Mainpage.md 2.53 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1 2
AMDiS    {#mainpage}
=====
3 4 5 6 7 8 9 10 11 12 13
The *Adaptive Multi-Dimensional Simulation Toolbox* (AMDiS) is a Finite-Element
discretization module allowing for fast prototyping of PDEs on adaptively refined
grids. It is implemented on top of the [Dune](https://dune-project.org) framework.

Tutorial
--------
As an example of usage, we want to solve an elliptic PDE, the Poisson equation,
\f$ -\Delta u = f \f$ in \f$ \Omega \f$ with \f$ u = g \f$ on a subset of the boundary
\f$ \Gamma\subset\partial\Omega \f$. For simplicity, we assume \f$ f(x) \equiv -1 \f$
and \f$ g(x) \equiv 0 \f$, the domain \f$ \Omega \f$ a square domain \f$ [0,1]^2 \f$
and \f$ \Gamma \f$ the lower and left edge of the boundary.
Praetorius, Simon's avatar
Praetorius, Simon committed
14 15

~~~~~~~~~~~~~~~{.cpp}
16 17 18 19
#include <amdis/AMDiS.hpp>
#include <amdis/AdaptInfo.hpp>
#include <amdis/AdaptStationary.hpp>
#include <amdis/ProblemStat.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
20

21
// The namespace all AMDiS classes and functions are defined in
Praetorius, Simon's avatar
Praetorius, Simon committed
22 23
using namespace AMDiS;

24
// A dune grid type
25
using Grid = Dune::AlbertaGrid<GRIDDIM, WORLDDIM>;
26 27 28

// A dune-functions globalBasis wrapped in a struct,
// here representing local polynomial shape functions of degree 1
Praetorius, Simon's avatar
Praetorius, Simon committed
29 30
using Traits = LagrangeBasis<Grid::LeafGridView, 1>;

31
int main(int argc, char* argv[])
Praetorius, Simon's avatar
Praetorius, Simon committed
32
{
33
  // Initialize linear-algebra backend and read parameters from file
Praetorius, Simon's avatar
Praetorius, Simon committed
34 35
  AMDiS::init(argc, argv);

36
  // Create a problem class containing all data for assembling
Praetorius, Simon's avatar
Praetorius, Simon committed
37
  ProblemStat<Traits> prob("poisson");
38
  // Initialize grid, globalBasis, solution vector and systenmatrix
Praetorius, Simon's avatar
Praetorius, Simon committed
39 40
  prob.initialize(INIT_ALL);

41
  // An operator representing the weak laplacian with coefficient = 1.0
Praetorius, Simon's avatar
Praetorius, Simon committed
42 43 44
  auto opL = makeOperator(tag::gradtest_gradtrial{}, 1.0);
  prob.addMatrixOperator(opL, 0, 0);

45
  // An rhs-operator representing an analytic function f(x) = -1
Praetorius, Simon's avatar
Praetorius, Simon committed
46 47 48
  auto opF = makeOperator(tag::test{}, [](auto const& x) { return -1.0; }, 0);
  prob.addVectorOperator(opF, 0);

49
  // Define the boundary Gamma
Praetorius, Simon's avatar
Praetorius, Simon committed
50
  auto predicate = [](auto const& x){ return x[0] < 1.e-8 || x[1] < 1.e-8; };
51
  // Set a value g(x) = 0 at this part of the boundary
Praetorius, Simon's avatar
Praetorius, Simon committed
52 53
  prob.addDirichletBC(predicate, 0, 0, 0.0);

54
  // assemble and solve the linear system
Praetorius, Simon's avatar
Praetorius, Simon committed
55 56 57 58 59 60 61
  AdaptInfo adaptInfo("adapt");
  AdaptStationary adapt("adapt", prob);
  adapt.adapt();

  AMDiS::finalize();
  return 0;
}
62 63 64 65 66 67 68 69 70
~~~~~~~~~~~~~~~

Notes
-----
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. A parameter file controlling several parts of the solution process, see \ref Initfile.