From da6b489f05fac57f0b0436ec779d770a586248b3 Mon Sep 17 00:00:00 2001
From: Timo Koch <timo.koch@iws.uni-stuttgart.de>
Date: Wed, 2 Aug 2017 13:37:22 +0200
Subject: [PATCH] [pmflow] No need for template specialization in the fluxvars.
 Greatly simplify.

---
 dumux/discretization/fluxvariablesbase.hh     |   2 +-
 .../implicit/fluxvariables.hh                 | 262 ++++--------------
 2 files changed, 57 insertions(+), 207 deletions(-)

diff --git a/dumux/discretization/fluxvariablesbase.hh b/dumux/discretization/fluxvariablesbase.hh
index b20dd3f15e..29a46c64de 100644
--- a/dumux/discretization/fluxvariablesbase.hh
+++ b/dumux/discretization/fluxvariablesbase.hh
@@ -101,7 +101,7 @@ public:
 
     //! Applies the upwind scheme to precalculated fluxes
     template<class UpwindTermFunction>
-    Scalar applyUpwindScheme(const UpwindTermFunction& upwindTerm, Scalar flux, int phaseIdx)
+    Scalar applyUpwindScheme(const UpwindTermFunction& upwindTerm, Scalar flux, int phaseIdx) const
     {
         //! Give the upwind scheme access to the cached variables
         return UpwindScheme::apply(*this, upwindTerm, flux, phaseIdx);
diff --git a/dumux/porousmediumflow/implicit/fluxvariables.hh b/dumux/porousmediumflow/implicit/fluxvariables.hh
index bb896cb687..54b21d7941 100644
--- a/dumux/porousmediumflow/implicit/fluxvariables.hh
+++ b/dumux/porousmediumflow/implicit/fluxvariables.hh
@@ -34,30 +34,20 @@ namespace Properties
 {
 // forward declaration
 NEW_PROP_TAG(NumPhases);
+NEW_PROP_TAG(NumComponents);
 NEW_PROP_TAG(EnableAdvection);
 NEW_PROP_TAG(EnableMolecularDiffusion);
 NEW_PROP_TAG(EnableEnergyBalance);
 }
 
-// forward declaration
-template<class TypeTag, bool enableAdvection, bool enableMolecularDiffusion, bool enableEnergyBalance>
-class PorousMediumFluxVariablesImpl;
-
 /*!
  * \ingroup ImplicitModel
- * \brief The flux variables class
- *        specializations are provided for combinations of physical processes
+ * \brief The porous medium flux variables class that computes advective / convective,
+ *        molecular diffusive and heat conduction fluxes.
  * \note  Not all specializations are currently implemented
  */
 template<class TypeTag>
-using PorousMediumFluxVariables = PorousMediumFluxVariablesImpl<TypeTag, GET_PROP_VALUE(TypeTag, EnableAdvection),
-                                                                         GET_PROP_VALUE(TypeTag, EnableMolecularDiffusion),
-                                                                         GET_PROP_VALUE(TypeTag, EnableEnergyBalance)>;
-
-
-// specialization for pure advective flow (e.g. 1p/2p/3p immiscible darcy flow)
-template<class TypeTag>
-class PorousMediumFluxVariablesImpl<TypeTag, true, false, false> : public FluxVariablesBase<TypeTag>
+class PorousMediumFluxVariables : public FluxVariablesBase<TypeTag>
 {
     using ParentType = FluxVariablesBase<TypeTag>;
     using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
@@ -67,238 +57,98 @@ class PorousMediumFluxVariablesImpl<TypeTag, true, false, false> : public FluxVa
     using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
     using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace);
     using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
-    using AdvectionType = typename GET_PROP_TYPE(TypeTag, AdvectionType);
     using VolumeVariables = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
     using ElementVolumeVariables = typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables);
     using ElementFluxVariablesCache = typename GET_PROP_TYPE(TypeTag, ElementFluxVariablesCache);
 
-    static const int dim = GridView::dimension;
-    static const int dimWorld = GridView::dimensionworld;
-
-public:
-
-    template<typename FunctionType>
-    Scalar advectiveFlux(const int phaseIdx, const FunctionType& upwindTerm)
-    {
-        Scalar flux = AdvectionType::flux(this->problem(),
-                                          this->element(),
-                                          this->fvGeometry(),
-                                          this->elemVolVars(),
-                                          this->scvFace(),
-                                          phaseIdx,
-                                          this->elemFluxVarsCache());
-
-        return this->applyUpwindScheme(upwindTerm, flux, phaseIdx);
-    }
-};
-
-
-// specialization for isothermal advection molecularDiffusion equations
-template<class TypeTag>
-class PorousMediumFluxVariablesImpl<TypeTag, true, true, false> : public FluxVariablesBase<TypeTag>
-{
-    using ParentType = FluxVariablesBase<TypeTag>;
-    using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
-    using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
-    using Element = typename GridView::template Codim<0>::Entity;
-    using IndexType = typename GridView::IndexSet::IndexType;
-    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
-    using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace);
-    using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
-    using ElementVolumeVariables = typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables);
-    using ElementFluxVariablesCache = typename GET_PROP_TYPE(TypeTag, ElementFluxVariablesCache);
-
     using AdvectionType = typename GET_PROP_TYPE(TypeTag, AdvectionType);
     using MolecularDiffusionType = typename GET_PROP_TYPE(TypeTag, MolecularDiffusionType);
+    using HeatConductionType = typename GET_PROP_TYPE(TypeTag, HeatConductionType);
 
     enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases),
            numComponents = GET_PROP_VALUE(TypeTag, NumComponents)
     };
 
+    static constexpr bool enableAdvection = GET_PROP_VALUE(TypeTag, EnableAdvection);
+    static constexpr bool enableMolecularDiffusion = GET_PROP_VALUE(TypeTag, EnableMolecularDiffusion);
+    static constexpr bool enableEnergyBalance = GET_PROP_VALUE(TypeTag, EnableEnergyBalance);
+
 public:
+
     //! The constructor
-    PorousMediumFluxVariablesImpl()
+    PorousMediumFluxVariables()
     {
-        advFluxCached_.reset();
-        advPreFlux_.fill(0.0);
+        advFluxIsCached_.reset();
+        advFluxBeforeUpwinding_.fill(0.0);
     }
 
     template<typename FunctionType>
-    Scalar advectiveFlux(const int phaseIdx, const FunctionType& upwindTerm)
+    Scalar advectiveFlux(const int phaseIdx, const FunctionType& upwindTerm) const
     {
-        if (!advFluxCached_[phaseIdx])
+        if (enableAdvection)
         {
-
-            advPreFlux_[phaseIdx] = AdvectionType::flux(this->problem(),
-                                                        this->element(),
-                                                        this->fvGeometry(),
-                                                        this->elemVolVars(),
-                                                        this->scvFace(),
-                                                        phaseIdx,
-                                                        this->elemFluxVarsCache());
-            advFluxCached_.set(phaseIdx, true);
+            if (!advFluxIsCached_[phaseIdx])
+            {
+
+                advFluxBeforeUpwinding_[phaseIdx] = AdvectionType::flux(this->problem(),
+                                                                        this->element(),
+                                                                        this->fvGeometry(),
+                                                                        this->elemVolVars(),
+                                                                        this->scvFace(),
+                                                                        phaseIdx,
+                                                                        this->elemFluxVarsCache());
+                advFluxIsCached_.set(phaseIdx, true);
+            }
+
+            return this->applyUpwindScheme(upwindTerm, advFluxBeforeUpwinding_[phaseIdx], phaseIdx);
         }
-
-        return this->applyUpwindScheme(upwindTerm, advPreFlux_[phaseIdx], phaseIdx);
-    }
-
-    Dune::FieldVector<Scalar, numComponents> molecularDiffusionFlux(const int phaseIdx)
-    {
-        return MolecularDiffusionType::flux(this->problem(),
-                                            this->element(),
-                                            this->fvGeometry(),
-                                            this->elemVolVars(),
-                                            this->scvFace(),
-                                            phaseIdx,
-                                            this->elemFluxVarsCache());
-    }
-
-private:
-    //! simple caching if advection flux is used twice with different upwind function
-    std::bitset<numPhases> advFluxCached_;
-    std::array<Scalar, numPhases> advPreFlux_;
-};
-
-// specialization for non-isothermal advective flow (e.g. non-isothermal one-phase darcy equation)
-template<class TypeTag>
-class PorousMediumFluxVariablesImpl<TypeTag, true, false, true> : public FluxVariablesBase<TypeTag>
-{
-    using ParentType = FluxVariablesBase<TypeTag>;
-    using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
-    using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
-    using Element = typename GridView::template Codim<0>::Entity;
-    using IndexType = typename GridView::IndexSet::IndexType;
-    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
-    using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace);
-    using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
-    using ElementVolumeVariables = typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables);
-    using ElementFluxVariablesCache = typename GET_PROP_TYPE(TypeTag, ElementFluxVariablesCache);
-
-    using AdvectionType = typename GET_PROP_TYPE(TypeTag, AdvectionType);
-    using HeatConductionType = typename GET_PROP_TYPE(TypeTag, HeatConductionType);
-
-    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
-
-public:
-    //! The constructor
-    PorousMediumFluxVariablesImpl()
-    {
-        advFluxCached_.reset();
-        advPreFlux_.fill(0.0);
-    }
-
-    template<typename FunctionType>
-    Scalar advectiveFlux(const int phaseIdx, const FunctionType& upwindTerm)
-    {
-        if (!advFluxCached_[phaseIdx])
+        else
         {
-
-            advPreFlux_[phaseIdx] = AdvectionType::flux(this->problem(),
-                                                        this->element(),
-                                                        this->fvGeometry(),
-                                                        this->elemVolVars(),
-                                                        this->scvFace(),
-                                                        phaseIdx,
-                                                        this->elemFluxVarsCache());
-            advFluxCached_.set(phaseIdx, true);
+            return 0.0;
         }
-
-        return this->applyUpwindScheme(upwindTerm, advPreFlux_[phaseIdx], phaseIdx);
     }
 
-    Scalar heatConductionFlux()
+    Dune::FieldVector<Scalar, numComponents> molecularDiffusionFlux(const int phaseIdx) const
     {
-        return HeatConductionType::flux(this->problem(),
-                                        this->element(),
-                                        this->fvGeometry(),
-                                        this->elemVolVars(),
-                                        this->scvFace(),
-                                        this->elemFluxVarsCache());
-    }
-
-private:
-    //! simple caching if advection flux is used twice with different upwind function
-    std::bitset<numPhases> advFluxCached_;
-    std::array<Scalar, numPhases> advPreFlux_;
-};
-
-// specialization for non-isothermal advection and difussion equations (e.g. non-isothermal three-phase three-component flow)
-template<class TypeTag>
-class PorousMediumFluxVariablesImpl<TypeTag, true, true, true> : public FluxVariablesBase<TypeTag>
-{
-    using ParentType = FluxVariablesBase<TypeTag>;
-    using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
-    using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
-    using Element = typename GridView::template Codim<0>::Entity;
-    using IndexType = typename GridView::IndexSet::IndexType;
-    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
-    using SubControlVolumeFace = typename GET_PROP_TYPE(TypeTag, SubControlVolumeFace);
-    using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
-    using ElementVolumeVariables = typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables);
-    using ElementFluxVariablesCache = typename GET_PROP_TYPE(TypeTag, ElementFluxVariablesCache);
-
-    using AdvectionType = typename GET_PROP_TYPE(TypeTag, AdvectionType);
-    using MolecularDiffusionType = typename GET_PROP_TYPE(TypeTag, MolecularDiffusionType);
-    using HeatConductionType = typename GET_PROP_TYPE(TypeTag, HeatConductionType);
-
-    enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases),
-           numComponents = GET_PROP_VALUE(TypeTag, NumComponents)
-    };
-
-public:
-    //! The constructor
-    PorousMediumFluxVariablesImpl()
-    {
-        advFluxCached_.reset();
-        advPreFlux_.fill(0.0);
-    }
-
-    template<typename FunctionType>
-    Scalar advectiveFlux(const int phaseIdx, const FunctionType& upwindTerm)
-    {
-        if (!advFluxCached_[phaseIdx])
+        if (enableMolecularDiffusion)
         {
-
-            advPreFlux_[phaseIdx] = AdvectionType::flux(this->problem(),
-                                                        this->element(),
-                                                        this->fvGeometry(),
-                                                        this->elemVolVars(),
-                                                        this->scvFace(),
-                                                        phaseIdx,
-                                                        this->elemFluxVarsCache());
-            advFluxCached_.set(phaseIdx, true);
+            return MolecularDiffusionType::flux(this->problem(),
+                                                this->element(),
+                                                this->fvGeometry(),
+                                                this->elemVolVars(),
+                                                this->scvFace(),
+                                                phaseIdx,
+                                                this->elemFluxVarsCache());
+        }
+        else
+        {
+            return Dune::FieldVector<Scalar, numComponents>(0.0);
         }
-
-        return this->applyUpwindScheme(upwindTerm, advPreFlux_[phaseIdx], phaseIdx);
     }
 
-    Dune::FieldVector<Scalar, numComponents> molecularDiffusionFlux(const int phaseIdx)
+    Scalar heatConductionFlux() const
     {
-        return MolecularDiffusionType::flux(this->problem(),
+        if (enableEnergyBalance)
+        {
+            return HeatConductionType::flux(this->problem(),
                                             this->element(),
                                             this->fvGeometry(),
                                             this->elemVolVars(),
                                             this->scvFace(),
-                                            phaseIdx,
                                             this->elemFluxVarsCache());
-    }
-
-    Scalar heatConductionFlux()
-    {
-        return HeatConductionType::flux(this->problem(),
-                                        this->element(),
-                                        this->fvGeometry(),
-                                        this->elemVolVars(),
-                                        this->scvFace(),
-                                        this->elemFluxVarsCache());
+        }
+        else
+        {
+            return 0.0;
+        }
     }
 
 private:
     //! simple caching if advection flux is used twice with different upwind function
-    std::bitset<numPhases> advFluxCached_;
-    std::array<Scalar, numPhases> advPreFlux_;
+    mutable std::bitset<numPhases> advFluxIsCached_;
+    mutable std::array<Scalar, numPhases> advFluxBeforeUpwinding_;
 };
 
-} // end namespace
+} // end namespace Dumux
 
 #endif
-- 
GitLab