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.58 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

#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>
Praetorius, Simon's avatar
Praetorius, Simon committed
20
#if HAVE_DUNE_UGGRID
Praetorius, Simon's avatar
Praetorius, Simon committed
21
#include <dune/grid/uggrid.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
22
#endif
Praetorius, Simon's avatar
Praetorius, Simon committed
23
24
#include <dune/grid/yaspgrid.hh>

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

using namespace Dune;
using namespace Dune::Functions;

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

template <class GridView>
void write (std::string prefix, GridView const& gridView)
Praetorius, Simon's avatar
Praetorius, Simon committed
42
43
44
45
{
  using namespace BasisFactory;
  auto basis = makeBasis(gridView, lagrange<1>());

Praetorius, Simon's avatar
Praetorius, Simon committed
46
47
48
49
  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
50

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

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

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

  for (auto const& test_case : test_cases) {
62
63
64
65
66
    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
67
    vtkWriter.write(prefix + "_" + std::to_string(GridView::dimensionworld) + "d_" + std::get<0>(test_case) + ".vtu");
68
  }
Praetorius, Simon's avatar
Praetorius, Simon committed
69
}
70
71
72
73
74
75
76
77

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

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

Praetorius, Simon's avatar
Praetorius, Simon committed
78
#if HAVE_DUNE_UGGRID
79
80
81
82
83
84
85
  // 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
86
      auto numElements = filledArray<dim.value,unsigned int>(8);
87
      auto gridPtr = StructuredGridFactory<GridType>::createSimplexGrid(lowerLeft, upperRight, numElements);
88
      gridPtr->loadBalance();
Praetorius, Simon's avatar
Praetorius, Simon committed
89
      write("vtkwriter_ug", gridPtr->leafGridView());
90
91
92
93
94
95
96
97
98
    }
  });
#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
99
100
    auto numElements = filledArray<dim.value,int>(8);
    GridType grid(upperRight, numElements, 0, 0);
Praetorius, Simon's avatar
Praetorius, Simon committed
101
    write("vtkwriter_yasp", grid.leafGridView());
102
  });
103
}