OGS
ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits > Class Template Referencefinal

Detailed Description

template<int DisplacementDim, typename ConstitutiveTraits>
class ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >

Global assembler for the monolithic scheme of the non-isothermal Richards flow coupled with mechanics.

Governing equations without vapor diffusion

The energy balance equation is given by

\[ (\rho c_p)^{eff}\dot T - \nabla (\mathbf{k}_T^{eff} \nabla T)+\rho^l c_p^l \nabla T \cdot \mathbf{v}^l = Q_T \]

with \(T\) the temperature, \((\rho c_p)^{eff}\) the effective volumetric heat capacity, \(\mathbf{k}_T^{eff} \) the effective thermal conductivity, \(\rho^l\) the density of liquid, \(c_p^l\) the specific heat capacity of liquid, \(\mathbf{v}^l\) the liquid velocity, and \(Q_T\) the point heat source. The effective volumetric heat can be considered as a composite of the contributions of solid phase and the liquid phase as

\[ (\rho c_p)^{eff} = (1-\phi) \rho^s c_p^s + S^l \phi \rho^l c_p^l \]

with \(\phi\) the porosity, \(S^l\) the liquid saturation, \(\rho^s \) the solid density, and \(c_p^s\) the specific heat capacity of solid. Similarly, the effective thermal conductivity is given by

\[ \mathbf{k}_T^{eff} = (1-\phi) \mathbf{k}_T^s + S^l \phi k_T^l \mathbf I \]

where \(\mathbf{k}_T^s\) is the thermal conductivity tensor of solid, \( k_T^l\) is the thermal conductivity of liquid, and \(\mathbf I\) is the identity tensor.

The mass balance equation is given by

\begin{eqnarray*} \left(S^l\beta - \phi\frac{\partial S}{\partial p_c}\right) \rho^l\dot p - S \left( \frac{\partial \rho^l}{\partial T} +\rho^l(\alpha_B -S) \alpha_T^s \right)\dot T\\ +\nabla (\rho^l \mathbf{v}^l) + S \alpha_B \rho^l \nabla \cdot \dot {\mathbf u}= Q_H \end{eqnarray*}

where \(p\) is the pore pressure, \(p_c\) is the capillary pressure, which is \(-p\) under the single phase assumption, \(\beta\) is a composite coefficient by the liquid compressibility and solid compressibility, \(\alpha_B\) is the Biot's constant, \(\alpha_T^s\) is the linear thermal expansivity of solid, \(Q_H\) is the point source or sink term, \( \mathbf u\) is the displacement, and \(H(S-1)\) is the Heaviside function. The liquid velocity \(\mathbf{v}^l\) is described by the Darcy's law as

\[ \mathbf{v}^l=-\frac{{\mathbf k} k_{ref}}{\mu} (\nabla p - \rho^l \mathbf g) \]

with \({\mathbf k}\) the intrinsic permeability, \(k_{ref}\) the relative permeability, \(\mathbf g\) the gravitational force.

The momentum balance equation takes the form of

\[ \nabla (\mathbf{\sigma}-b(S)\alpha_B p^l \mathbf I) +\mathbf f=0 \]

with \(\mathbf{\sigma}\) the effective stress tensor, \(b(S)\) the Bishop model, \(\mathbf f\) the body force, and \(\mathbf I\) the identity. The primary unknowns of the momentum balance equation are the displacement \(\mathbf u\), which is associated with the stress by the the generalized Hook's law as

\[ {\dot {\mathbf {\sigma}}} = C {\dot {\mathbf \epsilon}}^e = C ( {\dot {\mathbf \epsilon}} - {\dot {\mathbf \epsilon}}^T -{\dot {\mathbf \epsilon}}^p - {\dot {\mathbf \epsilon}}^{sw}-\cdots ) \]

with \(C\) the forth order elastic tensor, \({\dot {\mathbf \epsilon}}\) the total strain rate, \({\dot {\mathbf \epsilon}}^e\) the elastic strain rate, \({\dot {\mathbf \epsilon}}^T\) the thermal strain rate, \({\dot {\mathbf \epsilon}}^p\) the plastic strain rate, \({\dot {\mathbf \epsilon}}^{sw}\) the swelling strain rate.

The strain tensor is given by displacement vector as

\[ \mathbf \epsilon = \frac{1}{2} \left((\nabla \mathbf u)^{\text T}+\nabla \mathbf u\right) \]

where the superscript \({\text T}\) means transpose,

Definition at line 114 of file ThermoRichardsMechanicsProcess.h.

#include <ThermoRichardsMechanicsProcess.h>

Inheritance diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >:
[legend]
Collaboration diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >:
[legend]

Public Member Functions

 ThermoRichardsMechanicsProcess (std::string name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > &&process_variables, ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)
 
MathLib::MatrixSpecifications getMatrixSpecifications (const int process_id) const override
 
ODESystem interface
bool isLinear () const override
 
- Public Member Functions inherited from ProcessLib::Process
 Process (std::string name_, MeshLib::Mesh &mesh, std::unique_ptr< AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > &&process_variables, SecondaryVariableCollection &&secondary_variables, const bool use_monolithic_scheme=true)
 
void preTimestep (std::vector< GlobalVector * > const &x, const double t, const double delta_t, const int process_id)
 Preprocessing before starting assembly for new timestep.
 
void postTimestep (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, const double delta_t, int const process_id)
 Postprocessing after a complete timestep.
 
void postNonLinearSolver (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, double const dt, int const process_id)
 
void preIteration (const unsigned iter, GlobalVector const &x) final
 
void computeSecondaryVariable (double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
 compute secondary variables for the coupled equations or for output.
 
NumLib::IterationResult postIteration (GlobalVector const &x) final
 
void initialize (std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
 
void setInitialConditions (std::vector< GlobalVector * > &process_solutions, std::vector< GlobalVector * > const &process_solutions_prev, double const t, int const process_id)
 
MathLib::MatrixSpecifications getMatrixSpecifications (const int process_id) const override
 
void updateDeactivatedSubdomains (double const time, const int process_id)
 
virtual bool isMonolithicSchemeUsed () const
 
virtual void extrapolateIntegrationPointValuesToNodes (const double, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > &)
 
void preAssemble (const double t, double const dt, GlobalVector const &x) final
 
void assemble (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) final
 
void assembleWithJacobian (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) final
 
void preOutput (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id)
 
std::vector< NumLib::IndexValueVector< GlobalIndexType > > const * getKnownSolutions (double const t, GlobalVector const &x, int const process_id) const final
 
MeshLib::MeshgetMesh () const
 
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables (const int process_id) const
 
std::vector< std::size_t > const & getActiveElementIDs () const
 
SecondaryVariableCollection const & getSecondaryVariables () const
 
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > const & getIntegrationPointWriters () const
 
virtual Eigen::Vector3d getFlux (std::size_t, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &) const
 
virtual void solveReactionEquation (std::vector< GlobalVector * > &, std::vector< GlobalVector * > const &, double const, double const, NumLib::EquationSystem &, int const)
 
- Public Member Functions inherited from ProcessLib::SubmeshAssemblySupport
virtual ~SubmeshAssemblySupport ()=default
 

Private Types

using LocalAssemblerIF
 

Private Member Functions

void constructDofTable () override
 
void initializeConcreteProcess (NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
 Process specific initialization called by initialize().
 
void initializeBoundaryConditions (std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media) override
 
void setInitialConditionsConcreteProcess (std::vector< GlobalVector * > &x, double const t, int const process_id) override
 
void assembleConcreteProcess (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
 
void assembleWithJacobianConcreteProcess (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
 
void preTimestepConcreteProcess (std::vector< GlobalVector * > const &, const double, const double, const int) override
 
void postTimestepConcreteProcess (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, const int process_id) override
 
NumLib::LocalToGlobalIndexMap const & getDOFTable (const int process_id) const override
 
std::vector< std::string > initializeAssemblyOnSubmeshes (std::vector< std::reference_wrapper< MeshLib::Mesh > > const &meshes) override
 
void computeSecondaryVariableConcrete (double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id) override
 
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData () const override
 
bool hasMechanicalProcess (int const) const
 
- Private Member Functions inherited from ProcessLib::AssemblyMixin< ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits > >
void initializeAssemblyOnSubmeshes (const int process_id, std::vector< std::reference_wrapper< MeshLib::Mesh > > const &submeshes, std::vector< std::string > const &residuum_names)
 
void updateActiveElements (const int process_id)
 
void assemble (const double, double const, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, int const, GlobalMatrix &, GlobalMatrix &, GlobalVector &)
 
void assembleWithJacobian (const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)
 

Private Attributes

std::vector< MeshLib::Node * > base_nodes_
 
std::unique_ptr< MeshLib::MeshSubset const > mesh_subset_base_nodes_
 
ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > process_data_
 
std::vector< std::unique_ptr< LocalAssemblerIF > > local_assemblers_
 
std::unique_ptr< NumLib::LocalToGlobalIndexMaplocal_to_global_index_map_single_component_
 
std::unique_ptr< NumLib::LocalToGlobalIndexMaplocal_to_global_index_map_with_base_nodes_
 
GlobalSparsityPattern sparsity_pattern_with_linear_element_
 

Friends

class AssemblyMixin< ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits > >
 

Additional Inherited Members

- Public Attributes inherited from ProcessLib::Process
std::string const name
 
- Static Public Attributes inherited from ProcessLib::Process
static PROCESSLIB_EXPORT const std::string constant_one_parameter_name = "constant_one"
 
- Protected Member Functions inherited from ProcessLib::Process
NumLib::ExtrapolatorgetExtrapolator () const
 
NumLib::LocalToGlobalIndexMap const & getSingleComponentDOFTable () const
 
void initializeProcessBoundaryConditionsAndSourceTerms (const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
 
void constructMonolithicProcessDofTable ()
 
void constructDofTableOfSpecifiedProcessStaggeredScheme (const int specified_process_id)
 
std::vector< GlobalIndexTypegetIndicesOfResiduumWithoutInitialCompensation () const override
 
- Protected Attributes inherited from ProcessLib::Process
MeshLib::Mesh_mesh
 
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
 
std::unique_ptr< NumLib::LocalToGlobalIndexMap_local_to_global_index_map
 
SecondaryVariableCollection _secondary_variables
 
std::unique_ptr< ProcessLib::AbstractJacobianAssembler_jacobian_assembler
 
VectorMatrixAssembler _global_assembler
 
const bool _use_monolithic_scheme
 
unsigned const _integration_order
 
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
 
GlobalSparsityPattern _sparsity_pattern
 
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
 
std::vector< BoundaryConditionCollection_boundary_conditions
 

Member Typedef Documentation

◆ LocalAssemblerIF

template<int DisplacementDim, typename ConstitutiveTraits >
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::LocalAssemblerIF
private
Initial value:
LocalAssemblerInterface<DisplacementDim, ConstitutiveTraits>

Definition at line 160 of file ThermoRichardsMechanicsProcess.h.

Constructor & Destructor Documentation

◆ ThermoRichardsMechanicsProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::ThermoRichardsMechanicsProcess ( std::string name,
MeshLib::Mesh & mesh,
std::unique_ptr< ProcessLib::AbstractJacobianAssembler > && jacobian_assembler,
std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const & parameters,
unsigned const integration_order,
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > && process_variables,
ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > && process_data,
SecondaryVariableCollection && secondary_variables,
bool const use_monolithic_scheme )

Definition at line 31 of file ThermoRichardsMechanicsProcess.cpp.

46 : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
47 integration_order, std::move(process_variables),
48 std::move(secondary_variables), use_monolithic_scheme),
49 AssemblyMixin<
50 ThermoRichardsMechanicsProcess<DisplacementDim, ConstitutiveTraits>>{
52 process_data_(std::move(process_data))
53{
56 _integration_point_writer, integration_order,
58}
std::string const name
Definition Process.h:351
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
Definition Process.h:377
std::unique_ptr< ProcessLib::AbstractJacobianAssembler > _jacobian_assembler
Definition Process.h:363
ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > process_data_
void addReflectedIntegrationPointWriters(ReflData const &reflection_data, std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > &integration_point_writers, unsigned const integration_order, std::vector< std::unique_ptr< LocAsmIF > > const &local_assemblers)

References ProcessLib::Process::_jacobian_assembler.

Member Function Documentation

◆ assembleConcreteProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::assembleConcreteProcess ( const double t,
double const dt,
std::vector< GlobalVector * > const & x,
std::vector< GlobalVector * > const & x_prev,
int const process_id,
GlobalMatrix & M,
GlobalMatrix & K,
GlobalVector & b )
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 229 of file ThermoRichardsMechanicsProcess.cpp.

235{
236 OGS_FATAL(
237 "The Picard method or the Newton-Raphson method with numerical "
238 "Jacobian is not implemented for ThermoRichardsMechanics with the full "
239 "monolithic coupling scheme");
240}
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

◆ assembleWithJacobianConcreteProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::assembleWithJacobianConcreteProcess ( const double t,
double const dt,
std::vector< GlobalVector * > const & x,
std::vector< GlobalVector * > const & x_prev,
int const process_id,
GlobalMatrix & M,
GlobalMatrix & K,
GlobalVector & b,
GlobalMatrix & Jac )
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 243 of file ThermoRichardsMechanicsProcess.cpp.

248{
249 AssemblyMixin<ThermoRichardsMechanicsProcess<
250 DisplacementDim, ConstitutiveTraits>>::assembleWithJacobian(t, dt, x,
251 x_prev,
252 process_id,
253 M, K, b,
254 Jac);
255}
void assembleWithJacobian(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) final
Definition Process.cpp:273
ThermoRichardsMechanicsProcess(std::string name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > &&process_variables, ThermoRichardsMechanicsProcessData< DisplacementDim, ConstitutiveTraits > &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)

◆ computeSecondaryVariableConcrete()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::computeSecondaryVariableConcrete ( double const t,
double const dt,
std::vector< GlobalVector * > const & x,
GlobalVector const & x_prev,
int const process_id )
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 308 of file ThermoRichardsMechanicsProcess.cpp.

313{
314 DBUG("Compute the secondary variables for ThermoRichardsMechanicsProcess.");
315
316 auto get_a_dof_table_func = [this](const int processe_id) -> auto&
317 {
318 return getDOFTable(processe_id);
319 };
320 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
324 NumLib::getDOFTables(x.size(), get_a_dof_table_func), t, dt, x, x_prev,
325 process_id);
326}
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
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)
std::vector< std::size_t > const & getActiveElementIDs() const
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition Process.h:155
NumLib::LocalToGlobalIndexMap const & getDOFTable(const int process_id) const override
std::vector< NumLib::LocalToGlobalIndexMap const * > getDOFTables(int const number_of_processes, std::function< NumLib::LocalToGlobalIndexMap const &(const int)> get_single_dof_table)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)

References DBUG(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and NumLib::getDOFTables().

◆ constructDofTable()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::constructDofTable ( )
overrideprivatevirtual

This function is for general cases, in which all equations of the coupled processes have the same number of unknowns. For the general cases with the staggered scheme, all equations of the coupled processes share one DOF table hold by _local_to_global_index_map. Other cases can be considered by overloading this member function in the derived class.

Reimplemented from ProcessLib::Process.

Definition at line 79 of file ThermoRichardsMechanicsProcess.cpp.

80{
81 // Create single component dof in every of the mesh's nodes.
82 _mesh_subset_all_nodes = std::make_unique<MeshLib::MeshSubset>(
83 _mesh, _mesh.getNodes(), process_data_.use_TaylorHood_elements);
84 // Create single component dof in the mesh's base nodes.
86 mesh_subset_base_nodes_ = std::make_unique<MeshLib::MeshSubset>(
87 _mesh, base_nodes_, process_data_.use_TaylorHood_elements);
88
89 // TODO move the two data members somewhere else.
90 // for extrapolation of secondary variables of stress or strain
91 std::vector<MeshLib::MeshSubset> all__meshsubsets_single_component{
94 std::make_unique<NumLib::LocalToGlobalIndexMap>(
95 std::move(all__meshsubsets_single_component),
96 // by location order is needed for output
98
99 // For temperature, which is the first variable.
100 std::vector<MeshLib::MeshSubset> all__meshsubsets{*mesh_subset_base_nodes_};
101
102 // For pressure, which is the second variable
103 all__meshsubsets.push_back(*mesh_subset_base_nodes_);
104
105 // For displacement.
106 const int monolithic_process_id = 0;
107 std::generate_n(std::back_inserter(all__meshsubsets),
108 getProcessVariables(monolithic_process_id)[2]
109 .get()
110 .getNumberOfGlobalComponents(),
111 [&]() { return *_mesh_subset_all_nodes; });
112
113 std::vector<int> const vec_n_components{1, 1, DisplacementDim};
115 std::make_unique<NumLib::LocalToGlobalIndexMap>(
116 std::move(all__meshsubsets), vec_n_components,
119}
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:106
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:109
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition Process.h:355
MeshLib::Mesh & _mesh
Definition Process.h:354
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition Process.h:357
std::unique_ptr< NumLib::LocalToGlobalIndexMap > local_to_global_index_map_single_component_
std::vector< Node * > getBaseNodes(std::vector< Element * > const &elements)
Definition Utils.h:26
@ BY_LOCATION
Ordering data by spatial location.
auto & get(Tuples &... ts)
Definition Get.h:67

References NumLib::BY_LOCATION, and MeshLib::getBaseNodes().

◆ getDOFTable()

template<int DisplacementDim, typename ConstitutiveTraits >
NumLib::LocalToGlobalIndexMap const & ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::getDOFTable ( const int process_id) const
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 339 of file ThermoRichardsMechanicsProcess.cpp.

341{
343}

◆ getDOFTableForExtrapolatorData()

template<int DisplacementDim, typename ConstitutiveTraits >
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::getDOFTableForExtrapolatorData ( ) const
overrideprivatevirtual

Get the address of a LocalToGlobalIndexMap, and the status of its memory. If the LocalToGlobalIndexMap is created as new in this function, the function also returns a true boolean value to let Extrapolator manage the memory by the address returned by this function.

Returns
Address of a LocalToGlobalIndexMap and its memory status.

Reimplemented from ProcessLib::Process.

Definition at line 330 of file ThermoRichardsMechanicsProcess.cpp.

331{
332 const bool manage_storage = false;
333 return std::make_tuple(local_to_global_index_map_single_component_.get(),
334 manage_storage);
335}

◆ getMatrixSpecifications()

template<int DisplacementDim, typename ConstitutiveTraits >
MathLib::MatrixSpecifications ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::getMatrixSpecifications ( const int process_id) const
override

Get the size and the sparse pattern of the global matrix in order to create the global matrices and vectors for the system equations of this process.

Parameters
process_idProcess ID. If the monolithic scheme is applied, process_id = 0. For the staggered scheme, process_id = 0 represents the hydraulic (H) process, while process_id = 1 represents the mechanical (M) process.
Returns
Matrix specifications including size and sparse pattern.

Definition at line 70 of file ThermoRichardsMechanicsProcess.cpp.

71{
72 auto const& l = *_local_to_global_index_map;
73 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
74 &l.getGhostIndices(), &this->_sparsity_pattern};
75}
GlobalSparsityPattern _sparsity_pattern
Definition Process.h:379

◆ hasMechanicalProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
bool ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::hasMechanicalProcess ( int const ) const
inlineprivate

Check whether the process represented by process_id is/has mechanical process. In the present implementation, the mechanical process has process_id == 1 in the staggered scheme.

Definition at line 240 of file ThermoRichardsMechanicsProcess.h.

240{ return true; }

◆ initializeAssemblyOnSubmeshes()

template<int DisplacementDim, typename ConstitutiveTraits >
std::vector< std::string > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::initializeAssemblyOnSubmeshes ( std::vector< std::reference_wrapper< MeshLib::Mesh > > const & meshes)
overrideprivatevirtual

Initializes the assembly on submeshes

Parameters
meshesthe submeshes on whom the assembly shall proceed.
Attention
meshes must be a must be a non-overlapping cover of the entire simulation domain (bulk mesh)!
Returns
The names of the residuum vectors that will be assembled.

Reimplemented from ProcessLib::SubmeshAssemblySupport.

Definition at line 270 of file ThermoRichardsMechanicsProcess.cpp.

273{
274 INFO("TRM process initializeSubmeshOutput().");
275 const int process_id = 0;
276 std::vector<std::string> residuum_names{"HeatFlowRate", "MassFlowRate",
277 "NodalForces"};
278
279 AssemblyMixin<
280 ThermoRichardsMechanicsProcess<DisplacementDim, ConstitutiveTraits>>::
281 initializeAssemblyOnSubmeshes(process_id, meshes, residuum_names);
282
283 return residuum_names;
284}
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
std::vector< std::string > initializeAssemblyOnSubmeshes(std::vector< std::reference_wrapper< MeshLib::Mesh > > const &meshes) override

References INFO().

◆ initializeBoundaryConditions()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::initializeBoundaryConditions ( std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const & media)
overrideprivatevirtual

Member function to initialize the boundary conditions for all coupled processes. It is called by initialize().

Reimplemented from ProcessLib::Process.

Definition at line 201 of file ThermoRichardsMechanicsProcess.cpp.

205{
206 const int process_id = 0;
208 *_local_to_global_index_map, process_id, media);
209}
void initializeProcessBoundaryConditionsAndSourceTerms(const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
Definition Process.cpp:73

◆ initializeConcreteProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::initializeConcreteProcess ( NumLib::LocalToGlobalIndexMap const & dof_table,
MeshLib::Mesh const & mesh,
unsigned const integration_order )
overrideprivatevirtual

Process specific initialization called by initialize().

Implements ProcessLib::Process.

Definition at line 122 of file ThermoRichardsMechanicsProcess.cpp.

126{
127 createLocalAssemblers<DisplacementDim, ConstitutiveTraits>(
128 mesh.getElements(), dof_table, local_assemblers_,
129 NumLib::IntegrationOrder{integration_order}, mesh.isAxiallySymmetric(),
131
132 auto add_secondary_variable = [&](std::string const& name,
133 int const num_components,
134 auto get_ip_values_function)
135 {
137 name,
138 makeExtrapolator(num_components, getExtrapolator(),
140 std::move(get_ip_values_function)));
141 };
142
143 ProcessLib::Reflection::addReflectedSecondaryVariables<DisplacementDim>(
146
147 //
148 // enable output of internal variables defined by material models
149 //
151 LocalAssemblerIF>(process_data_.solid_materials,
152 add_secondary_variable);
153
156 process_data_.solid_materials, local_assemblers_,
157 _integration_point_writer, integration_order);
158
159 process_data_.element_saturation = MeshLib::getOrCreateMeshProperty<double>(
160 const_cast<MeshLib::Mesh&>(mesh), "saturation_avg",
162
163 process_data_.element_porosity = MeshLib::getOrCreateMeshProperty<double>(
164 const_cast<MeshLib::Mesh&>(mesh), "porosity_avg",
166
167 process_data_.element_liquid_density =
168 MeshLib::getOrCreateMeshProperty<double>(
169 const_cast<MeshLib::Mesh&>(mesh), "liquid_density_avg",
171
172 process_data_.element_viscosity = MeshLib::getOrCreateMeshProperty<double>(
173 const_cast<MeshLib::Mesh&>(mesh), "viscosity_avg",
175
176 process_data_.element_stresses = MeshLib::getOrCreateMeshProperty<double>(
177 const_cast<MeshLib::Mesh&>(mesh), "stress_avg",
180 DisplacementDim>::RowsAtCompileTime);
181
182 process_data_.pressure_interpolated =
183 MeshLib::getOrCreateMeshProperty<double>(
184 const_cast<MeshLib::Mesh&>(mesh), "pressure_interpolated",
186 process_data_.temperature_interpolated =
187 MeshLib::getOrCreateMeshProperty<double>(
188 const_cast<MeshLib::Mesh&>(mesh), "temperature_interpolated",
190
193
194 // Initialize local assemblers after all variables have been set.
198}
virtual void initialize(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
SecondaryVariableCollection _secondary_variables
Definition Process.h:359
NumLib::Extrapolator & getExtrapolator() const
Definition Process.h:196
void addSecondaryVariable(std::string const &internal_name, SecondaryVariableFunctions &&fcts)
LocalAssemblerInterface< DisplacementDim, ConstitutiveTraits > LocalAssemblerIF
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
void solidMaterialInternalVariablesToIntegrationPointWriter(std::map< int, std::unique_ptr< SolidMaterial > > const &solid_materials, std::vector< std::unique_ptr< LocalAssemblerInterface > > const &local_assemblers, std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > &integration_point_writer, int const integration_order)
void solidMaterialInternalToSecondaryVariables(std::map< int, std::unique_ptr< SolidMaterial > > const &solid_materials, AddSecondaryVariableCallback const &add_secondary_variable)
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
void setIPDataInitialConditions(std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > const &_integration_point_writer, MeshLib::Properties const &mesh_properties, LocalAssemblersVector &local_assemblers)
static void executeMemberOnDereferenced(Method method, Container const &container, Args &&... args)

References MeshLib::Cell, NumLib::SerialExecutor::executeMemberOnDereferenced(), MeshLib::Mesh::getElements(), MeshLib::Mesh::getProperties(), MeshLib::Mesh::isAxiallySymmetric(), ProcessLib::makeExtrapolator(), MeshLib::Node, ProcessLib::setIPDataInitialConditions(), ProcessLib::Deformation::solidMaterialInternalToSecondaryVariables(), and ProcessLib::Deformation::solidMaterialInternalVariablesToIntegrationPointWriter().

◆ isLinear()

template<int DisplacementDim, typename ConstitutiveTraits >
bool ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::isLinear ( ) const
override

Definition at line 62 of file ThermoRichardsMechanicsProcess.cpp.

63{
64 return false;
65}

◆ postTimestepConcreteProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::postTimestepConcreteProcess ( std::vector< GlobalVector * > const & x,
std::vector< GlobalVector * > const & x_prev,
double const t,
double const dt,
const int process_id )
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 287 of file ThermoRichardsMechanicsProcess.cpp.

292{
293 DBUG("PostTimestep ThermoRichardsMechanicsProcess.");
294
295 auto get_a_dof_table_func = [this](const int processe_id) -> auto&
296 {
297 return getDOFTable(processe_id);
298 };
299 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
303 NumLib::getDOFTables(x.size(), get_a_dof_table_func), x, x_prev, t, dt,
304 process_id);
305}
virtual void postTimestep(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)

References DBUG(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and NumLib::getDOFTables().

◆ preTimestepConcreteProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::preTimestepConcreteProcess ( std::vector< GlobalVector * > const & ,
const double ,
const double ,
const int process_id )
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 258 of file ThermoRichardsMechanicsProcess.cpp.

263{
264 AssemblyMixin<ThermoRichardsMechanicsProcess<
265 DisplacementDim, ConstitutiveTraits>>::updateActiveElements(process_id);
266}

◆ setInitialConditionsConcreteProcess()

template<int DisplacementDim, typename ConstitutiveTraits >
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::setInitialConditionsConcreteProcess ( std::vector< GlobalVector * > & x,
double const t,
int const process_id )
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 212 of file ThermoRichardsMechanicsProcess.cpp.

216{
217 DBUG("SetInitialConditions ThermoRichardsMechanicsProcess.");
218
219 auto get_a_dof_table_func = [this](const int process_id) -> auto&
220 {
221 return getDOFTable(process_id);
222 };
225 NumLib::getDOFTables(x.size(), get_a_dof_table_func), x, t, process_id);
226}
virtual void setInitialConditions(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, int const process_id)

References DBUG(), NumLib::SerialExecutor::executeMemberOnDereferenced(), and NumLib::getDOFTables().

Friends And Related Symbol Documentation

◆ AssemblyMixin< ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits > >

template<int DisplacementDim, typename ConstitutiveTraits >
friend class AssemblyMixin< ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits > >
friend

Definition at line 252 of file ThermoRichardsMechanicsProcess.h.

Member Data Documentation

◆ base_nodes_

template<int DisplacementDim, typename ConstitutiveTraits >
std::vector<MeshLib::Node*> ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::base_nodes_
private

Definition at line 208 of file ThermoRichardsMechanicsProcess.h.

◆ local_assemblers_

template<int DisplacementDim, typename ConstitutiveTraits >
std::vector<std::unique_ptr<LocalAssemblerIF> > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::local_assemblers_
private

Definition at line 213 of file ThermoRichardsMechanicsProcess.h.

◆ local_to_global_index_map_single_component_

template<int DisplacementDim, typename ConstitutiveTraits >
std::unique_ptr<NumLib::LocalToGlobalIndexMap> ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::local_to_global_index_map_single_component_
private

Definition at line 216 of file ThermoRichardsMechanicsProcess.h.

◆ local_to_global_index_map_with_base_nodes_

template<int DisplacementDim, typename ConstitutiveTraits >
std::unique_ptr<NumLib::LocalToGlobalIndexMap> ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::local_to_global_index_map_with_base_nodes_
private

Local to global index mapping for base nodes, which is used for linear interpolation for pressure in the staggered scheme.

Definition at line 221 of file ThermoRichardsMechanicsProcess.h.

◆ mesh_subset_base_nodes_

template<int DisplacementDim, typename ConstitutiveTraits >
std::unique_ptr<MeshLib::MeshSubset const> ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::mesh_subset_base_nodes_
private

Definition at line 209 of file ThermoRichardsMechanicsProcess.h.

◆ process_data_

template<int DisplacementDim, typename ConstitutiveTraits >
ThermoRichardsMechanicsProcessData<DisplacementDim, ConstitutiveTraits> ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::process_data_
private

Definition at line 211 of file ThermoRichardsMechanicsProcess.h.

◆ sparsity_pattern_with_linear_element_

template<int DisplacementDim, typename ConstitutiveTraits >
GlobalSparsityPattern ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim, ConstitutiveTraits >::sparsity_pattern_with_linear_element_
private

Sparsity pattern for the flow equation, and it is initialized only if the staggered scheme is used.

Definition at line 225 of file ThermoRichardsMechanicsProcess.h.


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