OGS
HTProcess.cpp
Go to the documentation of this file.
1
11#include "HTProcess.h"
12
13#include <cassert>
14
15#include "MonolithicHTFEM.h"
20#include "StaggeredHTFEM.h"
21
22namespace ProcessLib
23{
24namespace HT
25{
27 std::string name,
28 MeshLib::Mesh& mesh,
29 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
30 std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
31 unsigned const integration_order,
32 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
33 process_variables,
34 HTProcessData&& process_data,
35 SecondaryVariableCollection&& secondary_variables,
36 bool const use_monolithic_scheme,
37 std::unique_ptr<ProcessLib::SurfaceFluxData>&& surfaceflux)
38 : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
39 integration_order, std::move(process_variables),
40 std::move(secondary_variables), use_monolithic_scheme),
41 _process_data(std::move(process_data)),
42 _surfaceflux(std::move(surfaceflux))
43{
44}
45
47 NumLib::LocalToGlobalIndexMap const& dof_table,
48 MeshLib::Mesh const& mesh,
49 unsigned const integration_order)
50{
52 {
53 ProcessLib::createLocalAssemblers<MonolithicHTFEM>(
54 mesh.getDimension(), mesh.getElements(), dof_table,
55 _local_assemblers, mesh.isAxiallySymmetric(), integration_order,
57 }
58 else
59 {
60 ProcessLib::createLocalAssemblers<StaggeredHTFEM>(
61 mesh.getDimension(), mesh.getElements(), dof_table,
62 _local_assemblers, mesh.isAxiallySymmetric(), integration_order,
64 }
65
67 "darcy_velocity",
71}
72
73void HTProcess::assembleConcreteProcess(const double t, double const dt,
74 std::vector<GlobalVector*> const& x,
75 std::vector<GlobalVector*> const& xdot,
76 int const process_id, GlobalMatrix& M,
78{
79 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
80 dof_tables;
82 {
83 DBUG("Assemble HTProcess.");
84 dof_tables.emplace_back(*_local_to_global_index_map);
85 }
86 else
87 {
89 {
90 DBUG(
91 "Assemble the equations of heat transport process within "
92 "HTProcess.");
93 }
94 else
95 {
96 DBUG(
97 "Assemble the equations of single phase fully saturated "
98 "fluid flow process within HTProcess.");
99 }
100 dof_tables.emplace_back(*_local_to_global_index_map);
101 dof_tables.emplace_back(*_local_to_global_index_map);
102 }
103
104 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
105 // Call global assembler for each local assembly item.
108 pv.getActiveElementIDs(), dof_tables, t, dt, x, xdot, process_id, M, K,
109 b);
110}
111
113 const double t, double const dt, std::vector<GlobalVector*> const& x,
114 std::vector<GlobalVector*> const& xdot, int const process_id,
116{
117 DBUG("AssembleWithJacobian HTProcess.");
118
119 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
120 dof_tables;
122 {
123 dof_tables.emplace_back(std::ref(*_local_to_global_index_map));
124 }
125 else
126 {
127 dof_tables.emplace_back(std::ref(*_local_to_global_index_map));
128 dof_tables.emplace_back(std::ref(*_local_to_global_index_map));
129 }
130
131 // Call global assembler for each local assembly item.
132 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
135 _local_assemblers, pv.getActiveElementIDs(), dof_tables, t, dt, x, xdot,
136 process_id, M, K, b, Jac);
137}
138
140 int const process_id)
141{
142 DBUG("Set the coupled term for the staggered scheme to local assemblers.");
143
144 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
148}
149
150std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
152{
154 {
155 // For single-variable-single-component processes reuse the existing DOF
156 // table.
157 const bool manage_storage = false;
158 return std::make_tuple(_local_to_global_index_map.get(),
159 manage_storage);
160 }
161
162 // Otherwise construct a new DOF table.
163 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
165
166 const bool manage_storage = true;
167 return std::make_tuple(new NumLib::LocalToGlobalIndexMap(
168 std::move(all_mesh_subsets_single_component),
169 // by location order is needed for output
171 manage_storage);
172}
173
174Eigen::Vector3d HTProcess::getFlux(std::size_t element_id,
175 MathLib::Point3d const& p,
176 double const t,
177 std::vector<GlobalVector*> const& x) const
178{
179 // fetch local_x from primary variable
180 std::vector<GlobalIndexType> indices_cache;
181 auto const r_c_indices = NumLib::getRowColumnIndices(
182 element_id, *_local_to_global_index_map, indices_cache);
183 std::vector<std::vector<GlobalIndexType>> indices_of_all_coupled_processes{
184 x.size(), r_c_indices.rows};
185 auto const local_x =
186 getCoupledLocalSolutions(x, indices_of_all_coupled_processes);
187
188 return _local_assemblers[element_id]->getFlux(p, t, local_x);
189}
190
191// this is almost a copy of the implementation in the GroundwaterFlow
192void HTProcess::postTimestepConcreteProcess(std::vector<GlobalVector*> const& x,
193 const double t,
194 const double /*delta_t*/,
195 int const process_id)
196{
197 // For the monolithic scheme, process_id is always zero.
198 if (_use_monolithic_scheme && process_id != 0)
199 {
200 OGS_FATAL(
201 "The condition of process_id = 0 must be satisfied for monolithic "
202 "HTProcess, which is a single process.");
203 }
206 {
207 DBUG("This is the thermal part of the staggered HTProcess.");
208 return;
209 }
210 if (!_surfaceflux) // computing the surfaceflux is optional
211 {
212 return;
213 }
214
215 ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
216
217 _surfaceflux->integrate(x, t, *this, process_id, _integration_order, _mesh,
219}
220} // namespace HT
221} // namespace ProcessLib
#define OGS_FATAL(...)
Definition: Error.h:26
void DBUG(char const *fmt, Args const &... args)
Definition: Logging.h:29
Global vector based on Eigen vector.
Definition: EigenVector.h:28
bool isAxiallySymmetric() const
Definition: Mesh.h:131
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:103
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition: Mesh.h:76
void setStaggeredCoupledSolutions(std::size_t const, CoupledSolutionsForStaggeredScheme *const coupling_term)
virtual std::vector< double > const & getIntPtDarcyVelocity(const double, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) const =0
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, const double t, const double delta_t, int const process_id) override
Definition: HTProcess.cpp:192
void setCoupledTermForTheStaggeredSchemeToLocalAssemblers(int const process_id) override
Definition: HTProcess.cpp:139
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
Definition: HTProcess.cpp:46
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
Definition: HTProcess.cpp:73
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const override
Definition: HTProcess.cpp:151
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, bool const use_monolithic_scheme, std::unique_ptr< ProcessLib::SurfaceFluxData > &&surfaceflux)
Definition: HTProcess.cpp:26
std::unique_ptr< ProcessLib::SurfaceFluxData > _surfaceflux
Definition: HTProcess.h:117
Eigen::Vector3d getFlux(std::size_t element_id, MathLib::Point3d const &p, double const t, std::vector< GlobalVector * > const &x) const override
Definition: HTProcess.cpp:174
HTProcessData _process_data
Definition: HTProcess.h:113
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
Definition: HTProcess.cpp:112
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
Definition: HTProcess.h:115
std::vector< std::size_t > const & getActiveElementIDs() const
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:329
MeshLib::Mesh & _mesh
Definition: Process.h:328
SecondaryVariableCollection _secondary_variables
Definition: Process.h:333
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:341
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:146
VectorMatrixAssembler _global_assembler
Definition: Process.h:335
unsigned const _integration_order
Definition: Process.h:346
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:331
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:182
const bool _use_monolithic_scheme
Definition: Process.h:337
Handles configuration of several secondary variables from the project file.
void addSecondaryVariable(std::string const &internal_name, SecondaryVariableFunctions &&fcts)
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 &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)
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
NumLib::LocalToGlobalIndexMap::RowColumnIndices getRowColumnIndices(std::size_t const id, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< GlobalIndexType > &indices)
@ BY_LOCATION
Ordering data by spatial location.
std::vector< double > getCoupledLocalSolutions(std::vector< GlobalVector * > const &global_solutions, std::vector< std::vector< GlobalIndexType > > const &indices)
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)