OGS
HTProcess.cpp
Go to the documentation of this file.
1
11#include "HTProcess.h"
12
13#include <cassert>
14
15#include "MonolithicHTFEM.h"
21#include "StaggeredHTFEM.h"
22
23namespace ProcessLib
24{
25namespace HT
26{
28 std::string name,
29 MeshLib::Mesh& mesh,
30 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
31 std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
32 unsigned const integration_order,
33 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
34 process_variables,
35 HTProcessData&& process_data,
36 SecondaryVariableCollection&& secondary_variables,
37 bool const use_monolithic_scheme,
38 std::unique_ptr<ProcessLib::SurfaceFluxData>&& surfaceflux)
39 : Process(std::move(name), mesh, std::move(jacobian_assembler), parameters,
40 integration_order, std::move(process_variables),
41 std::move(secondary_variables), use_monolithic_scheme),
42 _process_data(std::move(process_data)),
43 _surfaceflux(std::move(surfaceflux))
44{
45}
46
48 NumLib::LocalToGlobalIndexMap const& dof_table,
49 MeshLib::Mesh const& mesh,
50 unsigned const integration_order)
51{
52 int const mesh_space_dimension = _process_data.mesh_space_dimension;
53
55 {
56 ProcessLib::createLocalAssemblers<MonolithicHTFEM>(
57 mesh_space_dimension, mesh.getElements(), dof_table,
60 }
61 else
62 {
63 ProcessLib::createLocalAssemblers<StaggeredHTFEM>(
64 mesh_space_dimension, mesh.getElements(), dof_table,
67 }
68
70 "darcy_velocity",
71 makeExtrapolator(mesh_space_dimension, getExtrapolator(),
74}
75
77 const double t, double const dt, std::vector<GlobalVector*> const& x,
78 std::vector<GlobalVector*> const& x_prev, int const process_id,
80{
81 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
83 {
84 DBUG("Assemble HTProcess.");
85 dof_tables.emplace_back(_local_to_global_index_map.get());
86 }
87 else
88 {
90 {
91 DBUG(
92 "Assemble the equations of heat transport process within "
93 "HTProcess.");
94 }
95 else
96 {
97 DBUG(
98 "Assemble the equations of single phase fully saturated "
99 "fluid flow process within HTProcess.");
100 }
101 dof_tables.emplace_back(_local_to_global_index_map.get());
102 dof_tables.emplace_back(_local_to_global_index_map.get());
103 }
104
105 // Call global assembler for each local assembly item.
108 getActiveElementIDs(), dof_tables, t, dt, x, x_prev, 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& x_prev, int const process_id,
115 GlobalVector& b, GlobalMatrix& Jac)
116{
117 DBUG("AssembleWithJacobian HTProcess.");
118
119 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
121 {
122 dof_tables.emplace_back(_local_to_global_index_map.get());
123 }
124 else
125 {
126 dof_tables.emplace_back(_local_to_global_index_map.get());
127 dof_tables.emplace_back(_local_to_global_index_map.get());
128 }
129
130 // Call global assembler for each local assembly item.
133 _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev,
134 process_id, &b, &Jac);
135}
136
137std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
139{
141 {
142 // For single-variable-single-component processes reuse the existing DOF
143 // table.
144 const bool manage_storage = false;
145 return std::make_tuple(_local_to_global_index_map.get(),
146 manage_storage);
147 }
148
149 // Otherwise construct a new DOF table.
150 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
152
153 const bool manage_storage = true;
154 return std::make_tuple(new NumLib::LocalToGlobalIndexMap(
155 std::move(all_mesh_subsets_single_component),
156 // by location order is needed for output
158 manage_storage);
159}
160
161Eigen::Vector3d HTProcess::getFlux(std::size_t element_id,
162 MathLib::Point3d const& p,
163 double const t,
164 std::vector<GlobalVector*> const& x) const
165{
166 // fetch local_x from primary variable
167 std::vector<GlobalIndexType> indices_cache;
168 auto const r_c_indices = NumLib::getRowColumnIndices(
169 element_id, *_local_to_global_index_map, indices_cache);
170 std::vector<std::vector<GlobalIndexType>> indices_of_all_coupled_processes{
171 x.size(), r_c_indices.rows};
172 auto const local_x =
173 getCoupledLocalSolutions(x, indices_of_all_coupled_processes);
174
175 return _local_assemblers[element_id]->getFlux(p, t, local_x);
176}
177
178// this is almost a copy of the implementation in the GroundwaterFlow
180 std::vector<GlobalVector*> const& x,
181 std::vector<GlobalVector*> const& /*x_prev*/,
182 const double t,
183 const double /*delta_t*/,
184 int const process_id)
185{
186 // For the monolithic scheme, process_id is always zero.
187 if (_use_monolithic_scheme && process_id != 0)
188 {
189 OGS_FATAL(
190 "The condition of process_id = 0 must be satisfied for monolithic "
191 "HTProcess, which is a single process.");
192 }
195 {
196 DBUG("This is the thermal part of the staggered HTProcess.");
197 return;
198 }
199 if (!_surfaceflux) // computing the surfaceflux is optional
200 {
201 return;
202 }
203
204 _surfaceflux->integrate(x, t, *this, process_id, _integration_order, _mesh,
206}
207} // namespace HT
208} // namespace ProcessLib
#define OGS_FATAL(...)
Definition Error.h:26
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
Global vector based on Eigen vector.
Definition EigenVector.h:25
bool isAxiallySymmetric() const
Definition Mesh.h:137
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:109
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 initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
Definition HTProcess.cpp:47
void assembleWithJacobianConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalVector &b, GlobalMatrix &Jac) override
std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const override
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:27
std::unique_ptr< ProcessLib::SurfaceFluxData > _surfaceflux
Definition HTProcess.h:113
Eigen::Vector3d getFlux(std::size_t element_id, MathLib::Point3d const &p, double const t, std::vector< GlobalVector * > const &x) const override
HTProcessData _process_data
Definition HTProcess.h:109
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
Definition HTProcess.cpp:76
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, const double delta_t, int const process_id) override
std::vector< std::unique_ptr< HTLocalAssemblerInterface > > _local_assemblers
Definition HTProcess.h:111
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition Process.h:366
MeshLib::Mesh & _mesh
Definition Process.h:365
std::vector< std::size_t > const & getActiveElementIDs() const
Definition Process.h:167
SecondaryVariableCollection _secondary_variables
Definition Process.h:370
VectorMatrixAssembler _global_assembler
Definition Process.h:377
unsigned const _integration_order
Definition Process.h:384
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition Process.h:368
NumLib::Extrapolator & getExtrapolator() const
Definition Process.h:208
const bool _use_monolithic_scheme
Definition Process.h:379
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< NumLib::LocalToGlobalIndexMap const * > 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)
void assembleWithJacobian(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id, GlobalVector *b, GlobalMatrix *Jac)
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 executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)