3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
staggeredvtkoutputmodule.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_STAGGERED_VTK_OUTPUT_MODULE_HH
25#define DUMUX_STAGGERED_VTK_OUTPUT_MODULE_HH
26
27#include <dune/common/fvector.hh>
28
33
34namespace Dumux {
35
36template<class Scalar, class GlobalPosition>
37class PointCloudVtkWriter;
38
47template<class GridVariables, class SolutionVector>
49: public VtkOutputModule<GridVariables, SolutionVector>
50{
52 using GridGeometry = typename GridVariables::GridGeometry;
53 using GridView = typename GridGeometry::GridView;
54 using Scalar = typename GridVariables::Scalar;
55 using FaceVariables = typename GridVariables::GridFaceVariables::FaceVariables;
56 using FVElementGeometry = typename GridGeometry::LocalView;
57 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
58
59 using Element = typename GridView::template Codim<0>::Entity;
60 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
61
62 struct FaceVarScalarDataInfo { std::function<Scalar(const FaceVariables&)> get; std::string name; };
63 struct FaceVarVectorDataInfo { std::function<GlobalPosition(const SubControlVolumeFace& scvf, const FaceVariables&)> get; std::string name; };
64
65 struct FaceFieldScalarDataInfo
66 {
67 FaceFieldScalarDataInfo(const std::vector<Scalar>& f, const std::string& n) : data(f), name(n) {}
68 const std::vector<Scalar>& data;
69 const std::string name;
70 };
71
72 struct FaceFieldVectorDataInfo
73 {
74 FaceFieldVectorDataInfo(const std::vector<GlobalPosition>& f, const std::string& n) : data(f), name(n) {}
75 const std::vector<GlobalPosition>& data;
76 const std::string name;
77 };
78
79public:
80
81 template<class Sol>
83 const Sol& sol,
84 const std::string& name,
85 const std::string& paramGroup = "",
86 Dune::VTK::DataMode dm = Dune::VTK::conforming,
87 bool verbose = true)
89 , faceWriter_(std::make_shared<PointCloudVtkWriter<Scalar, GlobalPosition>>(coordinates_))
90 , sequenceWriter_(faceWriter_, name + "-face", "","",
91 gridVariables.curGridVolVars().problem().gridGeometry().gridView().comm().rank(),
92 gridVariables.curGridVolVars().problem().gridGeometry().gridView().comm().size() )
93
94 {
95 static_assert(std::is_same<Sol, SolutionVector>::value, "Make sure that sol has the same type as SolutionVector."
96 "Use StaggeredVtkOutputModule<GridVariables, decltype(sol)> when calling the constructor.");
97
98 // enable velocity output per default
100 writeFaceVars_ = getParamFromGroup<bool>(paramGroup, "Vtk.WriteFaceData", false);
101 coordinatesInitialized_ = false;
102 }
103
108
112 void addFaceField(const std::vector<Scalar>& v, const std::string& name)
113 {
114 if (v.size() == this->gridGeometry().gridView().size(1))
115 faceFieldScalarDataInfo_.emplace_back(v, name);
116 else
117 DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
118 }
119
123 void addFaceField(const std::vector<GlobalPosition>& v, const std::string& name)
124 {
125 if (v.size() == this->gridGeometry().gridView().size(1))
126 faceFieldVectorDataInfo_.emplace_back(v, name);
127 else
128 DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
129 }
130
134 void addFaceVariable(std::function<Scalar(const FaceVariables&)>&& f, const std::string& name)
135 {
136 faceVarScalarDataInfo_.push_back(FaceVarScalarDataInfo{f, name});
137 }
138
142 void addFaceVariable(std::function<GlobalPosition(const SubControlVolumeFace& scvf, const FaceVariables&)>&& f, const std::string& name)
143 {
144 faceVarVectorDataInfo_.push_back(FaceVarVectorDataInfo{f, name});
145 }
146
150 void write(double time, Dune::VTK::OutputType type = Dune::VTK::ascii)
151 {
152 ParentType::write(time, type);
153 if(writeFaceVars_)
154 getFaceDataAndWrite_(time);
155 }
156
157
158private:
159
161 void updateCoordinates_()
162 {
163 coordinates_.resize(this->gridGeometry().numFaceDofs());
164 for(auto&& facet : facets(this->gridGeometry().gridView()))
165 {
166 const int dofIdxGlobal = this->gridGeometry().gridView().indexSet().index(facet);
167 coordinates_[dofIdxGlobal] = facet.geometry().center();
168 }
169 coordinatesInitialized_ = true;
170 }
171
174 void getFaceDataAndWrite_(const Scalar time)
175 {
176 const auto numPoints = this->gridGeometry().numFaceDofs();
177
178 // make sure not to iterate over the same dofs twice
179 std::vector<bool> dofVisited(numPoints, false);
180
181 // get fields for all primary coordinates and variables
182 if(!coordinatesInitialized_)
183 updateCoordinates_();
184
185 // prepare some containers to store the relevant data
186 std::vector<std::vector<Scalar>> faceVarScalarData;
187 std::vector<std::vector<GlobalPosition>> faceVarVectorData;
188
189 if(!faceVarScalarDataInfo_.empty())
190 faceVarScalarData.resize(faceVarScalarDataInfo_.size(), std::vector<Scalar>(numPoints));
191
192 if(!faceVarVectorDataInfo_.empty())
193 faceVarVectorData.resize(faceVarVectorDataInfo_.size(), std::vector<GlobalPosition>(numPoints));
194
195 for (const auto& element : elements(this->gridGeometry().gridView(), Dune::Partitions::interior))
196 {
197 auto fvGeometry = localView(this->gridGeometry());
198 auto elemFaceVars = localView(this->gridVariables().curGridFaceVars());
199
200 if (!faceVarScalarDataInfo_.empty() || !faceVarVectorDataInfo_.empty())
201 {
202 fvGeometry.bind(element);
203 elemFaceVars.bindElement(element, fvGeometry, this->sol());
204
205 for (auto&& scvf : scvfs(fvGeometry))
206 {
207 const auto dofIdxGlobal = scvf.dofIndex();
208 if(dofVisited[dofIdxGlobal])
209 continue;
210
211 dofVisited[dofIdxGlobal] = true;
212
213 const auto& faceVars = elemFaceVars[scvf];
214
215 // get the scalar-valued data
216 for (std::size_t i = 0; i < faceVarScalarDataInfo_.size(); ++i)
217 faceVarScalarData[i][dofIdxGlobal] = faceVarScalarDataInfo_[i].get(faceVars);
218
219 // get the vector-valued data
220 for (std::size_t i = 0; i < faceVarVectorDataInfo_.size(); ++i)
221 faceVarVectorData[i][dofIdxGlobal] = faceVarVectorDataInfo_[i].get(scvf, faceVars);
222 }
223 }
224 }
225
226 // transfer the data to the point writer
227 if(!faceVarScalarDataInfo_.empty())
228 for (std::size_t i = 0; i < faceVarScalarDataInfo_.size(); ++i)
229 faceWriter_->addPointData(faceVarScalarData[i], faceVarScalarDataInfo_[i].name);
230
231 if(!faceVarVectorDataInfo_.empty())
232 for (std::size_t i = 0; i < faceVarVectorDataInfo_.size(); ++i)
233 faceWriter_->addPointData(faceVarVectorData[i], faceVarVectorDataInfo_[i].name);
234
235 // account for the custom fields
236 for(const auto& field: faceFieldScalarDataInfo_)
237 faceWriter_->addPointData(field.data, field.name);
238
239 for(const auto& field: faceFieldVectorDataInfo_)
240 faceWriter_->addPointData(field.data, field.name);
241
242 // write for the current time step
243 sequenceWriter_.write(time);
244
245 // clear coordinates to save some memory
246 coordinates_.clear();
247 coordinates_.shrink_to_fit();
248 coordinatesInitialized_ = false;
249 }
250
251
252 std::shared_ptr<PointCloudVtkWriter<Scalar, GlobalPosition>> faceWriter_;
253
254 VTKSequenceWriter<PointCloudVtkWriter<Scalar, GlobalPosition>> sequenceWriter_;
255
256 bool writeFaceVars_;
257
258 std::vector<GlobalPosition> coordinates_;
259 bool coordinatesInitialized_;
260
261 std::vector<FaceVarScalarDataInfo> faceVarScalarDataInfo_;
262 std::vector<FaceVarVectorDataInfo> faceVarVectorDataInfo_;
263
264 std::vector<FaceFieldScalarDataInfo> faceFieldScalarDataInfo_;
265 std::vector<FaceFieldVectorDataInfo> faceFieldVectorDataInfo_;
266
267};
268
269} // end namespace Dumux
270
271#endif
Base class to write pvd-files which contains a list of all collected vtk-files. This is a modified ve...
A VTK writer specialized for staggered grid implementations with dofs on the faces.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
Definition: adapt.hh:29
Velocity output for staggered free-flow models.
Definition: discretization/staggered/freeflow/velocityoutput.hh:38
A VTK output module to simplify writing dumux simulation data to VTK format Specialization for stagge...
Definition: staggeredvtkoutputmodule.hh:50
void addFaceField(const std::vector< Scalar > &v, const std::string &name)
Definition: staggeredvtkoutputmodule.hh:112
void write(double time, Dune::VTK::OutputType type=Dune::VTK::ascii)
Definition: staggeredvtkoutputmodule.hh:150
StaggeredVtkOutputModule(const GridVariables &gridVariables, const Sol &sol, const std::string &name, const std::string &paramGroup="", Dune::VTK::DataMode dm=Dune::VTK::conforming, bool verbose=true)
Definition: staggeredvtkoutputmodule.hh:82
void addFaceField(const std::vector< GlobalPosition > &v, const std::string &name)
Definition: staggeredvtkoutputmodule.hh:123
void addFaceVariable(std::function< GlobalPosition(const SubControlVolumeFace &scvf, const FaceVariables &)> &&f, const std::string &name)
Definition: staggeredvtkoutputmodule.hh:142
void addFaceVariable(std::function< Scalar(const FaceVariables &)> &&f, const std::string &name)
Definition: staggeredvtkoutputmodule.hh:134
A VTK output module to simplify writing dumux simulation data to VTK format.
Definition: pointcloudvtkwriter.hh:51
const std::string & paramGroup() const
the parameter group for getting parameter from the parameter tree
Definition: io/vtkoutputmodule.hh:94
void write(double time, Dune::VTK::OutputType type=Dune::VTK::ascii)
Definition: io/vtkoutputmodule.hh:171
const std::string & name() const
Definition: io/vtkoutputmodule.hh:193
bool verbose() const
Definition: io/vtkoutputmodule.hh:192
A VTK output module to simplify writing dumux simulation data to VTK format.
Definition: io/vtkoutputmodule.hh:305
const auto & problem() const
Definition: io/vtkoutputmodule.hh:383
const GridVariables & gridVariables() const
Definition: io/vtkoutputmodule.hh:384
void addVelocityOutput(std::shared_ptr< VelocityOutputType > velocityOutput)
Add a velocity output policy.
Definition: io/vtkoutputmodule.hh:358
const SolutionVector & sol() const
Definition: io/vtkoutputmodule.hh:386
const GridGeometry & gridGeometry() const
Definition: io/vtkoutputmodule.hh:385
A VTK output module to simplify writing dumux simulation data to VTK format.