diff --git a/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.cpp b/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.cpp
index d6d1e4a8078c7ba6da0acd097730f7122b4d3065..046bb5eb2fb37cadf600c599be808eae2b7ad4a4 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.cpp
+++ b/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.cpp
@@ -39,7 +39,6 @@ void EffectiveVolumetricEnthalpyModel::dEval(
PorosityData const& porosity_data,
PorosityDerivativeData const& porosity_d_data,
SaturationData const& S_L_data,
- SaturationDataDeriv const& dS_L_dp_cap,
SolidDensityData const& solid_density_data,
SolidDensityDerivativeData const& solid_density_d_data,
SolidEnthalpyData const& solid_enthalpy_data,
@@ -51,11 +50,6 @@ void EffectiveVolumetricEnthalpyModel::dEval(
auto const phi_G = (1. - S_L_data.S_L) * porosity_data.phi;
double const phi_S = 1. - porosity_data.phi;
- // dphi_G_dp_GR = -ds_L_dp_GR * porosity_data.phi = 0;
- double const dphi_G_dp_cap = -dS_L_dp_cap() * porosity_data.phi;
- // dphi_L_dp_GR = ds_L_dp_GR * porosity_data.phi = 0;
- double const dphi_L_dp_cap = dS_L_dp_cap() * porosity_data.phi;
-
// From p_LR = p_GR - p_cap it follows for
// drho_LR/dp_GR = drho_LR/dp_LR * dp_LR/dp_GR
// = drho_LR/dp_LR * (dp_GR/dp_GR - dp_cap/dp_GR)
@@ -72,9 +66,11 @@ void EffectiveVolumetricEnthalpyModel::dEval(
fluid_enthalpy_data.h_L +*/
phi_L * drho_LR_dp_GR * fluid_enthalpy_data.h_L;
effective_volumetric_enthalpy_d_data.drho_h_eff_dp_cap =
- dphi_G_dp_cap * fluid_density_data.rho_GR * fluid_enthalpy_data.h_G +
+ porosity_d_data.dphi_L_dp_cap * fluid_density_data.rho_GR *
+ fluid_enthalpy_data.h_G +
/*phi_G * (drho_GR_dp_cap = 0) * fluid_enthalpy_data.h_G +*/
- dphi_L_dp_cap * fluid_density_data.rho_LR * fluid_enthalpy_data.h_L +
+ porosity_d_data.dphi_L_dp_cap * fluid_density_data.rho_LR *
+ fluid_enthalpy_data.h_L +
phi_L * drho_LR_dp_cap * fluid_enthalpy_data.h_L;
// TODO (naumov) Extend for temperature dependent porosities.
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.h b/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.h
index 7778ec5da52430734ea7bc6ab7825dbb3caefab0..4bff673ce9abeb66f8566e4bc089f066a584a2cb 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.h
+++ b/ProcessLib/TH2M/ConstitutiveRelations/Enthalpy.h
@@ -80,7 +80,6 @@ struct EffectiveVolumetricEnthalpyModel
PorosityData const& porosity_data,
PorosityDerivativeData const& porosity_d_data,
SaturationData const& S_L_data,
- SaturationDataDeriv const& dS_L_dp_cap,
SolidDensityData const& solid_density_data,
SolidDensityDerivativeData const& solid_density_d_data,
SolidEnthalpyData const& solid_enthalpy_data,
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.cpp b/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.cpp
index 36f740a8594314fbb226dc33218a7a76d418222a..2aa2b1d75526d4454a24d8d54939b946b9505b46 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.cpp
+++ b/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.cpp
@@ -39,7 +39,6 @@ void InternalEnergyModel::dEval(
PorosityData const& porosity_data,
PorosityDerivativeData const& porosity_d_data,
SaturationData const& S_L_data,
- SaturationDataDeriv const& dS_L_dp_cap,
SolidDensityData const& solid_density_data,
SolidDensityDerivativeData const& solid_density_d_data,
SolidEnthalpyData const& solid_enthalpy_data,
@@ -61,11 +60,6 @@ void InternalEnergyModel::dEval(
phi_S * solid_density_data.rho_SR * solid_heat_capacity_data() -
porosity_d_data.dphi_dT * solid_density_data.rho_SR * u_S;
- // dphi_G_dp_GR = -ds_L_dp_GR * porosity_data.phi = 0;
- double const dphi_G_dp_cap = -dS_L_dp_cap() * porosity_data.phi;
- // dphi_L_dp_GR = ds_L_dp_GR * porosity_data.phi = 0;
- double const dphi_L_dp_cap = dS_L_dp_cap() * porosity_data.phi;
-
// From p_LR = p_GR - p_cap it follows for
// drho_LR/dp_GR = drho_LR/dp_LR * dp_LR/dp_GR
// = drho_LR/dp_LR * (dp_GR/dp_GR - dp_cap/dp_GR)
@@ -85,9 +79,11 @@ void InternalEnergyModel::dEval(
phi_L * fluid_density_data.rho_LR * phase_transition_data.du_L_dp_GR;
effective_volumetric_internal_energy_d_data.drho_u_eff_dp_cap =
- dphi_G_dp_cap * fluid_density_data.rho_GR * phase_transition_data.uG +
+ -porosity_d_data.dphi_L_dp_cap * fluid_density_data.rho_GR *
+ phase_transition_data.uG +
/*phi_G * (drho_GR_dp_cap = 0) * phase_transition_data.uG +*/
- dphi_L_dp_cap * fluid_density_data.rho_LR * phase_transition_data.uL +
+ porosity_d_data.dphi_L_dp_cap * fluid_density_data.rho_LR *
+ phase_transition_data.uL +
phi_L * drho_LR_dp_cap * phase_transition_data.uL +
phi_L * fluid_density_data.rho_LR * phase_transition_data.du_L_dp_cap;
}
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.h b/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.h
index fcfd77f9c2e516b0dfa19727ca823f94e6ef2df4..d74dfd7167708c10cbe4a9beb146b30e7cabf3ec 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.h
+++ b/ProcessLib/TH2M/ConstitutiveRelations/InternalEnergy.h
@@ -48,7 +48,6 @@ struct InternalEnergyModel
PorosityData const& porosity_data,
PorosityDerivativeData const& porosity_d_data,
SaturationData const& S_L_data,
- SaturationDataDeriv const& dS_L_dp_cap,
SolidDensityData const& solid_density_data,
SolidDensityDerivativeData const& solid_density_d_data,
SolidEnthalpyData const& solid_enthalpy_data,
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/Porosity.cpp b/ProcessLib/TH2M/ConstitutiveRelations/Porosity.cpp
index f43cef3cf6e84708117413bc7f66f34c584c919e..2688dd02698f3e4072a3a465d36a12badb4ee7e7 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/Porosity.cpp
+++ b/ProcessLib/TH2M/ConstitutiveRelations/Porosity.cpp
@@ -27,8 +27,9 @@ void PorosityModel::eval(SpaceTimeData const& x_t,
mpl_porosity.template value<double>(variables, x_t.x, x_t.t, x_t.dt);
}
-void PorosityModel::dEval(SpaceTimeData const& x_t,
- MediaData const& media_data,
+void PorosityModel::dEval(SpaceTimeData const& x_t, MediaData const& media_data,
+ PorosityData const& porosity_data,
+ SaturationDataDeriv const& dS_L_dp_cap,
PorosityDerivativeData& porosity_d_data) const
{
MaterialPropertyLib::VariableArray variables;
@@ -39,6 +40,8 @@ void PorosityModel::dEval(SpaceTimeData const& x_t,
porosity_d_data.dphi_dT = mpl_porosity.template dValue<double>(
variables, MaterialPropertyLib::Variable::temperature, x_t.x, x_t.t,
x_t.dt);
+
+ porosity_d_data.dphi_L_dp_cap = dS_L_dp_cap() * porosity_data.phi;
}
template <int DisplacementDim>
@@ -74,9 +77,11 @@ template <int DisplacementDim>
void PorosityModelNonConstantSolidPhaseVolumeFraction<DisplacementDim>::dEval(
SpaceTimeData const& x_t,
MediaData const& media_data,
+ PorosityData const& porosity_data,
+ SaturationDataDeriv const& dS_L_dp_cap,
BiotData const& biot,
- StrainData<DisplacementDim> const& strain_data,
SolidThermalExpansionData<DisplacementDim> const& s_therm_exp_data,
+ StrainData<DisplacementDim> const& strain_data,
PorosityDerivativeData& porosity_d_data) const
{
MaterialPropertyLib::VariableArray variables;
@@ -100,6 +105,8 @@ void PorosityModelNonConstantSolidPhaseVolumeFraction<DisplacementDim>::dEval(
dphi_0_dT *
(1. + s_therm_exp_data.thermal_volume_strain - biot() * div_u) -
(1. - phi_0) * s_therm_exp_data.beta_T_SR;
+
+ porosity_d_data.dphi_L_dp_cap = dS_L_dp_cap() * porosity_data.phi;
}
template struct PorosityModelNonConstantSolidPhaseVolumeFraction<2>;
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/Porosity.h b/ProcessLib/TH2M/ConstitutiveRelations/Porosity.h
index aecbd3062625267c02562c1402ac23d1a5be7420..49884ed45eddd84dfe40308f5353bcda18538421 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/Porosity.h
+++ b/ProcessLib/TH2M/ConstitutiveRelations/Porosity.h
@@ -13,6 +13,7 @@
#include "Biot.h"
#include "ProcessLib/ConstitutiveRelations/StrainData.h"
#include "ProcessLib/Reflection/ReflectionData.h"
+#include "Saturation.h"
#include "SolidThermalExpansion.h"
namespace ProcessLib::TH2M
@@ -22,6 +23,10 @@ namespace ConstitutiveRelations
struct PorosityDerivativeData
{
double dphi_dT = nan;
+ // dphi_G_dp_GR = -ds_L_dp_GR * phi = 0;
+ // dphi_L_dp_GR = ds_L_dp_GR * phi = 0;
+ double dphi_L_dp_cap = nan;
+ // dphi_G_dp_cap = -dphi_L_dp_cap
};
struct PorosityData
@@ -45,6 +50,8 @@ struct PorosityModel
void dEval(SpaceTimeData const& x_t,
MediaData const& media_data,
+ PorosityData const& porosity_data,
+ SaturationDataDeriv const& dS_L_dp_cap,
PorosityDerivativeData& porosity_d_data) const;
};
@@ -62,9 +69,11 @@ struct PorosityModelNonConstantSolidPhaseVolumeFraction
void dEval(
SpaceTimeData const& x_t,
MediaData const& media_data,
+ PorosityData const& porosity_data,
+ SaturationDataDeriv const& dS_L_dp_cap,
BiotData const& biot,
- StrainData<DisplacementDim> const& strain_data,
SolidThermalExpansionData<DisplacementDim> const& s_therm_exp_data,
+ StrainData<DisplacementDim> const& strain_data,
PorosityDerivativeData& porosity_d_data) const;
};
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.cpp b/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.cpp
index 60002a217fb177f1579364fb4c5899463ffd3327..358ee797cc7005035ee0fdeeb6a1df98cca2f271 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.cpp
+++ b/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.cpp
@@ -40,7 +40,7 @@ void ThermalConductivityModel<DisplacementDim>::dEval(
SpaceTimeData const& x_t, MediaData const& media_data,
TemperatureData const& T_data, PorosityData const& porosity_data,
PorosityDerivativeData const& porosity_d_data,
- SaturationData const& S_L_data, SaturationDataDeriv const& dS_L_dp_cap,
+ SaturationData const& S_L_data,
ThermalConductivityDerivativeData<DisplacementDim>&
thermal_conductivity_d_data) const
{
@@ -98,13 +98,9 @@ void ThermalConductivityModel<DisplacementDim>::dEval(
x_t.t, x_t.dt))
: MPL::formEigenTensor<DisplacementDim>(0.);
- // dphi_G_dp_GR = -ds_L_dp_GR * phi = 0;
- double const dphi_G_dp_cap = -dS_L_dp_cap() * porosity_data.phi;
- // dphi_L_dp_GR = ds_L_dp_GR * phi = 0;
- double const dphi_L_dp_cap = dS_L_dp_cap() * porosity_data.phi;
-
thermal_conductivity_d_data.dlambda_dp_cap =
- dphi_G_dp_cap * lambdaGR + dphi_L_dp_cap * lambdaLR;
+ -porosity_d_data.dphi_L_dp_cap * lambdaGR +
+ porosity_d_data.dphi_L_dp_cap * lambdaLR;
double const phi_L = S_L_data.S_L * porosity_data.phi;
double const phi_G = (1. - S_L_data.S_L) * porosity_data.phi;
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h b/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h
index eefb94ef7077a4373d0675001770278538a1d5a7..35720adac3fee11f0e318e8111151868261f8cc8 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h
+++ b/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h
@@ -47,7 +47,6 @@ struct ThermalConductivityModel
TemperatureData const& T_data, PorosityData const& porosity_data,
PorosityDerivativeData const& porosity_d_data,
SaturationData const& S_L_data,
- SaturationDataDeriv const& dS_L_dp_cap,
ThermalConductivityDerivativeData<DisplacementDim>&
thermal_conductivity_d_data) const;
};
diff --git a/ProcessLib/TH2M/TH2MFEM-impl.h b/ProcessLib/TH2M/TH2MFEM-impl.h
index 5098124c862573ed2b3deac4f479957330a6d876..3e0beb96ac31e7d6465a447dd3307337e4d5166c 100644
--- a/ProcessLib/TH2M/TH2MFEM-impl.h
+++ b/ProcessLib/TH2M/TH2MFEM-impl.h
@@ -556,16 +556,16 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
ip_cv.phase_transition_data,
ip_dd.advection_d_data);
- models.porosity_model.dEval({pos, t, dt}, media_data,
+ models.porosity_model.dEval(
+ {pos, t, dt}, media_data, ip_out.porosity_data, ip_dd.dS_L_dp_cap,
#ifdef NON_CONSTANT_SOLID_PHASE_VOLUME_FRACTION
- ip_cv.biot_data, ip_out.eps_data,
- ip_cv.s_therm_exp_data,
+ ip_cv.biot_data, ip_out.eps_data, ip_cv.s_therm_exp_data,
#endif // NON_CONSTANT_SOLID_PHASE_VOLUME_FRACTION
- ip_dd.porosity_d_data);
+ ip_dd.porosity_d_data);
models.thermal_conductivity_model.dEval(
{pos, t, dt}, media_data, T_data, ip_out.porosity_data,
- ip_dd.porosity_d_data, current_state.S_L_data, ip_dd.dS_L_dp_cap,
+ ip_dd.porosity_d_data, current_state.S_L_data,
ip_dd.thermal_conductivity_d_data);
models.solid_density_model.dEval({pos, t, dt}, media_data, T_data,
@@ -581,7 +581,6 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
ip_out.porosity_data,
ip_dd.porosity_d_data,
current_state.S_L_data,
- ip_dd.dS_L_dp_cap,
ip_out.solid_density_data,
ip_dd.solid_density_d_data,
ip_out.solid_enthalpy_data,
@@ -595,7 +594,6 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
ip_out.porosity_data,
ip_dd.porosity_d_data,
current_state.S_L_data,
- ip_dd.dS_L_dp_cap,
ip_out.solid_density_data,
ip_dd.solid_density_d_data,
ip_out.solid_enthalpy_data,