DataTransferTest.hpp 4.93 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
31
32
33
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

#ifdef HAVE_CONFIG_H
# include "config.h"
#endif

#include <array>
#include <cmath>
#include <functional>
#include <memory>
#include <vector>

#include <dune/common/fvector.hh>
#include <dune/common/parallel/mpihelper.hh>

#ifdef HAVE_DUNE_UGGRID
#include <dune/grid/uggrid.hh>
#else
#include <dune/grid/yaspgrid.hh>
#endif

#include <amdis/ProblemStat.hpp>
#include <amdis/ProblemStatTraits.hpp>

#include "Tests.hpp"

using namespace AMDiS;

template <class BasisCreator>
auto makeGrid(bool simplex, int globalRefines = 0)
{
  using Grid = typename BasisCreator::GlobalBasis::GridView::Grid;
  static constexpr int d = Grid::dimension; // problem dimension
  using DomainType = typename Dune::FieldVector<double, d>;

  // constants
  DomainType lowerLeft; lowerLeft = 0.0;    // lower left grid corner
  DomainType upperRight; upperRight = 1.0;  // upper right grid corner
  std::array<unsigned int, d> s; s.fill(1); // number of elements on each axis

  // make grid
  std::unique_ptr<Grid> gridPtr;
  if (simplex)
  {
    gridPtr = std::unique_ptr<Grid>{
      Dune::StructuredGridFactory<Grid>::createSimplexGrid(lowerLeft, upperRight, s)};
  }
  else
  {
    gridPtr = std::unique_ptr<Grid>{
      Dune::StructuredGridFactory<Grid>::createCubeGrid(lowerLeft, upperRight, s)};
  }
  gridPtr->globalRefine(globalRefines);

  return gridPtr;
}

template <class BasisCreator, class Fcts>
auto makeProblem(typename BasisCreator::GlobalBasis::GridView::Grid& grid, Fcts const& funcs)
{
  using Problem = ProblemStat<BasisCreator>;

  // make problem
63
64
  auto prob = std::make_unique<Problem>("test", grid);
  prob->initialize(INIT_ALL);
65

66
  auto& globalBasis = *prob->globalBasis();
67
68
69
70
  auto localView = globalBasis.localView();

  // interpolate given function to initial grid
  int k = 0;
71
  for_each_leaf_node(localView.tree(), [&](auto const& node, auto tp)
72
  {
73
74
    auto gf = makeGridFunction(funcs[k], globalBasis.gridView());
    AMDiS::interpolate(globalBasis, prob->solution(tp).coefficients(), gf, tp);
75
76
77
    k++;
  });

78
  return std::move(prob);
79
80
81
82
83
}

template <class Problem, class Fcts>
double calcError(Problem const& prob, Fcts const& funcs)
{
84
  auto& globalBasis = *prob->globalBasis();
85
  auto localView = globalBasis.localView();
86
  auto const& sol = prob->solution().coefficients();
87
88

  auto ref = makeDOFVector(globalBasis, DataTransferOperation::NO_OPERATION);
89
90
91
  int k = 0;

  // interpolate given function onto reference vector
92
  for_each_leaf_node(localView.tree(), [&](auto const& node, auto tp)
93
  {
94
    auto gf = makeGridFunction(funcs[k], globalBasis.gridView());
95
    AMDiS::interpolate(globalBasis, ref, gf, tp);
96
97
98
99
    k++;
  });

  // compare the solution with the reference
100
101
  double maxError = 0;
  sol.forEach([&](auto dof, double coeff) {
102
    double error = std::abs(ref.at(dof) - coeff);
103
    maxError = std::max(maxError, error);
104
  });
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
  return maxError;
}

template <class BasisCreator>
  using Fcts = std::vector<std::function<double(
      Dune::FieldVector<double, BasisCreator::GlobalBasis::GridView::Grid::dimension>)> >;

/// Test data transfer for the case where no grid changes are made
template <class BasisCreator>
bool unchanged_test(Fcts<BasisCreator> const& funcs, bool simplex = true)
{
  using Grid = typename BasisCreator::GlobalBasis::GridView::Grid;
  static constexpr int d = Grid::dimension; // problem dimension

  auto gridPtr = makeGrid<BasisCreator>(simplex, (d > 2 ? 3 : 5));
  auto prob = makeProblem<BasisCreator, Fcts<BasisCreator>>(*gridPtr, funcs);
  // mark a single element for coarsening -> no adaptation is done
  auto e = *gridPtr->leafGridView().template begin<0>();
  gridPtr->mark(-1, e);
  AdaptInfo adaptInfo("adapt");
125
  prob->adaptGrid(adaptInfo);
126
127
128
129
130
131
132
133
134
135
136
137
138
139
  auto error = calcError(prob, funcs);

  return error < AMDIS_TEST_TOL;
}

/// Test data transfer for the case of global coarsening
template <class BasisCreator>
bool coarsen_test(Fcts<BasisCreator> const& funcs, bool simplex = true)
{
  using Grid = typename BasisCreator::GlobalBasis::GridView::Grid;
  static constexpr int d = Grid::dimension; // problem dimension

  auto gridPtr = makeGrid<BasisCreator>(simplex, (d > 2 ? 2 : 4));
  auto prob = makeProblem<BasisCreator, Fcts<BasisCreator>>(*gridPtr, funcs);
140
  prob->globalCoarsen(1);
141
142
143
144
145
146
147
148
149
150
151
152
153
154
  auto error = calcError(prob, funcs);

  return error < AMDIS_TEST_TOL;
}

/// Test data transfer for the case of global refinement
template <class BasisCreator>
bool refine_test(Fcts<BasisCreator> const& funcs, bool simplex = true)
{
  using Grid = typename BasisCreator::GlobalBasis::GridView::Grid;
  static constexpr int d = Grid::dimension; // problem dimension

  auto gridPtr = makeGrid<BasisCreator>(simplex, (d > 2 ? 1 : 3));
  auto prob = makeProblem<BasisCreator, Fcts<BasisCreator>>(*gridPtr, funcs);
155
  prob->globalRefine(1);
156
157
158
159
160
161
  auto error = calcError(prob, funcs);

  return error < AMDIS_TEST_TOL;
}

template <class Grid>
162
  using Lagrange3 = LagrangeBasis<Grid, 3>;
163
template <class Grid>
164
  using TaylorHood = TaylorHoodBasis<Grid>;
165
166

  constexpr double pi = 3.141592653589793238463;