OGS 6.2.1-76-gbb689931b
ProcessLib::HT::HTProcess Class Referencefinal

Detailed Description

HT process

The implementation uses a monolithic approach, i.e., both processes are assembled within one global system of equations.

Process Coupling

The advective term of the heat conduction equation is given by the confined groundwater flow process, i.e., the heat conduction depends on darcy velocity of the groundwater flow process. On the other hand the temperature dependencies of the viscosity and density in the groundwater flow couples the H process to the T process.

Note
At the moment there is not any coupling by source or sink terms, i.e., the coupling is implemented only by density changes due to temperature changes in the buoyance term in the groundwater flow. This coupling schema is referred to as the Boussinesq approximation.

Definition at line 49 of file HTProcess.h.

#include <HTProcess.h>

Inheritance diagram for ProcessLib::HT::HTProcess:
Collaboration diagram for ProcessLib::HT::HTProcess:

Public Member Functions

 HTProcess (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, HTProcessData &&process_data, SecondaryVariableCollection &&secondary_variables, NumLib::NamedFunctionCaller &&named_function_caller, bool const use_monolithic_scheme, std::unique_ptr< ProcessLib::SurfaceFluxData > &&surfaceflux, const int heat_transport_process_id, const int hydraulic_process_id)
 
Eigen::Vector3d getFlux (std::size_t element_id, MathLib::Point3d const &p, double const t, GlobalVector const &x) const override
 
void setCoupledTermForTheStaggeredSchemeToLocalAssemblers () override
 
void postTimestepConcreteProcess (GlobalVector const &x, const double t, const double delta_t, int const process_id) 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, NumLib::NamedFunctionCaller &&named_function_caller, const bool use_monolithic_scheme=true)
 
void preTimestep (GlobalVector const &x, const double t, const double delta_t, const int process_id)
 Preprocessing before starting assembly for new timestep. More...
 
void postTimestep (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, const double t, int const process_id)
 
void preIteration (const unsigned iter, GlobalVector const &x) final
 
void computeSecondaryVariable (const double t, GlobalVector const &x, 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 (const int process_id, const double t, GlobalVector &x)
 
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
 
void preAssemble (const double t, GlobalVector const &x) final
 
void assemble (const double t, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) final
 
void assembleWithJacobian (const double t, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) final
 
std::vector< NumLib::IndexValueVector< GlobalIndexType > > const * getKnownSolutions (double const t, GlobalVector const &x) 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
 
SecondaryVariableCollection const & getSecondaryVariables () const
 
std::vector< std::unique_ptr< IntegrationPointWriter > > const & getIntegrationPointWriter () const
 

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(). More...
 
void assembleConcreteProcess (const double t, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
 
void assembleWithJacobianConcreteProcess (const double t, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
 
void preTimestepConcreteProcess (GlobalVector const &x, double const t, double const dt, const int process_id) override
 
void setCoupledSolutionsOfPreviousTimeStepPerProcess (const int process_id)
 
void setCoupledSolutionsOfPreviousTimeStep ()
 
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData () const override
 

Private Attributes

HTProcessData _process_data
 
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
 
std::array< std::unique_ptr< GlobalVector >, 2 > _xs_previous_timestep
 Solutions of the previous time step. More...
 
std::unique_ptr< ProcessLib::SurfaceFluxData_surfaceflux
 
const int _heat_transport_process_id
 
const int _hydraulic_process_id
 

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
 
- 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)
 
virtual void constructDofTable ()
 
void constructMonolithicProcessDofTable ()
 
void constructDofTableOfSpecifiedProsessStaggerdScheme (const int specified_prosess_id)
 
- 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
 
NumLib::NamedFunctionCaller _named_function_caller
 
std::vector< std::unique_ptr< CachedSecondaryVariable > > _cached_secondary_variables
 
SecondaryVariableContext _secondary_variable_context
 
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
 

Constructor & Destructor Documentation

◆ HTProcess()

ProcessLib::HT::HTProcess::HTProcess ( 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,
HTProcessData &&  process_data,
SecondaryVariableCollection &&  secondary_variables,
NumLib::NamedFunctionCaller &&  named_function_caller,
bool const  use_monolithic_scheme,
std::unique_ptr< ProcessLib::SurfaceFluxData > &&  surfaceflux,
const int  heat_transport_process_id,
const int  hydraulic_process_id 
)

Definition at line 26 of file HTProcess.cpp.

41  : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
42  integration_order, std::move(process_variables),
43  std::move(secondary_variables), std::move(named_function_caller),
44  use_monolithic_scheme),
45  _process_data(std::move(process_data)),
46  _surfaceflux(std::move(surfaceflux)),
47  _heat_transport_process_id(heat_transport_process_id),
48  _hydraulic_process_id(hydraulic_process_id)
49 {
50 }
std::unique_ptr< ProcessLib::SurfaceFluxData > _surfaceflux
Definition: HTProcess.h:126
const int _hydraulic_process_id
Definition: HTProcess.h:129
HTProcessData _process_data
Definition: HTProcess.h:119
const int _heat_transport_process_id
Definition: HTProcess.h:128
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, NumLib::NamedFunctionCaller &&named_function_caller, const bool use_monolithic_scheme=true)
Definition: Process.cpp:24
std::string const name
Definition: Process.h:284

Member Function Documentation

◆ assembleConcreteProcess()

void ProcessLib::HT::HTProcess::assembleConcreteProcess ( const double  t,
GlobalVector const &  x,
GlobalMatrix &  M,
GlobalMatrix &  K,
GlobalVector &  b 
)
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 88 of file HTProcess.cpp.

References ProcessLib::Process::_coupled_solutions, ProcessLib::Process::_global_assembler, _heat_transport_process_id, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::Process::_use_monolithic_scheme, ProcessLib::VectorMatrixAssembler::assemble(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), ProcessLib::CoupledSolutionsForStaggeredScheme::process_id, and setCoupledSolutionsOfPreviousTimeStep().

93 {
94  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
95  dof_tables;
97  {
98  DBUG("Assemble HTProcess.");
99  dof_tables.emplace_back(*_local_to_global_index_map);
100  }
101  else
102  {
104  {
105  DBUG(
106  "Assemble the equations of heat transport process within "
107  "HTProcess.");
108  }
109  else
110  {
111  DBUG(
112  "Assemble the equations of single phase fully saturated "
113  "fluid flow process within HTProcess.");
114  }
116  dof_tables.emplace_back(*_local_to_global_index_map);
117  dof_tables.emplace_back(*_local_to_global_index_map);
118  }
119 
120  const int process_id =
122  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
123  // Call global assembler for each local assembly item.
126  pv.getActiveElementIDs(), dof_tables, t, x, M, K, b,
128 }
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:125
void setCoupledSolutionsOfPreviousTimeStep()
Definition: HTProcess.cpp:239
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, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:290
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:305
const int _heat_transport_process_id
Definition: HTProcess.h:128
const bool _use_monolithic_scheme
Definition: Process.h:301
VectorMatrixAssembler _global_assembler
Definition: Process.h:299
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
Definition: HTProcess.h:121
std::vector< std::size_t > & getActiveElementIDs() const

◆ assembleWithJacobianConcreteProcess()

void ProcessLib::HT::HTProcess::assembleWithJacobianConcreteProcess ( const double  t,
GlobalVector const &  x,
GlobalVector const &  xdot,
const double  dxdot_dx,
const double  dx_dx,
GlobalMatrix &  M,
GlobalMatrix &  K,
GlobalVector &  b,
GlobalMatrix &  Jac 
)
overrideprivatevirtual

Implements ProcessLib::Process.

Definition at line 130 of file HTProcess.cpp.

References ProcessLib::Process::_coupled_solutions, ProcessLib::Process::_global_assembler, _local_assemblers, ProcessLib::Process::_local_to_global_index_map, ProcessLib::Process::_use_monolithic_scheme, ProcessLib::VectorMatrixAssembler::assembleWithJacobian(), NumLib::SerialExecutor::executeSelectedMemberDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), ProcessLib::CoupledSolutionsForStaggeredScheme::process_id, and setCoupledSolutionsOfPreviousTimeStep().

134 {
135  DBUG("AssembleWithJacobian HTProcess.");
136 
137  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
138  dof_tables;
140  {
142  dof_tables.emplace_back(std::ref(*_local_to_global_index_map));
143  }
144  else
145  {
146  dof_tables.emplace_back(std::ref(*_local_to_global_index_map));
147  dof_tables.emplace_back(std::ref(*_local_to_global_index_map));
148  }
149 
150  // Call global assembler for each local assembly item.
151  const int process_id =
153  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
156  _local_assemblers, pv.getActiveElementIDs(), dof_tables, t, x, xdot,
157  dxdot_dx, dx_dx, M, K, b, Jac, _coupled_solutions);
158 }
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:125
void setCoupledSolutionsOfPreviousTimeStep()
Definition: HTProcess.cpp:239
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, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:290
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:305
const bool _use_monolithic_scheme
Definition: Process.h:301
VectorMatrixAssembler _global_assembler
Definition: Process.h:299
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
Definition: HTProcess.h:121
std::vector< std::size_t > & getActiveElementIDs() const

◆ getDOFTableForExtrapolatorData()

std::tuple< NumLib::LocalToGlobalIndexMap *, bool > ProcessLib::HT::HTProcess::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 200 of file HTProcess.cpp.

References ProcessLib::Process::_local_to_global_index_map, ProcessLib::Process::_mesh_subset_all_nodes, ProcessLib::Process::_use_monolithic_scheme, and NumLib::BY_LOCATION.

201 {
203  {
204  // For single-variable-single-component processes reuse the existing DOF
205  // table.
206  const bool manage_storage = false;
207  return std::make_tuple(_local_to_global_index_map.get(),
208  manage_storage);
209  }
210 
211  // Otherwise construct a new DOF table.
212  std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
214 
215  const bool manage_storage = true;
216  return std::make_tuple(new NumLib::LocalToGlobalIndexMap(
217  std::move(all_mesh_subsets_single_component),
218  // by location order is needed for output
220  manage_storage);
221 }
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:288
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:290
Ordering data by spatial location.
const bool _use_monolithic_scheme
Definition: Process.h:301

◆ getFlux()

Eigen::Vector3d ProcessLib::HT::HTProcess::getFlux ( std::size_t  element_id,
MathLib::Point3d const &  p,
double const  t,
GlobalVector const &  x 
) const
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 246 of file HTProcess.cpp.

References _local_assemblers, ProcessLib::Process::_local_to_global_index_map, and NumLib::getRowColumnIndices().

250 {
251  // fetch local_x from primary variable
252  std::vector<GlobalIndexType> indices_cache;
253  auto const r_c_indices = NumLib::getRowColumnIndices(
254  element_id, *_local_to_global_index_map, indices_cache);
255  std::vector<double> local_x(x.get(r_c_indices.rows));
256 
257  return _local_assemblers[element_id]->getFlux(p, t, local_x);
258 }
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:290
static const double p
NumLib::LocalToGlobalIndexMap::RowColumnIndices getRowColumnIndices(std::size_t const id, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< GlobalIndexType > &indices)
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
Definition: HTProcess.h:121

◆ initializeConcreteProcess()

void ProcessLib::HT::HTProcess::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 52 of file HTProcess.cpp.

References _heat_transport_process_id, _hydraulic_process_id, _local_assemblers, _process_data, ProcessLib::Process::_secondary_variables, ProcessLib::Process::_use_monolithic_scheme, ProcessLib::SecondaryVariableCollection::addSecondaryVariable(), MeshLib::Mesh::getDimension(), MeshLib::Mesh::getElements(), ProcessLib::Process::getExtrapolator(), ProcessLib::HT::HTLocalAssemblerInterface::getIntPtDarcyVelocity(), ProcessLib::Process::getProcessVariables(), ProcessLib::ProcessVariable::getShapeFunctionOrder(), MeshLib::Mesh::isAxiallySymmetric(), and ProcessLib::makeExtrapolator().

56 {
57  // For the staggered scheme, both processes are assumed to use the same
58  // element order. Therefore the order of shape function can be fetched from
59  // any set of the sets of process variables of the coupled processes. Here,
60  // we take the one from the first process by setting process_id = 0.
61  const int process_id = 0;
62  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
63 
65  {
66  ProcessLib::createLocalAssemblers<MonolithicHTFEM>(
67  mesh.getDimension(), mesh.getElements(), dof_table,
69  mesh.isAxiallySymmetric(), integration_order,
71  }
72  else
73  {
74  ProcessLib::createLocalAssemblers<StaggeredHTFEM>(
75  mesh.getDimension(), mesh.getElements(), dof_table,
77  mesh.isAxiallySymmetric(), integration_order, _process_data,
79  }
80 
82  "darcy_velocity",
83  makeExtrapolator(mesh.getDimension(), getExtrapolator(),
86 }
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:125
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
virtual std::vector< double > const & getIntPtDarcyVelocity(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &) const =0
const int _hydraulic_process_id
Definition: HTProcess.h:129
HTProcessData _process_data
Definition: HTProcess.h:119
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:152
SecondaryVariableCollection _secondary_variables
Definition: Process.h:292
unsigned getShapeFunctionOrder() const
const int _heat_transport_process_id
Definition: HTProcess.h:128
void addSecondaryVariable(std::string const &internal_name, SecondaryVariableFunctions &&fcts)
const bool _use_monolithic_scheme
Definition: Process.h:301
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
Definition: HTProcess.h:121

◆ isLinear()

bool ProcessLib::HT::HTProcess::isLinear ( ) const
inlineoverride

Definition at line 72 of file HTProcess.h.

72 { return false; }

◆ postTimestepConcreteProcess()

void ProcessLib::HT::HTProcess::postTimestepConcreteProcess ( GlobalVector const &  x,
const double  t,
const double  delta_t,
int const  process_id 
)
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 261 of file HTProcess.cpp.

References _hydraulic_process_id, ProcessLib::Process::_integration_order, ProcessLib::Process::_mesh, _surfaceflux, ProcessLib::Process::_use_monolithic_scheme, ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), and OGS_FATAL.

265 {
266  // For the monolithic scheme, process_id is always zero.
267  if (_use_monolithic_scheme && process_id != 0)
268  {
269  OGS_FATAL(
270  "The condition of process_id = 0 must be satisfied for "
271  "monolithic HTProcess, which is a single process.");
272  }
273  if (!_use_monolithic_scheme && process_id != _hydraulic_process_id)
274  {
275  DBUG("This is the thermal part of the staggered HTProcess.");
276  return;
277  }
278  if (!_surfaceflux) // computing the surfaceflux is optional
279  {
280  return;
281  }
282 
283  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
284 
285  _surfaceflux->integrate(x, t, *this, process_id, _integration_order, _mesh,
286  pv.getActiveElementIDs());
287  _surfaceflux->save(t);
288 }
MeshLib::Mesh & _mesh
Definition: Process.h:287
std::unique_ptr< ProcessLib::SurfaceFluxData > _surfaceflux
Definition: HTProcess.h:126
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:125
unsigned const _integration_order
Definition: Process.h:310
const int _hydraulic_process_id
Definition: HTProcess.h:129
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
const bool _use_monolithic_scheme
Definition: Process.h:301
std::vector< std::size_t > & getActiveElementIDs() const

◆ preTimestepConcreteProcess()

void ProcessLib::HT::HTProcess::preTimestepConcreteProcess ( GlobalVector const &  x,
double const  t,
double const  dt,
const int  process_id 
)
overrideprivatevirtual

Reimplemented from ProcessLib::Process.

Definition at line 160 of file HTProcess.cpp.

References ProcessLib::Process::_use_monolithic_scheme, _xs_previous_timestep, and MathLib::LinAlg::copy().

164 {
165  assert(process_id < 2);
166 
168  {
169  return;
170  }
171 
172  if (!_xs_previous_timestep[process_id])
173  {
174  _xs_previous_timestep[process_id] =
176  }
177  else
178  {
179  auto& x0 = *_xs_previous_timestep[process_id];
180  MathLib::LinAlg::copy(x, x0);
181  }
182 
183  auto& x0 = *_xs_previous_timestep[process_id];
184  MathLib::LinAlg::setLocalAccessibleVector(x0);
185 }
std::array< std::unique_ptr< GlobalVector >, 2 > _xs_previous_timestep
Solutions of the previous time step.
Definition: HTProcess.h:124
void copy(MatrixOrVector const &x, MatrixOrVector &y)
Copies x to y.
Definition: LinAlg.h:36
const bool _use_monolithic_scheme
Definition: Process.h:301

◆ setCoupledSolutionsOfPreviousTimeStep()

void ProcessLib::HT::HTProcess::setCoupledSolutionsOfPreviousTimeStep ( )
private

Set the solutions of the previous time step to the coupled term. It is only for the staggered scheme, and it must be called within the coupling loop because that the coupling term is only created there.

Definition at line 239 of file HTProcess.cpp.

References ProcessLib::Process::_coupled_solutions, _heat_transport_process_id, _hydraulic_process_id, ProcessLib::CoupledSolutionsForStaggeredScheme::coupled_xs_t0, and setCoupledSolutionsOfPreviousTimeStepPerProcess().

Referenced by assembleConcreteProcess(), and assembleWithJacobianConcreteProcess().

240 {
244 }
const int _hydraulic_process_id
Definition: HTProcess.h:129
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:305
const int _heat_transport_process_id
Definition: HTProcess.h:128
void setCoupledSolutionsOfPreviousTimeStepPerProcess(const int process_id)
Definition: HTProcess.cpp:223
std::vector< GlobalVector * > coupled_xs_t0
Pointers to the vector of the solutions of the previous time step.

◆ setCoupledSolutionsOfPreviousTimeStepPerProcess()

void ProcessLib::HT::HTProcess::setCoupledSolutionsOfPreviousTimeStepPerProcess ( const int  process_id)
private

Definition at line 223 of file HTProcess.cpp.

References ProcessLib::Process::_coupled_solutions, _xs_previous_timestep, ProcessLib::CoupledSolutionsForStaggeredScheme::coupled_xs_t0, and OGS_FATAL.

Referenced by setCoupledSolutionsOfPreviousTimeStep().

225 {
226  const auto& x_t0 = _xs_previous_timestep[process_id];
227  if (x_t0 == nullptr)
228  {
229  OGS_FATAL(
230  "Memory is not allocated for the global vector of the solution of "
231  "the previous time step for the staggered scheme.\n It can be done "
232  "by overriding Process::preTimestepConcreteProcess (ref. "
233  "HTProcess::preTimestepConcreteProcess) ");
234  }
235 
236  _coupled_solutions->coupled_xs_t0[process_id] = x_t0.get();
237 }
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:305
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
std::array< std::unique_ptr< GlobalVector >, 2 > _xs_previous_timestep
Solutions of the previous time step.
Definition: HTProcess.h:124
std::vector< GlobalVector * > coupled_xs_t0
Pointers to the vector of the solutions of the previous time step.

◆ setCoupledTermForTheStaggeredSchemeToLocalAssemblers()

void ProcessLib::HT::HTProcess::setCoupledTermForTheStaggeredSchemeToLocalAssemblers ( )
overridevirtual

Reimplemented from ProcessLib::Process.

Definition at line 187 of file HTProcess.cpp.

References ProcessLib::Process::_coupled_solutions, _local_assemblers, ProcessLib::Process::_use_monolithic_scheme, NumLib::SerialExecutor::executeSelectedMemberOnDereferenced(), ProcessLib::ProcessVariable::getActiveElementIDs(), ProcessLib::Process::getProcessVariables(), ProcessLib::CoupledSolutionsForStaggeredScheme::process_id, and ProcessLib::HT::HTLocalAssemblerInterface::setStaggeredCoupledSolutions().

188 {
189  DBUG("Set the coupled term for the staggered scheme to local assembers.");
190 
191  const int process_id =
193  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
197 }
void setStaggeredCoupledSolutions(std::size_t const, CoupledSolutionsForStaggeredScheme *const coupling_term)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:125
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:305
const bool _use_monolithic_scheme
Definition: Process.h:301
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
Definition: HTProcess.h:121
std::vector< std::size_t > & getActiveElementIDs() const

Member Data Documentation

◆ _heat_transport_process_id

const int ProcessLib::HT::HTProcess::_heat_transport_process_id
private

◆ _hydraulic_process_id

const int ProcessLib::HT::HTProcess::_hydraulic_process_id
private

◆ _local_assemblers

std::vector<std::unique_ptr<HTLocalAssemblerInterface> > ProcessLib::HT::HTProcess::_local_assemblers
private

◆ _process_data

HTProcessData ProcessLib::HT::HTProcess::_process_data
private

Definition at line 119 of file HTProcess.h.

Referenced by initializeConcreteProcess().

◆ _surfaceflux

std::unique_ptr<ProcessLib::SurfaceFluxData> ProcessLib::HT::HTProcess::_surfaceflux
private

Definition at line 126 of file HTProcess.h.

Referenced by postTimestepConcreteProcess().

◆ _xs_previous_timestep

std::array<std::unique_ptr<GlobalVector>, 2> ProcessLib::HT::HTProcess::_xs_previous_timestep
private

Solutions of the previous time step.

Definition at line 124 of file HTProcess.h.

Referenced by preTimestepConcreteProcess(), and setCoupledSolutionsOfPreviousTimeStepPerProcess().


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