OGS
ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess Class Referencefinal

Detailed Description

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

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, \(H(S-1)\) is the Heaviside function, and \( \mathbf u\) is the displacement. While this process does not contain a fully mechanical coupling, simplfied expressions can be given to approximate the latter term under certain stress conditions. 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.

Definition at line 81 of file ThermoRichardsFlowProcess.h.

#include <ThermoRichardsFlowProcess.h>

Inheritance diagram for ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess:
[legend]
Collaboration diagram for ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess:
[legend]

Public Member Functions

 ThermoRichardsFlowProcess (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, ThermoRichardsFlowProcessData &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)
 
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
 
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable (const int) const
 
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 std::vector< std::string > initializeAssemblyOnSubmeshes (std::vector< std::reference_wrapper< MeshLib::Mesh > > const &meshes)
 
virtual ~SubmeshAssemblySupport ()=default
 

Private Types

using LocalAssemblerIF = LocalAssemblerInterface
 

Private Member Functions

void initializeConcreteProcess (NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
 Process specific initialization called by initialize().
 
void setInitialConditionsConcreteProcess (std::vector< GlobalVector * > &x, double const t, int const) 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 postTimestepConcreteProcess (std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, const int process_id) override
 
void computeSecondaryVariableConcrete (double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id) override
 

Private Attributes

std::vector< MeshLib::Node * > _base_nodes
 
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_base_nodes
 
ThermoRichardsFlowProcessData _process_data
 
Assembly::ParallelVectorMatrixAssembler _pvma {*_jacobian_assembler}
 
std::vector< std::unique_ptr< LocalAssemblerIF > > _local_assemblers
 
GlobalSparsityPattern _sparsity_pattern_with_linear_element
 
MeshLib::PropertyVector< double > * _heat_flux = nullptr
 
MeshLib::PropertyVector< double > * _hydraulic_flow = nullptr
 

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)
 
virtual void constructDofTable ()
 
void constructMonolithicProcessDofTable ()
 
void constructDofTableOfSpecifiedProcessStaggeredScheme (const int specified_process_id)
 
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData () const
 
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

Constructor & Destructor Documentation

◆ ThermoRichardsFlowProcess()

ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::ThermoRichardsFlowProcess ( 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,
ThermoRichardsFlowProcessData && process_data,
SecondaryVariableCollection && secondary_variables,
bool const use_monolithic_scheme )

Definition at line 29 of file ThermoRichardsFlowProcess.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 _heat_flux = MeshLib::getOrCreateMeshProperty<double>(
46 mesh, "HeatFlowRate", MeshLib::MeshItemType::Node, 1);
47
48 _hydraulic_flow = MeshLib::getOrCreateMeshProperty<double>(
49 mesh, "MassFlowRate", MeshLib::MeshItemType::Node, 1);
50
51 // TODO (naumov) remove ip suffix. Probably needs modification of the mesh
52 // properties, s.t. there is no "overlapping" with cell/point data.
53 // See getOrCreateMeshProperty.
54 _integration_point_writer.emplace_back(
55 std::make_unique<MeshLib::IntegrationPointWriter>(
56 "saturation_ip", 1 /*n components*/, integration_order,
58
59 _integration_point_writer.emplace_back(
60 std::make_unique<MeshLib::IntegrationPointWriter>(
61 "porosity_ip", 1 /*n components*/, integration_order,
63}
std::string const name
Definition Process.h:351
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
Definition Process.h:377
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:28
std::vector< std::unique_ptr< LocalAssemblerIF > > _local_assemblers
virtual std::vector< double > getSaturation() const =0
virtual std::vector< double > getPorosity() const =0

References _heat_flux, _hydraulic_flow, ProcessLib::Process::_integration_point_writer, _local_assemblers, ProcessLib::ThermoRichardsFlow::LocalAssemblerInterface::getPorosity(), ProcessLib::ThermoRichardsFlow::LocalAssemblerInterface::getSaturation(), and MeshLib::Node.

Member Function Documentation

◆ assembleConcreteProcess()

void ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::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 132 of file ThermoRichardsFlowProcess.cpp.

136{
137 DBUG("Assemble the equations for ThermoRichardsFlowProcess.");
138
139 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
140 dof_table = {std::ref(*_local_to_global_index_map)};
141 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
142
143 // Call global assembler for each local assembly item.
146 pv.getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K,
147 b);
148}
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
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
VectorMatrixAssembler _global_assembler
Definition Process.h:364
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition Process.h:357
void assemble(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap > > const &dof_tables, double const 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)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)

References ProcessLib::Process::_global_assembler, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::VectorMatrixAssembler::assemble(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().

◆ assembleWithJacobianConcreteProcess()

void ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::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 150 of file ThermoRichardsFlowProcess.cpp.

154{
155 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
156 dof_tables;
157
158 DBUG(
159 "Assemble the Jacobian of ThermoRichardsFlow for the monolithic "
160 "scheme.");
161 dof_tables.emplace_back(*_local_to_global_index_map);
162
163 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
164
166 dof_tables, t, dt, x, x_prev, process_id, M, K,
167 b, Jac);
168
169 auto copyRhs = [&](int const variable_id, auto& output_vector)
170 {
171 transformVariableFromGlobalVector(b, variable_id, dof_tables[0],
172 output_vector, std::negate<double>());
173 };
174
175 copyRhs(0, *_heat_flux);
176 copyRhs(1, *_hydraulic_flow);
177}
void assembleWithJacobian(BaseLib::PolymorphicRandomAccessContainerView< LocalAssemblerInterface > const &local_assemblers, std::vector< std::size_t > const &active_elements, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap > > const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &xs, std::vector< GlobalVector * > const &x_prevs, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)
void transformVariableFromGlobalVector(GlobalVector const &input_vector, int const variable_id, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, MeshLib::PropertyVector< double > &output_vector, Functor map_function)

References _heat_flux, _hydraulic_flow, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, _pvma, ProcessLib::Assembly::ParallelVectorMatrixAssembler::assembleWithJacobian(), DBUG(), ProcessLib::ProcessVariable::getActiveElementIDs(), and ProcessLib::Process::getProcessVariables().

◆ computeSecondaryVariableConcrete()

void ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::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 203 of file ThermoRichardsFlowProcess.cpp.

206{
207 if (process_id != 0)
208 {
209 return;
210 }
211 DBUG(
212 "Compute the secondary variables for "
213 "ThermoRichardsFlowProcess.");
214
215 auto get_a_dof_table_func = [this](const int processe_id) -> auto&
216 {
217 return getDOFTable(processe_id);
218 };
219 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
220
224 NumLib::getDOFTables(x.size(), get_a_dof_table_func), t, dt, x, x_prev,
225 process_id);
226}
virtual void computeSecondaryVariable(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const
Definition Process.h:147
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 _local_assemblers, ProcessLib::LocalAssemblerInterface::computeSecondaryVariable(), DBUG(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getDOFTable(), NumLib::getDOFTables(), and ProcessLib::Process::getProcessVariables().

◆ initializeConcreteProcess()

void ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::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 65 of file ThermoRichardsFlowProcess.cpp.

69{
70 ProcessLib::createLocalAssemblers<ThermoRichardsFlowLocalAssembler>(
71 mesh.getDimension(), mesh.getElements(), dof_table, _local_assemblers,
72 NumLib::IntegrationOrder{integration_order}, mesh.isAxiallySymmetric(),
74
75 auto add_secondary_variable = [&](std::string const& name,
76 int const num_components,
77 auto get_ip_values_function)
78 {
80 name,
81 makeExtrapolator(num_components, getExtrapolator(),
83 std::move(get_ip_values_function)));
84 };
85
86 add_secondary_variable("velocity", mesh.getDimension(),
88
89 add_secondary_variable("saturation", 1,
91
92 add_secondary_variable("porosity", 1, &LocalAssemblerIF::getIntPtPorosity);
93
94 add_secondary_variable("dry_density_solid", 1,
96
97 _process_data.element_saturation = MeshLib::getOrCreateMeshProperty<double>(
98 const_cast<MeshLib::Mesh&>(mesh), "saturation_avg",
100
101 _process_data.element_porosity = MeshLib::getOrCreateMeshProperty<double>(
102 const_cast<MeshLib::Mesh&>(mesh), "porosity_avg",
104
107
108 // Initialize local assemblers after all variables have been set.
112}
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)
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)
virtual 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 =0
virtual 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 =0
virtual 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 =0
virtual 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 =0

References ProcessLib::Process::_integration_point_writer, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, _process_data, ProcessLib::Process::_secondary_variables, ProcessLib::SecondaryVariableCollection::addSecondaryVariable(), MeshLib::Cell, ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcessData::element_porosity, ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcessData::element_saturation, NumLib::SerialExecutor::executeMemberOnDereferenced(), MeshLib::Mesh::getDimension(), MeshLib::Mesh::getElements(), ProcessLib::Process::getExtrapolator(), ProcessLib::ThermoRichardsFlow::LocalAssemblerInterface::getIntPtDarcyVelocity(), ProcessLib::ThermoRichardsFlow::LocalAssemblerInterface::getIntPtDryDensitySolid(), ProcessLib::ThermoRichardsFlow::LocalAssemblerInterface::getIntPtPorosity(), ProcessLib::ThermoRichardsFlow::LocalAssemblerInterface::getIntPtSaturation(), MeshLib::Mesh::getProperties(), ProcessLib::LocalAssemblerInterface::initialize(), MeshLib::Mesh::isAxiallySymmetric(), ProcessLib::makeExtrapolator(), ProcessLib::Process::name, and ProcessLib::setIPDataInitialConditions().

◆ isLinear()

bool ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::isLinear ( ) const
inlineoverride

Definition at line 101 of file ThermoRichardsFlowProcess.h.

101{ return false; }

◆ postTimestepConcreteProcess()

void ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::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 179 of file ThermoRichardsFlowProcess.cpp.

183{
184 if (process_id != 0)
185 {
186 return;
187 }
188
189 DBUG("PostTimestep ThermoRichardsFlowProcess.");
190
191 auto get_a_dof_table_func = [this](const int processe_id) -> auto&
192 {
193 return getDOFTable(processe_id);
194 };
195 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
199 NumLib::getDOFTables(x.size(), get_a_dof_table_func), x, x_prev, t, dt,
200 process_id);
201}
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 _local_assemblers, DBUG(), NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getDOFTable(), NumLib::getDOFTables(), ProcessLib::Process::getProcessVariables(), and ProcessLib::LocalAssemblerInterface::postTimestep().

◆ setInitialConditionsConcreteProcess()

void ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::setInitialConditionsConcreteProcess ( std::vector< GlobalVector * > & x,
double const t,
int const process_id )
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 114 of file ThermoRichardsFlowProcess.cpp.

116{
117 if (process_id != 0)
118 {
119 return;
120 }
121 DBUG("SetInitialConditions ThermoRichardsFlowProcess.");
122
123 auto get_a_dof_table_func = [this](const int process_id) -> auto&
124 {
125 return getDOFTable(process_id);
126 };
129 NumLib::getDOFTables(x.size(), get_a_dof_table_func), x, t, process_id);
130}
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 _local_assemblers, DBUG(), NumLib::SerialExecutor::executeMemberOnDereferenced(), ProcessLib::Process::getDOFTable(), NumLib::getDOFTables(), and ProcessLib::LocalAssemblerInterface::setInitialConditions().

Member Data Documentation

◆ _base_nodes

std::vector<MeshLib::Node*> ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_base_nodes
private

Definition at line 134 of file ThermoRichardsFlowProcess.h.

◆ _heat_flux

MeshLib::PropertyVector<double>* ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_heat_flux = nullptr
private

◆ _hydraulic_flow

MeshLib::PropertyVector<double>* ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_hydraulic_flow = nullptr
private

◆ _local_assemblers

std::vector<std::unique_ptr<LocalAssemblerIF> > ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_local_assemblers
private

◆ _mesh_subset_base_nodes

std::unique_ptr<MeshLib::MeshSubset const> ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_mesh_subset_base_nodes
private

Definition at line 135 of file ThermoRichardsFlowProcess.h.

◆ _process_data

ThermoRichardsFlowProcessData ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_process_data
private

Definition at line 136 of file ThermoRichardsFlowProcess.h.

Referenced by initializeConcreteProcess().

◆ _pvma

Assembly::ParallelVectorMatrixAssembler ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_pvma {*_jacobian_assembler}
private

Definition at line 138 of file ThermoRichardsFlowProcess.h.

std::unique_ptr< ProcessLib::AbstractJacobianAssembler > _jacobian_assembler
Definition Process.h:363

Referenced by assembleWithJacobianConcreteProcess().

◆ _sparsity_pattern_with_linear_element

GlobalSparsityPattern ProcessLib::ThermoRichardsFlow::ThermoRichardsFlowProcess::_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 144 of file ThermoRichardsFlowProcess.h.


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