surface.cc 1.68 KB
Newer Older
1
#include "config.h"
Praetorius, Simon's avatar
Praetorius, Simon committed
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
#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>;
18
using Basis = LagrangeBasis<Grid, 1>;
Praetorius, Simon's avatar
Praetorius, Simon committed
19
20
21
22
23
24
25

struct UnitRadius
{
  template <class T>
  static double eval(T&&) { return 1.0; }
};

26
27
28
#if !HAVE_ALBERTA
#error "Need Alberta!!"
#endif
Praetorius, Simon's avatar
Praetorius, Simon committed
29
30
31
32

// solve the equation -laplace(u) - k^2 u = f on the sphere
int main(int argc, char** argv)
{
33
  Environment env(argc, argv);
Praetorius, Simon's avatar
Praetorius, Simon committed
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63

  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);

  return 0;
}