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 4.21 KB
Newer Older
1
2
3
4
5
6
7
#pragma once

#include <string>
#include <vector>
#include <array>
#include <memory>

8
#include <dune/grid/io/file/vtk/vtkwriter.hh>
9
10
11
12
13
14
15
16
17
#include <dune/grid/io/file/vtk/vtksequencewriter.hh>
#include <dune/geometry/genericgeometry/referenceelements.hh>

#include "Loops.hpp"
#include "Initfile.hpp"

namespace AMDiS
{

18
19
20
21
22
  template <class Traits>
  class FileWriter
  {
  private: // typedefs and static constants
    
23
24
25
26
27
28
29
30
31
32
33
    using Mesh     = typename Traits::Mesh;
    using MeshView = typename Mesh::LeafGridView;
    
    /// Dimension of the mesh
    static constexpr int dim = Traits::dim;
    
    /// Dimension of the world
    static constexpr int dow = Traits::dimworld;
    
    /// Number of problem components
    static constexpr int nComponents = Traits::nComponents;
34
35
36
      
      
  public:
37
    
38
    /// Constructor.
39
    FileWriter(std::string base, MeshView const& meshView, std::vector<std::string> const& names_)
40
41
      : meshView(meshView)
      , names(names_)
42
    {
43
44
45
46
47
48
49
      filename = "solution";
      Parameters::get(base + "->filename", filename);
      dir = "output";
      Parameters::get(base + "->output directory", dir);
      
      vtkWriter = std::make_shared<Dune::VTKWriter<MeshView>>(meshView);
      vtkSeqWriter = std::make_shared<Dune::VTKSequenceWriter<MeshView>>(vtkWriter, filename, dir, "");
50
51
    }

52
53
    
    /// default write method for time-depended data
54
55
56
    template <class SystemVectorType>
    void write(double time, SystemVectorType&& solutions)
    {
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
      vtkWriter->clear();
      // copy dofvector to vertex data
      for_<0, nComponents>([this, &solutions](const auto _i) 
      {
        this->dofVector2vertexVector(solutions[_i], std::get<_i>(data_vectors));
        vtkSeqWriter->addVertexData(std::get<_i>(data_vectors), names[_i]);
      });
      vtkSeqWriter->write(time/*, Dune::VTK::appendedraw*/);
    }
    

    /// default write method for stationary data
    template <class SystemVectorType>
    void write(SystemVectorType&& solutions)
    {
      vtkWriter->clear();
      // copy dofvector to vertex data
      for_<0, nComponents>([this, &solutions](const auto _i) 
      {
        this->dofVector2vertexVector(solutions[_i], std::get<_i>(data_vectors));
        vtkWriter->addVertexData(std::get<_i>(data_vectors), names[_i]);
      });
      vtkWriter->pwrite(filename, dir, "" /*, Dune::VTK::appendedraw*/);
80
81
    }
    
82
83
      
  protected:
84
85
    
    template <class DOFVector, class Vector>
86
    void dofVector2vertexVector(DOFVector const& dofvector, Vector& data)
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
119
120
121
122
123
      using Geometry = typename MeshView::template Codim<0>::Geometry;
      using RefElements = Dune::ReferenceElements<typename Geometry::ctype, Geometry::mydimension>;
      
      data.resize(meshView.size(dim));
      
      auto const& indexSet = meshView.indexSet();
      auto const& feSpace = dofvector.getFeSpace();
      
      auto localView = feSpace.localView();
      auto localIndexSet = feSpace.localIndexSet();
      
      // copy data to P1-vector
      for (auto const& element : elements(meshView)) {
        localView.bind(element);
        localIndexSet.bind(localView);
        auto const& localBasis = localView.tree().finiteElement().localBasis();
        
        auto const& refElement = RefElements::general(element.type());
        
        std::vector<Dune::FieldVector<double,1> > shapeValues;
        
        size_t nVertices = element.subEntities(dim);
        for (size_t i = 0; i < nVertices; ++i) {
          auto const& v = element.template subEntity<dim>(i);
          auto pos = refElement.position(i, dim);
          
          localBasis.evaluateFunction(pos, shapeValues);
          
          size_t idx = indexSet.index(v);
          data[idx] = 0.0;
          for (size_t j = 0; j < shapeValues.size(); ++j) {
            const auto global_idx = localIndexSet.index(j);
            data[idx] += dofvector[global_idx] * shapeValues[j];
          }
        }
      }
124
125
    }

126
  private:
127
128
129
130
131
132
    MeshView const& meshView;
    std::shared_ptr<Dune::VTKWriter<MeshView>> vtkWriter;
    std::shared_ptr<Dune::VTKSequenceWriter<MeshView>> vtkSeqWriter;
    
    std::vector<std::string> names;
    std::array<std::vector<double>, nComponents> data_vectors;
133
134
135
136
    
    std::string filename;
    std::string dir;
  };
137
138

} // end namespace AMDiS