diff --git a/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h b/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h
index 0965f8009c6c145b45c989ca9f5f97f00d7c1fed..f79ba86ea6ba9f48acf71b614a879cfdaecd3394 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h
+++ b/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h
@@ -28,6 +28,7 @@
 #include "Saturation.h"
 #include "SolidCompressibility.h"
 #include "SolidDensity.h"
+#include "SolidHeatCapacity.h"
 #include "SolidMechanics.h"
 #include "SolidThermalExpansion.h"
 #include "Swelling.h"
@@ -138,6 +139,7 @@ struct ConstitutiveTempData
     PhaseTransitionData phase_transition_data;
     PorosityDerivativeData porosity_d_data;
     SolidDensityDerivativeData solid_density_d_data;
+    SolidHeatCapacityData solid_heat_capacity_data;
 
     using DisplacementDimVector = Eigen::Matrix<double, DisplacementDim, 1>;
     using DisplacementDimMatrix =
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveModels.h b/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveModels.h
index b9526cdf6d34197e8049e328e2930dc1ec84772c..50fedd23598c0de90b6eb19efd0e3d9c24c153d4 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveModels.h
+++ b/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveModels.h
@@ -20,6 +20,7 @@
 #include "Saturation.h"
 #include "SolidCompressibility.h"
 #include "SolidDensity.h"
+#include "SolidHeatCapacity.h"
 #include "SolidMechanics.h"
 #include "SolidThermalExpansion.h"
 #include "Swelling.h"
@@ -69,6 +70,7 @@ struct ConstitutiveModels
     PorosityModel porosity_model;
     SolidDensityModel solid_density_model;
 #endif  // NON_CONSTANT_SOLID_PHASE_VOLUME_FRACTION
+    SolidHeatCapacityModel solid_heat_capacity_model;
 };
 }  // namespace ConstitutiveRelations
 }  // namespace ProcessLib::TH2M
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/SolidHeatCapacity.cpp b/ProcessLib/TH2M/ConstitutiveRelations/SolidHeatCapacity.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..7e1fcef50cd4f92202ad42ce1d7b93a14b4d08c3
--- /dev/null
+++ b/ProcessLib/TH2M/ConstitutiveRelations/SolidHeatCapacity.cpp
@@ -0,0 +1,36 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2024, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ */
+
+#include "SolidHeatCapacity.h"
+
+namespace ProcessLib::TH2M
+{
+namespace ConstitutiveRelations
+{
+
+void SolidHeatCapacityModel::eval(
+    SpaceTimeData const& x_t,
+    MediaData const& media_data,
+    TemperatureData const& T_data,
+    SolidHeatCapacityData& solid_heat_capacity) const
+{
+    namespace MPL = MaterialPropertyLib;
+
+    MPL::VariableArray variables;
+    variables.temperature = T_data.T;
+
+    auto const& mpl_cpS =
+        media_data.solid[MPL::PropertyType::specific_heat_capacity];
+
+    *solid_heat_capacity =
+        mpl_cpS.template value<double>(variables, x_t.x, x_t.t, x_t.dt);
+}
+
+}  // namespace ConstitutiveRelations
+}  // namespace ProcessLib::TH2M
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/SolidHeatCapacity.h b/ProcessLib/TH2M/ConstitutiveRelations/SolidHeatCapacity.h
new file mode 100644
index 0000000000000000000000000000000000000000..4f3b4fae0ad1903fac570986ca8e23fa34e746dc
--- /dev/null
+++ b/ProcessLib/TH2M/ConstitutiveRelations/SolidHeatCapacity.h
@@ -0,0 +1,32 @@
+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2024, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ */
+
+#pragma once
+
+#include "Base.h"
+#include "BaseLib/StrongType.h"
+
+namespace ProcessLib::TH2M
+{
+namespace ConstitutiveRelations
+{
+
+using SolidHeatCapacityData =
+    BaseLib::StrongType<double, struct SolidHeatCapacityTag>;
+
+struct SolidHeatCapacityModel
+{
+    void eval(SpaceTimeData const& x_t,
+              MediaData const& media_data,
+              TemperatureData const& T_data,
+              SolidHeatCapacityData& solid_heat_capacity) const;
+};
+
+}  // namespace ConstitutiveRelations
+}  // namespace ProcessLib::TH2M
diff --git a/ProcessLib/TH2M/TH2MFEM-impl.h b/ProcessLib/TH2M/TH2MFEM-impl.h
index 9203a3d5326e621bb29938ed33525e0f66ee59a6..df952164fb7695b12ca7af805ff2af798edac0fc 100644
--- a/ProcessLib/TH2M/TH2MFEM-impl.h
+++ b/ProcessLib/TH2M/TH2MFEM-impl.h
@@ -246,6 +246,9 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
 #endif  // NON_CONSTANT_SOLID_PHASE_VOLUME_FRACTION
             ip_out.solid_density_data, ip_cv.solid_density_d_data);
 
+        models.solid_heat_capacity_model.eval({pos, t, dt}, media_data, T_data,
+                                              ip_cv.solid_heat_capacity_data);
+
         MPL::VariableArray vars;
         MPL::VariableArray vars_prev;
         vars.temperature = T;
@@ -276,10 +279,7 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
                     MaterialPropertyLib::PropertyType::thermal_conductivity)
                 .value(vars, pos, t, dt));
 
-        auto const cpS =
-            solid_phase.property(MPL::PropertyType::specific_heat_capacity)
-                .template value<double>(vars, pos, t, dt);
-        ip_data.h_S = cpS * T;
+        ip_data.h_S = ip_cv.solid_heat_capacity_data() * T;
         auto const u_S = ip_data.h_S;
 
         ip_data.rho_u_eff = phi_G * ip_out.fluid_density_data.rho_GR * c.uG +
@@ -348,7 +348,8 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
             phi_L * drho_LR_dT * c.uL +
             phi_L * ip_out.fluid_density_data.rho_LR * c.du_L_dT +
             phi_S * ip_cv.solid_density_d_data.drho_SR_dT * u_S +
-            phi_S * ip_out.solid_density_data.rho_SR * cpS -
+            phi_S * ip_out.solid_density_data.rho_SR *
+                ip_cv.solid_heat_capacity_data() -
             ip_cv.porosity_d_data.dphi_dT * ip_out.solid_density_data.rho_SR *
                 u_S;
 
@@ -453,7 +454,8 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
             ip_cv.porosity_d_data.dphi_dT * ip_out.solid_density_data.rho_SR *
                 ip_data.h_S +
             phi_S * ip_cv.solid_density_d_data.drho_SR_dT * ip_data.h_S +
-            phi_S * ip_out.solid_density_data.rho_SR * cpS;
+            phi_S * ip_out.solid_density_data.rho_SR *
+                ip_cv.solid_heat_capacity_data();
 
         ip_cv.drho_u_eff_dp_GR =
             /*(dphi_G_dp_GR = 0) * ip_out.fluid_density_data.rho_GR * c.uG +*/