From ce23c725de5bf3ad33b51a3a6b92d171cb85a7d6 Mon Sep 17 00:00:00 2001
From: Dmitri Naumov <dmitri.naumov@ufz.de>
Date: Wed, 3 Apr 2024 19:31:10 +0200
Subject: [PATCH] [PL/TH2M] Collect thermal conductivities data
---
...fectiveThermalConductivityPorosityMixing.h | 6 ++--
.../ConstitutiveRelations/ConstitutiveData.h | 5 ++-
.../ThermalConductivity.h | 31 ++++++++++++++++++
ProcessLib/TH2M/IntegrationPointData.h | 1 -
ProcessLib/TH2M/TH2MFEM-impl.h | 32 +++++++++++--------
5 files changed, 55 insertions(+), 20 deletions(-)
create mode 100644 ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h
diff --git a/MaterialLib/MPL/Properties/EffectiveThermalConductivityPorosityMixing.h b/MaterialLib/MPL/Properties/EffectiveThermalConductivityPorosityMixing.h
index f5ef7fdbd78..01f80c95962 100644
--- a/MaterialLib/MPL/Properties/EffectiveThermalConductivityPorosityMixing.h
+++ b/MaterialLib/MPL/Properties/EffectiveThermalConductivityPorosityMixing.h
@@ -39,9 +39,9 @@ public:
double const dt) const override;
PropertyDataType dValue(VariableArray const& variable_array,
Variable const variable,
- ParameterLib::SpatialPosition const& /*pos*/,
- double const /*t*/,
- double const /*dt*/) const override;
+ ParameterLib::SpatialPosition const& pos,
+ double const t,
+ double const dt) const override;
private:
ParameterLib::CoordinateSystem const* const local_coordinate_system_;
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h b/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h
index f79ba86ea6b..a543613ec98 100644
--- a/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h
+++ b/ProcessLib/TH2M/ConstitutiveRelations/ConstitutiveData.h
@@ -32,6 +32,7 @@
#include "SolidMechanics.h"
#include "SolidThermalExpansion.h"
#include "Swelling.h"
+#include "ThermalConductivity.h"
#include "TotalStress.h"
#include "VapourPartialPressure.h"
#include "Viscosity.h"
@@ -140,14 +141,12 @@ struct ConstitutiveTempData
PorosityDerivativeData porosity_d_data;
SolidDensityDerivativeData solid_density_d_data;
SolidHeatCapacityData solid_heat_capacity_data;
+ ThermalConductivityData<DisplacementDim> thermal_conductivity_data;
using DisplacementDimVector = Eigen::Matrix<double, DisplacementDim, 1>;
using DisplacementDimMatrix =
Eigen::Matrix<double, DisplacementDim, DisplacementDim>;
- DisplacementDimMatrix dlambda_dp_GR;
- DisplacementDimMatrix dlambda_dp_cap;
- DisplacementDimMatrix dlambda_dT;
DisplacementDimVector drho_GR_h_w_eff_dp_GR_Npart;
DisplacementDimMatrix drho_GR_h_w_eff_dp_GR_gradNpart;
DisplacementDimVector drho_LR_h_w_eff_dp_cap_Npart;
diff --git a/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h b/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h
new file mode 100644
index 00000000000..2ec711957f0
--- /dev/null
+++ b/ProcessLib/TH2M/ConstitutiveRelations/ThermalConductivity.h
@@ -0,0 +1,31 @@
+/**
+ * \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"
+
+namespace ProcessLib::TH2M
+{
+namespace ConstitutiveRelations
+{
+
+template <int DisplacementDim>
+struct ThermalConductivityData
+{
+ GlobalDimMatrix<DisplacementDim> lambda;
+ // Currently unused, but there is a comment in TH2MFEM-impl.h referring to
+ // this matrix
+ // GlobalDimMatrix<DisplacementDim> dlambda_dp_GR;
+ GlobalDimMatrix<DisplacementDim> dlambda_dp_cap;
+ GlobalDimMatrix<DisplacementDim> dlambda_dT;
+};
+
+} // namespace ConstitutiveRelations
+} // namespace ProcessLib::TH2M
diff --git a/ProcessLib/TH2M/IntegrationPointData.h b/ProcessLib/TH2M/IntegrationPointData.h
index ec9de5fdb91..e33c29163ed 100644
--- a/ProcessLib/TH2M/IntegrationPointData.h
+++ b/ProcessLib/TH2M/IntegrationPointData.h
@@ -48,7 +48,6 @@ struct IntegrationPointData final
double rho_u_eff = std::numeric_limits<double>::quiet_NaN();
double rho_u_eff_prev = std::numeric_limits<double>::quiet_NaN();
- GlobalDimMatrixType lambda;
GlobalDimVectorType d_CG;
GlobalDimVectorType d_WG;
GlobalDimVectorType d_CL;
diff --git a/ProcessLib/TH2M/TH2MFEM-impl.h b/ProcessLib/TH2M/TH2MFEM-impl.h
index df952164fb7..380f15766ee 100644
--- a/ProcessLib/TH2M/TH2MFEM-impl.h
+++ b/ProcessLib/TH2M/TH2MFEM-impl.h
@@ -273,11 +273,12 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
double const phi_S = 1. - ip_out.porosity_data.phi;
// thermal conductivity
- ip_data.lambda = MaterialPropertyLib::formEigenTensor<DisplacementDim>(
- medium
- .property(
- MaterialPropertyLib::PropertyType::thermal_conductivity)
- .value(vars, pos, t, dt));
+ ip_cv.thermal_conductivity_data.lambda =
+ MaterialPropertyLib::formEigenTensor<DisplacementDim>(
+ medium
+ .property(
+ MaterialPropertyLib::PropertyType::thermal_conductivity)
+ .value(vars, pos, t, dt));
ip_data.h_S = ip_cv.solid_heat_capacity_data() * T;
auto const u_S = ip_data.h_S;
@@ -409,12 +410,12 @@ TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
.dValue(vars, MPL::Variable::temperature, pos, t, dt))
: MPL::formEigenTensor<DisplacementDim>(0.);
- ip_cv.dlambda_dp_cap =
+ ip_cv.thermal_conductivity_data.dlambda_dp_cap =
dphi_G_dp_cap * lambdaGR + dphi_L_dp_cap * lambdaLR;
- ip_cv.dlambda_dT = phi_G * dlambda_GR_dT + phi_L * dlambda_LR_dT +
- phi_S * dlambda_SR_dT -
- ip_cv.porosity_d_data.dphi_dT * lambdaSR;
+ ip_cv.thermal_conductivity_data.dlambda_dT =
+ phi_G * dlambda_GR_dT + phi_L * dlambda_LR_dT +
+ phi_S * dlambda_SR_dT - ip_cv.porosity_d_data.dphi_dT * lambdaSR;
// From p_LR = p_GR - p_cap it follows for
// drho_LR/dp_GR = drho_LR/dp_LR * dp_LR/dp_GR
@@ -1392,7 +1393,8 @@ void TH2MLocalAssembler<
MTu.noalias() += NTT * rho_h_eff * mT * Bu * w;
- KTT.noalias() += gradNTT * ip.lambda * gradNT * w;
+ KTT.noalias() +=
+ gradNTT * ip_cv.thermal_conductivity_data.lambda * gradNT * w;
fT.noalias() -= NTT * rho_u_eff_dot * w;
@@ -2053,7 +2055,8 @@ void TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
temperature_index, temperature_index)
.noalias() += NTT * ip_cv.drho_h_eff_dT * div_u_dot * NT * w;
- KTT.noalias() += gradNTT * ip.lambda * gradNT * w;
+ KTT.noalias() +=
+ gradNTT * ip_cv.thermal_conductivity_data.lambda * gradNT * w;
// d KTT/dp_GR * T
// TODO (naumov) always zero if lambda_xR have no derivatives wrt. p_GR.
@@ -2073,13 +2076,16 @@ void TH2MLocalAssembler<ShapeFunctionDisplacement, ShapeFunctionPressure,
local_Jac
.template block<temperature_size, W_size>(temperature_index,
W_index)
- .noalias() += gradNTT * ip_cv.dlambda_dp_cap * gradT * Np * w;
+ .noalias() += gradNTT *
+ ip_cv.thermal_conductivity_data.dlambda_dp_cap *
+ gradT * Np * w;
// d KTT/dT * T
local_Jac
.template block<temperature_size, temperature_size>(
temperature_index, temperature_index)
- .noalias() += gradNTT * ip_cv.dlambda_dT * gradT * NT * w;
+ .noalias() += gradNTT * ip_cv.thermal_conductivity_data.dlambda_dT *
+ gradT * NT * w;
// fT_1
fT.noalias() -= NTT * rho_u_eff_dot * w;
--
GitLab