From 55126561afe97c5a26afbae71d900f3480f407c6 Mon Sep 17 00:00:00 2001
From: Kilian Weishaupt <kilian.weishaupt@iws.uni-stuttgart.de>
Date: Thu, 29 Oct 2020 21:50:04 +0100
Subject: [PATCH] [test][3p3c] Adapt spatial params

---
 .../3p3c/implicit/columnxylol/params.input    |  9 +++
 .../implicit/columnxylol/spatialparams.hh     | 59 ++++++-------------
 .../3p3c/implicit/infiltration/params.input   |  8 +++
 .../3p3c/implicit/infiltration/problem.hh     |  9 ++-
 .../implicit/infiltration/spatialparams.hh    | 37 +++---------
 5 files changed, 47 insertions(+), 75 deletions(-)

diff --git a/test/porousmediumflow/3p3c/implicit/columnxylol/params.input b/test/porousmediumflow/3p3c/implicit/columnxylol/params.input
index f502ce89f7..0da3b34d3b 100644
--- a/test/porousmediumflow/3p3c/implicit/columnxylol/params.input
+++ b/test/porousmediumflow/3p3c/implicit/columnxylol/params.input
@@ -18,3 +18,12 @@ SolidDensity = 2650
 SolidThermalConductivity = 2.8
 SolidHeatCapacityFine = 850
 SolidHeatCapacityCoarse = 84000
+
+[SpatialParams]
+Swr = 0.12
+Snr = 0.10
+Sgr = 0.01
+ParkerVanGenuchtenN = 4
+ParkerVanGenuchtenRegardSnrForKrn = true
+Fine.ParkerVanGenuchtenAlpha = 0.0005
+Coarse.ParkerVanGenuchtenAlpha = 0.5
diff --git a/test/porousmediumflow/3p3c/implicit/columnxylol/spatialparams.hh b/test/porousmediumflow/3p3c/implicit/columnxylol/spatialparams.hh
index 4de10b8032..3fade2a52b 100644
--- a/test/porousmediumflow/3p3c/implicit/columnxylol/spatialparams.hh
+++ b/test/porousmediumflow/3p3c/implicit/columnxylol/spatialparams.hh
@@ -26,9 +26,7 @@
 
 #include <dumux/porousmediumflow/properties.hh>
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/3p/regularizedparkervangen3p.hh>
-#include <dumux/material/fluidmatrixinteractions/3p/regularizedparkervangen3pparams.hh>
-#include <dumux/material/fluidmatrixinteractions/3p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/3p/parkervangenuchten.hh>
 
 namespace Dumux {
 
@@ -50,15 +48,16 @@ class ColumnSpatialParams
                                        ColumnSpatialParams<GridGeometry, Scalar>>;
 
     using GlobalPosition = typename SubControlVolume::GlobalPosition;
-    using EffectiveLaw = RegularizedParkerVanGen3P<Scalar>;
+
+    using ThreePhasePcKrSw = FluidMatrix::ParkerVanGenuchten3PDefault<Scalar>;
 
 public:
-    using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-    using MaterialLawParams = typename MaterialLaw::Params;
     using PermeabilityType = Scalar;
 
     ColumnSpatialParams(std::shared_ptr<const GridGeometry> gridGeometry)
     : ParentType(gridGeometry)
+    , pcKrSwCurveFine_("SpatialParams.Fine")
+    , pcKrSwCurveCoarse_("SpatialParams.Coarse")
     {
         // intrinsic permeabilities
         fineK_ = 1.4e-11;
@@ -70,29 +69,6 @@ public:
         // specific heat capacities
         fineHeatCap_ = getParam<Scalar>("Component.SolidHeatCapacityFine", 850.0);
         coarseHeatCap_ = getParam<Scalar>("Component.SolidHeatCapacityCoarse", 84000.0);
-
-        // residual saturations
-        fineMaterialParams_.setSwr(0.12);
-        fineMaterialParams_.setSnr(0.10);
-        fineMaterialParams_.setSgr(0.01);
-        coarseMaterialParams_.setSwr(0.12);
-        coarseMaterialParams_.setSnr(0.10);
-        coarseMaterialParams_.setSgr(0.01);
-
-        // parameters for the 3phase van Genuchten law
-        fineMaterialParams_.setVgAlpha(0.0005);
-        coarseMaterialParams_.setVgAlpha(0.5);
-        fineMaterialParams_.setVgn(4.0);
-        coarseMaterialParams_.setVgn(4.0);
-
-        coarseMaterialParams_.setKrRegardsSnr(true);
-        fineMaterialParams_.setKrRegardsSnr(true);
-
-        // parameters for adsorption
-        coarseMaterialParams_.setKdNAPL(0.);
-        coarseMaterialParams_.setRhoBulk(1500.);
-        fineMaterialParams_.setKdNAPL(0.);
-        fineMaterialParams_.setRhoBulk(1500.);
     }
 
     /*!
@@ -116,32 +92,31 @@ public:
     }
 
     /*! \brief Defines the porosity in [-].
-   *
-   * \param globalPos The global position where we evaluate
-   */
+    *
+    * \param globalPos The global position where we evaluate
+    */
     Scalar porosityAtPos(const GlobalPosition& globalPos) const
     {
         return porosity_;
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the fluid-matrix interaction law.
      *
      * \param element The current element
      * \param scv The sub-control volume inside the element.
      * \param elemSol The solution at the dofs connected to the element.
-     * \return The material parameters object
      */
-    template<class ElementSolution>
-    const MaterialLawParams& materialLawParams(const Element& element,
-                                               const SubControlVolume& scv,
-                                               const ElementSolution& elemSol) const
+     template<class ElementSolution>
+     auto fluidMatrixInteraction(const Element& element,
+                                 const SubControlVolume& scv,
+                                 const ElementSolution& elemSol) const
     {
         const auto& globalPos = scv.dofPosition();
         if (isFineMaterial_(globalPos))
-            return fineMaterialParams_;
+            return makeFluidMatrixInteraction(pcKrSwCurveFine_);
         else
-            return coarseMaterialParams_;
+            return makeFluidMatrixInteraction(pcKrSwCurveCoarse_);
     }
 
     /*!
@@ -180,8 +155,8 @@ private:
     Scalar fineHeatCap_;
     Scalar coarseHeatCap_;
 
-    MaterialLawParams fineMaterialParams_;
-    MaterialLawParams coarseMaterialParams_;
+    const ThreePhasePcKrSw pcKrSwCurveFine_;
+    const ThreePhasePcKrSw pcKrSwCurveCoarse_;
 
     static constexpr Scalar eps_ = 1e-6;
 };
diff --git a/test/porousmediumflow/3p3c/implicit/infiltration/params.input b/test/porousmediumflow/3p3c/implicit/infiltration/params.input
index d910fd4511..9a76b58299 100644
--- a/test/porousmediumflow/3p3c/implicit/infiltration/params.input
+++ b/test/porousmediumflow/3p3c/implicit/infiltration/params.input
@@ -14,3 +14,11 @@ NumericDifferenceMethod = 0 # -1 backward differences, 0: central differences, +
 
 [Newton]
 MaxRelativeShift = 1e-6
+
+[SpatialParams]
+Swr = 0.12
+Snr = 0.07
+Sgr = 0.03
+ParkerVanGenuchtenN = 4
+ParkerVanGenuchtenAlpha = 0.0005
+ParkerVanGenuchtenRegardSnrForKrn = false
diff --git a/test/porousmediumflow/3p3c/implicit/infiltration/problem.hh b/test/porousmediumflow/3p3c/implicit/infiltration/problem.hh
index c710cf2b9b..ba40089f11 100644
--- a/test/porousmediumflow/3p3c/implicit/infiltration/problem.hh
+++ b/test/porousmediumflow/3p3c/implicit/infiltration/problem.hh
@@ -274,7 +274,7 @@ private:
             Scalar pc = 9.81 * 1000.0 * (y - (-5E-4*x+5));
             if (pc < 0.0) pc = 0.0;
 
-            sw = invertPcgw_(pc, this->spatialParams().materialLawParamsAtPos(globalPos));
+            sw = invertPcgw_(pc, this->spatialParams().fluidMatrixInteractionAtPos(globalPos));
             if (sw < swr) sw = swr;
             if (sw > 1.-sgr) sw = 1.-sgr;
 
@@ -289,10 +289,9 @@ private:
         return values;
     }
 
-    template<class MaterialLawParams>
-    static Scalar invertPcgw_(Scalar pcIn, const MaterialLawParams &pcParams)
+    template<class FluidMatrixInteraction>
+    static Scalar invertPcgw_(Scalar pcIn, const FluidMatrixInteraction& fluidmatrixinteraction)
     {
-        using MaterialLaw = typename ParentType::SpatialParams::MaterialLaw;
         Scalar lower,upper;
         int k;
         int maxIt = 50;
@@ -302,7 +301,7 @@ private:
         for (k=1; k<=25; k++)
         {
             sw = 0.5*(upper+lower);
-            pcgw = MaterialLaw::pcgw(pcParams, sw);
+            pcgw = fluidmatrixinteraction.pcgw(sw, 0.0/*sn*/);
             Scalar delta = pcgw-pcIn;
             if (delta<0.) delta*=-1.;
             if (delta<bisLimit)
diff --git a/test/porousmediumflow/3p3c/implicit/infiltration/spatialparams.hh b/test/porousmediumflow/3p3c/implicit/infiltration/spatialparams.hh
index 2b7f665c55..f47a1c37f3 100644
--- a/test/porousmediumflow/3p3c/implicit/infiltration/spatialparams.hh
+++ b/test/porousmediumflow/3p3c/implicit/infiltration/spatialparams.hh
@@ -28,10 +28,7 @@
 
 #include <dumux/porousmediumflow/properties.hh>
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/3p/regularizedparkervangen3p.hh>
-#include <dumux/material/fluidmatrixinteractions/3p/regularizedparkervangen3pparams.hh>
-#include <dumux/material/fluidmatrixinteractions/3p/efftoabslaw.hh>
-#include <dumux/io/plotmateriallaw3p.hh>
+#include <dumux/material/fluidmatrixinteractions/3p/parkervangenuchten.hh>
 
 namespace Dumux {
 
@@ -51,17 +48,16 @@ class InfiltrationThreePThreeCSpatialParams
     using ParentType = FVSpatialParams<GridGeometry, Scalar,
                                        InfiltrationThreePThreeCSpatialParams<GridGeometry, Scalar>>;
 
-    using EffectiveLaw = RegularizedParkerVanGen3P<Scalar>;
-
     using GlobalPosition = typename SubControlVolume::GlobalPosition;
 
+    using ThreePhasePcKrSw = FluidMatrix::ParkerVanGenuchten3PDefault<Scalar>;
+
 public:
-   using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-   using MaterialLawParams = typename MaterialLaw::Params;
    using PermeabilityType = Scalar;
 
     InfiltrationThreePThreeCSpatialParams(std::shared_ptr<const GridGeometry> gridGeometry)
     : ParentType(gridGeometry)
+    , pcKrSwCurve_("SpatialParams")
     {
         // intrinsic permeabilities
         fineK_ = 1.e-11;
@@ -69,20 +65,6 @@ public:
 
         // porosities
         porosity_ = 0.40;
-
-        // residual saturations
-        materialParams_.setSwr(0.12);
-        materialParams_.setSnr(0.07);
-        materialParams_.setSgr(0.03);
-
-        // parameters for the 3phase van Genuchten law
-        materialParams_.setVgAlpha(0.0005);
-        materialParams_.setVgn(4.);
-        materialParams_.setKrRegardsSnr(false);
-
-        // parameters for adsorption
-        materialParams_.setKdNAPL(0.);
-        materialParams_.setRhoBulk(1500.);
     }
 
     /*!
@@ -115,15 +97,14 @@ public:
         return porosity_;
     }
 
-
     /*!
-     * \brief Returns the parameter object for the material law which depends on the position
+     * \brief Returns the fluid-matrix interaction law at a given location
      *
-     * \param globalPos The position for which the material law should be evaluated
+     * \param globalPos The global coordinates for the given location
      */
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
     {
-        return materialParams_;
+        return makeFluidMatrixInteraction(pcKrSwCurve_);
     }
 
 private:
@@ -138,7 +119,7 @@ private:
 
     Scalar porosity_;
 
-    MaterialLawParams materialParams_;
+    const ThreePhasePcKrSw pcKrSwCurve_;
 };
 
 } // end namespace Dumux
-- 
GitLab