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 3.21 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
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
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
#pragma once

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

#include <dune/grid/io/file/vtk/vtksequencewriter.hh>
#include <dune/geometry/genericgeometry/referenceelements.hh>

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

namespace AMDiS
{

template <class Traits>
class FileWriter
{
    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;
    
    
public:
    FileWriter(std::string base, MeshView const& meshView, std::vector<std::string> const& names_)
	: meshView(meshView)
	, names(names_)
    {
	std::string filename = "solution";
	Parameters::get(base + "->filename", filename);
	std::string 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, "");
    }

    template <class SystemVectorType>
    void write(double time, SystemVectorType&& solutions)
    {
	vtkWriter->clear();
	// copy dofvector to vertex data
	For<0, nComponents>::loop([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*/);
    }
    
    
protected:
    template <class DOFVector, class Vector>
    void dofVector2vertexVector(DOFVector dofvector, Vector& data)
    {
	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];
		}
	    }
	}
    }

private:
    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;
};

} // end namespace AMDiS