From dfd7fb1c5151701627803ac64271c479be8ab5ed Mon Sep 17 00:00:00 2001
From: Roman Winter <rwinter@shrek.iws.uni-stuttgart.de>
Date: Wed, 28 Oct 2020 15:47:42 +0100
Subject: [PATCH] [feature] materiallaw-2p for course

---
 exercises/exercise-basic/params.input         | 10 +++-
 exercises/exercise-basic/spatialparams.hh     | 44 ++++++----------
 .../biominspatialparams.hh                    | 50 +++++++-----------
 .../exercise-biomineralization/params.input   |  9 ++++
 .../2pspatialparams.hh                        | 40 +++++---------
 .../models/params.input                       |  6 ++-
 .../turbulence/params.input                   |  6 ++-
 exercises/exercise-fluidsystem/aparams.input  | 10 ++++
 exercises/exercise-fluidsystem/bparams.input  | 10 ++++
 .../exercise-fluidsystem/spatialparams.hh     | 45 ++++++----------
 .../fracturespatialparams.hh                  | 49 ++++++++---------
 .../exercise-fractures/matrixspatialparams.hh | 34 ++++++------
 exercises/exercise-fractures/params.input     | 18 +++----
 exercises/exercise-grids/params.input         | 10 +++-
 exercises/exercise-grids/spatialparams.hh     | 43 ++++++---------
 exercises/exercise-properties/params.input    |  8 +++
 .../exercise-properties/spatialparams.hh      | 50 ++++++------------
 exercises/exercise-runtimeparams/params.input | 10 +++-
 .../exercise-runtimeparams/spatialparams.hh   | 44 ++++++----------
 .../solution/exercise-basic/params.input      | 10 +++-
 .../solution/exercise-basic/spatialparams.hh  | 43 ++++++---------
 .../biominspatialparams.hh                    | 52 +++++++------------
 .../exercise-biomineralization/params.input   |  8 +++
 .../2pspatialparams.hh                        | 40 +++++---------
 .../models/params.input                       |  6 ++-
 .../turbulence/params.input                   |  7 +--
 .../exercise-fluidsystem/aparams.input        | 10 ++++
 .../exercise-fluidsystem/bparams.input        | 10 ++++
 .../exercise-fluidsystem/spatialparams.hh     | 42 ++++++---------
 .../exercise_fractures_solution_a.input       | 18 +++----
 .../exercise_fractures_solution_b.input       | 18 +++----
 .../exercise_fractures_solution_c.input       | 18 +++----
 .../fracturespatialparams.hh                  | 52 +++++++++----------
 .../exercise-fractures/matrixspatialparams.hh | 33 ++++++------
 .../solution/exercise-grids/params.input      | 10 +++-
 .../solution/exercise-grids/spatialparams.hh  | 43 ++++++---------
 .../solution/exercise-properties/params.input |  8 +++
 .../exercise-properties/spatialparams.hh      | 50 ++++++------------
 .../exercise-runtimeparams/params.input       | 10 +++-
 .../exercise-runtimeparams/spatialparams.hh   | 44 ++++++----------
 40 files changed, 478 insertions(+), 550 deletions(-)

diff --git a/exercises/exercise-basic/params.input b/exercises/exercise-basic/params.input
index 61a16b57..b2e574e8 100644
--- a/exercises/exercise-basic/params.input
+++ b/exercises/exercise-basic/params.input
@@ -15,9 +15,15 @@ InjectionDuration = 2.628e6 # in seconds, i.e. one month
 
 [SpatialParams]
 PermeabilityAquitard = 1e-15 # m^2
-EntryPressureAquitard = 4.5e4 # Pa
+Aquitard.BrooksCoreyPcEntry = 4.5e4 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.2
+Aquitard.Snr = 0.0
 PermeabilityAquifer = 1e-12 # m^2
-EntryPressureAquifer = 1e4 # Pa
+Aquifer.BrooksCoreyPcEntry = 1e4 # Pa
+Aquifer.BrooksCoreyLambda = 2.0
+Aquifer.Swr = 0.2
+Aquifer.Snr = 0.0
 
 # these parameters are only used in the nonisothermal model. Uncomment them for that
 #[Component]
diff --git a/exercises/exercise-basic/spatialparams.hh b/exercises/exercise-basic/spatialparams.hh
index 715d02f2..354caeda 100644
--- a/exercises/exercise-basic/spatialparams.hh
+++ b/exercises/exercise-basic/spatialparams.hh
@@ -28,11 +28,10 @@
 #define DUMUX_EX_BASIC_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 
 #include <dumux/io/gnuplotinterface.hh>
-#include <dumux/io/plotmateriallaw.hh>
+#include <dumux/io/plotpckrsw.hh>
 
 namespace Dumux {
 
@@ -55,13 +54,12 @@ class InjectionSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     // export permeability type
     using PermeabilityType = Scalar;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -69,6 +67,8 @@ public:
      */
     InjectionSpatialParams(std::shared_ptr<const FVGridGeometry>& fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
+    , aquiferPcKrSwCurve_("SpatialParams.Aquifer")
     {
         aquiferHeightFromBottom_ = 30.0;
 
@@ -79,18 +79,6 @@ public:
         // porosities
         aquitardPorosity_ = 0.2;
         aquiferPorosity_ = 0.4;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.2);
-        aquitardMaterialParams_.setSnr(0.0);
-        aquiferMaterialParams_.setSwr(0.2);
-        aquiferMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquitard"));
-        aquiferMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquifer"));
-        aquitardMaterialParams_.setLambda(2.0);
-        aquiferMaterialParams_.setLambda(2.0);
     }
 
     /*!
@@ -120,17 +108,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInAquitard_(globalPos))
-            return aquitardMaterialParams_;
-        return aquiferMaterialParams_;
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(aquiferPcKrSwCurve_); 
     }
 
     /*!
@@ -158,12 +147,11 @@ private:
     Scalar aquiferK_;
     Scalar aquiferHeightFromBottom_;
 
-
     Scalar aquitardPorosity_;
     Scalar aquiferPorosity_;
 
-    MaterialLawParams aquitardMaterialParams_;
-    MaterialLawParams aquiferMaterialParams_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
+    const PcKrSwCurve aquiferPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/exercise-biomineralization/biominspatialparams.hh b/exercises/exercise-biomineralization/biominspatialparams.hh
index 2c901dc9..e147f8fe 100644
--- a/exercises/exercise-biomineralization/biominspatialparams.hh
+++ b/exercises/exercise-biomineralization/biominspatialparams.hh
@@ -26,8 +26,7 @@
 
 #include <dumux/material/spatialparams/fv.hh>
 #include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 #include <dumux/material/fluidmatrixinteractions/porosityprecipitation.hh>
 #include <dumux/material/fluidmatrixinteractions/permeabilitykozenycarman.hh>
 
@@ -46,7 +45,6 @@ class BioMinSpatialparams
     using FVElementGeometry = typename FVGridGeometry::LocalView;
     using SubControlVolume = typename FVElementGeometry::SubControlVolume;
     using ParentType = FVSpatialParams<FVGridGeometry, Scalar, BioMinSpatialparams<FVGridGeometry, Scalar, numFluidComps, numActiveSolidComps>>;
-    using EffectiveLaw = RegularizedBrooksCorey<Scalar>;
 
     using GridView = typename FVGridGeometry::GridView;
     using CoordScalar = typename GridView::ctype;
@@ -58,10 +56,10 @@ class BioMinSpatialparams
 
     using PoroLaw = PorosityPrecipitation<Scalar, numFluidComps, numActiveSolidComps>;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     using PermeabilityType = Tensor;
-    using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-    using MaterialLawParams = typename MaterialLaw::Params;
 
     /*!
     * \brief The constructor
@@ -70,6 +68,8 @@ public:
     */
     BioMinSpatialparams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , pcKrSwCurve_("SpatialParams")
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
     {
         //! set initial aquifer params
         initialPorosity_ = getParam<Scalar>("SpatialParams.InitialPorosity");
@@ -78,25 +78,9 @@ public:
         // set main diagonal entries of the permeability tensor to a value
         // setting to one value means: isotropic, homogeneous
 
-        // residual saturations
-        materialParams_.setSwr(0.2);
-        materialParams_.setSnr(0.05);
-
-        // parameters for the Brooks-Corey law
-        materialParams_.setPe(1e4);
-        materialParams_.setLambda(2.0);
-
         //! hard code specific params for aquitard layer
         aquitardPorosity_ = 0.1;
         aquitardPermeability_ = 1e-15;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.0);
-        aquitardMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(1e7);
-        aquitardMaterialParams_.setLambda(2.0);
     }
 
     template<class SolidSystem, class ElementSolution>
@@ -261,16 +245,18 @@ public:
     }
 
     /*!
-     * \brief return the parameter object for the Brooks-Corey material law which depends on the position
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
+     *
+     * \param globalPos The global coordinates for the given location
      */
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
-        if (isInAquitard_(globalPos) && !isFaultZone_(globalPos))
-            return aquitardMaterialParams_;
-        else
-            return materialParams_;
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
+        if (isInAquitard_(globalPos))
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(pcKrSwCurve_); 
     }
 
     // define which phase is to be considered as the wetting phase
@@ -298,11 +284,11 @@ private:
     Scalar initialPermeability_;
     std::vector< std::vector<PermeabilityType> > referencePermeability_;
 
-    MaterialLawParams materialParams_;
-
     Scalar aquitardPorosity_;
     Scalar aquitardPermeability_;
-    MaterialLawParams aquitardMaterialParams_;
+
+    const PcKrSwCurve pcKrSwCurve_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/exercise-biomineralization/params.input b/exercises/exercise-biomineralization/params.input
index fb2df919..5b3cd71f 100644
--- a/exercises/exercise-biomineralization/params.input
+++ b/exercises/exercise-biomineralization/params.input
@@ -26,6 +26,15 @@ ConcUrea = 60 # [kg/m³] injected urea concentration (max: 80)
 [SpatialParams]
 InitialPorosity = 0.35 # [-]
 InitialPermeability = 1e-13 # [m^2]
+BrooksCoreyPcEntry = 1e4 # Pa
+BrooksCoreyLambda = 2.0
+Swr = 0.2
+Snr = 0.05
+Aquitard.BrooksCoreyPcEntry = 1e7 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.0
+Aquitard.Snr = 0.0
+
 
 [BioCoefficients]
 RhoBiofilm = 6.9 # [kg/m³] density of biofilm
diff --git a/exercises/exercise-coupling-ff-pm/2pspatialparams.hh b/exercises/exercise-coupling-ff-pm/2pspatialparams.hh
index 23930847..7a6f7360 100644
--- a/exercises/exercise-coupling-ff-pm/2pspatialparams.hh
+++ b/exercises/exercise-coupling-ff-pm/2pspatialparams.hh
@@ -25,8 +25,7 @@
 #define DUMUX_TWOPHASE_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 #include <dumux/material/fluidmatrixinteractions/2p/thermalconductivity/somerton.hh>
 
 namespace Dumux {
@@ -47,28 +46,19 @@ class TwoPSpatialParams
     using ParentType = FVSpatialParams<FVGridGeometry, Scalar, TwoPSpatialParams<FVGridGeometry, Scalar>>;
 
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
 public:
-    using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-    using MaterialLawParams = typename MaterialLaw::Params;
     using PermeabilityType = Scalar;
 
     TwoPSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , pcKrSwCurve_("Darcy.SpatialParams")
     {
         permeability_ = getParam<Scalar>("Darcy.SpatialParams.Permeability");
         porosity_ = getParam<Scalar>("Darcy.SpatialParams.Porosity");
         alphaBJ_ = getParam<Scalar>("Darcy.SpatialParams.AlphaBeaversJoseph");
-
-        // residual saturations
-        params_.setSwr(getParam<Scalar>("Darcy.SpatialParams.Swr"));
-        params_.setSnr(getParam<Scalar>("Darcy.SpatialParams.Snr"));
-        // parameters for the vanGenuchten law
-        params_.setVgAlpha(getParam<Scalar>("Darcy.SpatialParams.VgAlpha"));
-        params_.setVgn(getParam<Scalar>("Darcy.SpatialParams.VgN"));
-        params_.setPcLowSw(params_.swr()*5.0);
-        params_.setPcHighSw(1.0-params_.snr()*5.0);
     }
 
     /*!
@@ -95,19 +85,17 @@ public:
     { return alphaBJ_; }
 
     /*!
-     * \brief Returns the parameter object for the Brooks-Corey material law.
-     *        In this test, we use element-wise distributed material parameters.
+     * \brief Returns the parameters for the material law at a given location
+     *
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \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
+     * \param globalPos The global coordinates for the given location
      */
-    template<class ElementSolutionVector>
-    const MaterialLawParams& materialLawParams(const Element& element,
-                                               const SubControlVolume& scv,
-                                               const ElementSolutionVector& elemSol) const
-    { return params_; }
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {
+        return makeFluidMatrixInteraction(pcKrSwCurve_);
+    }
 
     /*!
      * \brief Function for defining which phase is to be considered as the wetting phase.
@@ -123,7 +111,7 @@ private:
     Scalar permeability_;
     Scalar porosity_;
     Scalar alphaBJ_;
-    MaterialLawParams params_;
+    const PcKrSwCurve pcKrSwCurve_;
     static constexpr Scalar eps_ = 1.0e-7;
 };
 
diff --git a/exercises/exercise-coupling-ff-pm/models/params.input b/exercises/exercise-coupling-ff-pm/models/params.input
index 07bb4dd0..60720f80 100644
--- a/exercises/exercise-coupling-ff-pm/models/params.input
+++ b/exercises/exercise-coupling-ff-pm/models/params.input
@@ -32,10 +32,12 @@ Permeability = 2.65e-10 # m^2
 Porosity = 0.4 # -
 AlphaBeaversJoseph = 1.0 # -
 # EXNUMBER >= 1
+VanGenuchtenN = 8.0
+VanGenuchtenAlpha = 6.5e-4
 Swr = 0.005
 Snr = 0.01
-VgAlpha = 6.5e-4
-VgN = 8.0
+PcLowSw = Swr * 5.0
+PcHighSw = 1.0 - Snr * 5.0
 
 [Problem]
 Name = models_coupling
diff --git a/exercises/exercise-coupling-ff-pm/turbulence/params.input b/exercises/exercise-coupling-ff-pm/turbulence/params.input
index 00e02465..5be38764 100644
--- a/exercises/exercise-coupling-ff-pm/turbulence/params.input
+++ b/exercises/exercise-coupling-ff-pm/turbulence/params.input
@@ -40,10 +40,12 @@ InitPhasePresence = 3 # bothPhases
 Porosity = 0.41
 Permeability = 2.65e-10
 AlphaBeaversJoseph = 1.0
+VanGenuchtenN = 6.9
+VanGenuchtenAlpha = 6.371e-4
 Swr = 0.005
 Snr = 0.01
-VgAlpha = 6.371e-4
-VgN = 6.9
+PcLowSw = Swr * 5.0
+PcHighSw = 1.0 - Snr * 5.0
 
 [Problem]
 Name = ex_coupling_turbulence_ff-pm
diff --git a/exercises/exercise-fluidsystem/aparams.input b/exercises/exercise-fluidsystem/aparams.input
index 8b18ca89..60897262 100644
--- a/exercises/exercise-fluidsystem/aparams.input
+++ b/exercises/exercise-fluidsystem/aparams.input
@@ -5,6 +5,16 @@ DtInitial = 10 # initial time step size [s]
 [Problem]
 Name = exercise-fluidsystem_a # name will be given to e.g. to the vtk result files
 
+[SpatialParams]
+BrooksCoreyPcEntry = 5.0e2 # Pa
+BrooksCoreyLambda = 2.0
+Swr = 0.1
+Snr = 0.0
+Lens.BrooksCoreyPcEntry = 1e3 # Pa
+Lens.BrooksCoreyLambda = 2.0
+Lens.Swr = 0.1
+Lens.Snr = 0.0
+
 [Grid]
 UpperRight = 60 60 # x-/y-coordinates of the upper-right corner of the grid [m]
 Cells = 60 60 # x-/y-resolution of the grid
diff --git a/exercises/exercise-fluidsystem/bparams.input b/exercises/exercise-fluidsystem/bparams.input
index 3f2eea4c..2b6595f3 100644
--- a/exercises/exercise-fluidsystem/bparams.input
+++ b/exercises/exercise-fluidsystem/bparams.input
@@ -5,6 +5,16 @@ DtInitial = 10 # initial time step size [s]
 [Problem]
 Name = exercise-fluidsystem_b # name will be given to e.g. to the vtk result files
 
+[SpatialParams]
+BrooksCoreyPcEntry = 5.0e2 # Pa
+BrooksCoreyLambda = 2.0
+Swr = 0.1
+Snr = 0.0
+Lens.BrooksCoreyPcEntry = 1e3 # Pa
+Lens.BrooksCoreyLambda = 2.0
+Lens.Swr = 0.1
+Lens.Snr = 0.0
+
 [Grid]
 UpperRight = 60 60 # x-/y-coordinates of the upper-right corner of the grid [m]
 Cells = 60 60 # x-/y-resolution of the grid
diff --git a/exercises/exercise-fluidsystem/spatialparams.hh b/exercises/exercise-fluidsystem/spatialparams.hh
index b3caeaec..a9390531 100644
--- a/exercises/exercise-fluidsystem/spatialparams.hh
+++ b/exercises/exercise-fluidsystem/spatialparams.hh
@@ -29,8 +29,7 @@
 #include <dumux/material/spatialparams/fv.hh>
 
 // include material laws
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 #include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
 
 namespace Dumux {
@@ -52,15 +51,14 @@ class ExerciseFluidsystemSpatialParams
     static constexpr int dim = GridView::dimension;
     static constexpr int dimWorld = GridView::dimensionworld;
     using Element = typename GridView::template Codim<0>::Entity;
-    using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
+    using GlobalPosition = typename Element::Geometry::GlobalCoordinate;    
+
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
 
 public:
     // export permeability type
     using PermeabilityType = Dune::FieldMatrix<Scalar, dim, dim>;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -70,6 +68,8 @@ public:
     : ParentType(fvGridGeometry)
     , K_(0)
     , KLens_(0)
+    , pcKrSwCurve_("SpatialParams")
+    , lensPcKrSwCurve_("Lens.SpatialParams")
     {
         //set main diagonal entries of the permeability tensor to a value
         //setting to one value means: isotropic, homogeneous
@@ -78,18 +78,6 @@ public:
             K_[i][i] = 1e-7;
             KLens_[i][i] = 1e-10;
         }
-
-        //set residual saturations
-        materialParams_.setSwr(0.0);
-        materialParamsLens_.setSwr(0.1);
-        materialParams_.setSnr(0.0);
-        materialParamsLens_.setSnr(0.1);
-
-        //parameters of Brooks & Corey Law
-        materialParams_.setPe(500.0);
-        materialParamsLens_.setPe(1000.0);
-        materialParams_.setLambda(2);
-        materialParamsLens_.setLambda(2);
     }
 
 
@@ -119,17 +107,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInLens(globalPos))
-            return materialParamsLens_;
-        return materialParams_;
+            return makeFluidMatrixInteraction(pcKrSwCurve_); 
+        return makeFluidMatrixInteraction(lensPcKrSwCurve_); 
     }
 
     /*!
@@ -165,9 +154,9 @@ private:
 
     Dune::FieldMatrix<Scalar, dim, dim> K_;
     Dune::FieldMatrix<Scalar, dim, dim> KLens_;
-    // Object that holds the values/parameters of the selected material law.
-    MaterialLawParams materialParams_;
-    MaterialLawParams materialParamsLens_;
+    
+    const PcKrSwCurve pcKrSwCurve_;
+    const PcKrSwCurve lensPcKrSwCurve_;
 };
 } // end namespace Dumux
 #endif
diff --git a/exercises/exercise-fractures/fracturespatialparams.hh b/exercises/exercise-fractures/fracturespatialparams.hh
index 3431f351..85113fb0 100644
--- a/exercises/exercise-fractures/fracturespatialparams.hh
+++ b/exercises/exercise-fractures/fracturespatialparams.hh
@@ -30,8 +30,7 @@
 #include <dumux/io/grid/griddata.hh>
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 
 namespace Dumux {
 
@@ -54,8 +53,7 @@ class FractureSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
-    // use a regularized van-genuchten material law
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
     // we identify those fractures as barriers, that have a domain marker
     // of 2. This is what is set in the grid file (see grids/complex.geo)
@@ -65,31 +63,19 @@ public:
     //! export the type used for permeabilities
     using PermeabilityType = Scalar;
 
-    //! export the material law and parameters used
-    using MaterialLaw = EffToAbsLaw< EffectiveLaw >;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     //! the constructor
     FractureSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
                           std::shared_ptr<const Dumux::GridData<Grid>> gridData,
                           const std::string& paramGroup)
     : ParentType(fvGridGeometry)
     , gridDataPtr_(gridData)
+    , pcKrSwCurve_("Fracture.SpatialParams")
+    , barrierPcKrSwCurve_("Fracture.SpatialParams.Barrier")
     {
         porosity_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Porosity");
         porosityBarrier_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.PorosityBarrier");
         permeability_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Permeability");
         permeabilityBarrier_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.PermeabilityBarrier");
-
-        // set the material law parameters
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Snr"));
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Swr"));
-        materialLawParams_.setVgAlpha(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGAlpha"));
-        materialLawParams_.setVgn(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGN"));
-        materialLawParamsBarrier_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.SnrBarrier"));
-        materialLawParamsBarrier_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.SwrBarrier"));
-        materialLawParamsBarrier_.setVgAlpha(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGAlphaBarrier"));
-        materialLawParamsBarrier_.setVgn(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGNBarrier"));
     }
 
     //! Function for defining the (intrinsic) permeability \f$[m^2]\f$.
@@ -114,15 +100,24 @@ public:
         return porosity_;
     }
 
-    //! Return the material law parameters
-    template< class ElementSolution >
-    const MaterialLawParams& materialLawParams(const Element& element,
-                                               const SubControlVolume& scv,
-                                               const ElementSolution& elemSol) const
-    {
+    /*!
+     * \brief Returns the parameters for the material law for the sub-control volume
+     *
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
+     *
+     * \param element The current finite element
+     * \param scv The sub-control volume
+     * \param elemSol The current element solution
+     */
+    template<class ElementSolution>
+    auto fluidMatrixInteraction(const Element& element,
+                                const SubControlVolume& scv,
+                                const ElementSolution& elemSol) const
+    {  
         // TODO dumux-course-task B
         // Change fracture properties
-        return materialLawParams_;
+        return makeFluidMatrixInteraction(pcKrSwCurve_); 
     }
 
     //! Water is the wetting phase
@@ -146,8 +141,8 @@ private:
     Scalar porosityBarrier_;
     PermeabilityType permeability_;
     PermeabilityType permeabilityBarrier_;
-    MaterialLawParams materialLawParams_;
-    MaterialLawParams materialLawParamsBarrier_;
+    const PcKrSwCurve pcKrSwCurve_;
+    const PcKrSwCurve barrierPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/exercise-fractures/matrixspatialparams.hh b/exercises/exercise-fractures/matrixspatialparams.hh
index 4ca8c72e..b6c77c8e 100644
--- a/exercises/exercise-fractures/matrixspatialparams.hh
+++ b/exercises/exercise-fractures/matrixspatialparams.hh
@@ -30,8 +30,7 @@
 #include <dumux/io/grid/griddata.hh>
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 
 namespace Dumux {
 
@@ -53,32 +52,22 @@ class MatrixSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
-    // use a regularized van-genuchten material law
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
 public:
     //! export the type used for permeabilities
     using PermeabilityType = Scalar;
 
-    //! export the material law and parameters used
-    using MaterialLaw = EffToAbsLaw< EffectiveLaw >;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     //! the constructor
     MatrixSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
                         std::shared_ptr<const Dumux::GridData<Grid>> gridData,
                         const std::string& paramGroup)
     : ParentType(fvGridGeometry)
     , gridDataPtr_(gridData)
+    , pcKrSwCurve_("Matrix.SpatialParams")
     {
         porosity_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Porosity");
         permeability_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Permeability");
-
-        // set the material law parameters
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Snr"));
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Swr"));
-        materialLawParams_.setVgAlpha(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGAlpha"));
-        materialLawParams_.setVgn(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGN"));
     }
 
     //! Function for defining the (intrinsic) permeability \f$[m^2]\f$.
@@ -89,9 +78,18 @@ public:
     Scalar porosityAtPos(const GlobalPosition& globalPos) const
     { return porosity_; }
 
-    //! Return the material law parameters
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    { return materialLawParams_; }
+    /*!
+     * \brief Returns the parameters for the material law at a given location
+     *
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
+     *
+     * \param globalPos The global coordinates for the given location
+     */
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
+        return makeFluidMatrixInteraction(pcKrSwCurve_); 
+    }
 
     //! Water is the wetting phase
     template< class FluidSystem >
@@ -112,7 +110,7 @@ private:
 
     Scalar porosity_;
     PermeabilityType permeability_;
-    MaterialLawParams materialLawParams_;
+    const PcKrSwCurve pcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/exercise-fractures/params.input b/exercises/exercise-fractures/params.input
index 75a68cc2..2c3bd5db 100644
--- a/exercises/exercise-fractures/params.input
+++ b/exercises/exercise-fractures/params.input
@@ -16,8 +16,8 @@ Problem.BoundaryOverPressure = 5e5
 Problem.BoundarySaturation = 0.5
 SpatialParams.Permeability = 1e-12
 SpatialParams.Porosity = 0.1
-SpatialParams.VGAlpha = 1e-3
-SpatialParams.VGN = 3
+SpatialParams.VanGenuchtenAlpha = 1e-3
+SpatialParams.VanGenuchtenN = 3
 SpatialParams.Snr = 0.0
 SpatialParams.Swr = 0.0
 
@@ -28,11 +28,11 @@ SpatialParams.Permeability = 1e-7
 SpatialParams.PermeabilityBarrier = 1e-16
 SpatialParams.Porosity = 0.85
 SpatialParams.PorosityBarrier = 0.05
-SpatialParams.VGAlpha = 1e-1
-SpatialParams.VGAlphaBarrier = 1e-4
-SpatialParams.VGN = 2
-SpatialParams.VGNBarrier = 2.5
-SpatialParams.Snr = 0.0
-SpatialParams.SnrBarrier = 0.0
+SpatialParams.VanGenuchtenAlpha = 1e-1
+SpatialParams.VanGenuchtenN = 2
 SpatialParams.Swr = 0.0
-SpatialParams.SwrBarrier = 0.0
+SpatialParams.Snr = 0.0
+SpatialParams.Barrier.VanGenuchtenAlpha = 1e-4
+SpatialParams.Barrier.VanGenuchtenN = 2.5
+SpatialParams.Barrier.Snr = 0.0
+SpatialParams.Barrier.Swr = 0.0
diff --git a/exercises/exercise-grids/params.input b/exercises/exercise-grids/params.input
index 78155038..21d49997 100644
--- a/exercises/exercise-grids/params.input
+++ b/exercises/exercise-grids/params.input
@@ -19,6 +19,12 @@ InjectionDuration = 2.628e6 # in seconds, i.e. one month
 
 [SpatialParams]
 PermeabilityAquitard = 1e-15 # m^2
-EntryPressureAquitard = 4.5e4 # Pa
+Aquitard.BrooksCoreyPcEntry = 4.5e4 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.2
+Aquitard.Snr = 0.0
 PermeabilityAquifer = 1e-12 # m^2
-EntryPressureAquifer = 1e4 # Pa
+Aquifer.BrooksCoreyPcEntry = 1e4 # Pa
+Aquifer.BrooksCoreyLambda = 2.0
+Aquifer.Swr = 0.2
+Aquifer.Snr = 0.0
\ No newline at end of file
diff --git a/exercises/exercise-grids/spatialparams.hh b/exercises/exercise-grids/spatialparams.hh
index 0e72f8d5..ea7e3084 100644
--- a/exercises/exercise-grids/spatialparams.hh
+++ b/exercises/exercise-grids/spatialparams.hh
@@ -28,11 +28,10 @@
 #define DUMUX_EXGRIDS_INJECTION_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 
 #include <dumux/io/gnuplotinterface.hh>
-#include <dumux/io/plotmateriallaw.hh>
+#include <dumux/io/plotpckrsw.hh>
 
 namespace Dumux {
 
@@ -55,13 +54,12 @@ class InjectionSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     // export permeability type
     using PermeabilityType = Scalar;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -69,6 +67,8 @@ public:
      */
     InjectionSpatialParams(std::shared_ptr<const FVGridGeometry>& fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
+    , aquiferPcKrSwCurve_("SpatialParams.Aquifer")
     {
         // Aquifer Height, measured from the bottom
         aquiferHeightFromBottom_ = 30.0;
@@ -80,18 +80,6 @@ public:
         // porosities
         aquitardPorosity_ = 0.2;
         aquiferPorosity_ = 0.4;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.2);
-        aquitardMaterialParams_.setSnr(0.0);
-        aquiferMaterialParams_.setSwr(0.2);
-        aquiferMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquitard"));
-        aquiferMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquifer"));
-        aquitardMaterialParams_.setLambda(2.0);
-        aquiferMaterialParams_.setLambda(2.0);
     }
 
     /*!
@@ -121,17 +109,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-     const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInAquitard_(globalPos))
-            return aquitardMaterialParams_;
-        return aquiferMaterialParams_;
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(aquiferPcKrSwCurve_); 
     }
 
     /*!
@@ -163,8 +152,8 @@ private:
     Scalar aquitardPorosity_;
     Scalar aquiferPorosity_;
 
-    MaterialLawParams aquitardMaterialParams_;
-    MaterialLawParams aquiferMaterialParams_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
+    const PcKrSwCurve aquiferPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/exercise-properties/params.input b/exercises/exercise-properties/params.input
index 7ae1761e..8e25e5f8 100644
--- a/exercises/exercise-properties/params.input
+++ b/exercises/exercise-properties/params.input
@@ -10,6 +10,14 @@ Cells = 48 32
 [SpatialParams]
 LensLowerLeft = 1.0 2.0 # [m] coordinates of the lower left lens corner
 LensUpperRight = 4.0 3.0 # [m] coordinates of the upper right lens corner
+Lens.VanGenuchtenN = 7.3
+Lens.VanGenuchtenAlpha = 0.00045
+Lens.Swr = 0.18
+Lens.Snr = 0.0
+OuterDomain.VanGenuchtenN = 4.7
+OuterDomain.VanGenuchtenAlpha = 0.0037
+OuterDomain.Swr = 0.05
+OuterDomain.Snr = 0.0
 
 [Problem]
 Name = exercise_properties
diff --git a/exercises/exercise-properties/spatialparams.hh b/exercises/exercise-properties/spatialparams.hh
index 275984d8..fa106e87 100644
--- a/exercises/exercise-properties/spatialparams.hh
+++ b/exercises/exercise-properties/spatialparams.hh
@@ -24,8 +24,7 @@
 #define DUMUX_INCOMPRESSIBLE_TWOP_TEST_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 
 namespace Dumux {
 
@@ -47,32 +46,19 @@ class TwoPTestSpatialParams
     static constexpr int dimWorld = GridView::dimensionworld;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
 public:
-    using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-    using MaterialLawParams = typename MaterialLaw::Params;
     using PermeabilityType = Scalar;
 
     TwoPTestSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , lensPcKrSwCurve_("SpatialParams.Lens")
+    , outerPcKrSwCurve_("SpatialParams.OuterDomain")
     {
         lensLowerLeft_ = getParam<GlobalPosition>("SpatialParams.LensLowerLeft");
         lensUpperRight_ = getParam<GlobalPosition>("SpatialParams.LensUpperRight");
 
-        // residual saturations
-        lensMaterialParams_.setSwr(0.18);
-        lensMaterialParams_.setSnr(0.0);
-        outerMaterialParams_.setSwr(0.05);
-        outerMaterialParams_.setSnr(0.0);
-
-        // parameters for the Van Genuchten law
-        // alpha and n
-        lensMaterialParams_.setVgAlpha(0.00045);
-        lensMaterialParams_.setVgn(7.3);
-        outerMaterialParams_.setVgAlpha(0.0037);
-        outerMaterialParams_.setVgn(4.7);
-
         lensK_ = 9.05e-12;
         outerK_ = 4.6e-10;
     }
@@ -105,22 +91,19 @@ public:
     { return 0.4; }
 
     /*!
-     * \brief Returns the parameter object for the Brooks-Corey material law.
-     *        In this test, we use element-wise distributed material parameters.
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \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
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
+     *
+     * \param globalPos The global coordinates for the given location
      */
-    template<class ElementSolution>
-    const MaterialLawParams& materialLawParams(const Element& element,
-                                               const SubControlVolume& scv,
-                                               const ElementSolution& elemSol) const
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
     {
-        if (isInLens_(element.geometry().center()))
-            return lensMaterialParams_;
-        return outerMaterialParams_;
+        if (isInLens_(globalPos))
+            return makeFluidMatrixInteraction(lensPcKrSwCurve_);
+
+        return makeFluidMatrixInteraction(outerPcKrSwCurve_);
     }
 
     /*!
@@ -150,8 +133,9 @@ private:
 
     Scalar lensK_;
     Scalar outerK_;
-    MaterialLawParams lensMaterialParams_;
-    MaterialLawParams outerMaterialParams_;
+
+    const PcKrSwCurve lensPcKrSwCurve_;
+    const PcKrSwCurve outerPcKrSwCurve_;
 
     static constexpr Scalar eps_ = 1.5e-7;
 };
diff --git a/exercises/exercise-runtimeparams/params.input b/exercises/exercise-runtimeparams/params.input
index 33f3940b..7ec66fe2 100644
--- a/exercises/exercise-runtimeparams/params.input
+++ b/exercises/exercise-runtimeparams/params.input
@@ -16,6 +16,12 @@ InjectionDuration = 2.628e6 # in seconds, i.e. one month
 [SpatialParams]
 PermeabilityAquitard = 1e-15 # m^2
 # TODO: Task 1: Change the Aquitard's Entry Pressure
-EntryPressureAquitard = 4.5e4 # Pa
+Aquitard.BrooksCoreyPcEntry = 4.5e4 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.2
+Aquitard.Snr = 0.0
 PermeabilityAquifer = 1e-12 # m^2
-EntryPressureAquifer = 1e4 # Pa
+Aquifer.BrooksCoreyPcEntry = 1e4 # Pa
+Aquifer.BrooksCoreyLambda = 2.0
+Aquifer.Swr = 0.2
+Aquifer.Snr = 0.0
diff --git a/exercises/exercise-runtimeparams/spatialparams.hh b/exercises/exercise-runtimeparams/spatialparams.hh
index 634a1e9e..4a0327a6 100644
--- a/exercises/exercise-runtimeparams/spatialparams.hh
+++ b/exercises/exercise-runtimeparams/spatialparams.hh
@@ -28,11 +28,10 @@
 #define DUMUX_EXRUNTIMEPARAMS_INJECTION_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 
 #include <dumux/io/gnuplotinterface.hh>
-#include <dumux/io/plotmateriallaw.hh>
+#include <dumux/io/plotpckrsw.hh>
 
 namespace Dumux {
 
@@ -55,13 +54,12 @@ class InjectionSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     // export permeability type
     using PermeabilityType = Scalar;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -69,6 +67,8 @@ public:
      */
     InjectionSpatialParams(std::shared_ptr<const FVGridGeometry>& fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
+    , aquiferPcKrSwCurve_("SpatialParams.Aquifer")
     {
         // Aquifer Height, measured from the bottom
         aquiferHeightFromBottom_ = 30.0;
@@ -80,18 +80,6 @@ public:
         // porosities
         aquitardPorosity_ = 0.2;
         aquiferPorosity_ = 0.4;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.2);
-        aquitardMaterialParams_.setSnr(0.0);
-        aquiferMaterialParams_.setSwr(0.2);
-        aquiferMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquitard"));
-        aquiferMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquifer"));
-        aquitardMaterialParams_.setLambda(2.0);
-        aquiferMaterialParams_.setLambda(2.0);
     }
 
     /*!
@@ -121,17 +109,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-     const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInAquitard_(globalPos))
-            return aquitardMaterialParams_;
-        return aquiferMaterialParams_;
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(aquiferPcKrSwCurve_); 
     }
 
     /*!
@@ -159,12 +148,11 @@ private:
     Scalar aquiferK_;
     Scalar aquiferHeightFromBottom_;
 
-
     Scalar aquitardPorosity_;
     Scalar aquiferPorosity_;
 
-    MaterialLawParams aquitardMaterialParams_;
-    MaterialLawParams aquiferMaterialParams_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
+    const PcKrSwCurve aquiferPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/solution/exercise-basic/params.input b/exercises/solution/exercise-basic/params.input
index ada32b39..06e9e037 100644
--- a/exercises/solution/exercise-basic/params.input
+++ b/exercises/solution/exercise-basic/params.input
@@ -15,9 +15,15 @@ InjectionDuration = 2.628e6 # in seconds, i.e. one month
 
 [SpatialParams]
 PermeabilityAquitard = 1e-15 # m^2
-EntryPressureAquitard = 4.5e4 # Pa
+Aquitard.BrooksCoreyPcEntry = 4.5e4 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.2
+Aquitard.Snr = 0.0
 PermeabilityAquifer = 1e-12 # m^2
-EntryPressureAquifer = 1e4 # Pa
+Aquifer.BrooksCoreyPcEntry = 1e4 # Pa
+Aquifer.BrooksCoreyLambda = 2.0
+Aquifer.Swr = 0.2
+Aquifer.Snr = 0.0
 
 # these parameters are only used in the nonisothermal model. Uncomment them for that
 [Component]
diff --git a/exercises/solution/exercise-basic/spatialparams.hh b/exercises/solution/exercise-basic/spatialparams.hh
index 715d02f2..75d073b7 100644
--- a/exercises/solution/exercise-basic/spatialparams.hh
+++ b/exercises/solution/exercise-basic/spatialparams.hh
@@ -28,11 +28,10 @@
 #define DUMUX_EX_BASIC_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 
 #include <dumux/io/gnuplotinterface.hh>
-#include <dumux/io/plotmateriallaw.hh>
+#include <dumux/io/plotpckrsw.hh>
 
 namespace Dumux {
 
@@ -55,13 +54,12 @@ class InjectionSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     // export permeability type
     using PermeabilityType = Scalar;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -69,6 +67,8 @@ public:
      */
     InjectionSpatialParams(std::shared_ptr<const FVGridGeometry>& fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
+    , aquiferPcKrSwCurve_("SpatialParams.Aquifer")
     {
         aquiferHeightFromBottom_ = 30.0;
 
@@ -79,18 +79,6 @@ public:
         // porosities
         aquitardPorosity_ = 0.2;
         aquiferPorosity_ = 0.4;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.2);
-        aquitardMaterialParams_.setSnr(0.0);
-        aquiferMaterialParams_.setSwr(0.2);
-        aquiferMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquitard"));
-        aquiferMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquifer"));
-        aquitardMaterialParams_.setLambda(2.0);
-        aquiferMaterialParams_.setLambda(2.0);
     }
 
     /*!
@@ -120,17 +108,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInAquitard_(globalPos))
-            return aquitardMaterialParams_;
-        return aquiferMaterialParams_;
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(aquiferPcKrSwCurve_); 
     }
 
     /*!
@@ -162,8 +151,8 @@ private:
     Scalar aquitardPorosity_;
     Scalar aquiferPorosity_;
 
-    MaterialLawParams aquitardMaterialParams_;
-    MaterialLawParams aquiferMaterialParams_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
+    const PcKrSwCurve aquiferPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/solution/exercise-biomineralization/biominspatialparams.hh b/exercises/solution/exercise-biomineralization/biominspatialparams.hh
index 2c901dc9..1f70ef5d 100644
--- a/exercises/solution/exercise-biomineralization/biominspatialparams.hh
+++ b/exercises/solution/exercise-biomineralization/biominspatialparams.hh
@@ -26,8 +26,7 @@
 
 #include <dumux/material/spatialparams/fv.hh>
 #include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 #include <dumux/material/fluidmatrixinteractions/porosityprecipitation.hh>
 #include <dumux/material/fluidmatrixinteractions/permeabilitykozenycarman.hh>
 
@@ -46,7 +45,6 @@ class BioMinSpatialparams
     using FVElementGeometry = typename FVGridGeometry::LocalView;
     using SubControlVolume = typename FVElementGeometry::SubControlVolume;
     using ParentType = FVSpatialParams<FVGridGeometry, Scalar, BioMinSpatialparams<FVGridGeometry, Scalar, numFluidComps, numActiveSolidComps>>;
-    using EffectiveLaw = RegularizedBrooksCorey<Scalar>;
 
     using GridView = typename FVGridGeometry::GridView;
     using CoordScalar = typename GridView::ctype;
@@ -58,10 +56,10 @@ class BioMinSpatialparams
 
     using PoroLaw = PorosityPrecipitation<Scalar, numFluidComps, numActiveSolidComps>;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     using PermeabilityType = Tensor;
-    using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-    using MaterialLawParams = typename MaterialLaw::Params;
 
     /*!
     * \brief The constructor
@@ -70,6 +68,8 @@ public:
     */
     BioMinSpatialparams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , pcKrSwCurve_("SpatialParams")
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
     {
         //! set initial aquifer params
         initialPorosity_ = getParam<Scalar>("SpatialParams.InitialPorosity");
@@ -77,26 +77,10 @@ public:
 
         // set main diagonal entries of the permeability tensor to a value
         // setting to one value means: isotropic, homogeneous
-
-        // residual saturations
-        materialParams_.setSwr(0.2);
-        materialParams_.setSnr(0.05);
-
-        // parameters for the Brooks-Corey law
-        materialParams_.setPe(1e4);
-        materialParams_.setLambda(2.0);
-
+        
         //! hard code specific params for aquitard layer
         aquitardPorosity_ = 0.1;
         aquitardPermeability_ = 1e-15;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.0);
-        aquitardMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(1e7);
-        aquitardMaterialParams_.setLambda(2.0);
     }
 
     template<class SolidSystem, class ElementSolution>
@@ -261,16 +245,18 @@ public:
     }
 
     /*!
-     * \brief return the parameter object for the Brooks-Corey material law which depends on the position
+     * \brief Returns the parameters for the material law at a given location
+     *
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \param globalPos The global position
+     * \param globalPos The global coordinates for the given location
      */
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
-        if (isInAquitard_(globalPos) && !isFaultZone_(globalPos))
-            return aquitardMaterialParams_;
-        else
-            return materialParams_;
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
+        if (isInAquitard_(globalPos))
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(pcKrSwCurve_); 
     }
 
     // define which phase is to be considered as the wetting phase
@@ -298,11 +284,11 @@ private:
     Scalar initialPermeability_;
     std::vector< std::vector<PermeabilityType> > referencePermeability_;
 
-    MaterialLawParams materialParams_;
-
     Scalar aquitardPorosity_;
     Scalar aquitardPermeability_;
-    MaterialLawParams aquitardMaterialParams_;
+
+    const PcKrSwCurve pcKrSwCurve_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/solution/exercise-biomineralization/params.input b/exercises/solution/exercise-biomineralization/params.input
index fb2df919..f9b0a05f 100644
--- a/exercises/solution/exercise-biomineralization/params.input
+++ b/exercises/solution/exercise-biomineralization/params.input
@@ -26,6 +26,14 @@ ConcUrea = 60 # [kg/m³] injected urea concentration (max: 80)
 [SpatialParams]
 InitialPorosity = 0.35 # [-]
 InitialPermeability = 1e-13 # [m^2]
+BrooksCoreyPcEntry = 1e4 # Pa
+BrooksCoreyLambda = 2.0
+Swr = 0.2
+Snr = 0.05
+Aquitard.BrooksCoreyPcEntry = 1e7 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.0
+Aquitard.Snr = 0.0
 
 [BioCoefficients]
 RhoBiofilm = 6.9 # [kg/m³] density of biofilm
diff --git a/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh b/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh
index 23930847..7a6f7360 100644
--- a/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh
@@ -25,8 +25,7 @@
 #define DUMUX_TWOPHASE_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 #include <dumux/material/fluidmatrixinteractions/2p/thermalconductivity/somerton.hh>
 
 namespace Dumux {
@@ -47,28 +46,19 @@ class TwoPSpatialParams
     using ParentType = FVSpatialParams<FVGridGeometry, Scalar, TwoPSpatialParams<FVGridGeometry, Scalar>>;
 
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
 public:
-    using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-    using MaterialLawParams = typename MaterialLaw::Params;
     using PermeabilityType = Scalar;
 
     TwoPSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , pcKrSwCurve_("Darcy.SpatialParams")
     {
         permeability_ = getParam<Scalar>("Darcy.SpatialParams.Permeability");
         porosity_ = getParam<Scalar>("Darcy.SpatialParams.Porosity");
         alphaBJ_ = getParam<Scalar>("Darcy.SpatialParams.AlphaBeaversJoseph");
-
-        // residual saturations
-        params_.setSwr(getParam<Scalar>("Darcy.SpatialParams.Swr"));
-        params_.setSnr(getParam<Scalar>("Darcy.SpatialParams.Snr"));
-        // parameters for the vanGenuchten law
-        params_.setVgAlpha(getParam<Scalar>("Darcy.SpatialParams.VgAlpha"));
-        params_.setVgn(getParam<Scalar>("Darcy.SpatialParams.VgN"));
-        params_.setPcLowSw(params_.swr()*5.0);
-        params_.setPcHighSw(1.0-params_.snr()*5.0);
     }
 
     /*!
@@ -95,19 +85,17 @@ public:
     { return alphaBJ_; }
 
     /*!
-     * \brief Returns the parameter object for the Brooks-Corey material law.
-     *        In this test, we use element-wise distributed material parameters.
+     * \brief Returns the parameters for the material law at a given location
+     *
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \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
+     * \param globalPos The global coordinates for the given location
      */
-    template<class ElementSolutionVector>
-    const MaterialLawParams& materialLawParams(const Element& element,
-                                               const SubControlVolume& scv,
-                                               const ElementSolutionVector& elemSol) const
-    { return params_; }
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {
+        return makeFluidMatrixInteraction(pcKrSwCurve_);
+    }
 
     /*!
      * \brief Function for defining which phase is to be considered as the wetting phase.
@@ -123,7 +111,7 @@ private:
     Scalar permeability_;
     Scalar porosity_;
     Scalar alphaBJ_;
-    MaterialLawParams params_;
+    const PcKrSwCurve pcKrSwCurve_;
     static constexpr Scalar eps_ = 1.0e-7;
 };
 
diff --git a/exercises/solution/exercise-coupling-ff-pm/models/params.input b/exercises/solution/exercise-coupling-ff-pm/models/params.input
index 7dda6b1b..428750e4 100644
--- a/exercises/solution/exercise-coupling-ff-pm/models/params.input
+++ b/exercises/solution/exercise-coupling-ff-pm/models/params.input
@@ -32,10 +32,12 @@ Permeability = 2.65e-10 # m^2
 Porosity = 0.4 # -
 AlphaBeaversJoseph = 1.0 # -
 # EXNUMBER >= 1
+VanGenuchtenN = 8.0
+VanGenuchtenAlpha = 6.5e-4
 Swr = 0.005
 Snr = 0.01
-VgAlpha = 6.5e-4
-VgN = 8.0
+PcLowSw = Swr * 5.0
+PcHighSw = 1.0 - Snr * 5.0
 
 [Problem]
 Name = ex_models_coupling
diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input b/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input
index 00e02465..47d2331c 100644
--- a/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input
+++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input
@@ -40,11 +40,12 @@ InitPhasePresence = 3 # bothPhases
 Porosity = 0.41
 Permeability = 2.65e-10
 AlphaBeaversJoseph = 1.0
+VanGenuchtenN = 6.9
+VanGenuchtenAlpha = 6.371e-4
 Swr = 0.005
 Snr = 0.01
-VgAlpha = 6.371e-4
-VgN = 6.9
-
+PcLowSw = Swr * 5.0
+PcHighSw = 1.0 - Snr * 5.0
 [Problem]
 Name = ex_coupling_turbulence_ff-pm
 EnableGravity = true
diff --git a/exercises/solution/exercise-fluidsystem/aparams.input b/exercises/solution/exercise-fluidsystem/aparams.input
index ccfa1908..594ad15c 100644
--- a/exercises/solution/exercise-fluidsystem/aparams.input
+++ b/exercises/solution/exercise-fluidsystem/aparams.input
@@ -5,6 +5,16 @@ DtInitial = 10 # initial time step size [s]
 [Problem]
 Name = exercise-fluidsystem_a_solution # name will be given to e.g. to the vtk result files
 
+[SpatialParams]
+BrooksCoreyPcEntry = 5.0e2 # Pa
+BrooksCoreyLambda = 2.0
+Swr = 0.1
+Snr = 0.0
+Lens.BrooksCoreyPcEntry = 1e3 # Pa
+Lens.BrooksCoreyLambda = 2.0
+Lens.Swr = 0.1
+Lens.Snr = 0.0
+
 [Grid]
 UpperRight = 60 60 # x-/y-coordinates of the upper-right corner of the grid [m]
 Cells = 60 60 # x-/y-resolution of the grid
diff --git a/exercises/solution/exercise-fluidsystem/bparams.input b/exercises/solution/exercise-fluidsystem/bparams.input
index 139be637..a421776f 100644
--- a/exercises/solution/exercise-fluidsystem/bparams.input
+++ b/exercises/solution/exercise-fluidsystem/bparams.input
@@ -5,6 +5,16 @@ DtInitial = 10 # initial time step size [s]
 [Problem]
 Name = exercise-fluidsystem_b_solution # name will be given to e.g. to the vtk result files
 
+[SpatialParams]
+BrooksCoreyPcEntry = 5.0e2 # Pa
+BrooksCoreyLambda = 2.0
+Swr = 0.1
+Snr = 0.0
+Lens.BrooksCoreyPcEntry = 1e3 # Pa
+Lens.BrooksCoreyLambda = 2.0
+Lens.Swr = 0.1
+Lens.Snr = 0.0
+
 [Grid]
 UpperRight = 60 60 # x-/y-coordinates of the upper-right corner of the grid [m]
 Cells = 60 60 # x-/y-resolution of the grid
diff --git a/exercises/solution/exercise-fluidsystem/spatialparams.hh b/exercises/solution/exercise-fluidsystem/spatialparams.hh
index a627c0fb..edfc4958 100644
--- a/exercises/solution/exercise-fluidsystem/spatialparams.hh
+++ b/exercises/solution/exercise-fluidsystem/spatialparams.hh
@@ -29,8 +29,7 @@
 #include <dumux/material/spatialparams/fv.hh>
 
 // include material laws
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 #include <dumux/material/fluidmatrixinteractions/2p/linearmaterial.hh>
 
 namespace Dumux {
@@ -54,13 +53,12 @@ class ExerciseFluidsystemSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     // export permeability type
     using PermeabilityType = Dune::FieldMatrix<Scalar, dim, dim>;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -70,6 +68,8 @@ public:
     : ParentType(fvGridGeometry)
     , K_(0)
     , KLens_(0)
+    , pcKrSwCurve_("SpatialParams")
+    , lensPcKrSwCurve_("Lens.SpatialParams")
     {
         //set main diagonal entries of the permeability tensor to a value
         //setting to one value means: isotropic, homogeneous
@@ -78,21 +78,8 @@ public:
             K_[i][i] = 1e-7;
             KLens_[i][i] = 1e-10;
         }
-
-        //set residual saturations
-        materialParams_.setSwr(0.0);
-        materialParamsLens_.setSwr(0.1);
-        materialParams_.setSnr(0.0);
-        materialParamsLens_.setSnr(0.1);
-
-        //parameters of Brooks & Corey Law
-        materialParams_.setPe(500.0);
-        materialParamsLens_.setPe(1000.0);
-        materialParams_.setLambda(2);
-        materialParamsLens_.setLambda(2);
     }
 
-
     /*!
      * \brief Define the intrinsic permeability \f$\mathrm{[m^2]}\f$.
      *
@@ -119,17 +106,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInLens(globalPos))
-            return materialParamsLens_;
-        return materialParams_;
+            return makeFluidMatrixInteraction(pcKrSwCurve_); 
+        return makeFluidMatrixInteraction(lensPcKrSwCurve_); 
     }
 
     /*!
@@ -170,8 +158,8 @@ private:
     Dune::FieldMatrix<Scalar, dim, dim> KLens_;
 
     // Object that holds the values/parameters of the selected material law.
-    MaterialLawParams materialParams_;
-    MaterialLawParams materialParamsLens_;
+    const PcKrSwCurve pcKrSwCurve_;
+    const PcKrSwCurve lensPcKrSwCurve_;
 };
 } // end namespace Dumux
 #endif
diff --git a/exercises/solution/exercise-fractures/exercise_fractures_solution_a.input b/exercises/solution/exercise-fractures/exercise_fractures_solution_a.input
index e9294a23..e57f960d 100644
--- a/exercises/solution/exercise-fractures/exercise_fractures_solution_a.input
+++ b/exercises/solution/exercise-fractures/exercise_fractures_solution_a.input
@@ -19,8 +19,8 @@ Problem.BoundaryOverPressure = 5e5
 Problem.BoundarySaturation = 0.5
 SpatialParams.Permeability = 1e-12
 SpatialParams.Porosity = 0.1
-SpatialParams.VGAlpha = 1e-3
-SpatialParams.VGN = 3
+SpatialParams.VanGenuchtenAlpha = 1e-3
+SpatialParams.VanGenuchtenN = 3
 SpatialParams.Snr = 0.0
 SpatialParams.Swr = 0.0
 
@@ -31,11 +31,11 @@ SpatialParams.Permeability = 1e-7
 SpatialParams.PermeabilityBarrier = 1e-16
 SpatialParams.Porosity = 0.85
 SpatialParams.PorosityBarrier = 0.05
-SpatialParams.VGAlpha = 1e-1
-SpatialParams.VGAlphaBarrier = 1e-4
-SpatialParams.VGN = 2
-SpatialParams.VGNBarrier = 2.5
-SpatialParams.Snr = 0.0
-SpatialParams.SnrBarrier = 0.0
+SpatialParams.VanGenuchtenAlpha = 1e-1
+SpatialParams.VanGenuchtenN = 2
 SpatialParams.Swr = 0.0
-SpatialParams.SwrBarrier = 0.0
+SpatialParams.Snr = 0.0
+SpatialParams.Barrier.VanGenuchtenAlpha = 1e-4
+SpatialParams.Barrier.VanGenuchtenN = 2.5
+SpatialParams.Barrier.Snr = 0.0
+SpatialParams.Barrier.Swr = 0.0
diff --git a/exercises/solution/exercise-fractures/exercise_fractures_solution_b.input b/exercises/solution/exercise-fractures/exercise_fractures_solution_b.input
index 3797db3c..e8630924 100644
--- a/exercises/solution/exercise-fractures/exercise_fractures_solution_b.input
+++ b/exercises/solution/exercise-fractures/exercise_fractures_solution_b.input
@@ -19,8 +19,8 @@ Problem.BoundaryOverPressure = 5e5
 Problem.BoundarySaturation = 0.5
 SpatialParams.Permeability = 1e-12
 SpatialParams.Porosity = 0.1
-SpatialParams.VGAlpha = 1e-3
-SpatialParams.VGN = 3
+SpatialParams.VanGenuchtenAlpha = 1e-3
+SpatialParams.VanGenuchtenN = 3
 SpatialParams.Snr = 0.0
 SpatialParams.Swr = 0.0
 
@@ -31,11 +31,11 @@ SpatialParams.Permeability = 1e-7
 SpatialParams.PermeabilityBarrier = 1e-16
 SpatialParams.Porosity = 0.85
 SpatialParams.PorosityBarrier = 0.05
-SpatialParams.VGAlpha = 1e-1
-SpatialParams.VGAlphaBarrier = 1e-4
-SpatialParams.VGN = 2
-SpatialParams.VGNBarrier = 2.5
-SpatialParams.Snr = 0.0
-SpatialParams.SnrBarrier = 0.0
+SpatialParams.VanGenuchtenAlpha = 1e-1
+SpatialParams.VanGenuchtenN = 2
 SpatialParams.Swr = 0.0
-SpatialParams.SwrBarrier = 0.0
+SpatialParams.Snr = 0.0
+SpatialParams.Barrier.VanGenuchtenAlpha = 1e-4
+SpatialParams.Barrier.VanGenuchtenN = 2.5
+SpatialParams.Barrier.Snr = 0.0
+SpatialParams.Barrier.Swr = 0.0
diff --git a/exercises/solution/exercise-fractures/exercise_fractures_solution_c.input b/exercises/solution/exercise-fractures/exercise_fractures_solution_c.input
index 82dd9dac..51d92acf 100644
--- a/exercises/solution/exercise-fractures/exercise_fractures_solution_c.input
+++ b/exercises/solution/exercise-fractures/exercise_fractures_solution_c.input
@@ -19,8 +19,8 @@ Problem.BoundaryOverPressure = 5e5
 Problem.BoundarySaturation = 0.5
 SpatialParams.Permeability = 1e-12
 SpatialParams.Porosity = 0.1
-SpatialParams.VGAlpha = 1e-3
-SpatialParams.VGN = 3
+SpatialParams.VanGenuchtenAlpha = 1e-3
+SpatialParams.VanGenuchtenN = 3
 SpatialParams.Snr = 0.0
 SpatialParams.Swr = 0.0
 
@@ -31,11 +31,11 @@ SpatialParams.Permeability = 1e-7
 SpatialParams.PermeabilityBarrier = 1e-16
 SpatialParams.Porosity = 0.85
 SpatialParams.PorosityBarrier = 0.05
-SpatialParams.VGAlpha = 1e-1
-SpatialParams.VGAlphaBarrier = 1e-4
-SpatialParams.VGN = 2
-SpatialParams.VGNBarrier = 2.5
-SpatialParams.Snr = 0.0
-SpatialParams.SnrBarrier = 0.0
+SpatialParams.VanGenuchtenAlpha = 1e-1
+SpatialParams.VanGenuchtenN = 2
 SpatialParams.Swr = 0.0
-SpatialParams.SwrBarrier = 0.0
+SpatialParams.Snr = 0.0
+SpatialParams.Barrier.VanGenuchtenAlpha = 1e-4
+SpatialParams.Barrier.VanGenuchtenN = 2.5
+SpatialParams.Barrier.Snr = 0.0
+SpatialParams.Barrier.Swr = 0.0
diff --git a/exercises/solution/exercise-fractures/fracturespatialparams.hh b/exercises/solution/exercise-fractures/fracturespatialparams.hh
index e98e8aea..264c0b66 100644
--- a/exercises/solution/exercise-fractures/fracturespatialparams.hh
+++ b/exercises/solution/exercise-fractures/fracturespatialparams.hh
@@ -30,8 +30,7 @@
 #include <dumux/io/grid/griddata.hh>
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 
 namespace Dumux {
 
@@ -55,7 +54,7 @@ class FractureSpatialParams
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
     // use a regularized van-genuchten material law
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
     // we identify those fractures as barriers, that have a domain marker
     // of 2. This is what is set in the grid file (see grids/complex.geo)
@@ -65,16 +64,14 @@ public:
     //! export the type used for permeabilities
     using PermeabilityType = Scalar;
 
-    //! export the material law and parameters used
-    using MaterialLaw = EffToAbsLaw< EffectiveLaw >;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     //! the constructor
     FractureSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
                           std::shared_ptr<const Dumux::GridData<Grid>> gridData,
                           const std::string& paramGroup)
     : ParentType(fvGridGeometry)
     , gridDataPtr_(gridData)
+    , pcKrSwCurve_("Fracture.SpatialParams")
+    , barrierPcKrSwCurve_("Fracture.SpatialParams.Barrier")
     , isExercisePartA_(getParamFromGroup<bool>(paramGroup, "Problem.IsExercisePartA"))
     , isExercisePartB_(getParamFromGroup<bool>(paramGroup, "Problem.IsExercisePartB"))
     , isExercisePartC_(getParamFromGroup<bool>(paramGroup, "Problem.IsExercisePartC"))
@@ -83,16 +80,6 @@ public:
         porosityBarrier_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.PorosityBarrier");
         permeability_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Permeability");
         permeabilityBarrier_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.PermeabilityBarrier");
-
-        // set the material law parameters
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Snr"));
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Swr"));
-        materialLawParams_.setVgAlpha(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGAlpha"));
-        materialLawParams_.setVgn(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGN"));
-        materialLawParamsBarrier_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.SnrBarrier"));
-        materialLawParamsBarrier_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.SwrBarrier"));
-        materialLawParamsBarrier_.setVgAlpha(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGAlphaBarrier"));
-        materialLawParamsBarrier_.setVgn(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGNBarrier"));
     }
 
     //! Function for defining the (intrinsic) permeability \f$[m^2]\f$.
@@ -143,27 +130,36 @@ public:
         }
     }
 
-    //! Return the material law parameters
-    template< class ElementSolution >
-    const MaterialLawParams& materialLawParams(const Element& element,
-                                               const SubControlVolume& scv,
-                                               const ElementSolution& elemSol) const
+    /*!
+     * \brief Returns the parameters for the material law for the sub-control volume
+     *
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
+     *
+     * \param element The current finite element
+     * \param scv The sub-control volume
+     * \param elemSol The current element solution
+     */
+    template<class ElementSolution>
+    auto fluidMatrixInteraction(const Element& element,
+                                const SubControlVolume& scv,
+                                const ElementSolution& elemSol) const
     {
         // in the unmodified state and exercise part a return the non-barrier parameters
         if (!isExercisePartB_ && !isExercisePartC_)
-            return materialLawParams_;
+            return makeFluidMatrixInteraction(pcKrSwCurve_);
 
         // in exercise part b always return the barrier parameters
         else if (isExercisePartB_)
-            return materialLawParamsBarrier_;
+            return makeFluidMatrixInteraction(barrierPcKrSwCurve_);
 
         // in exercise part 3 return parameters depending on domain marker
         else
         {
             if (getElementDomainMarker(element) == barriersDomainMarker)
-                return materialLawParamsBarrier_;
+                return makeFluidMatrixInteraction(barrierPcKrSwCurve_);
             else
-                return materialLawParams_;
+                return makeFluidMatrixInteraction(pcKrSwCurve_);
         }
     }
 
@@ -192,8 +188,8 @@ private:
     Scalar porosityBarrier_;
     PermeabilityType permeability_;
     PermeabilityType permeabilityBarrier_;
-    MaterialLawParams materialLawParams_;
-    MaterialLawParams materialLawParamsBarrier_;
+    const PcKrSwCurve pcKrSwCurve_;
+    const PcKrSwCurve barrierPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/solution/exercise-fractures/matrixspatialparams.hh b/exercises/solution/exercise-fractures/matrixspatialparams.hh
index 4ca8c72e..e8f45aa8 100644
--- a/exercises/solution/exercise-fractures/matrixspatialparams.hh
+++ b/exercises/solution/exercise-fractures/matrixspatialparams.hh
@@ -30,8 +30,7 @@
 #include <dumux/io/grid/griddata.hh>
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 
 namespace Dumux {
 
@@ -54,31 +53,22 @@ class MatrixSpatialParams
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
     // use a regularized van-genuchten material law
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
 public:
     //! export the type used for permeabilities
     using PermeabilityType = Scalar;
 
-    //! export the material law and parameters used
-    using MaterialLaw = EffToAbsLaw< EffectiveLaw >;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     //! the constructor
     MatrixSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
                         std::shared_ptr<const Dumux::GridData<Grid>> gridData,
                         const std::string& paramGroup)
     : ParentType(fvGridGeometry)
     , gridDataPtr_(gridData)
+    , pcKrSwCurve_("Matrix.SpatialParams")
     {
         porosity_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Porosity");
         permeability_ = getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Permeability");
-
-        // set the material law parameters
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Snr"));
-        materialLawParams_.setSnr(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.Swr"));
-        materialLawParams_.setVgAlpha(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGAlpha"));
-        materialLawParams_.setVgn(getParamFromGroup<Scalar>(paramGroup, "SpatialParams.VGN"));
     }
 
     //! Function for defining the (intrinsic) permeability \f$[m^2]\f$.
@@ -89,9 +79,18 @@ public:
     Scalar porosityAtPos(const GlobalPosition& globalPos) const
     { return porosity_; }
 
-    //! Return the material law parameters
-    const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    { return materialLawParams_; }
+    /*!
+     * \brief Returns the parameters for the material law at a given location
+     *
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
+     *
+     * \param globalPos The global coordinates for the given location
+     */
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
+        return makeFluidMatrixInteraction(pcKrSwCurve_); 
+    }
 
     //! Water is the wetting phase
     template< class FluidSystem >
@@ -112,7 +111,7 @@ private:
 
     Scalar porosity_;
     PermeabilityType permeability_;
-    MaterialLawParams materialLawParams_;
+    const PcKrSwCurve pcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/solution/exercise-grids/params.input b/exercises/solution/exercise-grids/params.input
index 75af6c79..04d20959 100644
--- a/exercises/solution/exercise-grids/params.input
+++ b/exercises/solution/exercise-grids/params.input
@@ -40,6 +40,12 @@ InjectionDuration = 2.628e6 # in seconds, i.e. one month
 
 [SpatialParams]
 PermeabilityAquitard = 1e-15 # m^2
-EntryPressureAquitard = 4.5e4 # Pa
+Aquitard.BrooksCoreyPcEntry = 4.5e4 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.2
+Aquitard.Snr = 0.0
 PermeabilityAquifer = 1e-12 # m^2
-EntryPressureAquifer = 1e4 # Pa
+Aquifer.BrooksCoreyPcEntry = 1e4 # Pa
+Aquifer.BrooksCoreyLambda = 2.0
+Aquifer.Swr = 0.2
+Aquifer.Snr = 0.0
diff --git a/exercises/solution/exercise-grids/spatialparams.hh b/exercises/solution/exercise-grids/spatialparams.hh
index 0e72f8d5..ea7e3084 100644
--- a/exercises/solution/exercise-grids/spatialparams.hh
+++ b/exercises/solution/exercise-grids/spatialparams.hh
@@ -28,11 +28,10 @@
 #define DUMUX_EXGRIDS_INJECTION_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 
 #include <dumux/io/gnuplotinterface.hh>
-#include <dumux/io/plotmateriallaw.hh>
+#include <dumux/io/plotpckrsw.hh>
 
 namespace Dumux {
 
@@ -55,13 +54,12 @@ class InjectionSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     // export permeability type
     using PermeabilityType = Scalar;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -69,6 +67,8 @@ public:
      */
     InjectionSpatialParams(std::shared_ptr<const FVGridGeometry>& fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
+    , aquiferPcKrSwCurve_("SpatialParams.Aquifer")
     {
         // Aquifer Height, measured from the bottom
         aquiferHeightFromBottom_ = 30.0;
@@ -80,18 +80,6 @@ public:
         // porosities
         aquitardPorosity_ = 0.2;
         aquiferPorosity_ = 0.4;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.2);
-        aquitardMaterialParams_.setSnr(0.0);
-        aquiferMaterialParams_.setSwr(0.2);
-        aquiferMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquitard"));
-        aquiferMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquifer"));
-        aquitardMaterialParams_.setLambda(2.0);
-        aquiferMaterialParams_.setLambda(2.0);
     }
 
     /*!
@@ -121,17 +109,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-     const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInAquitard_(globalPos))
-            return aquitardMaterialParams_;
-        return aquiferMaterialParams_;
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(aquiferPcKrSwCurve_); 
     }
 
     /*!
@@ -163,8 +152,8 @@ private:
     Scalar aquitardPorosity_;
     Scalar aquiferPorosity_;
 
-    MaterialLawParams aquitardMaterialParams_;
-    MaterialLawParams aquiferMaterialParams_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
+    const PcKrSwCurve aquiferPcKrSwCurve_;
 };
 
 } // end namespace Dumux
diff --git a/exercises/solution/exercise-properties/params.input b/exercises/solution/exercise-properties/params.input
index 7ae1761e..8e25e5f8 100644
--- a/exercises/solution/exercise-properties/params.input
+++ b/exercises/solution/exercise-properties/params.input
@@ -10,6 +10,14 @@ Cells = 48 32
 [SpatialParams]
 LensLowerLeft = 1.0 2.0 # [m] coordinates of the lower left lens corner
 LensUpperRight = 4.0 3.0 # [m] coordinates of the upper right lens corner
+Lens.VanGenuchtenN = 7.3
+Lens.VanGenuchtenAlpha = 0.00045
+Lens.Swr = 0.18
+Lens.Snr = 0.0
+OuterDomain.VanGenuchtenN = 4.7
+OuterDomain.VanGenuchtenAlpha = 0.0037
+OuterDomain.Swr = 0.05
+OuterDomain.Snr = 0.0
 
 [Problem]
 Name = exercise_properties
diff --git a/exercises/solution/exercise-properties/spatialparams.hh b/exercises/solution/exercise-properties/spatialparams.hh
index 275984d8..fa106e87 100644
--- a/exercises/solution/exercise-properties/spatialparams.hh
+++ b/exercises/solution/exercise-properties/spatialparams.hh
@@ -24,8 +24,7 @@
 #define DUMUX_INCOMPRESSIBLE_TWOP_TEST_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedvangenuchten.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/vangenuchten.hh>
 
 namespace Dumux {
 
@@ -47,32 +46,19 @@ class TwoPTestSpatialParams
     static constexpr int dimWorld = GridView::dimensionworld;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
-    using EffectiveLaw = RegularizedVanGenuchten<Scalar>;
+    using PcKrSwCurve = FluidMatrix::VanGenuchtenDefault<Scalar>;
 
 public:
-    using MaterialLaw = EffToAbsLaw<EffectiveLaw>;
-    using MaterialLawParams = typename MaterialLaw::Params;
     using PermeabilityType = Scalar;
 
     TwoPTestSpatialParams(std::shared_ptr<const FVGridGeometry> fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , lensPcKrSwCurve_("SpatialParams.Lens")
+    , outerPcKrSwCurve_("SpatialParams.OuterDomain")
     {
         lensLowerLeft_ = getParam<GlobalPosition>("SpatialParams.LensLowerLeft");
         lensUpperRight_ = getParam<GlobalPosition>("SpatialParams.LensUpperRight");
 
-        // residual saturations
-        lensMaterialParams_.setSwr(0.18);
-        lensMaterialParams_.setSnr(0.0);
-        outerMaterialParams_.setSwr(0.05);
-        outerMaterialParams_.setSnr(0.0);
-
-        // parameters for the Van Genuchten law
-        // alpha and n
-        lensMaterialParams_.setVgAlpha(0.00045);
-        lensMaterialParams_.setVgn(7.3);
-        outerMaterialParams_.setVgAlpha(0.0037);
-        outerMaterialParams_.setVgn(4.7);
-
         lensK_ = 9.05e-12;
         outerK_ = 4.6e-10;
     }
@@ -105,22 +91,19 @@ public:
     { return 0.4; }
 
     /*!
-     * \brief Returns the parameter object for the Brooks-Corey material law.
-     *        In this test, we use element-wise distributed material parameters.
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \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
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
+     *
+     * \param globalPos The global coordinates for the given location
      */
-    template<class ElementSolution>
-    const MaterialLawParams& materialLawParams(const Element& element,
-                                               const SubControlVolume& scv,
-                                               const ElementSolution& elemSol) const
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
     {
-        if (isInLens_(element.geometry().center()))
-            return lensMaterialParams_;
-        return outerMaterialParams_;
+        if (isInLens_(globalPos))
+            return makeFluidMatrixInteraction(lensPcKrSwCurve_);
+
+        return makeFluidMatrixInteraction(outerPcKrSwCurve_);
     }
 
     /*!
@@ -150,8 +133,9 @@ private:
 
     Scalar lensK_;
     Scalar outerK_;
-    MaterialLawParams lensMaterialParams_;
-    MaterialLawParams outerMaterialParams_;
+
+    const PcKrSwCurve lensPcKrSwCurve_;
+    const PcKrSwCurve outerPcKrSwCurve_;
 
     static constexpr Scalar eps_ = 1.5e-7;
 };
diff --git a/exercises/solution/exercise-runtimeparams/params.input b/exercises/solution/exercise-runtimeparams/params.input
index 62cbcc86..542e29f6 100644
--- a/exercises/solution/exercise-runtimeparams/params.input
+++ b/exercises/solution/exercise-runtimeparams/params.input
@@ -17,6 +17,12 @@ TotalAreaSpecificInflow = -1e-4 # kg/(s m^2)
 [SpatialParams]
 PermeabilityAquitard = 1e-15 # m^2
 # TODO: Task 1: Change the Aquitard's Entry Pressure
-EntryPressureAquitard = 4.5e3 # Pa
+Aquitard.BrooksCoreyPcEntry = 4.5e3 # Pa
+Aquitard.BrooksCoreyLambda = 2.0
+Aquitard.Swr = 0.2
+Aquitard.Snr = 0.0
 PermeabilityAquifer = 1e-12 # m^2
-EntryPressureAquifer = 1e4 # Pa
+Aquifer.BrooksCoreyPcEntry = 1e4 # Pa
+Aquifer.BrooksCoreyLambda = 2.0
+Aquifer.Swr = 0.2
+Aquifer.Snr = 0.0
\ No newline at end of file
diff --git a/exercises/solution/exercise-runtimeparams/spatialparams.hh b/exercises/solution/exercise-runtimeparams/spatialparams.hh
index b42f8114..617cd87a 100644
--- a/exercises/solution/exercise-runtimeparams/spatialparams.hh
+++ b/exercises/solution/exercise-runtimeparams/spatialparams.hh
@@ -28,11 +28,10 @@
 #define DUMUX_RUNTIMEPARAMS_INJECTION_SPATIAL_PARAMS_HH
 
 #include <dumux/material/spatialparams/fv.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh>
-#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
+#include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh>
 
 #include <dumux/io/gnuplotinterface.hh>
-#include <dumux/io/plotmateriallaw.hh>
+#include <dumux/io/plotpckrsw.hh>
 
 namespace Dumux {
 
@@ -55,13 +54,12 @@ class InjectionSpatialParams
     using Element = typename GridView::template Codim<0>::Entity;
     using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
 
+    using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>;
+
 public:
     // export permeability type
     using PermeabilityType = Scalar;
 
-    using MaterialLaw = EffToAbsLaw<RegularizedBrooksCorey<Scalar>>;
-    using MaterialLawParams = typename MaterialLaw::Params;
-
     /*!
      * \brief The constructor
      *
@@ -69,6 +67,8 @@ public:
      */
     InjectionSpatialParams(std::shared_ptr<const FVGridGeometry>& fvGridGeometry)
     : ParentType(fvGridGeometry)
+    , aquitardPcKrSwCurve_("SpatialParams.Aquitard")
+    , aquiferPcKrSwCurve_("SpatialParams.Aquifer")
     {
         // Aquifer Height, measured from the bottom
         aquiferHeightFromBottom_ = 30.0;
@@ -80,18 +80,6 @@ public:
         // porosities
         aquitardPorosity_ = 0.2;
         aquiferPorosity_ = 0.4;
-
-        // residual saturations
-        aquitardMaterialParams_.setSwr(0.2);
-        aquitardMaterialParams_.setSnr(0.0);
-        aquiferMaterialParams_.setSwr(0.2);
-        aquiferMaterialParams_.setSnr(0.0);
-
-        // parameters for the Brooks-Corey law
-        aquitardMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquitard"));
-        aquiferMaterialParams_.setPe(getParam<Scalar>("SpatialParams.EntryPressureAquifer"));
-        aquitardMaterialParams_.setLambda(2.0);
-        aquiferMaterialParams_.setLambda(2.0);
     }
 
     /*!
@@ -121,17 +109,18 @@ public:
     }
 
     /*!
-     * \brief Function for defining the parameters needed by constitutive relationships (kr-sw, pc-sw, etc.).
+     * \brief Returns the parameters for the material law at a given location
      *
-     * \param globalPos The global position
+     * This method is not actually required by the Richards model, but provided
+     * for the convenience of the RichardsLensProblem
      *
-     * \return the material parameters object
+     * \param globalPos The global coordinates for the given location
      */
-     const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
-    {
+    auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const
+    {   
         if (isInAquitard_(globalPos))
-            return aquitardMaterialParams_;
-        return aquiferMaterialParams_;
+            return makeFluidMatrixInteraction(aquitardPcKrSwCurve_); 
+        return makeFluidMatrixInteraction(aquiferPcKrSwCurve_); 
     }
 
     /*!
@@ -159,12 +148,11 @@ private:
     Scalar aquiferK_;
     Scalar aquiferHeightFromBottom_;
 
-
     Scalar aquitardPorosity_;
     Scalar aquiferPorosity_;
 
-    MaterialLawParams aquitardMaterialParams_;
-    MaterialLawParams aquiferMaterialParams_;
+    const PcKrSwCurve aquitardPcKrSwCurve_;
+    const PcKrSwCurve aquiferPcKrSwCurve_;
 };
 
 } // end namespace Dumux
-- 
GitLab