periodic.cc 2.01 KB
Newer Older
1
2
3
4
5
#include "AMDiS.h"

using namespace std;
using namespace AMDiS;

6
/// Dirichlet boundary function
7
8
9
class G : public AbstractFunction<double, WorldVector<double> >
{
public:
10
11
12
  /// Implementation of AbstractFunction::operator().
  double operator()(const WorldVector<double>& x) const 
  {
13
14
    return 0.0;
  }
15
16
};

17
/// RHS function
18
19
20
class F : public AbstractFunction<double, WorldVector<double> >
{
public:
21
  F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {}
22

23
24
25
  /// Implementation of AbstractFunction::operator().
  double operator()(const WorldVector<double>& x) const 
  {
26
27
28
    int dim = x.getSize();
    double r2 = (x*x);
    double ux = exp(-10.0*r2);
29
    return -(400.0*r2 - 20.0*dim)*ux;
30
  }
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
};

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

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

  // ===== init parameters =====
  Parameters::init(false, argv[1]);

  // ===== create and init the scalar problem ===== 
  ProblemScal periodic("periodic");
  periodic.initialize(INIT_ALL);

  // === create adapt info ===
48
  AdaptInfo *adaptInfo = new AdaptInfo("periodic->adapt", 1);
49
50

  // === create adapt ===
51
  AdaptStationary *adapt = new AdaptStationary("periodic->adapt",
52
53
54
55
					       &periodic,
					       adaptInfo);
  
  // ===== create matrix operator =====
56
  Operator matrixOperator(Operator::MATRIX_OPERATOR, periodic.getFeSpace());
57
  matrixOperator.addSecondOrderTerm(new Laplace_SOT);
58
59
60
  periodic.addMatrixOperator(&matrixOperator);

  // ===== create rhs operator =====
61
62
  int degree = periodic.getFeSpace()->getBasisFcts()->getDegree();
  Operator rhsOperator(Operator::VECTOR_OPERATOR, periodic.getFeSpace());
63
  rhsOperator.addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree)));
64
65
  periodic.addVectorOperator(&rhsOperator);

66
67
68
69
70
71
  // ===== add boundary conditions =====
  periodic.addDirichletBC(1, new G);

  // ===== add periodic conditions =====
  periodic.addPeriodicBC(-1);

72
73
74
75
76
77
78
  // ===== start adaption loop =====
  adapt->adapt();

  periodic.writeFiles(adaptInfo, true);
}