[io][vtkoutput] let velocityOutput decide number of phase velocities

This commit will enable the geomechanics framework to not have to
state the number of phases in geomechanical model traits. Also, when
velocity output for the geomechanical models is implemented, we will
have to distinguish in the vtk output module anyway, as then we would
have to call numSolidPhases(). Outsourcing it directly to the velocity
output seems like an elegant solution.

dumux/discretization/staggered/freeflow/velocityoutput.hh

@@ -82,6 +82,10 @@ public:
     static std::string phaseName(int phaseIdx)
     { return GET_PROP_TYPE(TypeTag, FluidSystem)::phaseName(phaseIdx); }
 
+    // returns the number of phase velocities computed by this class
+    static constexpr int numPhaseVelocities()
+    { return GET_PROP_TYPE(TypeTag, ModelTraits)::numPhases(); }
+
     //! Return the problem boundary types
     auto problemBoundaryTypes(const Element& element, const SubControlVolumeFace& scvf) const
     { return problem_.boundaryTypes(element, scvf); }

dumux/io/vtkoutputmodule.hh

@@ -68,8 +68,7 @@ class VtkOutputModule
     using SolutionVector = typename GET_PROP_TYPE(TypeTag, SolutionVector);
     using FVGridGeometry = typename GET_PROP_TYPE(TypeTag, FVGridGeometry);
     using VelocityOutput = typename GET_PROP_TYPE(TypeTag, VelocityOutput);
-    using ModelTraits = typename GET_PROP_TYPE(TypeTag, ModelTraits);
-    static constexpr int numPhases = ModelTraits::numPhases();
+    static constexpr int numPhaseVelocities = VelocityOutput::numPhaseVelocities();
 
     using VV = typename GridVariables::VolumeVariables;
     using Scalar = typename GridVariables::Scalar;
@@ -263,7 +262,7 @@ private:
 
         // instatiate the velocity output
         VelocityOutput velocityOutput(problem_, gridGeom_, gridVariables_, sol_);
-        std::array<std::vector<VelocityVector>, numPhases> velocity;
+        std::array<std::vector<VelocityVector>, numPhaseVelocities> velocity;
 
         // process rank
         static bool addProcessRank = getParamFromGroup<bool>(paramGroup_, "Vtk.AddProcessRank");
@@ -291,7 +290,7 @@ private:
 
             if (velocityOutput.enableOutput())
             {
-                for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                 {
                     if(isBox && dim == 1)
                         velocity[phaseIdx].resize(numCells);
@@ -343,7 +342,7 @@ private:
 
                 // velocity output
                 if (velocityOutput.enableOutput())
-                    for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                    for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                         velocityOutput.calculateVelocity(velocity[phaseIdx], elemVolVars, fvGeometry, element, phaseIdx);
 
                 //! the rank
@@ -378,7 +377,7 @@ private:
             {
                 if (isBox && dim > 1)
                 {
-                    for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                    for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                         sequenceWriter_.addVertexData( Field(gridGeom_.gridView(), gridGeom_.vertexMapper(), velocity[phaseIdx],
                                                              "velocity_" + velocityOutput.phaseName(phaseIdx+phaseIdxOffset) + " (m/s)",
                                                              /*numComp*/dimWorld, /*codim*/dim).get() );
@@ -386,7 +385,7 @@ private:
                 // cell-centered models
                 else
                 {
-                    for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                    for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                         sequenceWriter_.addCellData( Field(gridGeom_.gridView(), gridGeom_.elementMapper(), velocity[phaseIdx],
                                                            "velocity_" + velocityOutput.phaseName(phaseIdx+phaseIdxOffset) + " (m/s)",
                                                            /*numComp*/dimWorld, /*codim*/0).get() );
@@ -432,7 +431,7 @@ private:
 
         // instatiate the velocity output
         VelocityOutput velocityOutput(problem_, gridGeom_, gridVariables_, sol_);
-        std::array<std::vector<VelocityVector>, numPhases> velocity;
+        std::array<std::vector<VelocityVector>, numPhaseVelocities> velocity;
 
         // process rank
         static bool addProcessRank = getParamFromGroup<bool>(paramGroup_, "Vtk.AddProcessRank");
@@ -462,7 +461,7 @@ private:
 
             if (velocityOutput.enableOutput())
             {
-                for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                 {
                     if(isBox && dim == 1)
                         velocity[phaseIdx].resize(numCells);
@@ -520,7 +519,7 @@ private:
 
                 // velocity output
                 if (velocityOutput.enableOutput())
-                    for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                    for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                         velocityOutput.calculateVelocity(velocity[phaseIdx], elemVolVars, fvGeometry, element, phaseIdx);
 
                 //! the rank
@@ -546,14 +545,14 @@ private:
             {
                 // node-wise velocities
                 if (dim > 1)
-                    for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                    for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                         sequenceWriter_.addVertexData( Field(gridGeom_.gridView(), gridGeom_.vertexMapper(), velocity[phaseIdx],
                                                              "velocity_" + velocityOutput.phaseName(phaseIdx+phaseIdxOffset) + " (m/s)",
                                                              /*numComp*/dimWorld, /*codim*/dim).get() );
 
                 // cell-wise velocities
                 else
-                    for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
+                    for (int phaseIdx = 0; phaseIdx < numPhaseVelocities; ++phaseIdx)
                         sequenceWriter_.addCellData( Field(gridGeom_.gridView(), gridGeom_.elementMapper(), velocity[phaseIdx],
                                                            "velocity_" + velocityOutput.phaseName(phaseIdx+phaseIdxOffset) + " (m/s)",
                                                            /*numComp*/dimWorld, /*codim*/0).get() );

dumux/porousmediumflow/velocityoutput.hh

@@ -111,6 +111,10 @@ public:
     static std::string phaseName(int phaseIdx)
     { return FluidSystem::phaseName(phaseIdx); }
 
+    // returns the number of phase velocities computed by this class
+    static constexpr int numPhaseVelocities()
+    { return GET_PROP_TYPE(TypeTag, ModelTraits)::numPhases(); }
+
     // The following SFINAE enable_if usage allows compilation, even if only a
     //
     // boundaryTypes(const Element&, const scv&)