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.

FileWriter.hpp 6.66 KB
Newer Older
1
2
3
#pragma once

#include <string>
4
5
//#include <vector>
//#include <array>
6
7
#include <memory>

8
#include <dune/functions/functionspacebases/lagrangebasis.hh>
9
#include <dune/grid/io/file/vtk/vtkwriter.hh>
10
#include <dune/grid/io/file/vtk/vtksequencewriter.hh>
11
12
//#include <dune/geometry/referenceelements.hh>
#include <dune/typetree/childextraction.hh>
13

14
#include <dune/amdis/Initfile.hpp>
15
#include <dune/amdis/io/FileWriterInterface.hpp>
16
#include <dune/amdis/utility/Filesystem.hpp>
17
18
19

namespace AMDiS
{
20
  template <class GlobalBasis, class TreePath>
21
  class FileWriter
22
      : public FileWriterInterface
23
24
  {
  private: // typedefs and static constants
25

26
    using GridView = typename GlobalBasis::GridView;
27

28
    /// Dimension of the mesh
29
    static constexpr int dim = GridView::dimension;
30

31
    /// Dimension of the world
32
    static constexpr int dow = GridView::dimensionworld;
33

34
    using Vector = DOFVector<GlobalBasis>;
35
36
37



38
  public:
39
    /// Constructor.
40
41
42
43
44
45
46
47
    FileWriter(std::string baseName,
               std::shared_ptr<GlobalBasis> const& basis,
               TreePath const& tp,
               std::shared_ptr<Vector> const& vector)
      : FileWriterInterface(baseName)
      , basis_(basis)
      , treePath_(tp)
      , vector_(vector)
48
    {
49
      //int subSampling = Parameters::get<int>(baseName + "->subsampling").value_or(0);
50

51
52
53
54
      //if (subSampling > 0)
      //  vtkWriter_ = std::make_shared<Dune::SubsamplingVTKWriter<GridView>>(basis_->gridView(), subSampling);
      //else
      vtkWriter_ = std::make_shared<Dune::VTKWriter<GridView>>(basis_->gridView());
55

56
57
58
      vtkSeqWriter_ = std::make_shared<Dune::VTKSequenceWriter<GridView>>(vtkWriter_, filename_, dir_, "");

      mode_ = Parameters::get<int>(baseName + "->ParaView mode").value_or(0);
59
60
    }

61

62
63
64
65
66
67
68
69
70
71
72
73
74
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
    /// Implements \ref FileWriterInterface::writeFiles
    virtual void writeFiles(AdaptInfo& adaptInfo, bool force) override
    {
      using Tree = typename GlobalBasis::LocalView::Tree;
      using Node = typename Dune::TypeTree::ChildForTreePath<Tree, TreePath>;

      writeVertexData(typename Node::NodeTag{}, index_<Node::CHILDREN>, [&,this]()
      {
        if (mode_ == 0)
          vtkSeqWriter_->write(adaptInfo.getTime(), Dune::VTK::ascii);
        else if (mode_ == 1)
          vtkSeqWriter_->write(adaptInfo.getTime(), Dune::VTK::appendedraw);
      });
    }

  protected:
    template <class W>
    void writeVertexData(Dune::TypeTree::LeafNodeTag, index_t<0>, W write)
    {
      using namespace Dune::Functions::BasisBuilder;
      auto fct = makeDiscreteFunction(basis_,treePath_,vector_);

      auto p1basis = makeBasis(basis_->gridView(), lagrange<1>());
      auto data = makeDOFVector(p1basis, name_);
      interpolate(p1basis, data, *fct);

      auto dataFct = Dune::Functions::makeDiscreteGlobalBasisFunction<double>(p1basis,data);

      vtkWriter_->addVertexData(dataFct, Dune::VTK::FieldInfo(name_, Dune::VTK::FieldInfo::Type::scalar, 1));
      write();
      vtkWriter_->clear();
    }

    template <std::size_t C, class W>
    void writeVertexData(Dune::TypeTree::PowerNodeTag, index_t<C>, W write)
    {
      using namespace Dune::Functions::BasisBuilder;
      assert( C == dow );
      auto fct = makeDiscreteFunction(basis_,treePath_,vector_);

      auto p1basis = makeBasis(basis_->gridView(), power<C>(lagrange<1>(), flatLexicographic()));
      auto data = makeDOFVector(p1basis, name_);
      interpolate(p1basis, data, *fct);

      using Range = Dune::FieldVector<double,C>;
      auto dataFct = Dune::Functions::makeDiscreteGlobalBasisFunction<Range>(p1basis,data);

      vtkWriter_->addVertexData(dataFct, Dune::VTK::FieldInfo(name_, Dune::VTK::FieldInfo::Type::vector, C));
      write();
      vtkWriter_->clear();
    }

    template <class NodeTag, std::size_t C, class W>
    void writeVertexData(NodeTag, index_t<C>, W) {}

#if 0
118
    /// default write method for time-depended data
119
    template <class SystemVectorType>
120
    void write(double time, SystemVectorType const& solutions)
121
    {
122
123
      vtkWriter->clear();
      // copy dofvector to vertex data
124
      forEach(range_<0, nComponents>, [this, &solutions](const auto _i)
125
      {
126
        this->dofVector2vertexVector(solutions, std::get<_i>(data_vectors));
127
128
129
130
        vtkSeqWriter->addVertexData(std::get<_i>(data_vectors), names[_i]);
      });
      vtkSeqWriter->write(time/*, Dune::VTK::appendedraw*/);
    }
131

132
133
134

    /// default write method for stationary data
    template <class SystemVectorType>
135
    void write(SystemVectorType const& solutions)
136
137
138
    {
      vtkWriter->clear();
      // copy dofvector to vertex data
139
      forEach(range_<0, nComponents>, [this, &solutions](const auto _i)
140
      {
141
        this->dofVector2vertexVector(solutions, std::get<_i>(data_vectors));
142
143
144
        vtkWriter->addVertexData(std::get<_i>(data_vectors), names[_i]);
      });
      vtkWriter->pwrite(filename, dir, "" /*, Dune::VTK::appendedraw*/);
145
    }
146

147
    template <class DOFVector, class Vector>
148
    void dofVector2vertexVector(DOFVector const& dofvector, Vector& data)
149
    {
150
151
      using Geometry = typename MeshView::template Codim<0>::Geometry;
      using RefElements = Dune::ReferenceElements<typename Geometry::ctype, Geometry::mydimension>;
152

153
      data.resize(meshView.size(dim));
154

155
156
      auto const& indexSet = meshView.indexSet();
      auto const& feSpace = dofvector.getFeSpace();
157

158
159
      auto localView = feSpace.localView();
      auto localIndexSet = feSpace.localIndexSet();
160

161
162
163
164
165
      // copy data to P1-vector
      for (auto const& element : elements(meshView)) {
        localView.bind(element);
        localIndexSet.bind(localView);
        auto const& localBasis = localView.tree().finiteElement().localBasis();
166

167
        auto const& refElement = RefElements::general(element.type());
168

169
        std::vector<Dune::FieldVector<double,1> > shapeValues;
170

171
172
        std::size_t nVertices = element.subEntities(dim);
        for (std::size_t i = 0; i < nVertices; ++i) {
173
174
          auto const& v = element.template subEntity<dim>(i);
          auto pos = refElement.position(i, dim);
175

176
          localBasis.evaluateFunction(pos, shapeValues);
177

178
          std::size_t idx = indexSet.index(v);
179
          data[idx] = 0.0;
180
          for (std::size_t j = 0; j < shapeValues.size(); ++j) {
181
182
183
184
185
            const auto global_idx = localIndexSet.index(j);
            data[idx] += dofvector[global_idx] * shapeValues[j];
          }
        }
      }
186
    }
187
#endif
188

189
  private:
190
191
192
193
194
195
    std::shared_ptr<GlobalBasis> basis_;
    TreePath treePath_;
    std::shared_ptr<Vector> vector_;

    std::shared_ptr<Dune::VTKWriter<GridView>> vtkWriter_;
    std::shared_ptr<Dune::VTKSequenceWriter<GridView>> vtkSeqWriter_;
196

197
    // std::vector<double> data_vector;
198

199
200
    // represents VTK::OutputType: 0...ascii, 1...appendedraw
    int mode_;
201
  };
202
203

} // end namespace AMDiS