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

Detailed Description

template<int DisplacementDim>
class ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >

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 109 of file ThermoRichardsMechanicsProcess.h.

#include <ThermoRichardsMechanicsProcess.h>

Inheritance diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >:
[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 > &&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. More...
 
void postTimestep (std::vector< GlobalVector * > const &x, const double t, const double delta_t, int const process_id)
 Postprocessing after a complete timestep. More...
 
void postNonLinearSolver (GlobalVector const &x, GlobalVector const &xdot, 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_dot, int const process_id)
 compute secondary variables for the coupled equations or for output. More...
 
NumLib::IterationResult postIteration (GlobalVector const &x) final
 
void initialize ()
 
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 setCoupledSolutionsForStaggeredScheme (CoupledSolutionsForStaggeredScheme *const coupled_solutions)
 
void updateDeactivatedSubdomains (double const time, const int process_id)
 
bool isMonolithicSchemeUsed () const
 
virtual void setCoupledTermForTheStaggeredSchemeToLocalAssemblers (int 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 &xdot, 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 &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) final
 
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
 
SecondaryVariableCollection const & getSecondaryVariables () const
 
std::vector< std::unique_ptr< 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)
 

Private Types

using LocalAssemblerIF = LocalAssemblerInterface< DisplacementDim >
 

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(). More...
 
void initializeBoundaryConditions () 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 &xdot, 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 &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
 
void postTimestepConcreteProcess (std::vector< GlobalVector * > const &x, double const t, double const dt, const int process_id) override
 
NumLib::LocalToGlobalIndexMap const & getDOFTable (const int process_id) const override
 
void computeSecondaryVariableConcrete (double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, int const process_id) override
 
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData () const override
 
std::vector< NumLib::LocalToGlobalIndexMap const * > getDOFTables (const int number_of_processes) const
 
bool hasMechanicalProcess (int const) const
 

Private Attributes

std::vector< MeshLib::Node * > base_nodes_
 
std::unique_ptr< MeshLib::MeshSubset const > mesh_subset_base_nodes_
 
ThermoRichardsMechanicsProcessData< DisplacementDim > 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_
 
MeshLib::PropertyVector< double > * nodal_forces_ = nullptr
 
MeshLib::PropertyVector< double > * hydraulic_flow_ = nullptr
 
MeshLib::PropertyVector< double > * heat_flux_ = nullptr
 

Additional Inherited Members

- Public Types inherited from ProcessLib::Process
using NonlinearSolver = NumLib::NonlinearSolverBase
 
using TimeDiscretization = NumLib::TimeDiscretization
 
- 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)
 
void constructMonolithicProcessDofTable ()
 
void constructDofTableOfSpecifiedProcessStaggeredScheme (const int specified_prosess_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
 
VectorMatrixAssembler _global_assembler
 
const bool _use_monolithic_scheme
 
CoupledSolutionsForStaggeredScheme_coupled_solutions
 
unsigned const _integration_order
 
std::vector< std::unique_ptr< 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>
using ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::LocalAssemblerIF = LocalAssemblerInterface<DisplacementDim>
private

Definition at line 148 of file ThermoRichardsMechanicsProcess.h.

Constructor & Destructor Documentation

◆ ThermoRichardsMechanicsProcess()

template<int DisplacementDim>
ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::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 > &&  process_data,
SecondaryVariableCollection &&  secondary_variables,
bool const  use_monolithic_scheme 
)

Definition at line 29 of file ThermoRichardsMechanicsProcess.cpp.

40  : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
41  integration_order, std::move(process_variables),
42  std::move(secondary_variables), use_monolithic_scheme),
43  process_data_(std::move(process_data))
44 {
45  nodal_forces_ = MeshLib::getOrCreateMeshProperty<double>(
46  mesh, "NodalForces", MeshLib::MeshItemType::Node, DisplacementDim);
47 
48  hydraulic_flow_ = MeshLib::getOrCreateMeshProperty<double>(
49  mesh, "MassFlowRate", MeshLib::MeshItemType::Node, 1);
50 
51  heat_flux_ = MeshLib::getOrCreateMeshProperty<double>(
52  mesh, "HeatFlowRate", MeshLib::MeshItemType::Node, 1);
53 
54  for (auto& [name, num_comp, ip_data_accessor] :
56  {
57  _integration_point_writer.emplace_back(
58  std::make_unique<IntegrationPointWriter>(
59  name, num_comp, integration_order, local_assemblers_,
60  ip_data_accessor));
61  }
62 }
std::string const name
Definition: Process.h:327
std::vector< std::unique_ptr< IntegrationPointWriter > > _integration_point_writer
Definition: Process.h:354
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)
Definition: Process.cpp:24
ThermoRichardsMechanicsProcessData< DisplacementDim > process_data_
static std::vector< IPDataAccessorForIPWriter > getIPDataAccessorsForIPWriter()

References ProcessLib::Process::_integration_point_writer, ProcessLib::ThermoRichardsMechanics::LocalAssemblerInterface< DisplacementDim >::getIPDataAccessorsForIPWriter(), ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::heat_flux_, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::hydraulic_flow_, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::local_assemblers_, ProcessLib::Process::name, ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::nodal_forces_, and MeshLib::Node.

Member Function Documentation

◆ assembleConcreteProcess()

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

Implements ProcessLib::Process.

Definition at line 224 of file ThermoRichardsMechanicsProcess.cpp.

229 {
230  OGS_FATAL(
231  "The Picard method or the Newton-Raphson method with numerical "
232  "Jacobian is not implemented for ThermoRichardsMechanics with the full "
233  "monolithic coupling scheme");
234 }
#define OGS_FATAL(...)
Definition: Error.h:26

References OGS_FATAL.

◆ assembleWithJacobianConcreteProcess()

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

Implements ProcessLib::Process.

Definition at line 237 of file ThermoRichardsMechanicsProcess.cpp.

244 {
245  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
246  dof_tables;
247 
248  DBUG(
249  "Assemble the Jacobian of ThermoRichardsMechanics for the monolithic "
250  "scheme.");
251  dof_tables.emplace_back(*_local_to_global_index_map);
252 
253  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
254 
257  local_assemblers_, pv.getActiveElementIDs(), dof_tables, t, dt, x, xdot,
258  process_id, M, K, b, Jac);
259 
260  auto copyRhs = [&](int const variable_id, auto& output_vector)
261  {
262  transformVariableFromGlobalVector(b, variable_id, dof_tables[0],
263  output_vector, std::negate<double>());
264  };
265 
266  copyRhs(0, *heat_flux_);
267  copyRhs(1, *hydraulic_flow_);
268  copyRhs(2, *nodal_forces_);
269 }
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition: Logging.h:29
std::vector< std::size_t > const & getActiveElementIDs() const
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:148
VectorMatrixAssembler _global_assembler
Definition: Process.h:337
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:333
void assembleWithJacobian(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)
static const double t
void transformVariableFromGlobalVector(GlobalVector const &input_vector, int const variable_id, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, MeshLib::PropertyVector< double > &output_vector, Functor mapFunction)
Definition: DOFTableUtil.h:59
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)

References ProcessLib::VectorMatrixAssembler::assembleWithJacobian(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), MathLib::t, and NumLib::transformVariableFromGlobalVector().

◆ computeSecondaryVariableConcrete()

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

Reimplemented from ProcessLib::Process.

Definition at line 288 of file ThermoRichardsMechanicsProcess.cpp.

293 {
294  DBUG("Compute the secondary variables for ThermoRichardsMechanicsProcess.");
295 
296  auto const dof_tables = getDOFTables(x.size());
297 
298  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
299 
302  pv.getActiveElementIDs(), dof_tables, t, dt, x, x_dot, process_id);
303 }
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)
std::vector< NumLib::LocalToGlobalIndexMap const * > getDOFTables(const int number_of_processes) const
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 MathLib::t.

◆ constructDofTable()

template<int DisplacementDim>
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::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 81 of file ThermoRichardsMechanicsProcess.cpp.

82 {
83  // Create single component dof in every of the mesh's nodes.
84  _mesh_subset_all_nodes = std::make_unique<MeshLib::MeshSubset>(
85  _mesh, _mesh.getNodes(), process_data_.use_TaylorHood_elements);
86  // Create single component dof in the mesh's base nodes.
88  mesh_subset_base_nodes_ = std::make_unique<MeshLib::MeshSubset>(
89  _mesh, base_nodes_, process_data_.use_TaylorHood_elements);
90 
91  // TODO move the two data members somewhere else.
92  // for extrapolation of secondary variables of stress or strain
93  std::vector<MeshLib::MeshSubset> all__meshsubsets_single_component{
96  std::make_unique<NumLib::LocalToGlobalIndexMap>(
97  std::move(all__meshsubsets_single_component),
98  // by location order is needed for output
100 
101  // For temperature, which is the first variable.
102  std::vector<MeshLib::MeshSubset> all__meshsubsets{*mesh_subset_base_nodes_};
103 
104  // For pressure, which is the second variable
105  all__meshsubsets.push_back(*mesh_subset_base_nodes_);
106 
107  // For displacement.
108  const int monolithic_process_id = 0;
109  std::generate_n(std::back_inserter(all__meshsubsets),
110  getProcessVariables(monolithic_process_id)[2]
111  .get()
112  .getNumberOfGlobalComponents(),
113  [&]() { return *_mesh_subset_all_nodes; });
114 
115  std::vector<int> const vec_n_components{1, 1, DisplacementDim};
117  std::make_unique<NumLib::LocalToGlobalIndexMap>(
118  std::move(all__meshsubsets), vec_n_components,
121 }
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:100
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:103
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:331
MeshLib::Mesh & _mesh
Definition: Process.h:330
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.

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

◆ getDOFTable()

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

Reimplemented from ProcessLib::Process.

Definition at line 316 of file ThermoRichardsMechanicsProcess.cpp.

318 {
320 }

◆ getDOFTableForExtrapolatorData()

template<int DisplacementDim>
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::getDOFTableForExtrapolatorData
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 307 of file ThermoRichardsMechanicsProcess.cpp.

308 {
309  const bool manage_storage = false;
310  return std::make_tuple(local_to_global_index_map_single_component_.get(),
311  manage_storage);
312 }

◆ getDOFTables()

template<int DisplacementDim>
std::vector< NumLib::LocalToGlobalIndexMap const * > ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::getDOFTables ( const int  number_of_processes) const
private

Definition at line 324 of file ThermoRichardsMechanicsProcess.cpp.

326 {
327  std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
328  dof_tables.reserve(number_of_processes);
329  std::generate_n(std::back_inserter(dof_tables), number_of_processes,
330  [&]() { return &getDOFTable(dof_tables.size()); });
331  return dof_tables;
332 }
NumLib::LocalToGlobalIndexMap const & getDOFTable(const int process_id) const override

◆ getMatrixSpecifications()

template<int DisplacementDim>
MathLib::MatrixSpecifications ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::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 72 of file ThermoRichardsMechanicsProcess.cpp.

74 {
75  auto const& l = *_local_to_global_index_map;
76  return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
77  &l.getGhostIndices(), &this->_sparsity_pattern};
78 }
GlobalSparsityPattern _sparsity_pattern
Definition: Process.h:356

◆ hasMechanicalProcess()

template<int DisplacementDim>
bool ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::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 217 of file ThermoRichardsMechanicsProcess.h.

217 { return true; }

◆ initializeBoundaryConditions()

template<int DisplacementDim>
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::initializeBoundaryConditions
overrideprivatevirtual

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

Reimplemented from ProcessLib::Process.

Definition at line 202 of file ThermoRichardsMechanicsProcess.cpp.

203 {
204  const int process_id = 0;
206  *_local_to_global_index_map, process_id);
207 }
void initializeProcessBoundaryConditionsAndSourceTerms(const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id)
Definition: Process.cpp:69

◆ initializeConcreteProcess()

template<int DisplacementDim>
void ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::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 124 of file ThermoRichardsMechanicsProcess.cpp.

128 {
129  createLocalAssemblersHM<DisplacementDim,
130  ThermoRichardsMechanicsLocalAssembler>(
131  mesh.getElements(), dof_table, local_assemblers_,
132  NumLib::IntegrationOrder{integration_order}, mesh.isAxiallySymmetric(),
133  process_data_);
134 
135  auto add_secondary_variable = [&](std::string const& name,
136  int const num_components,
137  auto get_ip_values_function)
138  {
140  name,
141  makeExtrapolator(num_components, getExtrapolator(),
143  std::move(get_ip_values_function)));
144  };
145 
146  for (auto& [name, num_comp, fct] : LocalAssemblerInterface<
147  DisplacementDim>::getIPDataAccessorsForExtrapolation())
148  {
149  add_secondary_variable(name, num_comp, fct);
150  };
151 
152  //
153  // enable output of internal variables defined by material models
154  //
156  LocalAssemblerIF>(process_data_.solid_materials,
157  add_secondary_variable);
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)
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:184
SecondaryVariableCollection _secondary_variables
Definition: Process.h:335
void addSecondaryVariable(std::string const &internal_name, SecondaryVariableFunctions &&fcts)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:48
void solidMaterialInternalToSecondaryVariables(std::map< int, std::unique_ptr< MaterialLib::Solids::MechanicsBase< DisplacementDim >>> const &solid_materials, AddSecondaryVariableCallback const &add_secondary_variable)
void setIPDataInitialConditions(std::vector< std::unique_ptr< IntegrationPointWriter >> const &_integration_point_writer, MeshLib::Properties const &mesh_properties, LocalAssemblersVector &local_assemblers)
void createLocalAssemblersHM(std::vector< MeshLib::Element * > const &mesh_elements, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< std::unique_ptr< LocalAssemblerInterface >> &local_assemblers, NumLib::IntegrationOrder const integration_order, ExtraCtorArgs &&... extra_ctor_args)
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
static void executeMemberOnDereferenced(Method method, Container const &container, Args &&... args)

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

◆ isLinear()

template<int DisplacementDim>
bool ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::isLinear
override

Definition at line 65 of file ThermoRichardsMechanicsProcess.cpp.

66 {
67  return false;
68 }

◆ postTimestepConcreteProcess()

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

Reimplemented from ProcessLib::Process.

Definition at line 272 of file ThermoRichardsMechanicsProcess.cpp.

276 {
277  DBUG("PostTimestep ThermoRichardsMechanicsProcess.");
278 
279  auto const dof_tables = getDOFTables(x.size());
280 
281  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
284  pv.getActiveElementIDs(), dof_tables, x, t, dt);
285 }
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)

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

◆ setInitialConditionsConcreteProcess()

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

Reimplemented from ProcessLib::Process.

Definition at line 210 of file ThermoRichardsMechanicsProcess.cpp.

214 {
215  DBUG("SetInitialConditions ThermoRichardsMechanicsProcess.");
216 
220  process_id);
221 }
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)
const bool _use_monolithic_scheme
Definition: Process.h:339

References DBUG(), NumLib::SerialExecutor::executeMemberOnDereferenced(), ProcessLib::setInitialConditions(), and MathLib::t.

Member Data Documentation

◆ base_nodes_

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

Definition at line 183 of file ThermoRichardsMechanicsProcess.h.

◆ heat_flux_

◆ hydraulic_flow_

◆ local_assemblers_

◆ local_to_global_index_map_single_component_

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

Definition at line 190 of file ThermoRichardsMechanicsProcess.h.

◆ local_to_global_index_map_with_base_nodes_

template<int DisplacementDim>
std::unique_ptr<NumLib::LocalToGlobalIndexMap> ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::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 195 of file ThermoRichardsMechanicsProcess.h.

◆ mesh_subset_base_nodes_

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

Definition at line 184 of file ThermoRichardsMechanicsProcess.h.

◆ nodal_forces_

◆ process_data_

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

Definition at line 185 of file ThermoRichardsMechanicsProcess.h.

◆ sparsity_pattern_with_linear_element_

template<int DisplacementDim>
GlobalSparsityPattern ProcessLib::ThermoRichardsMechanics::ThermoRichardsMechanicsProcess< DisplacementDim >::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 199 of file ThermoRichardsMechanicsProcess.h.


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