diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
index 04bd024fd98bbe6bc8306df8caa1d8d675c5a080..df26fa1c646638233a24be2564a998749562bb2d 100644
--- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
+++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h
@@ -372,10 +372,6 @@ void RichardsMechanicsLocalAssembler<
         auto const alpha =
             solid_phase.property(MPL::PropertyType::biot_coefficient)
                 .template value<double>(variables, x_position, t, dt);
-        auto const rho_SR =
-            solid_phase.property(MPL::PropertyType::density)
-                .template value<double>(variables, x_position, t, dt);
-
         auto const C_el = _ip_data[ip].computeElasticTangentStiffness(
             t, x_position, dt, temperature);
 
@@ -546,6 +542,13 @@ void RichardsMechanicsLocalAssembler<
         _ip_data[ip].updateConstitutiveRelation(variables, t, x_position, dt,
                                                 temperature);
 
+        // p_SR
+        variables[static_cast<int>(MPL::Variable::solid_grain_pressure)] =
+            p_FR - sigma_eff.dot(identity2) / (3 * (1 - phi));
+        auto const rho_SR =
+            solid_phase.property(MPL::PropertyType::density)
+                .template value<double>(variables, x_position, t, dt);
+
         //
         // displacement equation, displacement part
         //