Am Montag, 13. Mai 2022, finden Wartungsarbeiten am Gitlab-Server (Update auf neue Version statt). Der Dienst wird daher am Montag für einige Zeit nicht verfügbar sein.
On Monday, May 13th 2022, the Gitlab server will be updated. The service will therefore not be accessible for some time on Monday.

neumann.cc 2.14 KB
Newer Older
1
2
3
4
5
6
7
8
#include "AMDiS.h"

using namespace AMDiS;
using namespace std;

class N : public AbstractFunction<double, WorldVector<double> >
{
public:
9
10
  double operator()(const WorldVector<double>& x) const 
  {
11
12
    return 1.0;
  }
13
14
};

15
/// Dirichlet boundary function
16
17
18
class G : public AbstractFunction<double, WorldVector<double> >
{
public:
19
20
21
  /// Implementation of AbstractFunction::operator().
  double operator()(const WorldVector<double>& x) const 
  {
22
23
    return exp(-10.0 * (x * x));
  }
24
25
};

26
/// RHS function
27
28
29
class F : public AbstractFunction<double, WorldVector<double> >
{
public:
30
  F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {}
31

32
33
34
  /// Implementation of AbstractFunction::operator().
  double operator()(const WorldVector<double>& x) const 
  {
35
36
37
    int dow = x.getSize();
    double r2 = (x*x);
    double ux = exp(-10.0*r2);
38
39
    return -(400.0*r2 - 20.0*dow)*ux;
  }
40
41
42
43
44
45
46
47
48
49
50
};

// // ===== main program // 
int main(int argc, char* argv[])
{
  FUNCNAME("main");

  // ===== check for init file =====
  TEST_EXIT(argc == 2)("usage: neumann initfile\n");

  // ===== init parameters =====
51
  Parameters::init(argv[1]);
52
53

  // ===== create and init the scalar problem ===== 
54
  ProblemStat neumann("neumann");
55
56
57
  neumann.initialize(INIT_ALL);

  // === create adapt info ===
58
  AdaptInfo *adaptInfo = new AdaptInfo("neumann->adapt", neumann.getNumComponents());
59
60

  // === create adapt ===
61
  AdaptStationary *adapt = new AdaptStationary("neumann->adapt",
62
63
64
65
					       &neumann,
					       adaptInfo);
  
  // ===== create matrix operator =====
66
  Operator matrixOperator(neumann.getFeSpace());
67
68
  matrixOperator.addSecondOrderTerm(new Simple_SOT);
  neumann.addMatrixOperator(&matrixOperator, 0, 0);
69
70

  // ===== create rhs operator =====
71
  int degree = neumann.getFeSpace()->getBasisFcts()->getDegree();
72
  Operator rhsOperator(neumann.getFeSpace());
73
  rhsOperator.addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree)));
74
  neumann.addVectorOperator(&rhsOperator, 0);
75

76
  // ===== add boundary conditions =====
77
78
  neumann.addNeumannBC(1, 0, 0, new N);
  neumann.addDirichletBC(2, 0, 0, new G);
79

80
81
82
83
84
85
86
  // ===== start adaption loop =====
  adapt->adapt();

  neumann.writeFiles(adaptInfo, true);
}