From 15f060d11355b8fa0d1ba699a681586fb05fff75 Mon Sep 17 00:00:00 2001
From: Kilian Weishaupt <kilian.weishaupt@iws.uni-stuttgart.de>
Date: Wed, 11 Jan 2017 15:54:21 +0100
Subject: [PATCH] [staggeredGrid] Use standard ImplicitVolumeVariables

* Make use of [] operator of staggeredPriVars class
---
 .../staggered/globalvolumevariables.hh        |  14 +-
 .../staggered/volumevariables.hh              | 242 ------------------
 dumux/freeflow/staggered/volumevariables.hh   |  12 +-
 dumux/implicit/staggered/localjacobian.hh     |  13 +-
 dumux/implicit/staggered/primaryvariables.hh  |  51 +++-
 dumux/implicit/staggered/propertydefaults.hh  |   1 -
 6 files changed, 66 insertions(+), 267 deletions(-)
 delete mode 100644 dumux/discretization/staggered/volumevariables.hh

diff --git a/dumux/discretization/staggered/globalvolumevariables.hh b/dumux/discretization/staggered/globalvolumevariables.hh
index eaede1811d..310993f56c 100644
--- a/dumux/discretization/staggered/globalvolumevariables.hh
+++ b/dumux/discretization/staggered/globalvolumevariables.hh
@@ -55,6 +55,7 @@ class StaggeredGlobalVolumeVariables<TypeTag, /*enableGlobalVolVarsCache*/true>
     using VolumeVariables = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
     using ElementVolumeVariables = typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables);
     using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
+    using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
     using IndexType = typename GridView::IndexSet::IndexType;
 
     using DofTypeIndices = typename GET_PROP(TypeTag, DofTypeIndices);
@@ -64,6 +65,7 @@ class StaggeredGlobalVolumeVariables<TypeTag, /*enableGlobalVolVarsCache*/true>
     static const int dim = GridView::dimension;
     using Element = typename GridView::template Codim<0>::Entity;
     using CellCenterPrimaryVariables = typename GET_PROP_TYPE(TypeTag, CellCenterPrimaryVariables);
+    using PrimaryVariables = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
 
     enum { numEqCellCenter = GET_PROP_VALUE(TypeTag, NumEqCellCenter) };
 
@@ -82,7 +84,11 @@ public:
             fvGeometry.bindElement(element);
 
             for (auto&& scv : scvs(fvGeometry))
-                volumeVariables_[scv.index()].update(sol[cellCenterIdx][scv.dofIndex()], problem, element, scv);
+            {
+                PrimaryVariables priVars(0.0);
+                priVars[cellCenterIdx] = sol[cellCenterIdx][scv.dofIndex()];
+                volumeVariables_[scv.index()].update(priVars, problem, element, scv);
+            }
 
             // handle the boundary volume variables
             for (auto&& scvf : scvfs(fvGeometry))
@@ -95,14 +101,14 @@ public:
                 const auto insideScvIdx = scvf.insideScvIdx();
                 const auto& insideScv = fvGeometry.scv(insideScvIdx);
 
-                auto boundaryPriVars = CellCenterPrimaryVariables(0.0);
+                PrimaryVariables boundaryPriVars(0.0);
 
                 for(int eqIdx = 0; eqIdx < numEqCellCenter; ++eqIdx)
                 {
                     if(bcTypes.isDirichlet(eqIdx))
-                        boundaryPriVars[eqIdx] = problem.dirichlet(element, scvf)[eqIdx];
+                        boundaryPriVars[cellCenterIdx][eqIdx] = problem.dirichlet(element, scvf)[eqIdx];
                     else if(bcTypes.isNeumann(eqIdx))
-                        boundaryPriVars[eqIdx] = sol[cellCenterIdx][scvf.insideScvIdx()][eqIdx];
+                        boundaryPriVars[cellCenterIdx][eqIdx] = sol[cellCenterIdx][scvf.insideScvIdx()][eqIdx];
                     //TODO: this assumes a zero-gradient for e.g. the pressure on the boundary
                     // could be made more general by allowing a non-zero-gradient, provided in problem file
                     else
diff --git a/dumux/discretization/staggered/volumevariables.hh b/dumux/discretization/staggered/volumevariables.hh
deleted file mode 100644
index fbd8b4bc85..0000000000
--- a/dumux/discretization/staggered/volumevariables.hh
+++ /dev/null
@@ -1,242 +0,0 @@
-// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
-// vi: set et ts=4 sw=4 sts=4:
-/*****************************************************************************
- *   See the file COPYING for full copying permissions.                      *
- *                                                                           *
- *   This program is free software: you can redistribute it and/or modify    *
- *   it under the terms of the GNU General Public License as published by    *
- *   the Free Software Foundation, either version 2 of the License, or       *
- *   (at your option) any later version.                                     *
- *                                                                           *
- *   This program is distributed in the hope that it will be useful,         *
- *   but WITHOUT ANY WARRANTY; without even the implied warranty of          *
- *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the            *
- *   GNU General Public License for more details.                            *
- *                                                                           *
- *   You should have received a copy of the GNU General Public License       *
- *   along with this program.  If not, see <http://www.gnu.org/licenses/>.   *
- *****************************************************************************/
-/*!
- * \file
- * \brief Base class for the model specific class which provides
- *        access to all volume averaged quantities.
- */
-#ifndef DUMUX_DISCRETIZATION_STAGGERED_VOLUME_VARIABLES_HH
-#define DUMUX_DISCRETIZATION_STAGGERED_VOLUME_VARIABLES_HH
-
-#include <dumux/implicit/properties.hh>
-#include <dumux/common/valgrind.hh>
-
-namespace Dumux
-{
-
-namespace Properties
-{
-NEW_PROP_TAG(FluidSystem);
-NEW_PROP_TAG(Indices);
-NEW_PROP_TAG(EnableEnergyTransport);
-}
-
-// forward declaration
-template <class TypeTag, bool enableEnergyBalance>
-class StaggeredVolumeVariablesImplementation;
-
-/*!
- * \ingroup ImplicitVolumeVariables
- * \brief Base class for the model specific class which provides
- *        access to all volume averaged quantities. The volume variables base class
- *        is specialized for isothermal and non-isothermal models.
- */
-template <class TypeTag>
-using StaggeredVolumeVariables = StaggeredVolumeVariablesImplementation<TypeTag, GET_PROP_VALUE(TypeTag, EnableEnergyTransport)>;
-
-/*!
- * \ingroup StaggeredVolumeVariables
- * \brief The isothermal base class
- */
-template<class TypeTag>
-class StaggeredVolumeVariablesImplementation<TypeTag, false>
-{
-    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
-    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 SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
-    using CellCenterPrimaryVariables = typename GET_PROP_TYPE(TypeTag, CellCenterPrimaryVariables);
-    using FacePrimaryVariables = typename GET_PROP_TYPE(TypeTag, FacePrimaryVariables);
-
-public:
-
-    /*!
-     * \brief Update all quantities for a given control volume
-     *
-     * \param priVars A vector containing the primary variables for the control volume
-     * \param problem The object specifying the problem which ought to
-     *                be simulated
-     * \param element An element which contains part of the control volume
-     * \param fvGeometry The finite volume geometry for the element
-     * \param scvIdx Local index of the sub control volume which is inside the element
-     * \param isOldSol Specifies whether this is the previous solution or the current one
-     *
-     * \todo Eliminate the 'isOldSol' parameter. This implies that the
-     *       'pseudo-primary variables' must be somehow be stored
-     *       inside the PrimaryVariables. (e.g. we need to know the
-     *       phase state in the 2p2c model)
-     */
-    void update(const CellCenterPrimaryVariables &ccPriVars,
-                const Problem &problem,
-                const Element &element,
-                const SubControlVolume &scv)
-    {
-        extrusionFactor_ = problem.boxExtrusionFactor(element, scv);
-        ccPriVars_ = ccPriVars;
-    }
-
-
-
-    /*!
-     * \brief Return the vector of primary variables
-     */
-    CellCenterPrimaryVariables ccPriVars() const
-    { return ccPriVars_; }
-
-    /*!
-     * \brief Return a component of primary variable vector
-     *
-     * \param pvIdx The index of the primary variable of interest
-     */
-    Scalar ccPriVar(const int pvIdx) const
-    { return ccPriVars_[pvIdx]; }
-
-
-    /*!
-     * \brief Return how much the sub-control volume is extruded.
-     *
-     * This means the factor by which a lower-dimensional (1D or 2D)
-     * entity needs to be expanded to get a full dimensional cell. The
-     * default is 1.0 which means that 1D problems are actually
-     * thought as pipes with a cross section of 1 m^2 and 2D problems
-     * are assumed to extend 1 m to the back.
-     */
-    Scalar extrusionFactor() const
-    { return extrusionFactor_; }
-
-    //! The temperature is obtained from the problem as a constant for isothermal models
-    static Scalar temperature(const CellCenterPrimaryVariables &priVars,
-                              const Problem& problem,
-                              const Element &element,
-                              const SubControlVolume &scv)
-    {
-        return problem.temperatureAtPos(scv.dofPosition());
-    }
-
-    //! The phase enthalpy is zero for isothermal models
-    //! This is needed for completing the fluid state
-    template<class FluidState, class ParameterCache>
-    static Scalar enthalpy(const FluidState& fluidState,
-                           const ParameterCache& paramCache,
-                           const int phaseIdx)
-    {
-        return 0;
-    }
-
-private:
-    Scalar extrusionFactor_;
-
-    CellCenterPrimaryVariables ccPriVars_;
-};
-
-//! The non-isothermal implicit volume variables base class
-template <class TypeTag>
-class StaggeredVolumeVariablesImplementation<TypeTag, true>
-: public StaggeredVolumeVariablesImplementation<TypeTag, false>
-{
-    using ParentType = StaggeredVolumeVariablesImplementation<TypeTag, false>;
-    using Implementation = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
-    using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
-    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 SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
-    using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
-    using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
-    using CellCenterPrimaryVariables = typename GET_PROP_TYPE(TypeTag, CellCenterPrimaryVariables);
-
-    enum { temperatureIdx = Indices::temperatureIdx };
-
-public:
-
-    /*!
-     * \brief Update all quantities for a given control volume
-     *
-     * \param priVars A vector containing the primary variables for the control volume
-     * \param problem The object specifying the problem which ought to
-     *                be simulated
-     * \param element An element which contains part of the control volume
-     * \param fvGeometry The finite volume geometry for the element
-     * \param scvIdx Local index of the sub control volume which is inside the element
-     */
-    void update(const CellCenterPrimaryVariables &ccPriVars,
-                const Problem &problem,
-                const Element &element,
-                const SubControlVolume &scv)
-    {
-        ParentType::update(ccPriVars, problem, element, scv);
-    }
-
-    /*!
-     * \brief Returns the total internal energy of a phase in the
-     *        sub-control volume.
-     *
-     * \param phaseIdx The phase index
-     */
-    Scalar internalEnergy(const int phaseIdx) const
-    { return asImp_().fluidState().internalEnergy(phaseIdx); }
-
-    /*!
-     * \brief Returns the total enthalpy of a phase in the sub-control
-     *        volume.
-     *
-     * \param phaseIdx The phase index
-     */
-    Scalar enthalpy(const int phaseIdx) const
-    { return asImp_().fluidState().enthalpy(phaseIdx); }
-
-
-    /*!
-     * \brief Returns the thermal conductivity \f$\mathrm{[W/(m*K)]}\f$ of a fluid phase in
-     *        the sub-control volume.
-     */
-    Scalar fluidThermalConductivity(const int phaseIdx) const
-    { return FluidSystem::thermalConductivity(asImp_().fluidState(), phaseIdx); }
-
-
-    //! The temperature is a primary variable for non-isothermal models
-    static Scalar temperature(const CellCenterPrimaryVariables &priVars,
-                              const Problem& problem,
-                              const Element &element,
-                              const SubControlVolume &scv)
-    {
-        return priVars[temperatureIdx];
-    }
-
-    //! The phase enthalpy is zero for isothermal models
-    //! This is needed for completing the fluid state
-    template<class FluidState, class ParameterCache>
-    static Scalar enthalpy(const FluidState& fluidState,
-                           const ParameterCache& paramCache,
-                           const int phaseIdx)
-    {
-        return FluidSystem::enthalpy(fluidState, paramCache, phaseIdx);
-    }
-
-protected:
-    const Implementation &asImp_() const
-    { return *static_cast<const Implementation*>(this); }
-    Implementation &asImp_()
-    { return *static_cast<Implementation*>(this); }
-};
-
-} // end namespace Dumux
-
-#endif
diff --git a/dumux/freeflow/staggered/volumevariables.hh b/dumux/freeflow/staggered/volumevariables.hh
index 39b84dd995..10273615e4 100644
--- a/dumux/freeflow/staggered/volumevariables.hh
+++ b/dumux/freeflow/staggered/volumevariables.hh
@@ -25,7 +25,7 @@
 #define DUMUX_1P_VOLUME_VARIABLES_HH
 
 #include "properties.hh"
-#include <dumux/discretization/staggered/volumevariables.hh>
+#include <dumux/discretization/volumevariables.hh>
 
 #include <dumux/material/fluidstates/immiscible.hh>
 
@@ -39,15 +39,15 @@ namespace Dumux
  *        finite volume in the one-phase model.
  */
 template <class TypeTag>
-class NavierStokesVolumeVariables : public StaggeredVolumeVariables<TypeTag>
+class NavierStokesVolumeVariables : public ImplicitVolumeVariables<TypeTag>
 {
-    using ParentType = StaggeredVolumeVariables<TypeTag>;
+    using ParentType = ImplicitVolumeVariables<TypeTag>;
     using Implementation = typename GET_PROP_TYPE(TypeTag, VolumeVariables);
     using Scalar = typename GET_PROP_TYPE(TypeTag, Scalar);
     using Problem = typename GET_PROP_TYPE(TypeTag, Problem);
     using FVElementGeometry = typename GET_PROP_TYPE(TypeTag, FVElementGeometry);
     using SubControlVolume = typename GET_PROP_TYPE(TypeTag, SubControlVolume);
-    using CellCenterPrimaryVariables = typename GET_PROP_TYPE(TypeTag, CellCenterPrimaryVariables);
+    using PrimaryVariables = typename GET_PROP_TYPE(TypeTag, PrimaryVariables);
     using Indices = typename GET_PROP_TYPE(TypeTag, Indices);
     using FluidSystem = typename GET_PROP_TYPE(TypeTag, FluidSystem);
     using GridView = typename GET_PROP_TYPE(TypeTag, GridView);
@@ -60,7 +60,7 @@ public:
     /*!
      * \copydoc ImplicitVolumeVariables::update
      */
-    void update(const CellCenterPrimaryVariables &priVars,
+    void update(const PrimaryVariables &priVars,
                 const Problem &problem,
                 const Element &element,
                 const SubControlVolume& scv)
@@ -73,7 +73,7 @@ public:
     /*!
      * \copydoc ImplicitModel::completeFluidState
      */
-    static void completeFluidState(const CellCenterPrimaryVariables& priVars,
+    static void completeFluidState(const PrimaryVariables& priVars,
                                    const Problem& problem,
                                    const Element& element,
                                    const SubControlVolume& scv,
diff --git a/dumux/implicit/staggered/localjacobian.hh b/dumux/implicit/staggered/localjacobian.hh
index 94580cf6b5..7b3fd86e2b 100644
--- a/dumux/implicit/staggered/localjacobian.hh
+++ b/dumux/implicit/staggered/localjacobian.hh
@@ -255,9 +255,11 @@ private:
             for(int pvIdx = 0; pvIdx < numEqCellCenter; ++pvIdx)
             {
                 const Scalar eps = 1e-4; // TODO: do properly
-                CellCenterPrimaryVariables priVars(this->model_().curSol()[cellCenterIdx][globalJ]);
 
-                priVars[pvIdx] += eps;
+                PrimaryVariables priVars(CellCenterPrimaryVariables(this->model_().curSol()[cellCenterIdx][globalJ]),
+                                         FacePrimaryVariables(0.0));
+
+                priVars[cellCenterIdx][pvIdx] += eps;
                 curVolVars.update(priVars, this->problem_(), elementJ, scvJ);
 
                 this->localResidual().evalCellCenter(element, fvGeometry, scvI,
@@ -359,9 +361,12 @@ private:
                 for(int pvIdx = 0; pvIdx < numEqCellCenter; ++pvIdx)
                 {
                     const Scalar eps = 1e-4; // TODO: do properly
-                    CellCenterPrimaryVariables priVars(this->model_().curSol()[cellCenterIdx][globalJ]);
 
-                    priVars[pvIdx] += eps;
+                    PrimaryVariables priVars(CellCenterPrimaryVariables(this->model_().curSol()[cellCenterIdx][globalJ]),
+                                             FacePrimaryVariables(0.0));
+
+                    priVars[cellCenterIdx][pvIdx] += eps;
+
                     curVolVars.update(priVars, this->problem_(), elementJ, scvJ);
 
                     this->localResidual().evalFace(element, fvGeometry, scvf,
diff --git a/dumux/implicit/staggered/primaryvariables.hh b/dumux/implicit/staggered/primaryvariables.hh
index ed9ccb279d..f93d9f53e1 100644
--- a/dumux/implicit/staggered/primaryvariables.hh
+++ b/dumux/implicit/staggered/primaryvariables.hh
@@ -44,31 +44,62 @@ class StaggeredPrimaryVariables : public Dune::MultiTypeBlockVector<CellCenterPr
     typename DofTypeIndices::CellCenterIdx cellCenterIdx;
     typename DofTypeIndices::FaceIdx faceIdx;
 
-    using PrimaryVariables = Dune::MultiTypeBlockVector<CellCenterPrimaryVariables, FacePrimaryVariables>;
+    using ParentType = Dune::MultiTypeBlockVector<CellCenterPrimaryVariables, FacePrimaryVariables>;
 
 public:
     StaggeredPrimaryVariables() = default;
 
-    StaggeredPrimaryVariables(const Scalar value)
+     /*!
+     * \brief Constructor to initialize all entries with the same value
+     *
+     * \param value The value
+     */
+    StaggeredPrimaryVariables(const Scalar value) noexcept
     {
         (*this)[cellCenterIdx] = value;
         (*this)[faceIdx] = value;
     }
 
-    const auto& operator [](const unsigned int i) const
+     /*!
+     * \brief Constructor to initialize the cellcenter and face primary values with given values
+     *
+     * \param ccPriVars The cellcenter primary variables used for initialization
+     * \param facePriVars The face primary variables used for initialization
+     */
+    StaggeredPrimaryVariables(CellCenterPrimaryVariables&& ccPriVars, FacePrimaryVariables&& facePriVars) noexcept
     {
-        if(i < Indices::faceOffset)
-            return PrimaryVariables::operator[](cellCenterIdx)[i];
+        (*this)[cellCenterIdx] = std::move(ccPriVars);
+        (*this)[faceIdx] = std::move(facePriVars);
+    }
+
+     /*!
+     * \brief Operator overload which allows to automatically access the "right" priVars vector via pvIdx.
+     *        const version
+     * \note: the ParentType (DUNE multitypeblockvector) [] operator has to be visible (using ...)
+     *
+     * \param pvIdx The global index of the primary variable
+     */
+    using ParentType::operator [];
+    const Scalar& operator [](const unsigned int pvIdx) const
+    {
+        if(pvIdx < Indices::faceOffset)
+            return ParentType::operator[](cellCenterIdx)[pvIdx];
         else
-            return PrimaryVariables::operator[](faceIdx)[i - Indices::faceOffset];
+            return ParentType::operator[](faceIdx)[pvIdx - Indices::faceOffset];
     }
 
-    auto& operator [](const unsigned int i)
+     /*!
+     * \brief Operator overload which allows to automatically access the "right" priVars vector via pvIdx
+     *        non-const version
+     *
+     * \param pvIdx The global index of the primary variable
+     */
+    Scalar& operator [](const unsigned int pvIdx)
     {
-        if(i < Indices::faceOffset)
-            return PrimaryVariables::operator[](cellCenterIdx)[i];
+        if(pvIdx < Indices::faceOffset)
+            return ParentType::operator[](cellCenterIdx)[pvIdx];
         else
-            return PrimaryVariables::operator[](faceIdx)[i - Indices::faceOffset];
+            return ParentType::operator[](faceIdx)[pvIdx - Indices::faceOffset];
     }
 };
 
diff --git a/dumux/implicit/staggered/propertydefaults.hh b/dumux/implicit/staggered/propertydefaults.hh
index f21482a04a..f2739f0a57 100644
--- a/dumux/implicit/staggered/propertydefaults.hh
+++ b/dumux/implicit/staggered/propertydefaults.hh
@@ -39,7 +39,6 @@
 
 #include <dumux/discretization/staggered/elementvolumevariables.hh>
 #include <dumux/discretization/staggered/globalvolumevariables.hh>
-#include <dumux/discretization/staggered/volumevariables.hh>
 #include <dumux/discretization/staggered/globalfacevariables.hh>
 
 #include <dune/common/fvector.hh>
-- 
GitLab