diff --git a/dumux/io/vtkoutputmodule.hh b/dumux/io/vtkoutputmodule.hh
index 47ff05309160367b89e7799eaa7d8cb58113fa54..be8bab380eb4019bbcfad3c81af6fa32797d940e 100644
--- a/dumux/io/vtkoutputmodule.hh
+++ b/dumux/io/vtkoutputmodule.hh
@@ -242,6 +242,11 @@ class VtkOutputModule
public:
+ enum class FieldType : unsigned int
+ {
+ element, vertex, automatic
+ };
+
VtkOutputModule(const Problem& problem,
const FVGridGeometry& fvGridGeometry,
const GridVariables& gridVariables,
@@ -281,19 +286,51 @@ public:
volVarScalarDataInfo_.push_back(VolVarScalarDataInfo{f, name});
}
- //! Add a scalar or vector valued field
- //! \param v The scalar field to be added
- //! \param name The name of the vtk field
- //! \param nComp Number of components of the vector (maximum 3)
+ //! Add a scalar or vector valued vtk field
+ //! \param v The field to be added. Can be any indexable container. Its value type can be a number or itself an indexable container.
+ //! \param name The name of the field
+ //! \param fieldType The type of the field.
+ //! This determines whether the values are associated with vertices or elements.
+ //! By default, the method automatically deduces the correct type for the given input.
template<typename Vector>
- void addField(const Vector& v, const std::string& name, int nComp = 1)
+ void addField(const Vector& v, const std::string& name, FieldType fieldType = FieldType::automatic)
{
- if (v.size() == gridGeom_.gridView().size(0))
+ // Deduce the number of components from the given vector type
+ const auto nComp = getNumberOfComponents_(v);
+
+ const auto numElements = gridGeom_.gridView().size(0);
+ const auto numVertices = gridGeom_.gridView().size(dim);
+
+ // Automatically deduce the field type ...
+ if(fieldType == FieldType::automatic)
+ {
+ if(numElements == numVertices)
+ DUNE_THROW(Dune::InvalidStateException, "Automatic deduction of FieldType failed. Please explicitly specify FieldType::element or FieldType::vertex.");
+
+ if(v.size() == numElements)
+ fieldType = FieldType::element;
+ else if(v.size() == numVertices)
+ fieldType = FieldType::vertex;
+ else
+ DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
+ }
+ // ... or check if the user-specified type matches the size of v
+ else
+ {
+ if(fieldType == FieldType::element)
+ if(v.size() != numElements)
+ DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
+
+ if(fieldType == FieldType::vertex)
+ if(v.size() != numVertices)
+ DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
+ }
+
+ // add the appropriate field
+ if (fieldType == FieldType::element)
fields_.emplace_back(gridGeom_.gridView(), gridGeom_.elementMapper(), v, name, nComp, 0);
- else if (v.size() == gridGeom_.gridView().size(dim))
- fields_.emplace_back(gridGeom_.gridView(), gridGeom_.vertexMapper(), v, name, nComp, dim);
else
- DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
+ fields_.emplace_back(gridGeom_.gridView(), gridGeom_.vertexMapper(), v, name, nComp, dim);
}
//////////////////////////////////////////////////////////////////////////////////////////////
@@ -458,6 +495,15 @@ public:
}
private:
+ //! Deduces the number of components of the value type of a vector of values
+ template<class Vector, typename std::enable_if_t<Dune::is_indexable<decltype(std::declval<Vector>()[0])>::value, int> = 0>
+ std::size_t getNumberOfComponents_(const Vector& v)
+ { return v[0].size(); }
+
+ //! Deduces the number of components of the value type of a vector of values
+ template<class Vector, typename std::enable_if_t<!Dune::is_indexable<decltype(std::declval<Vector>()[0])>::value, int> = 0>
+ std::size_t getNumberOfComponents_(const Vector& v)
+ { return 1; }
template<typename Writer, typename... Args>
void addDofDataForWriter_(Writer& writer,
diff --git a/test/freeflow/staggered/test_angeli.cc b/test/freeflow/staggered/test_angeli.cc
index 7e53402957abacfbc5c9aa7eec866b6468713f83..6b7f620176f352d38c7ffe061fcab9d55a31e1c7 100644
--- a/test/freeflow/staggered/test_angeli.cc
+++ b/test/freeflow/staggered/test_angeli.cc
@@ -141,7 +141,6 @@ int main(int argc, char** argv) try
problem->setTimeLoop(timeLoop);
// the solution vector
- using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using SolutionVector = typename GET_PROP_TYPE(TypeTag, SolutionVector);
using DofTypeIndices = typename GET_PROP(TypeTag, DofTypeIndices);
typename DofTypeIndices::CellCenterIdx cellCenterIdx;
@@ -163,8 +162,8 @@ int main(int argc, char** argv) try
using VtkOutputFields = typename GET_PROP_TYPE(TypeTag, VtkOutputFields);
StaggeredVtkOutputModule<TypeTag> vtkWriter(*problem, *fvGridGeometry, *gridVariables, x, problem->name());
VtkOutputFields::init(vtkWriter); //! Add model specific output fields
- vtkWriter.addField(problem->getAnalyticalPressureSolution(), "pressureExact", 1);
- vtkWriter.addField(problem->getAnalyticalVelocitySolution(), "velocityExact", GridView::dimensionworld);
+ vtkWriter.addField(problem->getAnalyticalPressureSolution(), "pressureExact");
+ vtkWriter.addField(problem->getAnalyticalVelocitySolution(), "velocityExact");
vtkWriter.addFaceField(problem->getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact");
vtkWriter.write(0.0);
diff --git a/test/freeflow/staggered/test_donea.cc b/test/freeflow/staggered/test_donea.cc
index ca8684a1198d40aabb2ac2d25a434a236dcef8b2..173b09a2c3c3a2a9c61164fdf3fb77650b614042 100644
--- a/test/freeflow/staggered/test_donea.cc
+++ b/test/freeflow/staggered/test_donea.cc
@@ -125,7 +125,6 @@ int main(int argc, char** argv) try
auto problem = std::make_shared<Problem>(fvGridGeometry);
// the solution vector
- using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using SolutionVector = typename GET_PROP_TYPE(TypeTag, SolutionVector);
using DofTypeIndices = typename GET_PROP(TypeTag, DofTypeIndices);
typename DofTypeIndices::CellCenterIdx cellCenterIdx;
@@ -159,8 +158,8 @@ int main(int argc, char** argv) try
using VtkOutputFields = typename GET_PROP_TYPE(TypeTag, VtkOutputFields);
StaggeredVtkOutputModule<TypeTag> vtkWriter(*problem, *fvGridGeometry, *gridVariables, x, problem->name());
VtkOutputFields::init(vtkWriter); //! Add model specific output fields
- vtkWriter.addField(problem->getAnalyticalPressureSolution(), "pressureExact", 1);
- vtkWriter.addField(problem->getAnalyticalVelocitySolution(), "velocityExact", GridView::dimensionworld);
+ vtkWriter.addField(problem->getAnalyticalPressureSolution(), "pressureExact");
+ vtkWriter.addField(problem->getAnalyticalVelocitySolution(), "velocityExact");
vtkWriter.addFaceField(problem->getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact");
vtkWriter.write(0.0);
diff --git a/test/freeflow/staggered/test_kovasznay.cc b/test/freeflow/staggered/test_kovasznay.cc
index 16f5b81f7e7227f40419be6076a735d5e422c8df..b1c855a136d9605751c732f73eb91d507b87fd47 100644
--- a/test/freeflow/staggered/test_kovasznay.cc
+++ b/test/freeflow/staggered/test_kovasznay.cc
@@ -124,7 +124,6 @@ int main(int argc, char** argv) try
auto problem = std::make_shared<Problem>(fvGridGeometry);
// the solution vector
- using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
using SolutionVector = typename GET_PROP_TYPE(TypeTag, SolutionVector);
using DofTypeIndices = typename GET_PROP(TypeTag, DofTypeIndices);
typename DofTypeIndices::CellCenterIdx cellCenterIdx;
@@ -158,8 +157,8 @@ int main(int argc, char** argv) try
using VtkOutputFields = typename GET_PROP_TYPE(TypeTag, VtkOutputFields);
StaggeredVtkOutputModule<TypeTag> vtkWriter(*problem, *fvGridGeometry, *gridVariables, x, problem->name());
VtkOutputFields::init(vtkWriter); //! Add model specific output fields
- vtkWriter.addField(problem->getAnalyticalPressureSolution(), "pressureExact", 1);
- vtkWriter.addField(problem->getAnalyticalVelocitySolution(), "velocityExact", GridView::dimensionworld);
+ vtkWriter.addField(problem->getAnalyticalPressureSolution(), "pressureExact");
+ vtkWriter.addField(problem->getAnalyticalVelocitySolution(), "velocityExact");
vtkWriter.addFaceField(problem->getAnalyticalVelocitySolutionOnFace(), "faceVelocityExact");
vtkWriter.write(0.0);