OGS
ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits > Class Template Reference

Detailed Description

template<typename ShapeFunctionDisplacement, typename ShapeFunction, int DisplacementDim, typename ConstitutiveTraits>
class ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >

Definition at line 37 of file ThermoRichardsMechanicsFEM.h.

#include <ThermoRichardsMechanicsFEM.h>

Inheritance diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >:
[legend]
Collaboration diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >:
[legend]

Classes

class  LocalMatrices
 

Public Types

using ShapeMatricesTypeDisplacement = ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim >
 
using ShapeMatricesType = ShapeMatrixPolicyType< ShapeFunction, DisplacementDim >
 
using GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType
 
using GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType
 
using BMatricesType = BMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim >
 
using KelvinVectorType = typename BMatricesType::KelvinVectorType
 
using IpData = IntegrationPointData< ShapeMatricesTypeDisplacement, ShapeMatricesType, DisplacementDim, ShapeFunctionDisplacement::NPOINTS >
 
using Invariants = MathLib::KelvinVector::Invariants< KelvinVectorSize >
 
using SymmetricTensor = Eigen::Matrix< double, KelvinVectorSize, 1 >
 

Public Member Functions

 ThermoRichardsMechanicsLocalAssembler (ThermoRichardsMechanicsLocalAssembler const &)=delete
 
 ThermoRichardsMechanicsLocalAssembler (ThermoRichardsMechanicsLocalAssembler &&)=delete
 
 ThermoRichardsMechanicsLocalAssembler (MeshLib::Element const &e, std::size_t const, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > &process_data)
 
void setInitialConditionsConcrete (std::vector< double > const &local_x, double const t, bool const use_monolithic_scheme, int const process_id) 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 > &, std::vector< double > &, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
 
void initializeConcrete () override
 
void computeSecondaryVariableConcrete (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev) 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::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >
 LocalAssemblerInterface (MeshLib::Element const &e, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > &process_data)
 
std::size_t setIPDataInitialConditions (std::string_view name, double const *values, int const integration_order)
 
std::vector< double > getMaterialStateVariableInternalState (std::function< std::span< double >(typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables &)> const &get_values_span, int const &n_components) const
 
unsigned getNumberOfIntegrationPoints () const
 
int getMaterialID () const
 
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt (unsigned integration_point) const
 
void postTimestepConcrete (Eigen::VectorXd const &, Eigen::VectorXd const &, double const, double const, bool const, int const) override
 
- Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
virtual ~LocalAssemblerInterface ()=default
 
void setInitialConditions (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, bool const use_monolithic_scheme, 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 assemble (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_b_data)
 
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 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_b_data, std::vector< double > &local_Jac_data)
 
virtual 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)
 
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, bool const use_monolithic_scheme, 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, bool const use_monolithic_scheme, 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 Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix (const unsigned integration_point) const =0
 Provides the shape matrix at the given integration point.
 
virtual ~ExtrapolatableElement ()=default
 

Static Public Attributes

static int const KelvinVectorSize
 
static constexpr auto & N_u_op
 

Private Member Functions

void assembleWithJacobianSingleIP (double const t, double const dt, ParameterLib::SpatialPosition const &x_position, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, IpData const &ip_data, typename ConstitutiveTraits::ConstitutiveSetting &CS, MaterialPropertyLib::Medium &medium, LocalMatrices &out, typename ConstitutiveTraits::StatefulData &current_state, typename ConstitutiveTraits::StatefulDataPrev const &prev_state, MaterialStateData< DisplacementDim > &mat_state, typename ConstitutiveTraits::OutputData &output_data) const
 
void addToLocalMatrixData (double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, LocalMatrices const &loc_mat, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) const
 
void massLumping (LocalMatrices &loc_mat) const
 
void convertInitialStressType (unsigned const ip, double const t, ParameterLib::SpatialPosition const x_position, MaterialPropertyLib::Medium const &medium, MPL::VariableArray const &variables, double const p_at_ip)
 

Static Private Member Functions

static constexpr auto localDOF (std::vector< double > const &x)
 
static auto block_uu (auto &mat)
 
static auto block_up (auto &mat)
 
static auto block_uT (auto &mat)
 
static auto block_pu (auto &mat)
 
static auto block_pp (auto &mat)
 
static auto block_pT (auto &mat)
 
static auto block_Tp (auto &mat)
 
static auto block_TT (auto &mat)
 
static auto block_u (auto &vec)
 
static auto block_p (auto &vec)
 
static auto block_T (auto &vec)
 

Private Attributes

std::vector< IpDataip_data_
 

Static Private Attributes

static constexpr int temperature_index = 0
 
static constexpr int temperature_size = ShapeFunction::NPOINTS
 
static constexpr int pressure_index = temperature_size
 
static constexpr int pressure_size = ShapeFunction::NPOINTS
 
static constexpr int displacement_index = 2 * ShapeFunction::NPOINTS
 
static constexpr int displacement_size
 

Additional Inherited Members

- Static Public Member Functions inherited from ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >
static auto getReflectionDataForOutput ()
 
- Protected Attributes inherited from ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >
ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > & process_data_
 
std::vector< typename ConstitutiveTraits::StatefulData > current_states_
 
std::vector< typename ConstitutiveTraits::StatefulDataPrev > prev_states_
 
std::vector< MaterialStateData< DisplacementDim > > material_states_
 
NumLib::GenericIntegrationMethod const & integration_method_
 
MeshLib::Element const & element_
 
bool const is_axially_symmetric_
 
ConstitutiveTraits::SolidConstitutiveRelation const & solid_material_
 
std::vector< typename ConstitutiveTraits::OutputData > output_data_
 

Member Typedef Documentation

◆ BMatricesType

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::BMatricesType = BMatrixPolicyType<ShapeFunctionDisplacement, DisplacementDim>

Definition at line 59 of file ThermoRichardsMechanicsFEM.h.

◆ GlobalDimMatrixType

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType

Definition at line 56 of file ThermoRichardsMechanicsFEM.h.

◆ GlobalDimVectorType

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType

Definition at line 57 of file ThermoRichardsMechanicsFEM.h.

◆ Invariants

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::Invariants = MathLib::KelvinVector::Invariants<KelvinVectorSize>

Definition at line 69 of file ThermoRichardsMechanicsFEM.h.

◆ IpData

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::IpData = IntegrationPointData<ShapeMatricesTypeDisplacement, ShapeMatricesType, DisplacementDim, ShapeFunctionDisplacement::NPOINTS>

Definition at line 63 of file ThermoRichardsMechanicsFEM.h.

◆ KelvinVectorType

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::KelvinVectorType = typename BMatricesType::KelvinVectorType

Definition at line 61 of file ThermoRichardsMechanicsFEM.h.

◆ ShapeMatricesType

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ShapeMatricesType = ShapeMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 53 of file ThermoRichardsMechanicsFEM.h.

◆ ShapeMatricesTypeDisplacement

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ShapeMatricesTypeDisplacement = ShapeMatrixPolicyType<ShapeFunctionDisplacement, DisplacementDim>

Definition at line 49 of file ThermoRichardsMechanicsFEM.h.

◆ SymmetricTensor

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::SymmetricTensor = Eigen::Matrix<double, KelvinVectorSize, 1>

Definition at line 71 of file ThermoRichardsMechanicsFEM.h.

Constructor & Destructor Documentation

◆ ThermoRichardsMechanicsLocalAssembler() [1/3]

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ThermoRichardsMechanicsLocalAssembler ( ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits > const &  )
delete

◆ ThermoRichardsMechanicsLocalAssembler() [2/3]

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ThermoRichardsMechanicsLocalAssembler ( ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits > &&  )
delete

◆ ThermoRichardsMechanicsLocalAssembler() [3/3]

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ThermoRichardsMechanicsLocalAssembler ( MeshLib::Element const &  e,
std::size_t const  ,
NumLib::GenericIntegrationMethod const &  integration_method,
bool const  is_axially_symmetric,
ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > &  process_data 
)

Definition at line 41 of file ThermoRichardsMechanicsFEM-impl.h.

49 : LocalAssemblerInterface<DisplacementDim, ConstitutiveTraits>(
50 e, integration_method, is_axially_symmetric, process_data)
51{
52 unsigned const n_integration_points =
54
55 ip_data_.resize(n_integration_points);
56
57 auto const shape_matrices_u =
58 NumLib::initShapeMatrices<ShapeFunctionDisplacement,
60 DisplacementDim>(e, is_axially_symmetric,
61 integration_method);
62
63 auto const shape_matrices =
65 DisplacementDim>(e, is_axially_symmetric,
66 integration_method);
67
68 for (unsigned ip = 0; ip < n_integration_points; ip++)
69 {
70 auto& ip_data = ip_data_[ip];
71 auto const& sm_u = shape_matrices_u[ip];
72 ip_data_[ip].integration_weight =
73 integration_method.getWeightedPoint(ip).getWeight() *
74 sm_u.integralMeasure * sm_u.detJ;
75
76 ip_data.N_u = sm_u.N;
77 ip_data.dNdx_u = sm_u.dNdx;
78
79 // ip_data.N_p and ip_data.dNdx_p are used for both p and T variables
80 ip_data.N_p = shape_matrices[ip].N;
81 ip_data.dNdx_p = shape_matrices[ip].dNdx;
82 }
83}
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
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)

References NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::initShapeMatrices(), ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::integration_method_, and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ip_data_.

Member Function Documentation

◆ addToLocalMatrixData()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::addToLocalMatrixData ( double const  dt,
std::vector< double > const &  local_x,
std::vector< double > const &  local_x_prev,
LocalMatrices const &  loc_mat,
std::vector< double > &  local_rhs_data,
std::vector< double > &  local_Jac_data 
) const
private

Definition at line 294 of file ThermoRichardsMechanicsFEM-impl.h.

304{
305 constexpr auto local_matrix_dim =
307 assert(local_x.size() == local_matrix_dim);
308
309 auto local_Jac = MathLib::createZeroedMatrix<
310 typename ShapeMatricesTypeDisplacement::template MatrixType<
311 local_matrix_dim, local_matrix_dim>>(
312 local_Jac_data, local_matrix_dim, local_matrix_dim);
313
314 auto local_rhs =
315 MathLib::createZeroedVector<typename ShapeMatricesTypeDisplacement::
316 template VectorType<local_matrix_dim>>(
317 local_rhs_data, local_matrix_dim);
318
319 local_Jac.noalias() = loc_mat.Jac;
320 local_rhs.noalias() = -loc_mat.res;
321
322 //
323 // -- Jacobian
324 //
325 block_TT(local_Jac).noalias() += loc_mat.M_TT / dt + loc_mat.K_TT;
326 block_Tp(local_Jac).noalias() +=
327 loc_mat.K_Tp + loc_mat.dK_TT_dp + loc_mat.M_Tp / dt;
328
329 block_pT(local_Jac).noalias() += loc_mat.M_pT / dt + loc_mat.K_pT;
330 block_pp(local_Jac).noalias() +=
331 loc_mat.K_pp + loc_mat.storage_p_a_p / dt + loc_mat.storage_p_a_S_Jpp;
332 block_pu(local_Jac).noalias() = loc_mat.M_pu / dt;
333
334 //
335 // -- Residual
336 //
337 auto const [T, p_L, u] = localDOF(local_x);
338 auto const [T_prev, p_L_prev, u_prev] = localDOF(local_x_prev);
339
340 block_T(local_rhs).noalias() -= loc_mat.M_TT * (T - T_prev) / dt +
341 loc_mat.K_TT * T + loc_mat.K_Tp * p_L +
342 loc_mat.M_Tp * (p_L - p_L_prev) / dt;
343 block_p(local_rhs).noalias() -=
344 loc_mat.K_pp * p_L + loc_mat.K_pT * T +
345 (loc_mat.storage_p_a_p + loc_mat.storage_p_a_S) * (p_L - p_L_prev) /
346 dt +
347 loc_mat.M_pu * (u - u_prev) / dt + loc_mat.M_pT * (T - T_prev) / dt;
348}
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)
Definition: EigenMapTools.h:32

References MathLib::createZeroedMatrix(), and MathLib::createZeroedVector().

◆ assembleWithJacobian()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::assembleWithJacobian ( double const  t,
double const  dt,
std::vector< double > const &  local_x,
std::vector< double > const &  local_x_prev,
std::vector< double > &  ,
std::vector< double > &  ,
std::vector< double > &  local_rhs_data,
std::vector< double > &  local_Jac_data 
)
overridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 224 of file ThermoRichardsMechanicsFEM-impl.h.

232{
233 auto& medium =
234 *this->process_data_.media_map.getMedium(this->element_.getID());
235
236 LocalMatrices loc_mat;
237 loc_mat.setZero();
238 LocalMatrices loc_mat_current_ip;
239 loc_mat_current_ip.setZero(); // only to set the right matrix sizes
240
241 typename ConstitutiveTraits::ConstitutiveSetting constitutive_setting;
242
243 for (unsigned ip = 0; ip < this->integration_method_.getNumberOfPoints();
244 ++ip)
245 {
246 ParameterLib::SpatialPosition const x_position{
247 std::nullopt, this->element_.getID(), ip,
249 NumLib::interpolateCoordinates<ShapeFunctionDisplacement,
251 this->element_, ip_data_[ip].N_u))};
252
254 t, dt, x_position, //
255 local_x, local_x_prev, //
256 ip_data_[ip], constitutive_setting,
257 medium, //
258 loc_mat_current_ip, //
259 this->current_states_[ip], this->prev_states_[ip],
260 this->material_states_[ip], this->output_data_[ip]);
261 loc_mat += loc_mat_current_ip;
262 }
263
264 massLumping(loc_mat);
265
266 addToLocalMatrixData(dt, local_x, local_x_prev, loc_mat, local_rhs_data,
267 local_Jac_data);
268}
std::size_t getID() const
Returns the ID of the element.
Definition: Element.h:89
void addToLocalMatrixData(double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, LocalMatrices const &loc_mat, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) const
void assembleWithJacobianSingleIP(double const t, double const dt, ParameterLib::SpatialPosition const &x_position, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, IpData const &ip_data, typename ConstitutiveTraits::ConstitutiveSetting &CS, MaterialPropertyLib::Medium &medium, LocalMatrices &out, typename ConstitutiveTraits::StatefulData &current_state, typename ConstitutiveTraits::StatefulDataPrev const &prev_state, MaterialStateData< DisplacementDim > &mat_state, typename ConstitutiveTraits::OutputData &output_data) const
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
std::vector< typename ConstitutiveTraits::OutputData > output_data_
ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > & process_data_
std::vector< typename ConstitutiveTraits::StatefulDataPrev > prev_states_
std::vector< MaterialStateData< DisplacementDim > > material_states_
std::vector< typename ConstitutiveTraits::StatefulData > current_states_

References NumLib::interpolateCoordinates(), and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::LocalMatrices::setZero().

◆ assembleWithJacobianSingleIP()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::assembleWithJacobianSingleIP ( double const  t,
double const  dt,
ParameterLib::SpatialPosition const &  x_position,
std::vector< double > const &  local_x,
std::vector< double > const &  local_x_prev,
IpData const &  ip_data,
typename ConstitutiveTraits::ConstitutiveSetting &  CS,
MaterialPropertyLib::Medium medium,
LocalMatrices out,
typename ConstitutiveTraits::StatefulData &  current_state,
typename ConstitutiveTraits::StatefulDataPrev const &  prev_state,
MaterialStateData< DisplacementDim > &  mat_state,
typename ConstitutiveTraits::OutputData &  output_data 
) const
private

Definition at line 354 of file ThermoRichardsMechanicsFEM-impl.h.

372{
373 auto const& N_u = ip_data.N_u;
374 auto const& dNdx_u = ip_data.dNdx_u;
375
376 // N and dNdx are used for both p and T variables
377 auto const& N = ip_data.N_p;
378 auto const& dNdx = ip_data.dNdx_p;
379
380 auto const B =
381 LinearBMatrix::computeBMatrix<DisplacementDim,
382 ShapeFunctionDisplacement::NPOINTS,
384 dNdx_u, N_u, (*x_position.getCoordinates())[0],
386
387 auto const [T, p_L, u] = localDOF(local_x);
388 auto const [T_prev, p_L_prev, u_prev] = localDOF(local_x_prev);
389
390 GlobalDimVectorType const grad_T_ip = dNdx * T;
391
392 typename ConstitutiveTraits::ConstitutiveModels models(
393 this->process_data_, this->solid_material_);
394 typename ConstitutiveTraits::ConstitutiveTempData tmp;
395 typename ConstitutiveTraits::ConstitutiveData CD;
396
397 {
398 double const T_ip = N * T;
399 double const T_prev_ip = N * T_prev;
400
401 double const p_cap_ip = -N * p_L;
402 double const p_cap_prev_ip = -N * p_L_prev;
403 GlobalDimVectorType const grad_p_cap_ip = -dNdx * p_L;
404
405 KelvinVectorType eps = B * u;
406
407 CS.eval(models, t, dt, x_position, //
408 medium, //
409 {T_ip, T_prev_ip, grad_T_ip}, //
410 {p_cap_ip, p_cap_prev_ip, grad_p_cap_ip}, //
411 eps, current_state, prev_state, mat_state, tmp, output_data,
412 CD);
413 }
414
415 using NodalMatrix = typename ShapeMatricesType::NodalMatrixType;
416 NodalMatrix const NTN = N.transpose() * N;
417 NodalMatrix const dNTdN = dNdx.transpose() * dNdx;
418
419 // TODO is identity2.transpose() * B something like divergence?
420 auto const& identity2 = MathLib::KelvinVector::Invariants<
422 DisplacementDim)>::identity2;
423 typename ShapeMatricesTypeDisplacement::template MatrixType<
425 BTI2N = B.transpose() * identity2 * N;
426
427 /*
428 * Conventions:
429 *
430 * * use positive signs exclusively, any negative signs must be included in
431 * the coefficients coming from the constitutive setting
432 * * the used coefficients from the constitutive setting are named after the
433 * terms they appear in, e.g. K_TT_X_dNTdN means it is a scalar (X) that
434 * will be multiplied by dNdx.transpose() * dNdx. Placefolders for the
435 * coefficients are:
436 * * X -> scalar
437 * * V -> vector
438 * * K -> Kelvin vector
439 * * the Laplace terms have a special name, e.g., K_TT_Laplace
440 * * there shall be only one contribution to each of the LocalMatrices,
441 * assigned with = assignment; this point might be relaxed in the future
442 * * this method will overwrite the data in the passed LocalMatrices& out
443 * argument, not add to it
444 */
445
446 // residual, order T, p, u
447 block_p(out.res).noalias() =
448 dNdx.transpose() * std::get<EqPData<DisplacementDim>>(CD).rhs_p_dNT_V;
449 block_u(out.res).noalias() =
450 B.transpose() *
451 ProcessLib::Graph::get<TotalStressData<DisplacementDim>>(
452 CD, current_state)
453 .sigma_total -
454 static_cast<int>(this->process_data_.apply_body_force_for_deformation) *
455 N_u_op(N_u).transpose() *
456 std::get<GravityData<DisplacementDim>>(CD).volumetric_body_force;
457
458 // Storage matrices
459 out.storage_p_a_p.noalias() =
460 std::get<EqPData<DisplacementDim>>(CD).storage_p_a_p_X_NTN * NTN;
461 out.storage_p_a_S.noalias() =
462 std::get<TRMStorageData>(CD).storage_p_a_S_X_NTN * NTN;
463 out.storage_p_a_S_Jpp.noalias() =
464 std::get<TRMStorageData>(CD).storage_p_a_S_Jpp_X_NTN * NTN;
465
466 // M matrices, order T, p, u
467 out.M_TT.noalias() =
468 std::get<EqTData<DisplacementDim>>(CD).M_TT_X_NTN * NTN;
469 out.M_Tp.noalias() =
470 std::get<TRMVaporDiffusionData<DisplacementDim>>(CD).M_Tp_X_NTN * NTN;
471
472 out.M_pT.noalias() =
473 std::get<EqPData<DisplacementDim>>(CD).M_pT_X_NTN * NTN;
474 out.M_pu.noalias() =
475 std::get<EqPData<DisplacementDim>>(CD).M_pu_X_BTI2N * BTI2N.transpose();
476
477 // K matrices, order T, p, u
478 out.K_TT.noalias() =
479 dNdx.transpose() *
480 std::get<TRMHeatStorageAndFluxData<DisplacementDim>>(CD)
481 .K_TT_Laplace *
482 dNdx +
483 N.transpose() *
484 (std::get<EqTData<DisplacementDim>>(CD).K_TT_NT_V_dN.transpose() *
485 dNdx) +
486 std::get<TRMVaporDiffusionData<DisplacementDim>>(CD).K_TT_X_dNTdN *
487 dNTdN;
488
489 out.dK_TT_dp.noalias() =
490 N.transpose() *
491 (std::get<TRMHeatStorageAndFluxData<DisplacementDim>>(CD)
492 .K_Tp_NT_V_dN.transpose() *
493 dNdx) +
494 std::get<TRMHeatStorageAndFluxData<DisplacementDim>>(CD).K_Tp_X_NTN *
495 NTN;
496 out.K_Tp.noalias() =
497 dNdx.transpose() *
498 std::get<ThermoOsmosisData<DisplacementDim>>(CD).K_Tp_Laplace *
499 dNdx +
500 std::get<TRMVaporDiffusionData<DisplacementDim>>(CD).K_Tp_X_dNTdN *
501 dNTdN;
502
503 out.K_pp.noalias() =
504 dNdx.transpose() * std::get<EqPData<DisplacementDim>>(CD).K_pp_Laplace *
505 dNdx +
506 std::get<TRMVaporDiffusionData<DisplacementDim>>(CD).K_pp_X_dNTdN *
507 dNTdN;
508 out.K_pT.noalias() =
509 dNdx.transpose() *
510 std::get<ThermoOsmosisData<DisplacementDim>>(CD).K_pT_Laplace * dNdx;
511
512 // direct Jacobian contributions, order T, p, u
513 block_pT(out.Jac).noalias() =
514 std::get<TRMVaporDiffusionData<DisplacementDim>>(CD).J_pT_X_dNTdN *
515 dNTdN;
516 block_pp(out.Jac).noalias() =
517 std::get<TRMStorageData>(CD).J_pp_X_NTN * NTN +
518 std::get<EqPData<DisplacementDim>>(CD).J_pp_X_BTI2NT_u_dot_N *
519 BTI2N.transpose() * (u - u_prev) / dt *
520 N // TODO something with volumetric strain rate?
521 + dNdx.transpose() *
522 std::get<EqPData<DisplacementDim>>(CD).J_pp_dNT_V_N * N;
523
524 block_uT(out.Jac).noalias() =
525 B.transpose() *
526 std::get<SolidMechanicsDataStateless<DisplacementDim>>(CD).J_uT_BT_K_N *
527 N;
528 block_up(out.Jac).noalias() =
529 B.transpose() *
530 std::get<SolidMechanicsDataStateless<DisplacementDim>>(CD)
531 .J_up_BT_K_N *
532 N +
533 N_u_op(N_u).transpose() *
534 std::get<GravityData<DisplacementDim>>(CD).J_up_HT_V_N * N;
535 block_uu(out.Jac).noalias() =
536 B.transpose() *
537 std::get<SolidMechanicsDataStateless<DisplacementDim>>(CD)
538 .stiffness_tensor *
539 B;
540
541 out *= ip_data.integration_weight;
542}
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
Definition: BMatrixPolicy.h:55
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Definition: KelvinVector.h:24
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.
Definition: LinearBMatrix.h:42
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType
ConstitutiveTraits::SolidConstitutiveRelation const & solid_material_

References ProcessLib::LinearBMatrix::computeBMatrix(), ParameterLib::SpatialPosition::getCoordinates(), ProcessLib::ThermoRichardsMechanics::TRMVaporDiffusionData< DisplacementDim >::K_pp_X_dNTdN, ProcessLib::ThermoRichardsMechanics::TRMVaporDiffusionData< DisplacementDim >::K_Tp_X_dNTdN, ProcessLib::ThermoRichardsMechanics::TRMHeatStorageAndFluxData< DisplacementDim >::K_Tp_X_NTN, ProcessLib::ThermoRichardsMechanics::TRMVaporDiffusionData< DisplacementDim >::K_TT_X_dNTdN, and MathLib::KelvinVector::kelvin_vector_dimensions().

◆ block_p()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_p ( auto &  vec)
inlinestaticprivate

◆ block_pp()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_pp ( auto &  mat)
inlinestaticprivate

◆ block_pT()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_pT ( auto &  mat)
inlinestaticprivate

◆ block_pu()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_pu ( auto &  mat)
inlinestaticprivate

◆ block_T()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_T ( auto &  vec)
inlinestaticprivate

◆ block_Tp()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_Tp ( auto &  mat)
inlinestaticprivate

◆ block_TT()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_TT ( auto &  mat)
inlinestaticprivate

◆ block_u()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_u ( auto &  vec)
inlinestaticprivate

◆ block_up()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_up ( auto &  mat)
inlinestaticprivate

◆ block_uT()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_uT ( auto &  mat)
inlinestaticprivate

◆ block_uu()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::block_uu ( auto &  mat)
inlinestaticprivate

◆ computeSecondaryVariableConcrete()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::computeSecondaryVariableConcrete ( double const  t,
double const  dt,
Eigen::VectorXd const &  local_x,
Eigen::VectorXd const &  local_x_prev 
)
overridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 548 of file ThermoRichardsMechanicsFEM-impl.h.

552{
553 auto const T = block_T(local_x);
554 auto const p_L = block_p(local_x);
555 auto const u = block_u(local_x);
556
557 auto const T_prev = block_T(local_x_prev);
558 auto const p_L_prev = block_p(local_x_prev);
559
560 auto const e_id = this->element_.getID();
561 auto const& process_data = this->process_data_;
562 auto& medium = *process_data.media_map.getMedium(e_id);
563
564 unsigned const n_integration_points =
566
567 double saturation_avg = 0;
568 double porosity_avg = 0;
569 double liquid_density_avg = 0;
570 double viscosity_avg = 0;
571
573 KV sigma_avg = KV::Zero();
574
575 typename ConstitutiveTraits::ConstitutiveSetting constitutive_setting;
576
577 typename ConstitutiveTraits::ConstitutiveModels models(
578 process_data, this->solid_material_);
579 typename ConstitutiveTraits::ConstitutiveTempData tmp;
580 typename ConstitutiveTraits::ConstitutiveData CD;
581
582 for (unsigned ip = 0; ip < n_integration_points; ip++)
583 {
584 auto& current_state = this->current_states_[ip];
585 auto& output_data = this->output_data_[ip];
586
587 auto const& ip_data = ip_data_[ip];
588
589 // N is used for both p and T variables
590 auto const& N = ip_data.N_p;
591 auto const& N_u = ip_data.N_u;
592 auto const& dNdx_u = ip_data.dNdx_u;
593 auto const& dNdx = ip_data.dNdx_p;
594
595 ParameterLib::SpatialPosition const x_position{
596 std::nullopt, this->element_.getID(), ip,
598 NumLib::interpolateCoordinates<ShapeFunctionDisplacement,
600 this->element_, N_u))};
601 auto const x_coord =
602 NumLib::interpolateXCoordinate<ShapeFunctionDisplacement,
604 this->element_, N_u);
605 auto const B =
606 LinearBMatrix::computeBMatrix<DisplacementDim,
607 ShapeFunctionDisplacement::NPOINTS,
609 dNdx_u, N_u, x_coord, this->is_axially_symmetric_);
610
611 double const T_ip = N * T;
612 double const T_prev_ip = N * T_prev;
613 GlobalDimVectorType const grad_T_ip = dNdx * T;
614
615 double const p_cap_ip = -N * p_L;
616 double const p_cap_prev_ip = -N * p_L_prev;
617 GlobalDimVectorType const grad_p_cap_ip = -dNdx * p_L;
618
619 KelvinVectorType eps = B * u;
620
621 constitutive_setting.eval(models, //
622 t, dt, x_position, //
623 medium, //
624 {T_ip, T_prev_ip, grad_T_ip}, //
625 {p_cap_ip, p_cap_prev_ip, grad_p_cap_ip}, //
626 eps, current_state, this->prev_states_[ip],
627 this->material_states_[ip], tmp, output_data,
628 CD);
629
630 saturation_avg += std::get<SaturationData>(current_state).S_L;
631 porosity_avg += std::get<PorosityData>(current_state).phi;
632
633 liquid_density_avg += std::get<LiquidDensityData>(output_data).rho_LR;
634 viscosity_avg += std::get<LiquidViscosityData>(output_data).viscosity;
635 sigma_avg += ConstitutiveTraits::ConstitutiveSetting::statefulStress(
636 current_state);
637 }
638 saturation_avg /= n_integration_points;
639 porosity_avg /= n_integration_points;
640 viscosity_avg /= n_integration_points;
641 liquid_density_avg /= n_integration_points;
642 sigma_avg /= n_integration_points;
643
644 (*process_data.element_saturation)[e_id] = saturation_avg;
645 (*process_data.element_porosity)[e_id] = porosity_avg;
646 (*process_data.element_liquid_density)[e_id] = liquid_density_avg;
647 (*process_data.element_viscosity)[e_id] = viscosity_avg;
648
649 Eigen::Map<KV>(
650 &(*process_data.element_stresses)[e_id * KV::RowsAtCompileTime]) =
652
654 ShapeFunction, typename ShapeFunctionDisplacement::MeshElement,
655 DisplacementDim>(this->element_, this->is_axially_symmetric_, p_L,
656 *process_data.pressure_interpolated);
658 ShapeFunction, typename ShapeFunctionDisplacement::MeshElement,
659 DisplacementDim>(this->element_, this->is_axially_symmetric_, T,
660 *process_data.temperature_interpolated);
661}
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:49
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
void interpolateToHigherOrderNodes(MeshLib::Element const &element, bool const is_axially_symmetric, Eigen::MatrixBase< EigenMatrixType > const &node_values, MeshLib::PropertyVector< double > &interpolated_values_global_vector)

References ProcessLib::LinearBMatrix::computeBMatrix(), NumLib::interpolateCoordinates(), NumLib::interpolateToHigherOrderNodes(), NumLib::interpolateXCoordinate(), and MathLib::KelvinVector::kelvinVectorToSymmetricTensor().

◆ convertInitialStressType()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::convertInitialStressType ( unsigned const  ip,
double const  t,
ParameterLib::SpatialPosition const  x_position,
MaterialPropertyLib::Medium const &  medium,
MPL::VariableArray const &  variables,
double const  p_at_ip 
)
private

This function resets the initial stress type according to the input initial stress type, either total or effective. If subtype = StressSaturation_StrainPressureTemperature is bing used in the process setting, the initial effective stress is converted to total stress. Otherwise, the initial total stress is converted to effective stress.

Definition at line 182 of file ThermoRichardsMechanicsFEM-impl.h.

189{
190 bool constexpr is_strain_temperature_constitutive =
191 std::is_same<ConstitutiveStress_StrainTemperature::ConstitutiveTraits<
192 DisplacementDim>,
193 ConstitutiveTraits>::value;
194 if (is_strain_temperature_constitutive &&
195 this->process_data_.initial_stress.type ==
197 {
198 return;
199 }
200
201 if (!is_strain_temperature_constitutive &&
202 this->process_data_.initial_stress.type == InitialStress::Type::Total)
203 {
204 return;
205 }
206
207 double const alpha_b =
208 medium.property(MPL::PropertyType::biot_coefficient)
209 .template value<double>(variables, x_position, t, 0.0 /*dt*/);
210
211 double const bishop =
212 medium.property(MPL::PropertyType::bishops_effective_stress)
213 .template value<double>(variables, x_position, t, 0.0 /*dt*/);
214
215 ConstitutiveTraits::ConstitutiveSetting::convertInitialStressType(
216 this->current_states_[ip], this->prev_states_[ip],
217 bishop * alpha_b * p_at_ip * Invariants::identity2);
218}
static Eigen::Matrix< double, KelvinVectorSize, 1 > const identity2
Kelvin mapping of 2nd order identity tensor.
Definition: KelvinVector.h:107

References ProcessLib::InitialStress::Effective, MaterialPropertyLib::Medium::property(), and ProcessLib::InitialStress::Total.

◆ getShapeMatrix()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
Eigen::Map< const Eigen::RowVectorXd > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::getShapeMatrix ( const unsigned  integration_point) const
inlineoverridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 349 of file ThermoRichardsMechanicsFEM.h.

351 {
352 auto const& N_u = ip_data_[integration_point].N_u;
353
354 // assumes N is stored contiguously in memory
355 return Eigen::Map<const Eigen::RowVectorXd>(N_u.data(), N_u.size());
356 }

References ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ip_data_.

◆ initializeConcrete()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::initializeConcrete ( )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 276 of file ThermoRichardsMechanicsFEM.h.

277 {
278 unsigned const n_integration_points =
280 auto const time_independent = std::numeric_limits<double>::quiet_NaN();
281 auto const& medium =
282 *this->process_data_.media_map.getMedium(this->element_.getID());
283
284 for (unsigned ip = 0; ip < n_integration_points; ip++)
285 {
286 ParameterLib::SpatialPosition const x_position{
287 std::nullopt, this->element_.getID(), ip,
289 ShapeFunctionDisplacement,
291 this->element_, this->ip_data_[ip].N_u))};
292 auto& current_state = this->current_states_[ip];
293
294 // Set initial stress from parameter.
295 if (this->process_data_.initial_stress.value)
296 {
297 ConstitutiveTraits::ConstitutiveSetting::statefulStress(
298 current_state) =
300 DisplacementDim>(
301 (*this->process_data_.initial_stress.value)(
302 time_independent, x_position));
303 }
304
305 if (*this->process_data_.initialize_porosity_from_medium_property)
306 {
307 // Initial porosity. Could be read from integration point data
308 // or mesh.
309 std::get<PorosityData>(current_state).phi =
310 medium.property(MPL::porosity)
311 .template initialValue<double>(x_position,
312 time_independent);
313
314 if (medium.hasProperty(MPL::PropertyType::transport_porosity))
315 {
316 std::get<TransportPorosityData>(current_state).phi =
317 medium.property(MPL::transport_porosity)
318 .template initialValue<double>(x_position,
319 time_independent);
320 }
321 else
322 {
323 std::get<TransportPorosityData>(current_state).phi =
324 std::get<PorosityData>(current_state).phi;
325 }
326 }
327
328 double const t = 0; // TODO (naumov) pass t from top
329 this->solid_material_.initializeInternalStateVariables(
330 t, x_position,
331 *this->material_states_[ip].material_state_variables);
332 }
333
334 for (unsigned ip = 0; ip < n_integration_points; ip++)
335 {
336 this->material_states_[ip].pushBackState();
337 }
338
339 for (unsigned ip = 0; ip < n_integration_points; ip++)
340 {
341 this->prev_states_[ip] = this->current_states_[ip];
342 }
343 }
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
Definition: KelvinVector.h:201

References ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::current_states_, ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::element_, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::integration_method_, NumLib::interpolateCoordinates(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::ip_data_, ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::material_states_, MaterialPropertyLib::porosity, ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::prev_states_, ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::process_data_, ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits >::solid_material_, MathLib::KelvinVector::symmetricTensorToKelvinVector(), and MaterialPropertyLib::transport_porosity.

◆ localDOF()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
static constexpr auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::localDOF ( std::vector< double > const &  x)
inlinestaticconstexprprivate

Definition at line 361 of file ThermoRichardsMechanicsFEM.h.

362 {
363 return NumLib::localDOF<
364 ShapeFunction, ShapeFunction,
366 }
auto localDOF(ElementDOFVector const &x)
Definition: LocalDOF.h:64

References NumLib::localDOF().

◆ massLumping()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::massLumping ( LocalMatrices loc_mat) const
private

Definition at line 274 of file ThermoRichardsMechanicsFEM-impl.h.

278{
279 if (this->process_data_.apply_mass_lumping)
280 {
281 loc_mat.storage_p_a_p =
282 loc_mat.storage_p_a_p.colwise().sum().eval().asDiagonal();
283 loc_mat.storage_p_a_S =
284 loc_mat.storage_p_a_S.colwise().sum().eval().asDiagonal();
285 loc_mat.storage_p_a_S_Jpp =
286 loc_mat.storage_p_a_S_Jpp.colwise().sum().eval().asDiagonal();
287 }
288}

◆ setInitialConditionsConcrete()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::setInitialConditionsConcrete ( std::vector< double > const &  local_x,
double const  t,
bool const  use_monolithic_scheme,
int const  process_id 
)
overridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 89 of file ThermoRichardsMechanicsFEM-impl.h.

94{
95 assert(local_x.size() ==
97
98 auto const p_L = Eigen::Map<
99 typename ShapeMatricesType::template VectorType<pressure_size> const>(
100 local_x.data() + pressure_index, pressure_size);
101
102 auto const T = Eigen::Map<typename ShapeMatricesType::template VectorType<
103 temperature_size> const>(local_x.data() + temperature_index,
105
106 constexpr double dt = std::numeric_limits<double>::quiet_NaN();
107 auto const& medium =
108 this->process_data_.media_map.getMedium(this->element_.getID());
109 MPL::VariableArray variables;
110
111 auto const& solid_phase = medium->phase("Solid");
112
113 unsigned const n_integration_points =
115 for (unsigned ip = 0; ip < n_integration_points; ip++)
116 {
117 // N is used for both T and p variables.
118 auto const& N = ip_data_[ip].N_p;
119
120 ParameterLib::SpatialPosition const x_position{
121 std::nullopt, this->element_.getID(), ip,
123 NumLib::interpolateCoordinates<ShapeFunctionDisplacement,
125 this->element_, ip_data_[ip].N_u))};
126 double p_cap_ip;
127 NumLib::shapeFunctionInterpolate(-p_L, N, p_cap_ip);
128
129 variables.capillary_pressure = p_cap_ip;
130 variables.liquid_phase_pressure = -p_cap_ip;
131 // setting pG to 1 atm
132 // TODO : rewrite equations s.t. p_L = pG-p_cap
133 variables.gas_phase_pressure = 1.0e5;
134
135 double T_ip;
137 variables.temperature = T_ip;
138
139 double const S_L =
140 medium->property(MPL::PropertyType::saturation)
141 .template value<double>(variables, x_position, t, dt);
142
143 std::get<PrevState<SaturationData>>(this->prev_states_[ip])->S_L = S_L;
144
145 {
146 // Set eps_m_prev from potentially non-zero eps and sigma_sw from
147 // restart.
148 SpaceTimeData const x_t{x_position, t, dt};
149 ElasticTangentStiffnessData<DisplacementDim> C_el_data;
150 typename ConstitutiveTraits::ElasticTangentStiffnessModel{
151 this->solid_material_}
152 .eval(x_t, {T_ip, 0, {}}, C_el_data);
153
154 auto const& eps =
155 std::get<StrainData<DisplacementDim>>(this->current_states_[ip])
156 .eps;
157 auto const& sigma_sw =
158 std::get<SwellingDataStateful<DisplacementDim>>(
159 this->current_states_[ip])
160 .sigma_sw;
161 std::get<PrevState<MechanicalStrainData<DisplacementDim>>>(
162 this->prev_states_[ip])
163 ->eps_m.noalias() =
164 solid_phase.hasProperty(MPL::PropertyType::swelling_stress_rate)
165 ? eps + C_el_data.C_el.inverse() * sigma_sw
166 : eps;
167 }
168
169 if (this->process_data_.initial_stress.value)
170 {
171 variables.liquid_saturation = S_L;
172 convertInitialStressType(ip, t, x_position, *medium, variables,
173 -p_cap_ip);
174 }
175 }
176}
void convertInitialStressType(unsigned const ip, double const t, ParameterLib::SpatialPosition const x_position, MaterialPropertyLib::Medium const &medium, MPL::VariableArray const &variables, double const p_at_ip)
void shapeFunctionInterpolate(const NodalValues &, const ShapeMatrix &)
Definition: Interpolation.h:27

References ProcessLib::ThermoRichardsMechanics::ElasticTangentStiffnessData< DisplacementDim >::C_el, MaterialPropertyLib::VariableArray::capillary_pressure, MaterialPropertyLib::VariableArray::gas_phase_pressure, NumLib::interpolateCoordinates(), MaterialPropertyLib::VariableArray::liquid_phase_pressure, MaterialPropertyLib::VariableArray::liquid_saturation, NumLib::detail::shapeFunctionInterpolate(), and MaterialPropertyLib::VariableArray::temperature.

Member Data Documentation

◆ displacement_index

◆ displacement_size

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
constexpr int ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::displacement_size
staticconstexprprivate

◆ ip_data_

◆ KelvinVectorSize

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
int const ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::KelvinVectorSize
static
Initial value:

Definition at line 67 of file ThermoRichardsMechanicsFEM.h.

◆ N_u_op

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
constexpr auto& ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::N_u_op
staticconstexpr
Initial value:
DisplacementDim,
constexpr Eigen::CwiseNullaryOp< EigenBlockMatrixViewFunctor< D, M >, typename EigenBlockMatrixViewFunctor< D, M >::Matrix > eigenBlockMatrixView(const Eigen::MatrixBase< M > &matrix)
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType

Definition at line 73 of file ThermoRichardsMechanicsFEM.h.

◆ pressure_index

◆ pressure_size

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
constexpr int ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::pressure_size = ShapeFunction::NPOINTS
staticconstexprprivate

◆ temperature_index

◆ temperature_size

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim, typename ConstitutiveTraits >
constexpr int ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim, ConstitutiveTraits >::temperature_size = ShapeFunction::NPOINTS
staticconstexprprivate

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