Mainpage.md 4.26 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
AMDiS    {#mainpage}
=====
The *Adaptive Multi-Dimensional Simulation Toolbox* (AMDiS) is implemented as a
discretization module on top of the [Dune](https://dune-project.org) framework.

Example
-------

An AMDiS program consists of three main incredients:
1. A Problem class that holds all information necessary for assembling a linear
   system, see \ref ProblemStat.
2. Operators describing the (bi)linear-form of your PDE, see \ref operators.
3. Adaption-modules for the time- and space adaptive solution of the problem, see \ref Adaption.

**Poisson equation:**

The most simple elliptic PDE is the Poisson equation:
\f{eqnarray*}{
  -\Delta u &=& f(x),\quad\mbox{ in }\Omega \\
          u &=& g(x),\quad\mbox{ on }\partial\Omega
\f}
where \f$ f(x) \f$ and \f$ g(x) \f$ are parameters describing the volume and
boundary forces of the problem.

The corresponding weak form of the equation reads:
\f[
  \langle \nabla v, \nabla u\rangle_\Omega = (f(x)\,v)_\Omega,\quad\forall v\in V_0
\f]
with \f$ u\in V_g \f$.

Thus, we need to define a grid (discretization of \f$ \Omega \f$) and a finite
Praetorius, Simon's avatar
Praetorius, Simon committed
32
element space (discretization of \f$ V_0 \f$ and \f$ V_g \f$). This cab be
Praetorius, Simon's avatar
Praetorius, Simon committed
33
34
35
36
37
38
39
provided as `Traits` parameter in the ProblemStat class:

~~~~~~~~~~~~~~~{.cpp}
using Grid = Dune::AlbertaGrid<AMDIS_DIM, AMDIS_DOW>;
using Traits = LagrangeBasis<Grid::LeafGridView, 1>;
~~~~~~~~~~~~~~~

Praetorius, Simon's avatar
Praetorius, Simon committed
40
where `AlbertaGrid` defines a grid and `Traits` a Lagrange Finite-Element space with
Praetorius, Simon's avatar
Praetorius, Simon committed
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
local polynomial degree 1 on the elements.

All AMDiS programs start with initialization of the library, using `AMDiS::init`

~~~~~~~~~~~~~~~{.cpp}
using namespace AMDiS;
int main(int argc, char** argv)
{
  AMDiS::init(argc, argv);

  ProblemStat<Traits> prob("name");
~~~~~~~~~~~~~~~

The Problem class is initialized with a name, that is used as identifier in the
parameter files. In order to initialize the Finite-Element space, the grid and
all other parts of the problem class, call `initialize(Flag)` where `Flag`
specifies what to initialize. For now, we initialize everything: `INIT_ALL`:

~~~~~~~~~~~~~~~{.cpp}
prob.initialize(INIT_ALL);
~~~~~~~~~~~~~~~

Operators specify the (bi-)linear-form and the coefficient function in the term,
see \ref operators for a list of possible types. The bilinear-form in the Poisson
Praetorius, Simon's avatar
Praetorius, Simon committed
65
equation consists of a second-order term with coefficient = 1 and the linear-form
Praetorius, Simon's avatar
Praetorius, Simon committed
66
67
68
69
70
71
72
73
74
75
76
includes the function \f$ f(x)=-1 \f$:

~~~~~~~~~~~~~~~{.cpp}
auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0);
prob.addMatrixOperator(opL, 0, 0);

auto opF = makeOperator(tag::test{}, [](auto const& x) { return -1.0; }, 0);
prob.addVectorOperator(opF, 0);
~~~~~~~~~~~~~~~

Boundary conditions, in the example above a Dirichlet condition, is specified by
Praetorius, Simon's avatar
Praetorius, Simon committed
77
78
defining a predicate for the boundary \f$ \Gamma\subset\partial\Omega \f$ and the
values on the boundary \f$ g(x) = 0 \f$:
Praetorius, Simon's avatar
Praetorius, Simon committed
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95

~~~~~~~~~~~~~~~{.cpp}
auto predicate = [](auto const& x){ return x[0] < 1.e-8 || x[1] < 1.e-8; };
prob.addDirichletBC(predicate, 0, 0, 0.0);
~~~~~~~~~~~~~~~

The final step is the assembling and solution of the linear system. (Maybe including
grid adaption). This is realized using an \ref AdaptStationary class:

~~~~~~~~~~~~~~~{.cpp}
AdaptInfo adaptInfo("adapt");
AdaptStationary adapt("adapt", prob);
adapt.adapt(); // assemble and solve
~~~~~~~~~~~~~~~

Finally, finish the AMDiS program with `AMDiS::finish()`.

Praetorius, Simon's avatar
Praetorius, Simon committed
96
The complete program then reads:
Praetorius, Simon's avatar
Praetorius, Simon committed
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
~~~~~~~~~~~~~~~{.cpp}
#include <dune/amdis/AMDiS.hpp>
#include <dune/amdis/AdaptInfo.hpp>
#include <dune/amdis/AdaptStationary.hpp>
#include <dune/amdis/ProblemStat.hpp>

using namespace AMDiS;

using Grid = Dune::AlbertaGrid<AMDIS_DIM, AMDIS_DOW>;
using Traits = LagrangeBasis<Grid::LeafGridView, 1>;

int main(int argc, char** argv)
{
  AMDiS::init(argc, argv);

  ProblemStat<Traits> prob("poisson");
  prob.initialize(INIT_ALL);

  auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0);
  prob.addMatrixOperator(opL, 0, 0);

  auto opF = makeOperator(tag::test{}, [](auto const& x) { return -1.0; }, 0);
  prob.addVectorOperator(opF, 0);

  // set boundary condition
  auto predicate = [](auto const& x){ return x[0] < 1.e-8 || x[1] < 1.e-8; };
  prob.addDirichletBC(predicate, 0, 0, 0.0);

  // assemble and solve
  AdaptInfo adaptInfo("adapt");
  AdaptStationary adapt("adapt", prob);
  adapt.adapt();

  AMDiS::finalize();
  return 0;
}
~~~~~~~~~~~~~~~