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.

vtkwriter.cc 3.51 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1
2
3
4
5
6
7
8
9
10
11
12
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:

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

#include <iostream>
#include <vector>

#include <dune/common/parallel/mpihelper.hh> // An initializer of MPI
#include <dune/common/exceptions.hh> // We use exceptions
13
#include <dune/common/version.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
14
15
16
17
18
19
20
21
22

#include <dune/functions/functionspacebases/defaultglobalbasis.hh>
#include <dune/functions/functionspacebases/lagrangebasis.hh>
#include <dune/functions/functionspacebases/interpolate.hh>
#include <dune/functions/gridfunctions/analyticgridviewfunction.hh>
#include <dune/functions/gridfunctions/discreteglobalbasisfunction.hh>
#include <dune/grid/uggrid.hh>
#include <dune/grid/yaspgrid.hh>

23
#include <dune/vtk/vtkwriter.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
24
25
26
27

using namespace Dune;
using namespace Dune::Functions;

28
29
30
31
32
33
34
35
36
37
38
39
using TestCases = std::set<std::tuple<std::string,Vtk::FormatTypes,Vtk::DataTypes>>;
static TestCases test_cases = {
  {"ascii32", Vtk::ASCII, Vtk::FLOAT32},
  {"bin32", Vtk::BINARY, Vtk::FLOAT32},
  {"zlib32", Vtk::COMPRESSED, Vtk::FLOAT32},
  {"ascii64", Vtk::ASCII, Vtk::FLOAT64},
  {"bin64", Vtk::BINARY, Vtk::FLOAT64},
  {"zlib64", Vtk::COMPRESSED, Vtk::FLOAT64},
};

template <class GridView>
void write (std::string prefix, GridView const& gridView)
Praetorius, Simon's avatar
Praetorius, Simon committed
40
{
41
#if DUNE_VERSION_LT(DUNE_FUNCTIONS, 2, 6)
42
43
  using namespace BasisBuilder;
#else
Praetorius, Simon's avatar
Praetorius, Simon committed
44
  using namespace BasisFactory;
45
#endif
Praetorius, Simon's avatar
Praetorius, Simon committed
46
47
  auto basis = makeBasis(gridView, lagrange<1>());

Praetorius, Simon's avatar
Praetorius, Simon committed
48
49
50
51
  FieldVector<double,GridView::dimensionworld> c;
  if (GridView::dimensionworld > 0) c[0] = 11.0;
  if (GridView::dimensionworld > 1) c[1] = 7.0;
  if (GridView::dimensionworld > 2) c[2] = 3.0;
Praetorius, Simon's avatar
Praetorius, Simon committed
52

Praetorius, Simon's avatar
Praetorius, Simon committed
53
  assert(basis.dimension() > 0);
54
55
  std::vector<double> vec(basis.dimension());
  interpolate(basis, vec, [&c](auto const& x) { return c.dot(x); });
Praetorius, Simon's avatar
Praetorius, Simon committed
56
57

  // write discrete global-basis function
58
  auto p1Interpol = makeDiscreteGlobalBasisFunction<double>(basis, vec);
Praetorius, Simon's avatar
Praetorius, Simon committed
59
60

  // write analytic function
61
62
63
  auto p1Analytic = makeAnalyticGridViewFunction([&c](auto const& x) { return c.dot(x); }, gridView);

  for (auto const& test_case : test_cases) {
64
65
66
67
68
    VtkWriter<GridView> vtkWriter(gridView, std::get<1>(test_case), std::get<2>(test_case));
    vtkWriter.addPointData(p1Interpol, "p1");
    vtkWriter.addCellData(p1Interpol, "p0");
    vtkWriter.addPointData(p1Analytic, "q1");
    vtkWriter.addCellData(p1Analytic, "q0");
Praetorius, Simon's avatar
Praetorius, Simon committed
69
    vtkWriter.write(prefix + "_" + std::to_string(GridView::dimensionworld) + "d_" + std::get<0>(test_case) + ".vtu");
70
  }
Praetorius, Simon's avatar
Praetorius, Simon committed
71
}
72
73
74
75
76
77
78
79

template <int I>
using int_ = std::integral_constant<int,I>;

int main (int argc, char** argv)
{
  Dune::MPIHelper::instance(argc, argv);

80
#if HAVE_UG
81
82
83
84
85
86
87
  // Test VtkWriter for UGGrid
  Hybrid::forEach(std::make_tuple(int_<2>{}, int_<3>{}), [](auto dim)
  {
    using GridType = UGGrid<dim.value>;
    {
      FieldVector<double,dim.value> lowerLeft; lowerLeft = 0.0;
      FieldVector<double,dim.value> upperRight; upperRight = 1.0;
Praetorius, Simon's avatar
Praetorius, Simon committed
88
      auto numElements = filledArray<dim.value,unsigned int>(8);
89
      auto gridPtr = StructuredGridFactory<GridType>::createSimplexGrid(lowerLeft, upperRight, numElements);
90
      gridPtr->loadBalance();
Praetorius, Simon's avatar
Praetorius, Simon committed
91
      write("vtkwriter_ug", gridPtr->leafGridView());
92
93
94
95
96
97
98
99
100
    }
  });
#endif

  // Test VtkWriter for YaspGrid
  Hybrid::forEach(std::make_tuple(int_<1>{}, int_<2>{}, int_<3>{}), [](auto dim)
  {
    using GridType = YaspGrid<dim.value>;
    FieldVector<double,dim.value> upperRight; upperRight = 1.0;
Praetorius, Simon's avatar
Praetorius, Simon committed
101
102
    auto numElements = filledArray<dim.value,int>(8);
    GridType grid(upperRight, numElements, 0, 0);
Praetorius, Simon's avatar
Praetorius, Simon committed
103
    write("vtkwriter_yasp", grid.leafGridView());
104
  });
105
}