BallProject.cpp 2.46 KB
Newer Older
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
#include "BallProject.h"
#include "AdaptStationary.h"
/// Dirichlet boundary function
class G : public AbstractFunction<double, WorldVector<double> >
{
public:
  
  /// Implementation of AbstractFunction::operator().
  double operator()(const WorldVector<double>& x) const 
  {
    return exp(-10.0 * (x * x));
  }
};

/// RHS function
class F : public AbstractFunction<double, WorldVector<double> >
{
public:
  F(int degree) : AbstractFunction<double, WorldVector<double> >(degree) {}

  /// Implementation of AbstractFunction::operator().
  double operator()(const WorldVector<double>& x) const 
  {
    int dim = x.getSize();
    double r2 = x * x;
    double ux = exp(-10.0*r2);
    return -(400.0 * r2 - 20.0 * dim) * ux;
  }
};

31
Balldemo::Balldemo():
32
33
34
35
36
37
38
39
40
  ball("ball"),
  ballCenter(NULL),
  adaptInfo(NULL),
  adapt(NULL),
  matrixOperator(NULL),
  rhsOperator(NULL)
{
}

41
Balldemo::~Balldemo() {
42
43
44
45
46
47
48
49
50
51
52
53
  if(matrixOperator != NULL)
    delete matrixOperator;
  if(rhsOperator != NULL)
    delete rhsOperator;
  if(ballCenter != NULL)
    delete ballCenter;
  if(adapt != NULL)
    delete adapt;
  if(adaptInfo != NULL)
    delete adaptInfo;
}

54
void Balldemo::create(std::string& filename)  {
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
  // ===== init parameters =====
  Parameters::init(false, filename);
  ballCenter = new WorldVector< double >();
  ballCenter->set(0.0);

  // ===== create projection =====
  new AMDiS::BallProject(1, 
		  BOUNDARY_PROJECTION, 
		  *ballCenter, 
		  1.0);

  // ===== create and init the scalar problem ===== 
  ball.initialize(INIT_ALL);

  // === create adapt info ===
  adaptInfo = new AdaptInfo("ball->adapt");

  // === create adapt ===
  adapt = new AdaptStationary("ball->adapt", &ball, adaptInfo);
  
  // ===== create matrix operator =====
  matrixOperator=new Operator(ball.getFeSpace());
  matrixOperator->addSecondOrderTerm(new Laplace_SOT);
  ball.addMatrixOperator(matrixOperator);

  // ===== create rhs operator =====
  int degree = ball.getFeSpace()->getBasisFcts()->getDegree();
  rhsOperator = new Operator(ball.getFeSpace());
  rhsOperator->addZeroOrderTerm(new CoordsAtQP_ZOT(new F(degree)));
  ball.addVectorOperator(rhsOperator);

  // ===== add boundary conditions =====
  ball.addDirichletBC(1, new G);


}

92
int Balldemo::solve(SolutionInformation& solinfo) {
93
94
95
96
97
98
99
  assert(adaptInfo!=NULL);
  assert(adapt!=NULL);
  assert(matrixOperator!=NULL);
  assert(rhsOperator!=NULL);

  int retCode = adapt->adapt();
  solinfo.dofVec = ball.getSolution();
100
  //ball.writeFiles(adaptInfo,true);
101
102
  return retCode;
}
103