OGS
ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, int DisplacementDim>
class ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >

Local assembler of ThermoMechanics process.

In this assembler, the solid density can be expressed as an exponential function of temperature, if a temperature dependent solid density is assumed. The theory behind this is given below.

During the temperature variation, the solid mass balance is

\[ \rho_s^{n+1} V^{n+1} = \rho_s^n V^n, \]

with \(\rho_s\) the density, \(V\), the volume, and \(n\) the index of the time step. Under pure thermo-mechanics condition, the volume change along with the temperature change is given by

\[ V^{n+1} = V^{n} + \alpha_T dT V^{n} = (1+\alpha_T dT)V^{n}, \]

where \(\alpha_T\) is the volumetric thermal expansivity of the solid. This gives

\[ \rho_s^{n+1} + \alpha_T dT \rho_s^{n+1} = \rho_s^n. \]

Therefore, we obtain the differential expression of the temperature dependent solid density as

\[ \frac{d \rho_s}{d T} = -\alpha_T \rho_s. \]

The solution of the above ODE, i.e the density expression, is given by

\[ \rho_s = {\rho_0} \mathrm{exp} (- \alpha_T (T-T0)), \]

with reference density \(\rho_0\) at a reference temperature of \(T_0\).

An MPL property with the type of Exponential (see MaterialPropertyLib::Exponential) can be used for the temperature dependent solid property.

Definition at line 127 of file ThermoMechanicsFEM.h.

#include <ThermoMechanicsFEM.h>

Inheritance diagram for ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >:
[legend]

Public Types

using ShapeMatricesType
 
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
 
using GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType
 
using BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>
 
using Invariants = MathLib::KelvinVector::Invariants<KelvinVectorSize>
 
using NodalForceVectorType = typename BMatricesType::NodalForceVectorType
 
using RhsVector
 
using JacobianMatrix
 
using IpData
 

Public Member Functions

 ThermoMechanicsLocalAssembler (ThermoMechanicsLocalAssembler const &)=delete
 
 ThermoMechanicsLocalAssembler (ThermoMechanicsLocalAssembler &&)=delete
 
 ThermoMechanicsLocalAssembler (MeshLib::Element const &e, std::size_t const, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, ThermoMechanicsProcessData< DisplacementDim > &process_data)
 
std::size_t setIPDataInitialConditions (std::string_view const name, double const *values, int const integration_order) override
 Returns number of read integration points.
 
void assemble (double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
 
void assembleWithJacobianForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
 
void assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
 
void initializeConcrete () override
 
void postTimestepConcrete (Eigen::VectorXd const &, Eigen::VectorXd const &, double const, double const, int const) override
 
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix (const unsigned integration_point) const override
 Provides the shape matrix at the given integration point.
 
- Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
virtual ~LocalAssemblerInterface ()=default
 
virtual void setInitialConditions (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, int const process_id)
 
virtual void initialize (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
 
virtual void preAssemble (double const, double const, std::vector< double > const &)
 
virtual void assembleForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
 
virtual void computeSecondaryVariable (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
 
virtual void preTimestep (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, double const delta_t)
 
virtual void postTimestep (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
 
void postNonLinearSolver (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
 
virtual Eigen::Vector3d getFlux (MathLib::Point3d const &, double const, std::vector< double > const &) const
 
virtual Eigen::Vector3d getFlux (MathLib::Point3d const &, double const, std::vector< std::vector< double > > const &) const
 Fits to staggered scheme.
 
- Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default
 

Static Public Attributes

static int const KelvinVectorSize
 

Private Member Functions

void assembleWithJacobianForDeformationEquations (const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
 
void assembleWithJacobianForHeatConductionEquations (const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
 
std::size_t setSigma (double const *values)
 
std::vector< double > getSigma () const override
 
std::vector< double > const & getIntPtSigma (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
 
std::size_t setEpsilon (double const *values)
 
std::vector< double > getEpsilon () const override
 
std::vector< double > const & getIntPtEpsilon (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilonMechanical (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
 
std::size_t setEpsilonMechanical (double const *values)
 
std::vector< double > getEpsilonMechanical () const override
 
unsigned getNumberOfIntegrationPoints () const override
 
int getMaterialID () const override
 
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt (unsigned integration_point) const override
 

Private Attributes

ThermoMechanicsProcessData< DisplacementDim > & _process_data
 
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
 
NumLib::GenericIntegrationMethod const & _integration_method
 
MeshLib::Element const & _element
 
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
 
bool const _is_axially_symmetric
 

Static Private Attributes

static const int temperature_index = 0
 
static const int temperature_size = ShapeFunction::NPOINTS
 
static const int displacement_index = ShapeFunction::NPOINTS
 
static const int displacement_size
 

Member Typedef Documentation

◆ BMatricesType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 138 of file ThermoMechanicsFEM.h.

◆ GlobalDimMatrixType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType

Definition at line 137 of file ThermoMechanicsFEM.h.

◆ Invariants

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::Invariants = MathLib::KelvinVector::Invariants<KelvinVectorSize>

Definition at line 142 of file ThermoMechanicsFEM.h.

◆ IpData

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::IpData
Initial value:
IntegrationPointData<BMatricesType, ShapeMatricesType, DisplacementDim>

Definition at line 150 of file ThermoMechanicsFEM.h.

◆ JacobianMatrix

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::JacobianMatrix
Initial value:
typename ShapeMatricesType::template MatrixType<
ShapeFunction::NPOINTS + ShapeFunction::NPOINTS * DisplacementDim,
ShapeFunction::NPOINTS + ShapeFunction::NPOINTS * DisplacementDim>

Definition at line 147 of file ThermoMechanicsFEM.h.

◆ NodalForceVectorType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::NodalForceVectorType = typename BMatricesType::NodalForceVectorType

Definition at line 144 of file ThermoMechanicsFEM.h.

◆ RhsVector

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::RhsVector
Initial value:
typename ShapeMatricesType::template VectorType<
ShapeFunction::NPOINTS + ShapeFunction::NPOINTS * DisplacementDim>

Definition at line 145 of file ThermoMechanicsFEM.h.

◆ ShapeMatrices

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices

Definition at line 136 of file ThermoMechanicsFEM.h.

◆ ShapeMatricesType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::ShapeMatricesType

Constructor & Destructor Documentation

◆ ThermoMechanicsLocalAssembler() [1/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::ThermoMechanicsLocalAssembler ( ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim > const & )
delete

◆ ThermoMechanicsLocalAssembler() [2/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::ThermoMechanicsLocalAssembler ( ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim > && )
delete

◆ ThermoMechanicsLocalAssembler() [3/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::ThermoMechanicsLocalAssembler ( MeshLib::Element const & e,
std::size_t const ,
NumLib::GenericIntegrationMethod const & integration_method,
bool const is_axially_symmetric,
ThermoMechanicsProcessData< DisplacementDim > & process_data )

Definition at line 27 of file ThermoMechanicsFEM-impl.h.

34 : _process_data(process_data),
35 _integration_method(integration_method),
36 _element(e),
37 _is_axially_symmetric(is_axially_symmetric)
38{
39 unsigned const n_integration_points =
41
42 _ip_data.reserve(n_integration_points);
43 _secondary_data.N.resize(n_integration_points);
44
45 auto const shape_matrices =
47 DisplacementDim>(e, is_axially_symmetric,
49
51 _process_data.solid_materials, _process_data.material_ids, e.getID());
52
53 for (unsigned ip = 0; ip < n_integration_points; ip++)
54 {
55 _ip_data.emplace_back(solid_material);
56 auto& ip_data = _ip_data[ip];
57 ip_data.integration_weight =
59 shape_matrices[ip].integralMeasure * shape_matrices[ip].detJ;
60
61 static const int kelvin_vector_size =
63 // Initialize current time step values
64 ip_data.sigma.setZero(kelvin_vector_size);
65 ip_data.eps.setZero(kelvin_vector_size);
66
67 // Previous time step values are not initialized and are set later.
68 ip_data.sigma_prev.resize(kelvin_vector_size);
69 ip_data.eps_prev.resize(kelvin_vector_size);
70
71 ip_data.eps_m.setZero(kelvin_vector_size);
72 ip_data.eps_m_prev.setZero(kelvin_vector_size);
74 x_position.setElementID(_element.getID());
75 ip_data.N = shape_matrices[ip].N;
76 ip_data.dNdx = shape_matrices[ip].dNdx;
77
78 _secondary_data.N[ip] = shape_matrices[ip].N;
79 }
80}
double getWeight() const
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:89
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
void setElementID(std::size_t element_id)
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
ThermoMechanicsProcessData< DisplacementDim > & _process_data
NumLib::GenericIntegrationMethod const & _integration_method
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
auto & selectSolidConstitutiveRelation(SolidMaterialsMap const &constitutive_relations, MeshLib::PropertyVector< int > const *const material_ids, std::size_t const element_id)
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
std::vector< typename ShapeMatricesType::ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > initShapeMatrices(MeshLib::Element const &e, bool const is_axially_symmetric, IntegrationMethod const &integration_method)
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N

References ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_element, ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_integration_method, ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_process_data, ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_secondary_data, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::initShapeMatrices(), MathLib::KelvinVector::kelvin_vector_dimensions(), ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N, MaterialLib::Solids::selectSolidConstitutiveRelation(), and ParameterLib::SpatialPosition::setElementID().

Member Function Documentation

◆ assemble()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::assemble ( double const ,
double const ,
std::vector< double > const & ,
std::vector< double > const & ,
std::vector< double > & ,
std::vector< double > & ,
std::vector< double > &  )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 170 of file ThermoMechanicsFEM.h.

176 {
177 OGS_FATAL(
178 "ThermoMechanicsLocalAssembler: assembly without jacobian is not "
179 "implemented.");
180 }
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

◆ assembleWithJacobian()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobian ( double const t,
double const dt,
std::vector< double > const & local_x,
std::vector< double > const & local_x_prev,
std::vector< double > & local_M_data,
std::vector< double > & local_K_data,
std::vector< double > & local_rhs_data,
std::vector< double > & local_Jac_data )
overridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 115 of file ThermoMechanicsFEM-impl.h.

123{
124 auto const local_matrix_size = local_x.size();
125 assert(local_matrix_size == temperature_size + displacement_size);
126
127 auto T = Eigen::Map<typename ShapeMatricesType::template VectorType<
128 temperature_size> const>(local_x.data() + temperature_index,
130
131 auto u = Eigen::Map<typename ShapeMatricesType::template VectorType<
132 displacement_size> const>(local_x.data() + displacement_index,
134 bool const is_u_non_zero = u.norm() > 0.0;
135
136 auto T_prev = Eigen::Map<typename ShapeMatricesType::template VectorType<
137 temperature_size> const>(local_x_prev.data() + temperature_index,
139
140 auto local_Jac = MathLib::createZeroedMatrix<JacobianMatrix>(
141 local_Jac_data, local_matrix_size, local_matrix_size);
142
143 auto local_rhs = MathLib::createZeroedVector<RhsVector>(local_rhs_data,
144 local_matrix_size);
145
146 typename ShapeMatricesType::template MatrixType<displacement_size,
148 KuT;
150
153
156
157 unsigned const n_integration_points =
159
160 MPL::VariableArray variables;
161 MPL::VariableArray variables_prev;
163 x_position.setElementID(_element.getID());
164
165 auto const& medium = _process_data.media_map.getMedium(_element.getID());
166 auto const& solid_phase = medium->phase("Solid");
167
168 for (unsigned ip = 0; ip < n_integration_points; ip++)
169 {
170 x_position.setIntegrationPoint(ip);
171 auto const& w = _ip_data[ip].integration_weight;
172 auto const& N = _ip_data[ip].N;
173 auto const& dNdx = _ip_data[ip].dNdx;
174
175 auto const x_coord =
176 NumLib::interpolateXCoordinate<ShapeFunction, ShapeMatricesType>(
177 _element, N);
178 auto const& B =
179 LinearBMatrix::computeBMatrix<DisplacementDim,
180 ShapeFunction::NPOINTS,
182 dNdx, N, x_coord, _is_axially_symmetric);
183
184 auto& sigma = _ip_data[ip].sigma;
185 auto const& sigma_prev = _ip_data[ip].sigma_prev;
186
187 auto& eps = _ip_data[ip].eps;
188 auto const& eps_prev = _ip_data[ip].eps_prev;
189
190 auto& eps_m = _ip_data[ip].eps_m;
191 auto const& eps_m_prev = _ip_data[ip].eps_m_prev;
192
193 auto& state = _ip_data[ip].material_state_variables;
194
195 const double T_ip = N.dot(T); // T at integration point
196 double const T_prev_ip = N.dot(T_prev);
197
198 // Consider also anisotropic thermal expansion.
199 auto const solid_linear_thermal_expansivity_vector =
200 MPL::formKelvinVector<DisplacementDim>(
201 solid_phase
202 .property(
204 .value(variables, x_position, t, dt));
205
207 dthermal_strain =
208 solid_linear_thermal_expansivity_vector * (T_ip - T_prev_ip);
209
210 //
211 // displacement equation, displacement part
212 //
213 // For the restart computation, the displacement may not be
214 // reloaded but the initial strains are always available. For such case,
215 // the following computation is skipped.
216 if (is_u_non_zero)
217 {
218 eps.noalias() = B * u;
219 }
220
221 eps_m.noalias() = eps_m_prev + eps - eps_prev - dthermal_strain;
222
223 variables_prev.stress
225 sigma_prev);
226 variables_prev.mechanical_strain
228 eps_m_prev);
229 variables_prev.temperature = T_ip;
230 variables.mechanical_strain
232 eps_m);
233 variables.temperature = T_ip;
234
235 auto&& solution = _ip_data[ip].solid_material.integrateStress(
236 variables_prev, variables, t, x_position, dt, *state);
237
238 if (!solution)
239 {
240 OGS_FATAL("Computation of local constitutive relation failed.");
241 }
242
244 std::tie(sigma, state, C) = std::move(*solution);
245
246 local_Jac
247 .template block<displacement_size, displacement_size>(
249 .noalias() += B.transpose() * C * B * w;
250
251 typename ShapeMatricesType::template MatrixType<DisplacementDim,
253 N_u = ShapeMatricesType::template MatrixType<
254 DisplacementDim, displacement_size>::Zero(DisplacementDim,
256
257 for (int i = 0; i < DisplacementDim; ++i)
258 {
259 N_u.template block<1, displacement_size / DisplacementDim>(
260 i, i * displacement_size / DisplacementDim)
261 .noalias() = N;
262 }
263
264 auto const rho_s =
265 solid_phase.property(MPL::PropertyType::density)
266 .template value<double>(variables, x_position, t, dt);
267
268 auto const& b = _process_data.specific_body_force;
269 local_rhs.template block<displacement_size, 1>(displacement_index, 0)
270 .noalias() -=
271 (B.transpose() * sigma - N_u.transpose() * rho_s * b) * w;
272
273 //
274 // displacement equation, temperature part
275 // The computation of KuT can be ignored.
276 auto const alpha_T_tensor =
278 solid_linear_thermal_expansivity_vector);
279 KuT.noalias() += B.transpose() * (C * alpha_T_tensor) * N * w;
280
281 if (_ip_data[ip].solid_material.getConstitutiveModel() ==
283 {
284 auto const s = Invariants::deviatoric_projection * sigma;
285 double const norm_s = Invariants::FrobeniusNorm(s);
286 const double creep_coefficient =
287 _ip_data[ip].solid_material.getTemperatureRelatedCoefficient(
288 t, dt, x_position, T_ip, norm_s);
289 KuT.noalias() += creep_coefficient * B.transpose() * s * N * w;
290 }
291
292 //
293 // temperature equation, temperature part;
294 //
295 auto const lambda =
296 solid_phase
297 .property(
299 .value(variables, x_position, t, dt);
300
302 MaterialPropertyLib::formEigenTensor<DisplacementDim>(lambda);
303
304 KTT.noalias() += dNdx.transpose() * thermal_conductivity * dNdx * w;
305
306 auto const c =
307 solid_phase
308 .property(
310 .template value<double>(variables, x_position, t, dt);
311 DTT.noalias() += N.transpose() * rho_s * c * N * w;
312 }
313
314 // temperature equation, temperature part
315 local_Jac
316 .template block<temperature_size, temperature_size>(temperature_index,
318 .noalias() += KTT + DTT / dt;
319
320 // displacement equation, temperature part
321 local_Jac
322 .template block<displacement_size, temperature_size>(displacement_index,
324 .noalias() -= KuT;
325
326 local_rhs.template block<temperature_size, 1>(temperature_index, 0)
327 .noalias() -= KTT * T + DTT * (T - T_prev) / dt;
328}
std::variant< std::monostate, Eigen::Matrix< double, 4, 1 >, Eigen::Matrix< double, 6, 1 > > mechanical_strain
std::variant< std::monostate, Eigen::Matrix< double, 4, 1 >, Eigen::Matrix< double, 6, 1 > > stress
void setIntegrationPoint(unsigned integration_point)
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
BMatrixType computeBMatrix(DNDX_Type const &dNdx, N_Type const &N, const double radius, const bool is_axially_symmetric)
Fills a B-matrix based on given shape function dN/dx values.
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType
static double FrobeniusNorm(Eigen::Matrix< double, KelvinVectorSize, 1 > const &deviatoric_v)
Get the norm of the deviatoric stress.
static Eigen::Matrix< double, KelvinVectorSize, KelvinVectorSize > const deviatoric_projection

References ProcessLib::LinearBMatrix::computeBMatrix(), MaterialLib::Solids::CreepBGRa, MathLib::KelvinVector::kelvinVectorToSymmetricTensor(), MaterialPropertyLib::VariableArray::mechanical_strain, OGS_FATAL, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::specific_heat_capacity, MaterialPropertyLib::VariableArray::stress, MaterialPropertyLib::VariableArray::temperature, MaterialPropertyLib::thermal_conductivity, and MaterialPropertyLib::thermal_expansivity.

◆ assembleWithJacobianForDeformationEquations()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobianForDeformationEquations ( const double t,
double const dt,
Eigen::VectorXd const & local_x,
Eigen::VectorXd const & local_x_prev,
std::vector< double > & local_b_data,
std::vector< double > & local_Jac_data )
private

Assemble local matrices and vectors arise from the linearized discretized weak form of the residual of the momentum balance equation,

\[ \nabla (\sigma - \mathbf{D} \alpha_T (T-T_0) \mathrm{I}) = f \]

where \( \sigma\) is the effective stress tensor, \(\mathbf{D}\) is the tangential operator, \(T\) is the temperature, \(T_0\) is the initial temperature, \(\alpha_T\) is the linear thermal expansion, \(\mathrm{I}\) is the identity tensor, and \(f\) is the body force.

Parameters
tTime
dtTime increment
local_xNodal values of \(x\) of all processes of an element.
local_x_prevNodal values of \(x_{prev}\) of all processes of an element.
local_b_dataRight hand side vector of an element.
local_Jac_dataElement Jacobian matrix for the Newton-Raphson method.

Definition at line 353 of file ThermoMechanicsFEM-impl.h.

358{
359 auto const local_T =
360 local_x.template segment<temperature_size>(temperature_index);
361
362 auto const local_T_prev =
363 local_x_prev.template segment<temperature_size>(temperature_index);
364
365 auto const local_u =
366 local_x.template segment<displacement_size>(displacement_index);
367 bool const is_u_non_zero = local_u.norm() > 0.0;
368
369 auto local_Jac = MathLib::createZeroedMatrix<
370 typename ShapeMatricesType::template MatrixType<displacement_size,
372 local_Jac_data, displacement_size, displacement_size);
373
374 auto local_rhs = MathLib::createZeroedVector<
375 typename ShapeMatricesType::template VectorType<displacement_size>>(
376 local_b_data, displacement_size);
377
378 unsigned const n_integration_points =
380
381 MPL::VariableArray variables;
382 MPL::VariableArray variables_prev;
384 x_position.setElementID(_element.getID());
385 auto const& medium = _process_data.media_map.getMedium(_element.getID());
386 auto const& solid_phase = medium->phase("Solid");
387
388 for (unsigned ip = 0; ip < n_integration_points; ip++)
389 {
390 x_position.setIntegrationPoint(ip);
391 auto const& w = _ip_data[ip].integration_weight;
392 auto const& N = _ip_data[ip].N;
393 auto const& dNdx = _ip_data[ip].dNdx;
394
395 auto const x_coord =
396 NumLib::interpolateXCoordinate<ShapeFunction, ShapeMatricesType>(
397 _element, N);
398 auto const& B =
399 LinearBMatrix::computeBMatrix<DisplacementDim,
400 ShapeFunction::NPOINTS,
402 dNdx, N, x_coord, _is_axially_symmetric);
403
404 auto& sigma = _ip_data[ip].sigma;
405 auto const& sigma_prev = _ip_data[ip].sigma_prev;
406
407 auto& eps = _ip_data[ip].eps;
408 auto const& eps_prev = _ip_data[ip].eps_prev;
409
410 auto& eps_m = _ip_data[ip].eps_m;
411 auto const& eps_m_prev = _ip_data[ip].eps_m_prev;
412
413 auto& state = _ip_data[ip].material_state_variables;
414
415 const double T_ip = N.dot(local_T); // T at integration point
416 variables.temperature = T_ip;
417 double const dT_ip = T_ip - N.dot(local_T_prev);
418
419 //
420 // displacement equation, displacement part
421 //
422 // For the restart computation, the displacement may not be
423 // reloaded but the initial strains are always available. For such case,
424 // the following computation is skipped.
425 if (is_u_non_zero)
426 {
427 eps.noalias() = B * local_u;
428 }
429
430 // Consider also anisotropic thermal expansion.
431 auto const solid_linear_thermal_expansivity_vector =
432 MPL::formKelvinVector<DisplacementDim>(
433 solid_phase
434 .property(
436 .value(variables, x_position, t, dt));
437
439 dthermal_strain = solid_linear_thermal_expansivity_vector * dT_ip;
440
441 eps_m.noalias() = eps_m_prev + eps - eps_prev - dthermal_strain;
442
443 variables_prev.stress
445 sigma_prev);
446 variables_prev.mechanical_strain
448 eps_m_prev);
449 variables_prev.temperature = T_ip;
450 variables.mechanical_strain
452 eps_m);
453
454 auto&& solution = _ip_data[ip].solid_material.integrateStress(
455 variables_prev, variables, t, x_position, dt, *state);
456
457 if (!solution)
458 {
459 OGS_FATAL("Computation of local constitutive relation failed.");
460 }
461
463 std::tie(sigma, state, C) = std::move(*solution);
464
465 local_Jac.noalias() += B.transpose() * C * B * w;
466
467 typename ShapeMatricesType::template MatrixType<DisplacementDim,
469 N_u = ShapeMatricesType::template MatrixType<
470 DisplacementDim, displacement_size>::Zero(DisplacementDim,
472
473 for (int i = 0; i < DisplacementDim; ++i)
474 {
475 N_u.template block<1, displacement_size / DisplacementDim>(
476 i, i * displacement_size / DisplacementDim)
477 .noalias() = N;
478 }
479
480 auto const rho_s =
481 solid_phase.property(MPL::PropertyType::density)
482 .template value<double>(variables, x_position, t, dt);
483
484 auto const& b = _process_data.specific_body_force;
485 local_rhs.noalias() -=
486 (B.transpose() * sigma - N_u.transpose() * rho_s * b) * w;
487 }
488}
Eigen::Map< Vector > createZeroedVector(std::vector< double > &data, Eigen::VectorXd::Index size)
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)

References ProcessLib::LinearBMatrix::computeBMatrix(), MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), MaterialPropertyLib::VariableArray::mechanical_strain, OGS_FATAL, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::VariableArray::stress, MaterialPropertyLib::VariableArray::temperature, and MaterialPropertyLib::thermal_expansivity.

◆ assembleWithJacobianForHeatConductionEquations()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobianForHeatConductionEquations ( const double t,
double const dt,
Eigen::VectorXd const & local_x,
Eigen::VectorXd const & local_x_prev,
std::vector< double > & local_b_data,
std::vector< double > & local_Jac_data )
private

Assemble local matrices and vectors arise from the linearized discretized weak form of the residual of the energy balance equation,

\[ \rho c_p \cdot{T} - \nabla (\mathbf{K} (\nabla T) = Q_T \]

where \( rho\) is the solid density, \( c_p\) is the specific heat capacity, \( \mathbf{K} \) is the thermal conductivity, and \( Q_T\) is the source/sink term.

Parameters
tTime
dtTime increment
local_xNodal values of \(x\) of all processes of an element.
local_x_prevNodal values of \(x_{prev}\) of all processes of an element.
local_b_dataRight hand side vector of an element.
local_Jac_dataElement Jacobian matrix for the Newton-Raphson method.

Definition at line 491 of file ThermoMechanicsFEM-impl.h.

496{
497 auto const local_T =
498 local_x.template segment<temperature_size>(temperature_index);
499
500 auto const local_T_prev =
501 local_x_prev.template segment<temperature_size>(temperature_index);
502
503 auto local_Jac = MathLib::createZeroedMatrix<
504 typename ShapeMatricesType::template MatrixType<temperature_size,
506 local_Jac_data, temperature_size, temperature_size);
507
508 auto local_rhs = MathLib::createZeroedVector<
509 typename ShapeMatricesType::template VectorType<temperature_size>>(
510 local_b_data, temperature_size);
511
513 mass.setZero(temperature_size, temperature_size);
514
515 typename ShapeMatricesType::NodalMatrixType laplace;
516 laplace.setZero(temperature_size, temperature_size);
517
518 unsigned const n_integration_points =
520
522 x_position.setElementID(_element.getID());
523 auto const& medium = _process_data.media_map.getMedium(_element.getID());
524 auto const& solid_phase = medium->phase("Solid");
525 MPL::VariableArray variables;
526
527 for (unsigned ip = 0; ip < n_integration_points; ip++)
528 {
529 x_position.setIntegrationPoint(ip);
530 auto const& w = _ip_data[ip].integration_weight;
531 auto const& N = _ip_data[ip].N;
532 auto const& dNdx = _ip_data[ip].dNdx;
533
534 const double T_ip = N.dot(local_T); // T at integration point
535 variables.temperature = T_ip;
536
537 auto const rho_s =
538 solid_phase.property(MPL::PropertyType::density)
539 .template value<double>(variables, x_position, t, dt);
540 auto const c_p =
541 solid_phase.property(MPL::PropertyType::specific_heat_capacity)
542 .template value<double>(variables, x_position, t, dt);
543
544 mass.noalias() += N.transpose() * rho_s * c_p * N * w;
545
546 auto const lambda =
547 solid_phase
548 .property(
550 .value(variables, x_position, t, dt);
551
553 MaterialPropertyLib::formEigenTensor<DisplacementDim>(lambda);
554
555 laplace.noalias() += dNdx.transpose() * thermal_conductivity * dNdx * w;
556 }
557 local_Jac.noalias() += laplace + mass / dt;
558
559 local_rhs.noalias() -=
560 laplace * local_T + mass * (local_T - local_T_prev) / dt;
561}

References MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::VariableArray::temperature, and MaterialPropertyLib::thermal_conductivity.

◆ assembleWithJacobianForStaggeredScheme()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobianForStaggeredScheme ( double const t,
double const dt,
Eigen::VectorXd const & local_x,
Eigen::VectorXd const & local_x_prev,
int const process_id,
std::vector< double > & local_M_data,
std::vector< double > & local_K_data,
std::vector< double > & local_b_data,
std::vector< double > & local_Jac_data )
overridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 331 of file ThermoMechanicsFEM-impl.h.

338{
339 // For the equations with pressure
340 if (process_id == _process_data.heat_conduction_process_id)
341 {
343 t, dt, local_x, local_x_prev, local_b_data, local_Jac_data);
344 return;
345 }
346
347 // For the equations with deformation
348 assembleWithJacobianForDeformationEquations(t, dt, local_x, local_x_prev,
349 local_b_data, local_Jac_data);
350}
void assembleWithJacobianForHeatConductionEquations(const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
void assembleWithJacobianForDeformationEquations(const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)

◆ getEpsilon()

template<typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getEpsilon ( ) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 607 of file ThermoMechanicsFEM-impl.h.

608{
609 constexpr int kelvin_vector_size =
611
612 return transposeInPlace<kelvin_vector_size>(
613 [this](std::vector<double>& values)
614 { return getIntPtEpsilon(0, {}, {}, values); });
615}
std::vector< double > const & getIntPtEpsilon(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override

References MathLib::KelvinVector::kelvin_vector_dimensions().

◆ getEpsilonMechanical()

template<typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getEpsilonMechanical ( ) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 651 of file ThermoMechanicsFEM-impl.h.

652{
653 constexpr int kelvin_vector_size =
655
656 return transposeInPlace<kelvin_vector_size>(
657 [this](std::vector<double>& values)
658 { return getIntPtEpsilonMechanical(0, {}, {}, values); });
659}
std::vector< double > const & getIntPtEpsilonMechanical(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override

References MathLib::KelvinVector::kelvin_vector_dimensions().

◆ getIntPtEpsilon()

template<typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getIntPtEpsilon ( const double t,
std::vector< GlobalVector * > const & x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & dof_table,
std::vector< double > & cache ) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 619 of file ThermoMechanicsFEM-impl.h.

624{
625 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
626 _ip_data, &IpData::eps, cache);
627}

◆ getIntPtEpsilonMechanical()

template<typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getIntPtEpsilonMechanical ( const double t,
std::vector< GlobalVector * > const & x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & dof_table,
std::vector< double > & cache ) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 631 of file ThermoMechanicsFEM-impl.h.

637{
638 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
639 _ip_data, &IpData::eps_m, cache);
640}

◆ getIntPtSigma()

template<typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getIntPtSigma ( const double t,
std::vector< GlobalVector * > const & x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & dof_table,
std::vector< double > & cache ) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 586 of file ThermoMechanicsFEM-impl.h.

591{
592 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
593 _ip_data, &IpData::sigma, cache);
594}

◆ getMaterialID()

template<typename ShapeFunction , int DisplacementDim>
int ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getMaterialID ( ) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 670 of file ThermoMechanicsFEM-impl.h.

671{
672 return _process_data.material_ids == nullptr
673 ? 0
674 : (*_process_data.material_ids)[_element.getID()];
675}

◆ getMaterialStateVariablesAt()

template<typename ShapeFunction , int DisplacementDim>
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getMaterialStateVariablesAt ( unsigned integration_point) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 680 of file ThermoMechanicsFEM-impl.h.

682{
683 return *_ip_data[integration_point].material_state_variables;
684}

◆ getNumberOfIntegrationPoints()

template<typename ShapeFunction , int DisplacementDim>
unsigned ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getNumberOfIntegrationPoints ( ) const
overrideprivatevirtual

◆ getShapeMatrix()

template<typename ShapeFunction , int DisplacementDim>
Eigen::Map< const Eigen::RowVectorXd > ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getShapeMatrix ( const unsigned integration_point) const
inlineoverridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 245 of file ThermoMechanicsFEM.h.

247 {
248 auto const& N = _secondary_data.N[integration_point];
249
250 // assumes N is stored contiguously in memory
251 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
252 }

References ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_secondary_data, and ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N.

◆ getSigma()

template<typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::getSigma ( ) const
overrideprivatevirtual

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 574 of file ThermoMechanicsFEM-impl.h.

575{
576 constexpr int kelvin_vector_size =
578
579 return transposeInPlace<kelvin_vector_size>(
580 [this](std::vector<double>& values)
581 { return getIntPtSigma(0, {}, {}, values); });
582}
std::vector< double > const & getIntPtSigma(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override

References MathLib::KelvinVector::kelvin_vector_dimensions().

◆ initializeConcrete()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::initializeConcrete ( )
inlineoverridevirtual

Set initial stress from parameter.

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 197 of file ThermoMechanicsFEM.h.

198 {
199 unsigned const n_integration_points =
201 for (unsigned ip = 0; ip < n_integration_points; ip++)
202 {
203 auto& ip_data = _ip_data[ip];
204
205 ParameterLib::SpatialPosition const x_position{
206 std::nullopt, _element.getID(), ip,
208 NumLib::interpolateCoordinates<ShapeFunction,
210 _element, ip_data.N))};
211
213 if (_process_data.initial_stress != nullptr)
214 {
215 ip_data.sigma =
217 DisplacementDim>((*_process_data.initial_stress)(
218 std::numeric_limits<
219 double>::quiet_NaN() /* time independent */,
220 x_position));
221 }
222
223 double const t = 0; // TODO (naumov) pass t from top
224 ip_data.solid_material.initializeInternalStateVariables(
225 t, x_position, *ip_data.material_state_variables);
226
227 ip_data.pushBackState();
228 }
229 }
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)

References ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_element, ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_integration_method, ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_process_data, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), NumLib::interpolateCoordinates(), and MathLib::KelvinVector::symmetricTensorToKelvinVector().

◆ postTimestepConcrete()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::postTimestepConcrete ( Eigen::VectorXd const & ,
Eigen::VectorXd const & ,
double const ,
double const ,
int const  )
inlineoverridevirtual

◆ setEpsilon()

template<typename ShapeFunction , int DisplacementDim>
std::size_t ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::setEpsilon ( double const * values)
private

Definition at line 598 of file ThermoMechanicsFEM-impl.h.

600{
601 return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
602 values, _ip_data, &IpData::eps);
603}

◆ setEpsilonMechanical()

template<typename ShapeFunction , int DisplacementDim>
std::size_t ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::setEpsilonMechanical ( double const * values)
private

Definition at line 644 of file ThermoMechanicsFEM-impl.h.

645{
646 return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
647 values, _ip_data, &IpData::eps_m);
648}

◆ setIPDataInitialConditions()

template<typename ShapeFunction , int DisplacementDim>
std::size_t ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::setIPDataInitialConditions ( std::string_view const name,
double const * values,
int const integration_order )
overridevirtual

Returns number of read integration points.

Implements ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssemblerInterface< DisplacementDim >.

Definition at line 83 of file ThermoMechanicsFEM-impl.h.

87{
88 if (integration_order !=
89 static_cast<int>(_integration_method.getIntegrationOrder()))
90 {
92 "Setting integration point initial conditions; The integration "
93 "order of the local assembler for element {:d} is different from "
94 "the integration order in the initial condition.",
95 _element.getID());
96 }
97
98 if (name == "sigma")
99 {
100 return setSigma(values);
101 }
102 if (name == "epsilon")
103 {
104 return setEpsilon(values);
105 }
106 if (name == "epsilon_m")
107 {
108 return setEpsilonMechanical(values);
109 }
110
111 return 0;
112}

References OGS_FATAL.

◆ setSigma()

template<typename ShapeFunction , int DisplacementDim>
std::size_t ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::setSigma ( double const * values)
private

Definition at line 565 of file ThermoMechanicsFEM-impl.h.

567{
568 return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
569 values, _ip_data, &IpData::sigma);
570}

Member Data Documentation

◆ _element

◆ _integration_method

◆ _ip_data

◆ _is_axially_symmetric

template<typename ShapeFunction , int DisplacementDim>
bool const ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::_is_axially_symmetric
private

Definition at line 356 of file ThermoMechanicsFEM.h.

◆ _process_data

◆ _secondary_data

◆ displacement_index

template<typename ShapeFunction , int DisplacementDim>
const int ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::displacement_index = ShapeFunction::NPOINTS
staticprivate

Definition at line 360 of file ThermoMechanicsFEM.h.

◆ displacement_size

template<typename ShapeFunction , int DisplacementDim>
const int ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::displacement_size
staticprivate
Initial value:
=
ShapeFunction::NPOINTS * DisplacementDim

Definition at line 361 of file ThermoMechanicsFEM.h.

◆ KelvinVectorSize

template<typename ShapeFunction , int DisplacementDim>
int const ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::KelvinVectorSize
static
Initial value:

Definition at line 140 of file ThermoMechanicsFEM.h.

◆ temperature_index

template<typename ShapeFunction , int DisplacementDim>
const int ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::temperature_index = 0
staticprivate

Definition at line 358 of file ThermoMechanicsFEM.h.

◆ temperature_size

template<typename ShapeFunction , int DisplacementDim>
const int ProcessLib::ThermoMechanics::ThermoMechanicsLocalAssembler< ShapeFunction, DisplacementDim >::temperature_size = ShapeFunction::NPOINTS
staticprivate

Definition at line 359 of file ThermoMechanicsFEM.h.


The documentation for this class was generated from the following files: