ISTLCommTest.cpp 2.67 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1
#include <config.h>
2 3 4 5 6 7 8 9 10 11 12 13 14 15

#include <array>
#include <bitset>
#include <cmath>
#include <string>
#include <vector>

#include <dune/common/fvector.hh>

#include <dune/functions/functionspacebases/interpolate.hh>

#include <dune/grid/yaspgrid.hh>

#include <amdis/linearalgebra/istl/Communication.hpp>
16
#include <amdis/AdaptiveGrid.hpp>
17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34
#include <amdis/Environment.hpp>
#include <amdis/ProblemStatTraits.hpp>

#include "Tests.hpp"

using namespace AMDiS;

template <class Basis>
bool test(Basis& basis, std::string const& prefix)
{
  using GV = typename Basis::GridView;
  using Grid = typename GV::Grid;
  using Coord = typename Grid::template Codim<0>::Geometry::GlobalCoordinate;
  using T = double;
  static constexpr int d = Grid::dimension;

  bool passed = true;

35
#if HAVE_MPI
36 37 38 39 40 41
  // Interpolate test function f onto dofs and reference
  auto f = [&](const Coord& x) -> T { return std::pow(x[0],3) + x[1] + (d == 3 ? x[2] : 0.0); };
  std::vector<T> ref(basis.size());
  Dune::Functions::interpolate(basis, ref, f);
  std::vector<T> dofs = ref;

42 43
  // Exchange dofs
  basis.communicator().get().copyOwnerToAll(dofs, dofs);
44 45 46 47 48

  // Compare post-exchange data with reference
  for (std::size_t i = 0; i < dofs.size(); ++i)
    if (std::abs(dofs[i]-ref[i]) > AMDIS_TEST_TOL)
      passed = false;
49
#endif
50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77

  return passed;
}

int main(int argc, char** argv)
{
  using namespace Dune::Functions::BasisFactory;

  Environment env(argc, argv);

  if (Environment::mpiSize() == 1)
    return 0;

  using YaspGrid2d = Dune::YaspGrid< 2, Dune::EquidistantOffsetCoordinates<double, 2> >;
  using YaspGrid3d = Dune::YaspGrid< 3, Dune::EquidistantOffsetCoordinates<double, 3> >;

  // constants
  using Coord2d = Dune::FieldVector<double, 2>;
  using Coord3d = Dune::FieldVector<double, 3>;
  Coord2d lower2d(1.0);
  Coord2d upper2d(2.0);
  std::array<int, 2> nElements2d{16,3};
  Coord3d lower3d(1.0);
  Coord3d upper3d(2.0);
  std::array<int, 3> nElements3d{16,3,3};
  int ovlp = 1;

  {
Praetorius, Simon's avatar
Praetorius, Simon committed
78
    YaspGrid2d hostGrid(lower2d, upper2d, nElements2d, std::bitset<2>(), ovlp);
Praetorius, Simon's avatar
Praetorius, Simon committed
79
    AdaptiveGrid<YaspGrid2d> grid(hostGrid);
80

Praetorius, Simon's avatar
Praetorius, Simon committed
81
    auto l1 = LagrangeBasis<YaspGrid2d, 1>::create(grid.leafGridView());
82
    AMDIS_TEST(test(l1, "Yasp2d_L1_Ovlp"));
83

Praetorius, Simon's avatar
Praetorius, Simon committed
84
    auto th = TaylorHoodBasis<YaspGrid2d, 1>::create(grid.leafGridView());
85 86 87
    AMDIS_TEST(test(th, "Yasp2d_TH_Ovlp"));
  }
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
88
    YaspGrid3d hostGrid(lower3d, upper3d, nElements3d, std::bitset<3>(), ovlp);
Praetorius, Simon's avatar
Praetorius, Simon committed
89
    AdaptiveGrid<YaspGrid3d> grid(hostGrid);
90

Praetorius, Simon's avatar
Praetorius, Simon committed
91
    auto l1 = LagrangeBasis<YaspGrid3d, 1>::create(grid.leafGridView());
92
    AMDIS_TEST(test(l1, "Yasp3d_L1_Ovlp"));
93

Praetorius, Simon's avatar
Praetorius, Simon committed
94
    auto th = TaylorHoodBasis<YaspGrid3d, 1>::create(grid.leafGridView());
95 96 97 98 99
    AMDIS_TEST(test(th, "Yasp3d_TH_Ovlp"));
  }

  return report_errors();
}