Skip to content
GitLab
Menu
Projects
Groups
Snippets
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
amdis
amdis-core
Commits
107bcb89
Commit
107bcb89
authored
Apr 26, 2019
by
Praetorius, Simon
Browse files
Example/surface grid
parent
344e1aa3
Changes
5
Hide whitespace changes
Inline
Side-by-side
examples/CMakeLists.txt
View file @
107bcb89
...
...
@@ -27,3 +27,8 @@ add_dependencies(examples
navier_stokes.2d
convection_diffusion.2d
cahn_hilliard.2d
)
if
(
dune-foamgrid_FOUND
)
add_amdis_executable
(
NAME surface.2d SOURCES surface.cc DIM 2 DOW 3
)
add_dependencies
(
examples surface.2d
)
endif
()
\ No newline at end of file
examples/FactoryParametrization.hpp
0 → 100644
View file @
107bcb89
#pragma once
#include
<type_traits>
#include
<vector>
#include
<dune/grid/common/gridfactory.hh>
#include
<amdis/Output.hpp>
namespace
AMDiS
{
template
<
class
GridType
,
class
Projection
>
struct
Parametrization
{
using
Grid
=
GridType
;
using
ctype
=
typename
Grid
::
ctype
;
static
constexpr
int
dimension
=
Grid
::
dimension
;
static
constexpr
int
dimensionworld
=
Grid
::
dimensionworld
;
using
LeafIntersection
=
typename
Grid
::
LeafIntersection
;
template
<
int
cd
>
using
Codim
=
typename
Grid
::
template
Codim
<
cd
>;
};
}
// end namespace AMDiS
namespace
Dune
{
/// \brief A factory class for \ref SimpleGrid.
template
<
class
Grid
,
class
Projection
>
class
GridFactory
<
AMDiS
::
Parametrization
<
Grid
,
Projection
>>
:
public
GridFactoryInterface
<
AMDiS
::
Parametrization
<
Grid
,
Projection
>>
{
using
GridWrapper
=
AMDiS
::
Parametrization
<
Grid
,
Projection
>
;
/// Type used by the grid for coordinates
using
ctype
=
typename
Grid
::
ctype
;
static
constexpr
int
dim
=
Grid
::
dimension
;
static
constexpr
int
dow
=
Grid
::
dimensionworld
;
using
LocalCoordinate
=
FieldVector
<
ctype
,
dim
>
;
using
GlobalCoordinate
=
FieldVector
<
ctype
,
dow
>
;
using
ProjectionBase
=
VirtualFunction
<
LocalCoordinate
,
GlobalCoordinate
>
;
public:
/// Default constructor
GridFactory
()
{}
/// Insert a vertex into the coarse grid
virtual
void
insertVertex
(
GlobalCoordinate
const
&
pos
)
override
{
GlobalCoordinate
new_pos
(
pos
);
Projection
::
project
(
new_pos
);
coordinates
.
push_back
(
new_pos
);
factory
.
insertVertex
(
new_pos
);
}
/// Insert an element into the coarse grid.
virtual
void
insertElement
(
GeometryType
const
&
type
,
std
::
vector
<
unsigned
int
>
const
&
vertices
)
override
{
assert
(
type
.
isSimplex
()
);
std
::
vector
<
GlobalCoordinate
>
corners
(
vertices
.
size
());
for
(
std
::
size_t
i
=
0
;
i
<
vertices
.
size
();
++
i
)
corners
[
i
]
=
coordinates
[
vertices
[
i
]];
factory
.
insertElement
(
type
,
vertices
,
std
::
shared_ptr
<
ProjectionBase
>
(
new
Projection
(
std
::
move
(
corners
)))
);
}
virtual
void
insertBoundarySegment
(
std
::
vector
<
unsigned
int
>
const
&
/*vertices*/
)
override
{
AMDiS
::
warning
(
"Not implemented!"
);
}
/// Finalize grid creation and hand over the grid.
Grid
*
create
()
{
return
factory
.
createGrid
();
}
virtual
GridWrapper
*
createGrid
()
override
{
AMDiS
::
warning
(
"Should not be created. Use non-virtual method `create()` instead, to create the underlying grid!"
);
return
nullptr
;
}
private:
// buffers for the mesh data
std
::
vector
<
GlobalCoordinate
>
coordinates
;
GridFactory
<
Grid
>
factory
;
};
}
// end namespace Dune
examples/SphereMapping.hpp
0 → 100644
View file @
107bcb89
#pragma once
#include
<algorithm>
#include
<array>
#include
<dune/common/function.hh>
#include
<dune/common/fvector.hh>
namespace
AMDiS
{
/**
* \brief Mapping class mapping from a triangle with points on the unit sphere onto the sphere
* with theta in [0, pi] and phi in [0, 2*pi)
*/
template
<
int
dim
,
int
dow
,
class
Radius
,
class
ct
=
double
>
class
SphereMapping
:
public
Dune
::
VirtualFunction
<
Dune
::
FieldVector
<
ct
,
dim
>
,
Dune
::
FieldVector
<
ct
,
dow
>
>
{
using
LocalCoordinate
=
Dune
::
FieldVector
<
ct
,
dim
>
;
using
GlobalCoordinate
=
Dune
::
FieldVector
<
ct
,
dow
>
;
public:
template
<
class
VertexContainer
>
SphereMapping
(
VertexContainer
&&
vertices
)
{
std
::
copy_n
(
vertices
.
begin
(),
dim
+
1
,
vertices_
.
begin
());
}
void
evaluate
(
LocalCoordinate
const
&
x
,
GlobalCoordinate
&
y
)
const
{
// calculate global coordinate
auto
shapeFunctions
=
evaluateShapeFunctions
(
x
);
y
=
0.0
;
for
(
size_t
i
=
0
;
i
<
dim
+
1
;
i
++
)
for
(
size_t
j
=
0
;
j
<
dow
;
j
++
)
y
[
j
]
+=
vertices_
[
i
][
j
]
*
shapeFunctions
[
i
];
project
(
y
);
}
GlobalCoordinate
evaluateShapeFunctions
(
LocalCoordinate
const
&
x
)
const
{
GlobalCoordinate
out
;
out
[
0
]
=
1.0
;
for
(
size_t
i
=
0
;
i
<
2
;
i
++
)
{
out
[
0
]
-=
x
[
i
];
out
[
i
+
1
]
=
x
[
i
];
}
return
out
;
}
static
void
project
(
GlobalCoordinate
&
y
)
{
// project it on the unit sphere
y
/=
y
.
two_norm
()
/
Radius
::
eval
(
y
);
}
private:
std
::
array
<
GlobalCoordinate
,
dim
+
1
>
vertices_
;
};
}
// end namespace AMDiS
examples/init/surface.dat.2d
0 → 100644
View file @
107bcb89
surfaceMesh->macro file name: ./macro/sphere_flat.3d
surfaceMesh->global refinements: 1
surface->mesh: surfaceMesh
surface->solver->name: bicgstab
surface->solver->max iteration: 10000
surface->solver->relative tolerance: 1.e-6
surface->solver->info: 10
surface->solver->left precon: ilu
surface->output[0]->filename: surface.2d
surface->output[0]->name: u
surface->output[0]->output directory: ./output
examples/surface.cc
0 → 100644
View file @
107bcb89
#include
<amdis/AMDiS.hpp>
#include
<amdis/Integrate.hpp>
#include
<amdis/LocalOperators.hpp>
#include
<amdis/ProblemStat.hpp>
#include
<amdis/common/Literals.hpp>
#include
<dune/grid/albertagrid/albertareader.hh>
#include
<dune/foamgrid/foamgrid.hh>
#include
"FactoryParametrization.hpp"
#include
"SphereMapping.hpp"
using
namespace
AMDiS
;
using
namespace
Dune
::
Indices
;
using
Grid
=
Dune
::
FoamGrid
<
2
,
3
>
;
using
Basis
=
LagrangeBasis
<
typename
Grid
::
LeafGridView
,
1
>
;
struct
UnitRadius
{
template
<
class
T
>
static
double
eval
(
T
&&
)
{
return
1.0
;
}
};
// solve the equation -laplace(u) - k^2 u = f on the sphere
int
main
(
int
argc
,
char
**
argv
)
{
AMDiS
::
init
(
argc
,
argv
);
std
::
string
gridName
=
Parameters
::
get
<
std
::
string
>
(
"surface->mesh"
).
value
();
std
::
string
gridFileName
=
Parameters
::
get
<
std
::
string
>
(
gridName
+
"->macro file name"
).
value
();
using
Param
=
Parametrization
<
Grid
,
SphereMapping
<
2
,
3
,
UnitRadius
>>
;
Dune
::
GridFactory
<
Param
>
gridFactory
;
Dune
::
AlbertaReader
<
Param
>
().
readGrid
(
gridFileName
,
gridFactory
);
std
::
unique_ptr
<
Grid
>
grid
(
gridFactory
.
create
());
ProblemStat
<
Basis
>
prob
(
"surface"
,
*
grid
);
prob
.
initialize
(
INIT_ALL
);
auto
opL
=
makeOperator
(
tag
::
gradtest_gradtrial
{},
1.0
);
prob
.
addMatrixOperator
(
opL
,
0
,
0
);
double
kappa
=
Parameters
::
get
<
double
>
(
"surface->kappa"
).
value_or
(
1.0
);
auto
opM
=
makeOperator
(
tag
::
test_trial
{},
-
kappa
*
kappa
);
prob
.
addMatrixOperator
(
opM
,
0
,
0
);
auto
opForce
=
makeOperator
(
tag
::
test
{},
X
(
0
)
+
X
(
1
)
+
X
(
2
));
prob
.
addVectorOperator
(
opForce
,
0
);
AdaptInfo
adaptInfo
(
"adapt"
);
prob
.
assemble
(
adaptInfo
);
prob
.
solve
(
adaptInfo
);
prob
.
writeFiles
(
adaptInfo
);
AMDiS
::
finalize
();
return
0
;
}
Write
Preview
Supports
Markdown
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment