test-vtkwriter.cc 3.97 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:

#if HAVE_CONFIG_H
#include "config.h" // autoconf defines, needed by the dune headers
#endif

#include <algorithm>
#include <iostream>
#include <ostream>
#include <sstream>
#include <string>
#include <vector>

#include <dune/common/fvector.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/functions/gridfunctions/analyticgridviewfunction.hh>
#include <dune/grid/yaspgrid.hh>
Praetorius, Simon's avatar
Praetorius, Simon committed
19
#if HAVE_DUNE_UGGRID
20
  #include <dune/grid/uggrid.hh>
21
  #include <dune/grid/utility/structuredgridfactory.hh>
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
#endif
#include <dune/vtk/vtkwriter.hh>

#include "checkvtkfile.hh"

// accumulate exit status
void acc (int &accresult, int result)
{
  if (accresult == 0 || (accresult == 77 && result != 0))
    accresult = result;
}

struct Acc
{
  int operator() (int v1, int v2) const
  {
    acc(v1, v2);
    return v1;
  }
};

template <class GridView>
int testGridView (std::string prefix, Dune::VTKChecker& vtkChecker, GridView const& gridView)
{
  enum { dim = GridView :: dimension };

  Dune::VtkWriter<GridView> vtk(gridView);

  auto f1 = Dune::Functions::makeAnalyticGridViewFunction([](const auto& x) { return std::sin(x.two_norm()); },gridView);
  vtk.addCellData(f1, "scalar-valued lambda");

  auto f2 = Dune::Functions::makeAnalyticGridViewFunction([](const auto& x) { return x; },gridView);
  vtk.addPointData(f2, "vector-valued lambda");

  int result = 0;

  // ASCII files
  {
60
    vtk.setFormat(Dune::Vtk::FormatTypes::ASCII);
61
62
63
64
65
66
    std::string name = vtk.write(prefix + "_ascii");
    if (gridView.comm().rank() == 0) vtkChecker.push(name);
  }

  // BINARY files
  {
67
    vtk.setFormat(Dune::Vtk::FormatTypes::BINARY);
68
69
70
71
72
73
    std::string name = vtk.write(prefix + "_binary");
    if (gridView.comm().rank() == 0) vtkChecker.push(name);
  }

  // COMPRESSED files
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
74
    vtk.setFormat(Dune::Vtk::FormatTypes::COMPRESSED);
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
    std::string name = vtk.write(prefix + "_compressed");
    if (gridView.comm().rank() == 0) vtkChecker.push(name);
  }

  return result;
}

template <class Grid>
int testGrid (std::string prefix, Dune::VTKChecker& vtkChecker, Grid& grid)
{
  if (grid.comm().rank() == 0)
    std::cout << "vtkCheck(" << prefix << ")" << std::endl;

  grid.globalRefine(1);

  int result = 0;

  acc(result, testGridView( prefix + "_leaf",   vtkChecker, grid.leafGridView() ));
  acc(result, testGridView( prefix + "_level0", vtkChecker, grid.levelGridView(0) ));
  acc(result, testGridView( prefix + "_level1", vtkChecker, grid.levelGridView(grid.maxLevel()) ));

  return result;
}

int main (int argc, char** argv)
{
  using namespace Dune;
  const MPIHelper& mpiHelper = MPIHelper::instance(argc, argv);

  if (mpiHelper.rank() == 0)
    std::cout << "vtktest: MPI_Comm_size == " << mpiHelper.size() << std::endl;

  int result = 0; // pass by default

  VTKChecker vtkChecker;

  YaspGrid<1> yaspGrid1({1.0}, {8});
  YaspGrid<2> yaspGrid2({1.0,2.0}, {8,4});
  YaspGrid<3> yaspGrid3({1.0,2.0,3.0}, {8,4,4});

  acc(result, testGrid("yaspgrid_1d", vtkChecker, yaspGrid1));
  acc(result, testGrid("yaspgrid_2d", vtkChecker, yaspGrid2));
  acc(result, testGrid("yaspgrid_3d", vtkChecker, yaspGrid3));

Praetorius, Simon's avatar
Praetorius, Simon committed
119
#if HAVE_DUNE_UGGRID
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
  auto ugGrid2a = StructuredGridFactory<UGGrid<2>>::createSimplexGrid({0.0,0.0}, {1.0,2.0}, {2u,4u});
  auto ugGrid3a = StructuredGridFactory<UGGrid<3>>::createSimplexGrid({0.0,0.0,0.0}, {1.0,2.0,3.0}, {2u,4u,6u});

  acc(result, testGrid("uggrid_simplex_2d", vtkChecker, *ugGrid2a));
  acc(result, testGrid("uggrid_simplex_3d", vtkChecker, *ugGrid3a));

  auto ugGrid2b = StructuredGridFactory<UGGrid<2>>::createCubeGrid({0.0,0.0}, {1.0,2.0}, {2u,4u});
  auto ugGrid3b = StructuredGridFactory<UGGrid<3>>::createCubeGrid({0.0,0.0,0.0}, {1.0,2.0,3.0}, {2u,4u,6u});

  acc(result, testGrid("uggrid_cube_2d", vtkChecker, *ugGrid2b));
  acc(result, testGrid("uggrid_cube_3d", vtkChecker, *ugGrid3b));
#endif

  acc(result, vtkChecker.check());

  mpiHelper.getCollectiveCommunication().allreduce<Acc>(&result, 1);

  return result;
}