diff --git a/Applications/ApplicationsLib/ProjectData.cpp b/Applications/ApplicationsLib/ProjectData.cpp
index 17f63414460d7c0ce0e0ec55e490d32300aac4ac..c95940db3276f6828e5ab43da5bff32614a45ea4 100644
--- a/Applications/ApplicationsLib/ProjectData.cpp
+++ b/Applications/ApplicationsLib/ProjectData.cpp
@@ -903,7 +903,7 @@ void ProjectData::parseProcesses(BaseLib::ConfigTree const& processes_config,
                             name, *_mesh_vec[0], std::move(jacobian_assembler),
                             _process_variables, _parameters,
                             _local_coordinate_system, integration_order,
-                            process_config);
+                            process_config, _media);
                     break;
                 case 3:
                     process = ProcessLib::RichardsMechanics::
@@ -911,7 +911,7 @@ void ProjectData::parseProcesses(BaseLib::ConfigTree const& processes_config,
                             name, *_mesh_vec[0], std::move(jacobian_assembler),
                             _process_variables, _parameters,
                             _local_coordinate_system, integration_order,
-                            process_config);
+                            process_config, _media);
                     break;
             }
         }
diff --git a/ProcessLib/RichardsMechanics/CMakeLists.txt b/ProcessLib/RichardsMechanics/CMakeLists.txt
index 78346233e41f233c2471db93d0166d79433d4cd9..ca248293a8e7936e226451d7e7e8fd6d34ef5057 100644
--- a/ProcessLib/RichardsMechanics/CMakeLists.txt
+++ b/ProcessLib/RichardsMechanics/CMakeLists.txt
@@ -5,8 +5,9 @@ if(BUILD_SHARED_LIBS)
     install(TARGETS RichardsMechanics
             LIBRARY DESTINATION ${CMAKE_INSTALL_LIBDIR})
 endif()
-target_link_libraries(RichardsMechanics
-                      PUBLIC ProcessLib
-                      PRIVATE RichardsFlow ParameterLib)
+target_link_libraries(
+    RichardsMechanics
+    PUBLIC ProcessLib
+    PRIVATE ParameterLib)
 
 include(Tests.cmake)
diff --git a/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.cpp b/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.cpp
index 4daa4d0bb3d6150051e03340870a440d20186551..612874870f05d5f0b207517e588109a7c45309b7 100644
--- a/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.cpp
+++ b/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.cpp
@@ -12,13 +12,15 @@
 
 #include <cassert>
 
+#include "MaterialLib/MPL/CreateMaterialSpatialDistributionMap.h"
+#include "MaterialLib/MPL/MaterialSpatialDistributionMap.h"
+#include "MaterialLib/MPL/Medium.h"
 #include "MaterialLib/SolidModels/CreateConstitutiveRelation.h"
 #include "MaterialLib/SolidModels/MechanicsBase.h"
 #include "ParameterLib/Utils.h"
 #include "ProcessLib/Output/CreateSecondaryVariables.h"
 #include "ProcessLib/RichardsFlow/CreateRichardsFlowMaterialProperties.h"
 #include "ProcessLib/Utils/ProcessUtils.h"
-
 #include "RichardsMechanicsProcess.h"
 #include "RichardsMechanicsProcessData.h"
 
@@ -36,7 +38,8 @@ std::unique_ptr<Process> createRichardsMechanicsProcess(
     boost::optional<ParameterLib::CoordinateSystem> const&
         local_coordinate_system,
     unsigned const integration_order,
-    BaseLib::ConfigTree const& config)
+    BaseLib::ConfigTree const& config,
+    std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const& media)
 {
     //! \ogs_file_param{prj__processes__process__type}
     config.checkConfigParameter("type", "RICHARDS_MECHANICS");
@@ -110,37 +113,6 @@ std::unique_ptr<Process> createRichardsMechanicsProcess(
         MaterialLib::Solids::createConstitutiveRelations<DisplacementDim>(
             parameters, local_coordinate_system, config);
 
-    // Fluid bulk modulus
-    auto& fluid_bulk_modulus = ParameterLib::findParameter<double>(
-        config,
-        //! \ogs_file_param_special{prj__processes__process__RICHARDS_MECHANICS__fluid_bulk_modulus}
-        "fluid_bulk_modulus", parameters, 1, &mesh);
-    DBUG("Use '%s' as fluid bulk modulus parameter.",
-         fluid_bulk_modulus.name.c_str());
-
-    // Biot coefficient
-    auto& biot_coefficient = ParameterLib::findParameter<double>(
-        config,
-        //! \ogs_file_param_special{prj__processes__process__RICHARDS_MECHANICS__biot_coefficient}
-        "biot_coefficient", parameters, 1, &mesh);
-    DBUG("Use '%s' as Biot coefficient parameter.",
-         biot_coefficient.name.c_str());
-
-    // Solid density
-    auto& solid_density = ParameterLib::findParameter<double>(
-        config,
-        //! \ogs_file_param_special{prj__processes__process__RICHARDS_MECHANICS__solid_density}
-        "solid_density", parameters, 1, &mesh);
-    DBUG("Use '%s' as solid density parameter.", solid_density.name.c_str());
-
-    // Solid bulk modulus
-    auto& solid_bulk_modulus = ParameterLib::findParameter<double>(
-        config,
-        //! \ogs_file_param_special{prj__processes__process__RICHARDS_MECHANICS__solid_bulk_modulus}
-        "solid_bulk_modulus", parameters, 1, &mesh);
-    DBUG("Use '%s' as solid bulk modulus parameter.",
-         solid_bulk_modulus.name.c_str());
-
     // Specific body force
     Eigen::Matrix<double, DisplacementDim, 1> specific_body_force;
     {
@@ -160,6 +132,21 @@ std::unique_ptr<Process> createRichardsMechanicsProcess(
         std::copy_n(b.data(), b.size(), specific_body_force.data());
     }
 
+    auto media_map =
+        MaterialPropertyLib::createMaterialSpatialDistributionMap(media, mesh);
+
+    std::array const requiredGasProperties = {MaterialPropertyLib::viscosity};
+    std::array const requiredSolidProperties = {
+        MaterialPropertyLib::porosity, MaterialPropertyLib::biot_coefficient,
+        MaterialPropertyLib::density};
+    for (auto const& m : media)
+    {
+        checkRequiredProperties(m.second->phase("AqueousLiquid"),
+                                requiredGasProperties);
+        checkRequiredProperties(m.second->phase("Solid"),
+                                requiredSolidProperties);
+    }
+
     // Initial stress conditions
     auto const initial_stress = ParameterLib::findOptionalTagParameter<double>(
         //! \ogs_file_param_special{prj__processes__process__RICHARDS_MECHANICS__initial_stress}
@@ -168,49 +155,15 @@ std::unique_ptr<Process> createRichardsMechanicsProcess(
         MathLib::KelvinVector::KelvinVectorDimensions<DisplacementDim>::value,
         &mesh);
 
-    auto& temperature = ParameterLib::findParameter<double>(
-        config,
-        //! \ogs_file_param_special{prj__processes__process__RICHARDS_MECHANICS__temperature}
-        "temperature", parameters, 1, &mesh);
-
-    auto const& flow_material_config =
-        //! \ogs_file_param{prj__processes__process__RICHARDS_MECHANICS__material_property}
-        config.getConfigSubtree("material_property");
-
-    auto const material_ids =
-        mesh.getProperties().existsPropertyVector<int>("MaterialIDs")
-            ? mesh.getProperties().getPropertyVector<int>("MaterialIDs")
-            : nullptr;
-    if (material_ids != nullptr)
-    {
-        INFO(
-            "MaterialIDs vector found; the Richards flow is in heterogeneous "
-            "porous media.");
-    }
-    else
-    {
-        INFO(
-            "MaterialIDs vector not found; the Richards flow is in homogeneous "
-            "porous media.");
-    }
-    auto flow_material =
-        ProcessLib::RichardsFlow::createRichardsFlowMaterialProperties(
-            flow_material_config, material_ids, parameters);
-
     auto const mass_lumping =
         //! \ogs_file_param{prj__processes__process__RICHARDS_MECHANICS__mass_lumping}
         config.getConfigParameter<bool>("mass_lumping", false);
 
     RichardsMechanicsProcessData<DisplacementDim> process_data{
-        material_ids,
-        std::move(flow_material),
+        materialIDs(mesh),
+        std::move(media_map),
         std::move(solid_constitutive_relations),
         initial_stress,
-        fluid_bulk_modulus,
-        biot_coefficient,
-        solid_density,
-        solid_bulk_modulus,
-        temperature,
         specific_body_force,
         mass_lumping};
 
@@ -234,7 +187,8 @@ template std::unique_ptr<Process> createRichardsMechanicsProcess<2>(
     boost::optional<ParameterLib::CoordinateSystem> const&
         local_coordinate_system,
     unsigned const integration_order,
-    BaseLib::ConfigTree const& config);
+    BaseLib::ConfigTree const& config,
+    std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const& media);
 
 template std::unique_ptr<Process> createRichardsMechanicsProcess<3>(
     std::string name,
@@ -245,7 +199,8 @@ template std::unique_ptr<Process> createRichardsMechanicsProcess<3>(
     boost::optional<ParameterLib::CoordinateSystem> const&
         local_coordinate_system,
     unsigned const integration_order,
-    BaseLib::ConfigTree const& config);
+    BaseLib::ConfigTree const& config,
+    std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const& media);
 
 }  // namespace RichardsMechanics
 }  // namespace ProcessLib
diff --git a/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.h b/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.h
index 986cd0b61e46346d436e63fdbbcc2e9ba6269469..77f7a2196beb39dce7b37e6b063cf1b5abd9d4d5 100644
--- a/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.h
+++ b/ProcessLib/RichardsMechanics/CreateRichardsMechanicsProcess.h
@@ -27,6 +27,10 @@ namespace MathLib
 {
 class PiecewiseLinearInterpolation;
 }
+namespace MaterialPropertyLib
+{
+class Medium;
+}
 namespace ParameterLib
 {
 struct CoordinateSystem;
@@ -53,7 +57,8 @@ std::unique_ptr<Process> createRichardsMechanicsProcess(
     boost::optional<ParameterLib::CoordinateSystem> const&
         local_coordinate_system,
     unsigned const integration_order,
-    BaseLib::ConfigTree const& config);
+    BaseLib::ConfigTree const& config,
+    std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const& media);
 
 extern template std::unique_ptr<Process> createRichardsMechanicsProcess<2>(
     std::string name,
@@ -64,7 +69,8 @@ extern template std::unique_ptr<Process> createRichardsMechanicsProcess<2>(
     boost::optional<ParameterLib::CoordinateSystem> const&
         local_coordinate_system,
     unsigned const integration_order,
-    BaseLib::ConfigTree const& config);
+    BaseLib::ConfigTree const& config,
+    std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const& media);
 
 extern template std::unique_ptr<Process> createRichardsMechanicsProcess<3>(
     std::string name,
@@ -75,6 +81,7 @@ extern template std::unique_ptr<Process> createRichardsMechanicsProcess<3>(
     boost::optional<ParameterLib::CoordinateSystem> const&
         local_coordinate_system,
     unsigned const integration_order,
-    BaseLib::ConfigTree const& config);
+    BaseLib::ConfigTree const& config,
+    std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const& media);
 }  // namespace RichardsMechanics
 }  // namespace ProcessLib
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
index 4176bf6bc40d4b3fb9579c0933a997f47f3afd98..51e1b0eb14686ada9358027db0ce69a2d70d8304 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
@@ -15,6 +15,8 @@
 
 #include <cassert>
 
+#include "MaterialLib/MPL/Medium.h"
+#include "MaterialLib/MPL/Utils/FormEigenTensor.h"
 #include "MaterialLib/SolidModels/SelectSolidConstitutiveRelation.h"
 #include "MathLib/KelvinVector.h"
 #include "NumLib/Function/Interpolation.h"
@@ -157,8 +159,10 @@ void RichardsMechanicsLocalAssembler<
             displacement_size + pressure_size>>(
         local_rhs_data, displacement_size + pressure_size);
 
-    auto const material_id =
-        _process_data.flow_material->getMaterialID(_element.getID());
+    auto const& medium = _process_data.media_map->getMedium(_element.getID());
+    auto const& liquid_phase = medium->phase("AqueousLiquid");
+    auto const& solid_phase = medium->phase("Solid");
+    MPL::VariableArray variables;
 
     ParameterLib::SpatialPosition x_position;
     x_position.setElementID(_element.getID());
@@ -356,8 +360,10 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
             pressure_size, displacement_size>::Zero(pressure_size,
                                                     displacement_size);
 
-    auto const material_id =
-        _process_data.flow_material->getMaterialID(_element.getID());
+    auto const& medium = _process_data.media_map->getMedium(_element.getID());
+    auto const& liquid_phase = medium->phase("AqueousLiquid");
+    auto const& solid_phase = medium->phase("Solid");
+    MPL::VariableArray variables;
 
     ParameterLib::SpatialPosition x_position;
     x_position.setElementID(_element.getID());
@@ -673,8 +679,10 @@ std::vector<double> const& RichardsMechanicsLocalAssembler<
             pressure_size> const>(local_x.data() + pressure_index,
                                   pressure_size);
 
-    auto const material_id =
-        _process_data.flow_material->getMaterialID(_element.getID());
+    auto const& medium = _process_data.media_map->getMedium(_element.getID());
+    auto const& liquid_phase = medium->phase("AqueousLiquid");
+    auto const& solid_phase = medium->phase("Solid");
+    MPL::VariableArray variables;
 
     unsigned const n_integration_points =
         _integration_method.getNumberOfPoints();
@@ -812,6 +820,8 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
         Eigen::Map<typename ShapeMatricesTypeDisplacement::template VectorType<
             displacement_size> const>(local_x.data() + displacement_offset,
                                       displacement_size);
+    auto const& medium = _process_data.media_map->getMedium(_element.getID());
+    MPL::VariableArray variables;
     ParameterLib::SpatialPosition x_position;
     x_position.setElementID(_element.getID());
 
@@ -853,8 +863,9 @@ void RichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
         Eigen::Map<typename ShapeMatricesTypePressure::template VectorType<
             pressure_size> const>(local_x.data() + pressure_index,
                                   pressure_size);
-    auto const material_id =
-        _process_data.flow_material->getMaterialID(_element.getID());
+
+    auto const& medium = _process_data.media_map->getMedium(_element.getID());
+    MPL::VariableArray variables;
 
     ParameterLib::SpatialPosition x_position;
     x_position.setElementID(_element.getID());
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
index cb3e26fb1f474e4f62d72ce29660147459faa635..fd9da55cb0ef74922f6aab0129612e935a0d01e4 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h
@@ -32,6 +32,8 @@ namespace ProcessLib
 {
 namespace RichardsMechanics
 {
+namespace MPL = MaterialPropertyLib;
+
 /// Used for the extrapolation of the integration point values. It is ordered
 /// (and stored) by integration points.
 template <typename ShapeMatrixType>
diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h b/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h
index bf8537f66aebe51584d430ead9646134bb9dc055..97dfcfb101e0d6503fb654863cb32e719284d27d 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcessData.h
@@ -17,7 +17,7 @@
 
 #include <Eigen/Dense>
 
-#include "ProcessLib/RichardsFlow/RichardsFlowMaterialProperties.h"
+#include "MaterialLib/MPL/MaterialSpatialDistributionMap.h"
 
 namespace MaterialLib
 {
@@ -34,10 +34,11 @@ namespace RichardsMechanics
 template <int DisplacementDim>
 struct RichardsMechanicsProcessData
 {
+
     MeshLib::PropertyVector<int> const* const material_ids = nullptr;
 
-    std::unique_ptr<ProcessLib::RichardsFlow::RichardsFlowMaterialProperties>
-        flow_material;
+    std::unique_ptr<MaterialPropertyLib::MaterialSpatialDistributionMap>
+        media_map = nullptr;
 
     /// The constitutive relation for the mechanical part.
     std::map<
@@ -49,20 +50,6 @@ struct RichardsMechanicsProcessData
     /// representation of length 4 or 6, ParameterLib::Parameter<double>.
     ParameterLib::Parameter<double> const* const initial_stress;
 
-    /// Fluid's bulk modulus. A scalar quantity,
-    /// ParameterLib::Parameter<double>.
-    ParameterLib::Parameter<double> const& fluid_bulk_modulus;
-    /// Biot coefficient. A scalar quantity, ParameterLib::Parameter<double>.
-    ParameterLib::Parameter<double> const& biot_coefficient;
-    /// Density of the solid. A scalar quantity,
-    /// ParameterLib::Parameter<double>.
-    ParameterLib::Parameter<double> const& solid_density;
-    /// Solid's bulk modulus. A scalar quantity,
-    /// ParameterLib::Parameter<double>.
-    ParameterLib::Parameter<double> const& solid_bulk_modulus;
-    /// Reference temperature for material properties. A scalar quantity,
-    /// ParameterLib::Parameter<double>.
-    ParameterLib::Parameter<double> const& temperature;
     /// Specific body forces applied to solid and fluid.
     /// It is usually used to apply gravitational forces.
     /// A vector of displacement dimension's length.