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

Detailed Description

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

Definition at line 48 of file ThermoRichardsMechanicsFEM.h.

#include <ThermoRichardsMechanicsFEM.h>

Inheritance diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >:
[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 > &process_data)
 
std::size_t setIPDataInitialConditions (std::string const &name, double const *values, int const integration_order) override
 
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_xdot, std::vector< double > &, std::vector< double > &, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
 
void initializeConcrete () override
 
void postTimestepConcrete (Eigen::VectorXd const &, double const, double const) override
 
void computeSecondaryVariableConcrete (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_dot) override
 
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix (const unsigned integration_point) const override
 Provides the shape matrix at the given integration point. More...
 
std::vector< double > getSigma () const override
 
std::vector< double > const & getIntPtDarcyVelocity (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 > getSaturation () const override
 
std::vector< double > const & getIntPtSaturation (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 > getPorosity () const override
 
std::vector< double > const & getIntPtPorosity (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 > getTransportPorosity () const override
 
std::vector< double > const & getIntPtTransportPorosity (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 & getIntPtSigma (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 > getSwellingStress () const override
 
std::vector< double > const & getIntPtSwellingStress (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 > 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 & getIntPtDryDensitySolid (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 & getIntPtLiquidDensity (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 & getIntPtViscosity (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) 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_xdot, 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_xdot, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
 
virtual void assembleWithJacobianForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, 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_dot, 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, double const t, double const dt)
 
void postNonLinearSolver (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, GlobalVector const &xdot, 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. More...
 
- Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default
 

Static Public Attributes

static int const KelvinVectorSize
 

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_xdot, IpData const &ip_data, ConstitutiveSetting< DisplacementDim > &CS, MaterialPropertyLib::Medium &medium, LocalMatrices &out, StatefulData< DisplacementDim > &current_state, StatefulData< DisplacementDim > const &prev_state, MaterialStateData< DisplacementDim > &mat_state, OutputData< DisplacementDim > &output_data) const
 
void addToLocalMatrixData (double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, LocalMatrices const &loc_mat, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) const
 
void massLumping (LocalMatrices &loc_mat) const
 
auto localDOF (std::vector< double > const &local_dof_data) const
 
unsigned getNumberOfIntegrationPoints () const override
 
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt (unsigned integration_point) const override
 

Static Private Member Functions

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

ThermoRichardsMechanicsProcessData< DisplacementDim > & process_data_
 
std::vector< StatefulData< DisplacementDim > > current_states_
 
std::vector< StatefulData< DisplacementDim > > prev_states_
 
std::vector< MaterialStateData< DisplacementDim > > material_states_
 
std::vector< IpDataip_data_
 
NumLib::GenericIntegrationMethod const & integration_method_
 
MeshLib::Element const & element_
 
bool const is_axially_symmetric_
 
SecondaryData< typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeTypesecondary_data_
 
MaterialLib::Solids::MechanicsBase< DisplacementDim > const & solid_material_
 
std::vector< OutputData< DisplacementDim > > output_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
 

Member Typedef Documentation

◆ BMatricesType

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

Definition at line 70 of file ThermoRichardsMechanicsFEM.h.

◆ GlobalDimMatrixType

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

Definition at line 67 of file ThermoRichardsMechanicsFEM.h.

◆ GlobalDimVectorType

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

Definition at line 68 of file ThermoRichardsMechanicsFEM.h.

◆ Invariants

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

Definition at line 80 of file ThermoRichardsMechanicsFEM.h.

◆ IpData

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

Definition at line 74 of file ThermoRichardsMechanicsFEM.h.

◆ KelvinVectorType

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

Definition at line 72 of file ThermoRichardsMechanicsFEM.h.

◆ ShapeMatricesType

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

Definition at line 64 of file ThermoRichardsMechanicsFEM.h.

◆ ShapeMatricesTypeDisplacement

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

Definition at line 60 of file ThermoRichardsMechanicsFEM.h.

◆ SymmetricTensor

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

Definition at line 82 of file ThermoRichardsMechanicsFEM.h.

Constructor & Destructor Documentation

◆ ThermoRichardsMechanicsLocalAssembler() [1/3]

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

◆ ThermoRichardsMechanicsLocalAssembler() [2/3]

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

◆ ThermoRichardsMechanicsLocalAssembler() [3/3]

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

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

45  : process_data_(process_data),
46  integration_method_(integration_method),
47  element_(e),
48  is_axially_symmetric_(is_axially_symmetric),
50  process_data_.solid_materials, process_data_.material_ids, e.getID()))
51 {
52  unsigned const n_integration_points =
54 
55  current_states_.resize(n_integration_points);
56  prev_states_.resize(n_integration_points);
57  ip_data_.resize(n_integration_points);
58  secondary_data_.N_u.resize(n_integration_points);
59  output_data_.resize(n_integration_points);
60 
61  material_states_.reserve(n_integration_points);
62  for (unsigned ip = 0; ip < n_integration_points; ++ip)
63  {
64  material_states_.emplace_back(solid_material_);
65  }
66 
67  auto const shape_matrices_u =
68  NumLib::initShapeMatrices<ShapeFunctionDisplacement,
70  DisplacementDim>(e, is_axially_symmetric,
72 
73  auto const shape_matrices =
75  DisplacementDim>(e, is_axially_symmetric,
77 
78  auto const& medium = process_data_.media_map->getMedium(element_.getID());
79 
81  x_position.setElementID(element_.getID());
82  for (unsigned ip = 0; ip < n_integration_points; ip++)
83  {
84  x_position.setIntegrationPoint(ip);
85  auto& ip_data = ip_data_[ip];
86  auto const& sm_u = shape_matrices_u[ip];
87  ip_data_[ip].integration_weight =
89  sm_u.integralMeasure * sm_u.detJ;
90 
91  ip_data.N_u_op = ShapeMatricesTypeDisplacement::template MatrixType<
92  DisplacementDim, displacement_size>::Zero(DisplacementDim,
94  for (int i = 0; i < DisplacementDim; ++i)
95  {
96  ip_data.N_u_op
97  .template block<1, displacement_size / DisplacementDim>(
98  i, i * displacement_size / DisplacementDim)
99  .noalias() = sm_u.N;
100  }
101 
102  ip_data.N_u = sm_u.N;
103  ip_data.dNdx_u = sm_u.dNdx;
104 
105  // ip_data.N_p and ip_data.dNdx_p are used for both p and T variables
106  ip_data.N_p = shape_matrices[ip].N;
107  ip_data.dNdx_p = shape_matrices[ip].dNdx;
108 
109  auto& current_state = current_states_[ip];
110 
111  // Initial porosity. Could be read from integration point data or mesh.
112  current_state.poro_data.phi =
113  medium->property(MPL::porosity)
114  .template initialValue<double>(
115  x_position,
116  std::numeric_limits<
117  double>::quiet_NaN() /* t independent */);
118 
119  if (medium->hasProperty(MPL::PropertyType::transport_porosity))
120  {
121  current_state.transport_poro_data.phi =
122  medium->property(MPL::transport_porosity)
123  .template initialValue<double>(
124  x_position,
125  std::numeric_limits<
126  double>::quiet_NaN() /* t independent */);
127  }
128  else
129  {
130  current_state.transport_poro_data.phi = current_state.poro_data.phi;
131  }
132 
133  secondary_data_.N_u[ip] = shape_matrices_u[ip].N;
134  }
135 }
double getWeight() const
Definition: WeightedPoint.h:80
virtual std::size_t getID() const final
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)
void setIntegrationPoint(unsigned integration_point)
SecondaryData< typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeType > secondary_data_
ThermoRichardsMechanicsProcessData< DisplacementDim > & process_data_
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
MaterialLib::Solids::MechanicsBase< DisplacementDim > const & solid_material_
MechanicsBase< DisplacementDim > & selectSolidConstitutiveRelation(std::map< int, std::unique_ptr< MechanicsBase< DisplacementDim >>> const &constitutive_relations, MeshLib::PropertyVector< int > const *const material_ids, std::size_t const element_id)
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_u

References ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::current_states_, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::displacement_size, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::element_, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::initShapeMatrices(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::integration_method_, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::ip_data_, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::material_states_, ProcessLib::HydroMechanics::SecondaryData< ShapeMatrixType >::N_u, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::output_data_, MaterialPropertyLib::porosity, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::prev_states_, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::process_data_, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::secondary_data_, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::solid_material_, and MaterialPropertyLib::transport_porosity.

Member Function Documentation

◆ addToLocalMatrixData()

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

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

348 {
349  constexpr auto local_matrix_dim =
351  assert(local_x.size() == local_matrix_dim);
352 
353  auto local_Jac = MathLib::createZeroedMatrix<
354  typename ShapeMatricesTypeDisplacement::template MatrixType<
355  local_matrix_dim, local_matrix_dim>>(
356  local_Jac_data, local_matrix_dim, local_matrix_dim);
357 
358  auto local_rhs =
359  MathLib::createZeroedVector<typename ShapeMatricesTypeDisplacement::
360  template VectorType<local_matrix_dim>>(
361  local_rhs_data, local_matrix_dim);
362 
363  local_Jac.noalias() = loc_mat.Jac;
364  local_rhs.noalias() = -loc_mat.res;
365 
366  //
367  // -- Jacobian
368  //
369  block_TT(local_Jac).noalias() += loc_mat.M_TT / dt + loc_mat.K_TT;
370  block_Tp(local_Jac).noalias() +=
371  loc_mat.K_Tp + loc_mat.dK_TT_dp + loc_mat.M_Tp / dt;
372 
373  block_pT(local_Jac).noalias() += loc_mat.M_pT / dt + loc_mat.K_pT;
374  block_pp(local_Jac).noalias() +=
375  loc_mat.K_pp + loc_mat.storage_p_a_p / dt + loc_mat.storage_p_a_S_Jpp;
376  block_pu(local_Jac).noalias() = loc_mat.M_pu / dt;
377 
378  //
379  // -- Residual
380  //
381  auto const [T, p_L, u] = localDOF(local_x);
382  auto const [T_dot, p_L_dot, u_dot] = localDOF(local_xdot);
383 
384  block_T(local_rhs).noalias() -= loc_mat.M_TT * T_dot + loc_mat.K_TT * T +
385  loc_mat.K_Tp * p_L + loc_mat.M_Tp * p_L_dot;
386  block_p(local_rhs).noalias() -=
387  loc_mat.K_pp * p_L + loc_mat.K_pT * T +
388  (loc_mat.storage_p_a_p + loc_mat.storage_p_a_S) * p_L_dot +
389  loc_mat.M_pu * u_dot + loc_mat.M_pT * T_dot;
390 }
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
Definition: EigenMapTools.h:32
static const double u
Eigen::Map< Vector > createZeroedVector(std::vector< double > &data, Eigen::VectorXd::Index size)

References MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::dK_TT_dp, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::Jac, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_pp, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_pT, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_Tp, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_TT, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_pT, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_pu, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_Tp, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_TT, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::res, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_p, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_S, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_S_Jpp, and MathLib::u.

◆ assembleWithJacobian()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::assembleWithJacobian ( double const  t,
double const  dt,
std::vector< double > const &  local_x,
std::vector< double > const &  local_xdot,
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 275 of file ThermoRichardsMechanicsFEM-impl.h.

283 {
284  auto& medium = *process_data_.media_map->getMedium(element_.getID());
285 
287  x_position.setElementID(element_.getID());
288 
289  LocalMatrices loc_mat;
290  loc_mat.setZero();
291  LocalMatrices loc_mat_current_ip;
292  loc_mat_current_ip.setZero(); // only to set the right matrix sizes
293 
294  ConstitutiveSetting<DisplacementDim> constitutive_setting(solid_material_,
295  process_data_);
296 
297  for (unsigned ip = 0; ip < integration_method_.getNumberOfPoints(); ++ip)
298  {
299  x_position.setIntegrationPoint(ip);
300 
301  assembleWithJacobianSingleIP(t, dt, x_position, //
302  local_x, local_xdot, //
303  ip_data_[ip], constitutive_setting,
304  medium, //
305  loc_mat_current_ip, //
306  current_states_[ip], prev_states_[ip],
307  material_states_[ip], output_data_[ip]);
308  loc_mat += loc_mat_current_ip;
309  }
310 
311  massLumping(loc_mat);
312 
313  addToLocalMatrixData(dt, local_x, local_xdot, loc_mat, local_rhs_data,
314  local_Jac_data);
315 }
void addToLocalMatrixData(double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, 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_xdot, IpData const &ip_data, ConstitutiveSetting< DisplacementDim > &CS, MaterialPropertyLib::Medium &medium, LocalMatrices &out, StatefulData< DisplacementDim > &current_state, StatefulData< DisplacementDim > const &prev_state, MaterialStateData< DisplacementDim > &mat_state, OutputData< DisplacementDim > &output_data) const
static const double t

References ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::setZero(), and MathLib::t.

◆ assembleWithJacobianSingleIP()

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

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

413 {
414  auto const& N_u = ip_data.N_u;
415  auto const& N_u_op = ip_data.N_u_op;
416  auto const& dNdx_u = ip_data.dNdx_u;
417 
418  // N and dNdx are used for both p and T variables
419  auto const& N = ip_data.N_p;
420  auto const& dNdx = ip_data.dNdx_p;
421 
422  auto const x_coord =
423  NumLib::interpolateXCoordinate<ShapeFunctionDisplacement,
425  N_u);
426  auto const B =
427  LinearBMatrix::computeBMatrix<DisplacementDim,
428  ShapeFunctionDisplacement::NPOINTS,
429  typename BMatricesType::BMatrixType>(
430  dNdx_u, N_u, x_coord, is_axially_symmetric_);
431 
432  auto const [T, p_L, u] = localDOF(local_x);
433  auto const [T_dot, p_L_dot, u_dot] = localDOF(local_xdot);
434 
435  GlobalDimVectorType const grad_T_ip = dNdx * T;
436 
437  ConstitutiveModels<DisplacementDim> models(process_data_, solid_material_);
438  ConstitutiveTempData<DisplacementDim> tmp;
439  ConstitutiveData<DisplacementDim> CD;
440 
441  {
442  double const T_ip = N * T;
443  double const T_dot_ip = N * T_dot;
444 
445  double const p_cap_ip = -N * p_L;
446  double const p_cap_dot_ip = -N * p_L_dot;
447  GlobalDimVectorType const grad_p_cap_ip = -dNdx * p_L;
448 
449  KelvinVectorType eps = B * u;
450  // TODO conceptual consistency check. introduced for volumetric strain
451  // rate computation
452  KelvinVectorType eps_prev = eps - B * (u_dot * dt);
453 
454  CS.eval(models, t, dt, x_position, //
455  medium, //
456  {T_ip, T_dot_ip, grad_T_ip}, //
457  {p_cap_ip, p_cap_dot_ip, grad_p_cap_ip}, //
458  eps, eps_prev, current_state, prev_state, mat_state, tmp,
459  output_data, CD);
460  }
461 
462  using NodalMatrix = typename ShapeMatricesType::NodalMatrixType;
463  NodalMatrix const NTN = N.transpose() * N;
464  NodalMatrix const dNTdN = dNdx.transpose() * dNdx;
465 
466  // TODO is identity2.transpose() * B something like divergence?
467  auto const& identity2 = MathLib::KelvinVector::Invariants<
469  DisplacementDim)>::identity2;
470  typename ShapeMatricesTypeDisplacement::template MatrixType<
472  BTI2N = B.transpose() * identity2 * N;
473 
474  /*
475  * Conventions:
476  *
477  * * use positive signs exclusively, any negative signs must be included in
478  * the coefficients coming from the constitutive setting
479  * * the used coefficients from the constitutive setting are named after the
480  * terms they appear in, e.g. K_TT_X_dNTdN means it is a scalar (X) that
481  * will be multiplied by dNdx.transpose() * dNdx. Placefolders for the
482  * coefficients are:
483  * * X -> scalar
484  * * V -> vector
485  * * K -> Kelvin vector
486  * * the Laplace terms have a special name, e.g., K_TT_Laplace
487  * * there shall be only one contribution to each of the LocalMatrices,
488  * assigned with = assignment; this point might be relaxed in the future
489  * * this method will overwrite the data in the passed LocalMatrices& out
490  * argument, not add to it
491  */
492 
493  // residual, order T, p, u
494  block_p(out.res).noalias() = dNdx.transpose() * CD.eq_p_data.rhs_p_dNT_V;
495  block_u(out.res).noalias() =
496  B.transpose() * CD.s_mech_data.sigma_total -
497  N_u_op.transpose() * CD.grav_data.volumetric_body_force;
498 
499  // Storage matrices
500  out.storage_p_a_p.noalias() = CD.eq_p_data.storage_p_a_p_X_NTN * NTN;
501  out.storage_p_a_S.noalias() = CD.storage_data.storage_p_a_S_X_NTN * NTN;
502  out.storage_p_a_S_Jpp.noalias() =
503  CD.storage_data.storage_p_a_S_Jpp_X_NTN * NTN;
504 
505  // M matrices, order T, p, u
506  out.M_TT.noalias() = CD.eq_T_data.M_TT_X_NTN * NTN;
507  out.M_Tp.noalias() = CD.vap_data.M_Tp_X_NTN * NTN;
508 
509  out.M_pT.noalias() = CD.eq_p_data.M_pT_X_NTN * NTN;
510  out.M_pu.noalias() = CD.eq_p_data.M_pu_X_BTI2N * BTI2N.transpose();
511 
512  // K matrices, order T, p, u
513  out.K_TT.noalias() =
514  dNdx.transpose() * CD.heat_data.K_TT_Laplace * dNdx +
515  N.transpose() * (CD.eq_T_data.K_TT_NT_V_dN.transpose() * dNdx) +
516  CD.vap_data.K_TT_X_dNTdN * dNTdN;
517 
518  out.dK_TT_dp.noalias() =
519  N.transpose() * (CD.heat_data.K_Tp_NT_V_dN.transpose() * dNdx) +
520  CD.heat_data.K_Tp_X_NTN * NTN + CD.vap_data.K_Tp_X_dNTdN * dNTdN;
521  out.K_Tp.noalias() =
522  dNdx.transpose() * CD.th_osmosis_data.K_Tp_Laplace * dNdx;
523 
524  out.K_pp.noalias() = dNdx.transpose() * CD.eq_p_data.K_pp_Laplace * dNdx +
525  CD.vap_data.K_pp_X_dNTdN * dNTdN;
526  out.K_pT.noalias() =
527  dNdx.transpose() * CD.th_osmosis_data.K_pT_Laplace * dNdx;
528 
529  // direct Jacobian contributions, order T, p, u
530  block_pT(out.Jac).noalias() = CD.vap_data.J_pT_X_dNTdN * dNTdN;
531  block_pp(out.Jac).noalias() =
532  CD.storage_data.J_pp_X_NTN * NTN +
533  CD.eq_p_data.J_pp_X_BTI2NT_u_dot_N * BTI2N.transpose() * u_dot *
534  N // TODO something with volumetric strain rate?
535  + dNdx.transpose() * CD.eq_p_data.J_pp_dNT_V_N * N;
536 
537  block_uT(out.Jac).noalias() =
538  B.transpose() * CD.s_mech_data.J_uT_BT_K_N * N;
539  block_up(out.Jac).noalias() =
540  B.transpose() * CD.swelling_data.J_up_BT_K_N * N +
541  CD.eq_u_data.J_up_X_BTI2N * BTI2N +
542  N_u_op.transpose() * CD.eq_u_data.J_up_HT_V_N * N;
543  block_uu(out.Jac).noalias() =
544  B.transpose() * CD.s_mech_data.stiffness_tensor * B;
545 
546  out *= ip_data.integration_weight;
547 }
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:23
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
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

References ProcessLib::LinearBMatrix::computeBMatrix(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::dK_TT_dp, ProcessLib::ThermoRichardsMechanics::IntegrationPointData< ShapeMatrixTypeDisplacement, ShapeMatricesTypePressure, DisplacementDim, NPoints >::dNdx_p, ProcessLib::ThermoRichardsMechanics::IntegrationPointData< ShapeMatrixTypeDisplacement, ShapeMatricesTypePressure, DisplacementDim, NPoints >::dNdx_u, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::eq_p_data, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::eq_T_data, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::eq_u_data, ProcessLib::ThermoRichardsMechanics::ConstitutiveSetting< DisplacementDim >::eval(), ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::grav_data, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::heat_data, ProcessLib::ThermoRichardsMechanics::IntegrationPointData< ShapeMatrixTypeDisplacement, ShapeMatricesTypePressure, DisplacementDim, NPoints >::integration_weight, NumLib::interpolateXCoordinate(), ProcessLib::ThermoRichardsMechanics::TRMStorageData::J_pp_X_NTN, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::Jac, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_pp, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_pT, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_Tp, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::K_TT, MathLib::KelvinVector::kelvin_vector_dimensions(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_pT, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_pu, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_Tp, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::M_TT, ProcessLib::ThermoRichardsMechanics::IntegrationPointData< ShapeMatrixTypeDisplacement, ShapeMatricesTypePressure, DisplacementDim, NPoints >::N_p, ProcessLib::ThermoRichardsMechanics::IntegrationPointData< ShapeMatrixTypeDisplacement, ShapeMatricesTypePressure, DisplacementDim, NPoints >::N_u, ProcessLib::ThermoRichardsMechanics::IntegrationPointData< ShapeMatrixTypeDisplacement, ShapeMatricesTypePressure, DisplacementDim, NPoints >::N_u_op, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::res, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::s_mech_data, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::storage_data, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_p, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_S, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_S_Jpp, ProcessLib::ThermoRichardsMechanics::TRMStorageData::storage_p_a_S_Jpp_X_NTN, ProcessLib::ThermoRichardsMechanics::TRMStorageData::storage_p_a_S_X_NTN, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::swelling_data, MathLib::t, ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::th_osmosis_data, MathLib::u, and ProcessLib::ThermoRichardsMechanics::ConstitutiveData< DisplacementDim >::vap_data.

◆ block_p()

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

◆ block_pp()

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

◆ block_pT()

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

◆ block_pu()

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

◆ block_T()

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

◆ block_Tp()

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

◆ block_TT()

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

◆ block_u()

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

◆ block_up()

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

◆ block_uT()

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

◆ block_uu()

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

◆ computeSecondaryVariableConcrete()

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

Reimplemented from ProcessLib::LocalAssemblerInterface.

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

815 {
816  auto const T = block_T(local_x);
817  auto const p_L = block_p(local_x);
818  auto const u = block_u(local_x);
819 
820  auto const T_dot = block_T(local_x_dot);
821  auto const p_L_dot = block_p(local_x_dot);
822  auto const u_dot = block_u(local_x_dot);
823 
824  auto const e_id = element_.getID();
825  auto& medium = *process_data_.media_map->getMedium(e_id);
826 
828  x_position.setElementID(e_id);
829 
830  unsigned const n_integration_points =
832 
833  double saturation_avg = 0;
834  double porosity_avg = 0;
835  double liquid_density_avg = 0;
836  double viscosity_avg = 0;
837 
839  KV sigma_avg = KV::Zero();
840 
841  ConstitutiveSetting<DisplacementDim> constitutive_setting(solid_material_,
842  process_data_);
843 
844  ConstitutiveModels<DisplacementDim> models(process_data_, solid_material_);
845  ConstitutiveTempData<DisplacementDim> tmp;
846  ConstitutiveData<DisplacementDim> CD;
847 
848  for (unsigned ip = 0; ip < n_integration_points; ip++)
849  {
850  x_position.setIntegrationPoint(ip);
851 
852  auto const& ip_data = ip_data_[ip];
853 
854  // N is used for both p and T variables
855  auto const& N = ip_data.N_p;
856  auto const& N_u = ip_data.N_u;
857  auto const& dNdx_u = ip_data.dNdx_u;
858  auto const& dNdx = ip_data.dNdx_p;
859 
860  auto const x_coord =
861  NumLib::interpolateXCoordinate<ShapeFunctionDisplacement,
863  element_, N_u);
864  auto const B =
865  LinearBMatrix::computeBMatrix<DisplacementDim,
866  ShapeFunctionDisplacement::NPOINTS,
867  typename BMatricesType::BMatrixType>(
868  dNdx_u, N_u, x_coord, is_axially_symmetric_);
869 
870  double const T_ip = N * T;
871  double const T_dot_ip = N * T_dot;
872  GlobalDimVectorType const grad_T_ip = dNdx * T;
873 
874  double const p_cap_ip = -N * p_L;
875  double const p_cap_dot_ip = -N * p_L_dot;
876  GlobalDimVectorType const grad_p_cap_ip = -dNdx * p_L;
877 
878  KelvinVectorType eps = B * u;
879  // TODO conceptual consistency check. introduced for volumetric strain
880  // rate computation
881  KelvinVectorType eps_prev = eps - B * (u_dot * dt);
882 
883  constitutive_setting.eval(models, //
884  t, dt, x_position, //
885  medium, //
886  {T_ip, T_dot_ip, grad_T_ip}, //
887  {p_cap_ip, p_cap_dot_ip, grad_p_cap_ip}, //
888  eps, eps_prev, current_states_[ip],
889  prev_states_[ip], material_states_[ip], tmp,
890  output_data_[ip], CD);
891 
892  saturation_avg += current_states_[ip].S_L_data.S_L;
893  porosity_avg += current_states_[ip].poro_data.phi;
894 
895  liquid_density_avg += output_data_[ip].rho_L_data.rho_LR;
896  viscosity_avg += output_data_[ip].mu_L_data.viscosity;
897  sigma_avg += current_states_[ip].s_mech_data.sigma_eff;
898  }
899  saturation_avg /= n_integration_points;
900  porosity_avg /= n_integration_points;
901  viscosity_avg /= n_integration_points;
902  liquid_density_avg /= n_integration_points;
903  sigma_avg /= n_integration_points;
904 
905  (*process_data_.element_saturation)[e_id] = saturation_avg;
906  (*process_data_.element_porosity)[e_id] = porosity_avg;
907  (*process_data_.element_liquid_density)[e_id] = liquid_density_avg;
908  (*process_data_.element_viscosity)[e_id] = viscosity_avg;
909 
910  Eigen::Map<KV>(
911  &(*process_data_.element_stresses)[e_id * KV::RowsAtCompileTime]) =
913 
915  ShapeFunction, typename ShapeFunctionDisplacement::MeshElement,
916  DisplacementDim>(element_, is_axially_symmetric_, p_L,
917  *process_data_.pressure_interpolated);
919  ShapeFunction, typename ShapeFunctionDisplacement::MeshElement,
920  DisplacementDim>(element_, is_axially_symmetric_, T,
921  *process_data_.temperature_interpolated);
922 }
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
Definition: KelvinVector.h:48
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(), ProcessLib::ThermoRichardsMechanics::ConstitutiveSetting< DisplacementDim >::eval(), NumLib::interpolateToHigherOrderNodes(), NumLib::interpolateXCoordinate(), MathLib::KelvinVector::kelvinVectorToSymmetricTensor(), ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MathLib::t, and MathLib::u.

◆ getEpsilon()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getEpsilon
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

626 {
627  constexpr int kelvin_vector_size =
629 
630  return transposeInPlace<kelvin_vector_size>(
631  [this](std::vector<double>& values)
632  { return getIntPtEpsilon(0, {}, {}, values); });
633 }
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().

◆ getIntPtDarcyVelocity()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtDarcyVelocity ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

659 {
660  unsigned const n_integration_points =
662 
663  cache.clear();
664  auto cache_matrix = MathLib::createZeroedMatrix<Eigen::Matrix<
665  double, DisplacementDim, Eigen::Dynamic, Eigen::RowMajor>>(
666  cache, DisplacementDim, n_integration_points);
667 
668  for (unsigned ip = 0; ip < n_integration_points; ip++)
669  {
670  cache_matrix.col(ip).noalias() = output_data_[ip].darcy_data.v_darcy;
671  }
672 
673  return cache;
674 }

References MathLib::createZeroedMatrix().

◆ getIntPtDryDensitySolid()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtDryDensitySolid ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

767 {
769  output_data_,
770  [](auto const& out) -> auto const& {
771  return out.rho_S_data.dry_density_solid;
772  },
773  cache);
774 }
std::vector< double > const & getIntegrationPointScalarData(IntegrationPointDataVector const &ip_data_vector, MemberType IpData::*const member, std::vector< double > &cache)

References ProcessLib::getIntegrationPointScalarData().

◆ getIntPtEpsilon()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtEpsilon ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

644 {
645  return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
647  [](auto const& cs) -> auto const& { return cs.eps_data.eps; }, cache);
648 }

◆ getIntPtLiquidDensity()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtLiquidDensity ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

785 {
787  output_data_,
788  [](auto const& out) -> auto const& { return out.rho_L_data.rho_LR; },
789  cache);
790 }

References ProcessLib::getIntegrationPointScalarData().

◆ getIntPtPorosity()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtPorosity ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

723 {
726  [](auto const& state) -> auto const& { return state.poro_data.phi; },
727  cache);
728 }

References ProcessLib::getIntegrationPointScalarData().

◆ getIntPtSaturation()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtSaturation ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

696 {
699  [](auto const& state) -> auto const& { return state.S_L_data.S_L; },
700  cache);
701 }

References ProcessLib::getIntegrationPointScalarData().

◆ getIntPtSigma()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtSigma ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

571 {
572  return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
574  [](auto const& cs) -> auto const& { return cs.s_mech_data.sigma_eff; },
575  cache);
576 }

◆ getIntPtSwellingStress()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtSwellingStress ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

601 {
602  constexpr int kelvin_vector_size =
604  auto const n_integration_points = ip_data_.size();
605 
606  cache.clear();
607  auto cache_mat = MathLib::createZeroedMatrix<Eigen::Matrix<
608  double, kelvin_vector_size, Eigen::Dynamic, Eigen::RowMajor>>(
609  cache, kelvin_vector_size, n_integration_points);
610 
611  for (unsigned ip = 0; ip < n_integration_points; ++ip)
612  {
613  auto const& sigma_sw = current_states_[ip].swelling_data.sigma_sw;
614  cache_mat.col(ip) =
616  }
617 
618  return cache;
619 }

References MathLib::createZeroedMatrix(), MathLib::KelvinVector::kelvin_vector_dimensions(), and MathLib::KelvinVector::kelvinVectorToSymmetricTensor().

◆ getIntPtTransportPorosity()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtTransportPorosity ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

750 {
752  current_states_, [](auto const& state) -> auto const& {
753  return state.transport_poro_data.phi;
754  },
755  cache);
756 }

References ProcessLib::getIntegrationPointScalarData().

◆ getIntPtViscosity()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getIntPtViscosity ( const double  t,
std::vector< GlobalVector * > const &  x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table,
std::vector< double > &  cache 
) const
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

801 {
803  output_data_,
804  [](auto const& out) -> auto const& { return out.mu_L_data.viscosity; },
805  cache);
806 }

References ProcessLib::getIntegrationPointScalarData().

◆ getMaterialStateVariablesAt()

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

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

940 {
941  return *material_states_[integration_point].material_state_variables;
942 }

◆ getNumberOfIntegrationPoints()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
unsigned ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getNumberOfIntegrationPoints
overrideprivatevirtual

◆ getPorosity()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getPorosity
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

708 {
709  std::vector<double> result;
710  getIntPtPorosity(0, {}, {}, result);
711  return result;
712 }
std::vector< double > const & getIntPtPorosity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override

◆ getSaturation()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getSaturation
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

681 {
682  std::vector<double> result;
683  getIntPtSaturation(0, {}, {}, result);
684  return result;
685 }
std::vector< double > const & getIntPtSaturation(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override

◆ getShapeMatrix()

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

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 325 of file ThermoRichardsMechanicsFEM.h.

327  {
328  auto const& N_u = secondary_data_.N_u[integration_point];
329 
330  // assumes N is stored contiguously in memory
331  return Eigen::Map<const Eigen::RowVectorXd>(N_u.data(), N_u.size());
332  }

References ProcessLib::HydroMechanics::SecondaryData< ShapeMatrixType >::N_u, and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::secondary_data_.

◆ getSigma()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getSigma
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

553 {
554  constexpr int kelvin_vector_size =
556 
557  return transposeInPlace<kelvin_vector_size>(
558  [this](std::vector<double>& values)
559  { return getIntPtSigma(0, {}, {}, values); });
560 }
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().

◆ getSwellingStress()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getSwellingStress
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

583 {
584  constexpr int kelvin_vector_size =
586 
587  return transposeInPlace<kelvin_vector_size>(
588  [this](std::vector<double>& values)
589  { return getIntPtSwellingStress(0, {}, {}, values); });
590 }
std::vector< double > const & getIntPtSwellingStress(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().

◆ getTransportPorosity()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector< double > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::getTransportPorosity
overridevirtual

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

735 {
736  std::vector<double> result;
737  getIntPtTransportPorosity(0, {}, {}, result);
738  return result;
739 }
std::vector< double > const & getIntPtTransportPorosity(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override

◆ initializeConcrete()

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

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 274 of file ThermoRichardsMechanicsFEM.h.

275  {
276  unsigned const n_integration_points =
278 
279  for (unsigned ip = 0; ip < n_integration_points; ip++)
280  {
281  // Set initial stress from parameter.
282  if (process_data_.initial_stress != nullptr)
283  {
284  ParameterLib::SpatialPosition const x_position{
285  std::nullopt, element_.getID(), ip,
287  ShapeFunctionDisplacement,
289  element_, ip_data_[ip].N_u))};
290 
291  current_states_[ip].s_mech_data.sigma_eff =
293  DisplacementDim>((*process_data_.initial_stress)(
294  std::numeric_limits<
295  double>::quiet_NaN() /* time independent */,
296  x_position));
297  }
298 
299  material_states_[ip].pushBackState();
300  }
301 
303  }
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
Definition: KelvinVector.h:170
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)

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

◆ localDOF()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
auto ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::localDOF ( std::vector< double > const &  local_dof_data) const
inlineprivate

Makes local d.o.f.s more accessible. Order of the returned Eigen vectors: T, p_L, u.

Definition at line 255 of file ThermoRichardsMechanicsFEM.h.

256  {
257  static_assert(temperature_size == pressure_size);
258 
259  using NodalTOrPVec =
260  typename ShapeMatricesType::template VectorType<temperature_size>;
261  using NodalDispVec =
262  typename ShapeMatricesTypeDisplacement::template VectorType<
264 
265  return std::tuple<Eigen::Map<NodalTOrPVec const>,
266  Eigen::Map<NodalTOrPVec const>,
267  Eigen::Map<NodalDispVec const>>(
268  {local_dof_data.data() + temperature_index, temperature_size},
269  {local_dof_data.data() + pressure_index, pressure_size},
270  {local_dof_data.data() + displacement_index, displacement_size});
271  };

References ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::displacement_index, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::displacement_size, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::pressure_index, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::pressure_size, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::temperature_index, and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::temperature_size.

◆ massLumping()

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

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

324 {
325  if (process_data_.apply_mass_lumping)
326  {
327  loc_mat.storage_p_a_p =
328  loc_mat.storage_p_a_p.colwise().sum().eval().asDiagonal();
329  loc_mat.storage_p_a_S =
330  loc_mat.storage_p_a_S.colwise().sum().eval().asDiagonal();
331  loc_mat.storage_p_a_S_Jpp =
332  loc_mat.storage_p_a_S_Jpp.colwise().sum().eval().asDiagonal();
333  }
334 }

References ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_p, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_S, and ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::LocalMatrices::storage_p_a_S_Jpp.

◆ postTimestepConcrete()

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

◆ setInitialConditionsConcrete()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::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 206 of file ThermoRichardsMechanicsFEM-impl.h.

211 {
212  assert(local_x.size() ==
214 
215  auto const p_L = Eigen::Map<
216  typename ShapeMatricesType::template VectorType<pressure_size> const>(
217  local_x.data() + pressure_index, pressure_size);
218 
219  auto const T = Eigen::Map<typename ShapeMatricesType::template VectorType<
220  temperature_size> const>(local_x.data() + temperature_index,
222 
223  constexpr double dt = std::numeric_limits<double>::quiet_NaN();
224  auto const& medium = process_data_.media_map->getMedium(element_.getID());
225  MPL::VariableArray variables;
226 
228  x_position.setElementID(element_.getID());
229 
230  auto const& solid_phase = medium->phase("Solid");
231 
232  unsigned const n_integration_points =
234  for (unsigned ip = 0; ip < n_integration_points; ip++)
235  {
236  x_position.setIntegrationPoint(ip);
237 
238  // N is used for both T and p variables.
239  auto const& N = ip_data_[ip].N_p;
240 
241  double p_cap_ip;
242  NumLib::shapeFunctionInterpolate(-p_L, N, p_cap_ip);
243 
244  variables[static_cast<int>(MPL::Variable::capillary_pressure)] =
245  p_cap_ip;
246  variables[static_cast<int>(MPL::Variable::phase_pressure)] = -p_cap_ip;
247 
248  double T_ip;
250  variables[static_cast<int>(MPL::Variable::temperature)] = T_ip;
251 
252  prev_states_[ip].S_L_data.S_L =
253  medium->property(MPL::PropertyType::saturation)
254  .template value<double>(variables, x_position, t, dt);
255 
256  // Set eps_m_prev from potentially non-zero eps and sigma_sw from
257  // restart.
258  SpaceTimeData const x_t{x_position, t, dt};
259  ElasticTangentStiffnessData<DisplacementDim> C_el_data;
260  ElasticTangentStiffnessModel<DisplacementDim>{solid_material_}.eval(
261  x_t, {T_ip, 0, {}}, C_el_data);
262 
263  auto const& eps = current_states_[ip].eps_data.eps;
264  auto const& sigma_sw = current_states_[ip].swelling_data.sigma_sw;
265  prev_states_[ip].s_mech_data.eps_m.noalias() =
266  solid_phase.hasProperty(MPL::PropertyType::swelling_stress_rate)
267  ? eps + C_el_data.C_el.inverse() * sigma_sw
268  : eps;
269  }
270 }
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
Definition: VariableType.h:110
void shapeFunctionInterpolate(const NodalValues &nodal_values, const ShapeMatrix &shape_matrix_N, double &interpolated_value, ScalarTypes &... interpolated_values)
Definition: Interpolation.h:80

References ProcessLib::ThermoRichardsMechanics::ElasticTangentStiffnessData< DisplacementDim >::C_el, MaterialPropertyLib::capillary_pressure, ProcessLib::ThermoRichardsMechanics::ElasticTangentStiffnessModel< DisplacementDim >::eval(), MaterialPropertyLib::saturation, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), NumLib::shapeFunctionInterpolate(), MaterialPropertyLib::swelling_stress_rate, MathLib::t, and MaterialPropertyLib::temperature.

◆ setIPDataInitialConditions()

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::size_t ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::setIPDataInitialConditions ( std::string const &  name,
double const *  values,
int const  integration_order 
)
overridevirtual
Returns
the number of read integration points.

Implements ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >.

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

144 {
145  if (integration_order !=
146  static_cast<int>(integration_method_.getIntegrationOrder()))
147  {
148  OGS_FATAL(
149  "Setting integration point initial conditions; The integration "
150  "order of the local assembler for element {:d} is different "
151  "from the integration order in the initial condition.",
152  element_.getID());
153  }
154 
155  if (name == "sigma_ip")
156  {
157  if (process_data_.initial_stress != nullptr)
158  {
159  OGS_FATAL(
160  "Setting initial conditions for stress from integration "
161  "point data and from a parameter '{:s}' is not possible "
162  "simultaneously.",
163  process_data_.initial_stress->name);
164  }
165  return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
166  values, current_states_,
167  [](auto& cs) -> auto& { return cs.s_mech_data.sigma_eff; });
168  }
169 
170  if (name == "saturation_ip")
171  {
173  values, current_states_,
174  [](auto& state) -> auto& { return state.S_L_data.S_L; });
175  }
176  if (name == "porosity_ip")
177  {
179  values, current_states_,
180  [](auto& state) -> auto& { return state.poro_data.phi; });
181  }
182  if (name == "transport_porosity_ip")
183  {
185  values, current_states_,
186  [](auto& state) -> auto& { return state.transport_poro_data.phi; });
187  }
188  if (name == "swelling_stress_ip")
189  {
190  return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
191  values, current_states_,
192  [](auto& cs) -> auto& { return cs.swelling_data.sigma_sw; });
193  }
194  if (name == "epsilon_ip")
195  {
196  return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
197  values, current_states_,
198  [](auto& cs) -> auto& { return cs.eps_data.eps; });
199  }
200  return 0;
201 }
#define OGS_FATAL(...)
Definition: Error.h:26
std::size_t setIntegrationPointScalarData(double const *values, IntegrationPointDataVector &ip_data_vector, MemberType IpData::*const member)

References MaterialPropertyLib::name, OGS_FATAL, and ProcessLib::setIntegrationPointScalarData().

Member Data Documentation

◆ current_states_

◆ displacement_index

◆ displacement_size

◆ element_

◆ integration_method_

◆ ip_data_

◆ is_axially_symmetric_

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
bool const ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::is_axially_symmetric_
private

Definition at line 424 of file ThermoRichardsMechanicsFEM.h.

◆ KelvinVectorSize

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

Definition at line 78 of file ThermoRichardsMechanicsFEM.h.

◆ material_states_

◆ output_data_

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
std::vector<OutputData<DisplacementDim> > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::output_data_
private

◆ pressure_index

◆ pressure_size

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

◆ prev_states_

◆ process_data_

◆ secondary_data_

◆ solid_material_

template<typename ShapeFunctionDisplacement , typename ShapeFunction , int DisplacementDim>
MaterialLib::Solids::MechanicsBase<DisplacementDim> const& ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsLocalAssembler< ShapeFunctionDisplacement, ShapeFunction, DisplacementDim >::solid_material_
private

◆ temperature_index

◆ temperature_size

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

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