[mm][heatpipe] adapt to new material laws

lecture/mm/heatpipe/heatpipespatialparams.hh

@@ -38,47 +38,40 @@ class HeatPipeSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::KrPcHeatPipeDefault<Scalar>;
+
 public:
     using PermeabilityType = Scalar;
-    using MaterialLaw = KrPcHeatpipe<Scalar>;
-    using MaterialLawParams = typename MaterialLaw::Params;
 
     HeatPipeSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
     : ParentType(fvGridGeometry)
     {
         permeability_ = getParam<Scalar>("Problem.Permeability");
         porosity_ = 0.4;
+        Scalar p0 = std::pow((porosity_/permeability_), 0.5);
 
-        materialParams_.setSwr(0.15);
-        materialParams_.setSnr(0.0);
-        materialParams_.setP0(std::pow((porosity_/permeability_), 0.5));
+        typename PcKrSwCurve::BasicParams params(0.15, 0.0, p0);
+        pcKrSwCurve_ = std::make_unique<PcKrSwCurve>(params);
     }
 
     PermeabilityType permeabilityAtPos(const GlobalPosition& globalPos) const
-    {
-        return permeability_;
-    }
+    { return permeability_; }
 
     Scalar porosityAtPos(const GlobalPosition& globalPos) const
-    {
-        return porosity_;
-    }
+    { return porosity_; }
+
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    { return makeFluidMatrixInteraction(*pcKrSwCurve_); }
 
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
-        return materialParams_;
-    }
 
     template<class FluidSystem>
     int wettingPhaseAtPos(const GlobalPosition& globalPos) const
-    {
-        return FluidSystem::phase0Idx;
-    }
+    { return FluidSystem::phase0Idx; }
 
 private:
     PermeabilityType permeability_;
     Scalar porosity_;
-    MaterialLawParams materialParams_;
+    std::unique_ptr<const PcKrSwCurve> pcKrSwCurve_;
 };
 
 }

lecture/mm/heatpipe/krpcheatpipe.hh

@@ -19,13 +19,15 @@
 #ifndef DUMUX_KR_PC_HEATPIPE_HH
 #define DUMUX_KR_PC_HEATPIPE_HH
 
-#include <dumux/common/parameters.hh>
-#include "krpcheatpipeparams.hh"
-
 #include <algorithm>
+#include <cmath>
 
-namespace Dumux
-{
+#include <dumux/common/parameters.hh>
+#include <dumux/common/spline.hh>
+#include <dumux/common/optionalscalar.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/materiallaw.hh>
+
+namespace Dumux::FluidMatrix {
 /*!
  *
  * \brief Implementation of the capillary pressure <-> saturation
@@ -34,18 +36,57 @@ namespace Dumux
  *        cap-press based on the function of Leverett.
  *
  */
-template <class ScalarT, class ParamsT = KrPcHeatpipeParams<ScalarT> >
-class KrPcHeatpipe
+class KrPcHeatPipe
 {
 public:
-    using Params = ParamsT;
-    using Scalar = ScalarT;
+    template<class Scalar>
+    struct Params
+    {
+        Params(Scalar swr, Scalar snr, Scalar p0)
+        : swr_(swr)
+        , snr_(snr)
+        , p0_(p0)
+        {}
+
+        Scalar swr() const{ return swr_; }
+        void setSwr(Scalar swr){ swr_ = swr; }
+
+        Scalar snr() const { return snr_; }
+        void setSnr(Scalar snr) { snr_ = snr; }
+
+        Scalar p0() const { return p0_; }
+        void setp0(Scalar p0) { p0_ = p0; }
+
+        bool operator== (const Params& p) const
+        {
+            return Dune::FloatCmp::eq(snr(), p.snr(), 1e-6)
+                   && Dune::FloatCmp::eq(snr(), p.snr(), 1e-6)
+                   && Dune::FloatCmp::eq(p0(), p.p0(), 1e-6);
+        }
+
+    private:
+        Scalar swr_, snr_, p0_;
+    };
+
+    /*!
+     * \brief Construct from a subgroup from the global parameter tree
+     * \note This will give you nice error messages if a mandatory parameter is missing
+     */
+    template<class Scalar = double>
+    static Params<Scalar> makeParams(const std::string& paramGroup)
+    {
+        const auto swr = getParamFromGroup<Scalar>(paramGroup, "Swr");
+        const auto snr = getParamFromGroup<Scalar>(paramGroup, "Snr");
+        const auto p0 = getParamFromGroup<Scalar>(paramGroup, "P0");
+        return {swr, snr, p0};
+    }
 
     /*!
      * \brief The capillary pressure-saturation curve according to Leverett.
      *
      */
-    static Scalar pc(const Params &params, Scalar sw)
+    template<class Scalar>
+    static Scalar pc(Scalar sw, const Params<Scalar>& params)
     {
         if(sw<0.) sw=0.;
          /* effective values */
@@ -72,7 +113,8 @@ public:
      *        the medium according to the Fatt-Klikoff
      *        parameterization.
      */
-    static Scalar krw(const Params &params, Scalar sw)
+    template<class Scalar>
+    static Scalar krw(Scalar sw, const Params<Scalar>& params)
     {
          /*** effective saturation ***/
          Scalar Se = (sw-params.swr())/(1-params.snr()-params.swr());
@@ -88,7 +130,8 @@ public:
      *        the medium according to the Fatt-Klikoff
      *        parameterization.
      */
-    static Scalar krn(const Params &params, Scalar sw)
+    template<class Scalar>
+    static Scalar krn(Scalar sw, const Params<Scalar>& params)
     {
          /*** effective saturation ***/
          Scalar Se = (sw-params.swr())/(1-params.snr()-params.swr());
@@ -98,9 +141,15 @@ public:
          /* compute and return value */
          return std::pow(1.-Se,3);
     };
+};
 
+/*!
+ * \ingroup Fluidmatrixinteractions
+ * \brief A default configuration for using the VanGenuchten material law
+ */
+template<typename Scalar = double>
+using KrPcHeatPipeDefault = TwoPMaterialLaw<Scalar, KrPcHeatPipe, NoRegularization, TwoPEffToAbsDefaultPolicy>;
 
-};
 }
 
 #endif

lecture/mm/heatpipe/krpcheatpipeparams.hh

-// -*- 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 Specification of the material parameters
- *       for the Brooks Corey constitutive relations.
- */
-#ifndef DUMUX_KR_PC_HEATPIPE_PARAMS_HH
-#define DUMUX_KR_PC_HEATPIPE_PARAMS_HH
-
-namespace Dumux {
-
-/*!
- * \brief Specification of the material parameters
- *       for the kr-pc constitutive relations for the heatpipe problem.
- *
- */
-template <class ScalarT>
-class KrPcHeatpipeParams
-{
-public:
-    typedef ScalarT Scalar;
-
-    KrPcHeatpipeParams() = default;
-
-    KrPcHeatpipeParams(Scalar swr, Scalar snr)
-        : swr_(swr), snr_(snr)
-    {
-    }
-
-    /*!
-     * \brief Returns the entry pressure [Pa]
-     */
-    Scalar swr() const
-    { return swr_; }
-
-    /*!
-     * \brief Set the entry pressure [Pa]
-     */
-    void setSwr(Scalar v)
-    { swr_ = v; }
-
-
-    /*!
-     * \brief Returns the lambda shape parameter
-     */
-    Scalar snr() const
-    { return snr_; }
-
-    /*!
-     * \brief Set the lambda shape parameter
-     */
-    void setSnr(Scalar v)
-    { snr_ = v; }
-
-    /*!
-     * \brief Returns the lambda shape parameter
-     */
-    Scalar p0() const
-    { return p0_; }
-
-    /*!
-     * \brief Set the lambda shape parameter
-     */
-    void setP0(Scalar v)
-    { p0_ = v; }
-
-private:
-    Scalar swr_;
-    Scalar snr_;
-    Scalar p0_;
-};
-} // namespace Dumux
-
-#endif