OGS
Process.h
Go to the documentation of this file.
1
11#pragma once
12
13#include <string>
14#include <tuple>
15
26#include "ProcessVariable.h"
29#include "processlib_export.h"
30
31namespace MeshLib
32{
33class Mesh;
34}
35
36namespace ProcessLib
37{
38
40 : public NumLib::ODESystem< // TODO: later on use a simpler ODE system
41 NumLib::ODESystemTag::FirstOrderImplicitQuasilinear,
42 NumLib::NonlinearSolverTag::Newton>,
44{
45public:
47
48 Process(std::string name_,
49 MeshLib::Mesh& mesh,
50 std::unique_ptr<AbstractJacobianAssembler>&& jacobian_assembler,
51 std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const&
52 parameters,
53 unsigned const integration_order,
54 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
55 process_variables,
56 SecondaryVariableCollection&& secondary_variables,
57 const bool use_monolithic_scheme = true);
58
60 void preTimestep(std::vector<GlobalVector*> const& x, const double t,
61 const double delta_t, const int process_id);
62
64 void postTimestep(std::vector<GlobalVector*> const& x,
65 std::vector<GlobalVector*> const& x_prev, const double t,
66 const double delta_t, int const process_id);
67
70 void postNonLinearSolver(std::vector<GlobalVector*> const& x,
71 std::vector<GlobalVector*> const& x_prev,
72 const double t, double const dt,
73 int const process_id);
74
75 void preIteration(const unsigned iter, GlobalVector const& x) final;
76
78 void computeSecondaryVariable(double const t,
79 double const dt,
80 std::vector<GlobalVector*> const& x,
81 GlobalVector const& x_prev,
82 int const process_id);
83
85
86 void initialize(
87 std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const&
88 media);
89
91 std::vector<GlobalVector*>& process_solutions,
92 std::vector<GlobalVector*> const& process_solutions_prev,
93 double const t,
94 int const process_id);
95
97 const int process_id) const override;
98
99 void updateDeactivatedSubdomains(double const time, const int process_id);
100
101 virtual bool isMonolithicSchemeUsed() const
102 {
104 }
105
107 const double /*t*/,
108 std::vector<GlobalVector*> const& /*integration_point_values_vectors*/,
109 std::vector<GlobalVector*>& /*nodal_values_vectors*/)
110 {
111 }
112
113 void preAssemble(const double t, double const dt,
114 GlobalVector const& x) final;
115
116 void assemble(const double t, double const dt,
117 std::vector<GlobalVector*> const& x,
118 std::vector<GlobalVector*> const& x_prev,
119 int const process_id, GlobalMatrix& M, GlobalMatrix& K,
120 GlobalVector& b) final;
121
122 void assembleWithJacobian(const double t, double const dt,
123 std::vector<GlobalVector*> const& x,
124 std::vector<GlobalVector*> const& x_prev,
125 int const process_id, GlobalMatrix& M,
127 GlobalMatrix& Jac) final;
128
129 /* Computes data necessary for output.
130 *
131 * \pre Must be called before postTimestep() since the latter will overwrite
132 * previous states with current ones possibly affecting the data computed by
133 * this method.
134 */
135 void preOutput(const double t, double const dt,
136 std::vector<GlobalVector*> const& x,
137 std::vector<GlobalVector*> const& x_prev,
138 int const process_id);
139
140 std::vector<NumLib::IndexValueVector<GlobalIndexType>> const*
141 getKnownSolutions(double const t, GlobalVector const& x,
142 int const process_id) const final
143 {
144 return _boundary_conditions[process_id].getKnownSolutions(t, x);
145 }
146
148 const int /*process_id*/) const
149 {
151 }
152
153 MeshLib::Mesh& getMesh() const { return _mesh; }
154 std::vector<std::reference_wrapper<ProcessVariable>> const&
155 getProcessVariables(const int process_id) const
156 {
157 return _process_variables[process_id];
158 }
159
161 {
163 }
164
165 std::vector<std::unique_ptr<MeshLib::IntegrationPointWriter>> const&
167 {
169 }
170
171 // Used as a call back for CalculateSurfaceFlux process.
172
173 virtual Eigen::Vector3d getFlux(
174 std::size_t /*element_id*/,
175 MathLib::Point3d const& /*p*/,
176 double const /*t*/,
177 std::vector<GlobalVector*> const& /*x*/) const
178 {
179 return Eigen::Vector3d{};
180 }
181
183 std::vector<GlobalVector*>& /*x*/,
184 std::vector<GlobalVector*> const& /*x_prev*/, double const /*t*/,
185 double const /*dt*/, NumLib::EquationSystem& /*ode_sys*/,
186 int const /*process_id*/)
187 {
188 }
189
190protected:
192 {
194 }
195
197 {
199 }
200
207 const NumLib::LocalToGlobalIndexMap& dof_table, const int process_id,
208 std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const&
209 media);
210
211private:
214 NumLib::LocalToGlobalIndexMap const& dof_table,
215 MeshLib::Mesh const& mesh,
216 unsigned const integration_order) = 0;
217
220 virtual void initializeBoundaryConditions(
221 std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const&
222 media);
223
225 std::vector<GlobalVector*>& /*x*/,
226 double const /*t*/,
227 int const /*process_id*/)
228 {
229 }
230
231 virtual void preAssembleConcreteProcess(const double /*t*/,
232 double const /*dt*/,
233 GlobalVector const& /*x*/)
234 {
235 }
236
238 const double t, double const dt, std::vector<GlobalVector*> const& x,
239 std::vector<GlobalVector*> const& x_prev, int const process_id,
240 GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b) = 0;
241
243 const double t, double const dt, std::vector<GlobalVector*> const& x,
244 std::vector<GlobalVector*> const& x_prev, int const process_id,
246 GlobalMatrix& Jac) = 0;
247
249 std::vector<GlobalVector*> const& /*x*/,
250 const double /*t*/,
251 const double /*dt*/,
252 const int /*process_id*/)
253 {
254 }
255
257 std::vector<GlobalVector*> const& /*x*/,
258 std::vector<GlobalVector*> const& /*x_prev*/,
259 const double /*t*/,
260 const double /*dt*/,
261 int const /*process_id*/)
262 {
263 }
264
266 std::vector<GlobalVector*> const& /*x*/,
267 std::vector<GlobalVector*> const& /*x_prev*/, const double /*t*/,
268 double const /*dt*/, int const /*process_id*/)
269 {
270 }
271
272 virtual void preIterationConcreteProcess(const unsigned /*iter*/,
273 GlobalVector const& /*x*/)
274 {
275 }
276
278 double const /*t*/,
279 double const /*dt*/,
280 std::vector<GlobalVector*> const& /*x*/,
281 GlobalVector const& /*x_prev*/,
282 int const /*process_id*/)
283 {
284 }
285
287 GlobalVector const& /*x*/)
288 {
290 }
291
293 const double /*t*/, double const /*dt*/,
294 std::vector<GlobalVector*> const& /*x*/,
295 std::vector<GlobalVector*> const& /*x_prev*/, int const /*process_id*/)
296 {
297 }
298
299protected:
306 virtual void constructDofTable();
307
313
320 const int specified_process_id);
321
330 virtual std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
332
335 std::vector<GlobalIndexType>
337
338private:
340
344
345public:
346 std::string const name;
347
348protected:
350 std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_all_nodes;
351
352 std::unique_ptr<NumLib::LocalToGlobalIndexMap> _local_to_global_index_map;
353
355
356 // TODO (CL) switch to the parallel vector matrix assembler here, once
357 // feature complete, or use the assembly mixin and remove these two members.
358 std::unique_ptr<ProcessLib::AbstractJacobianAssembler> _jacobian_assembler;
360
362
366 unsigned const _integration_order;
367
371 std::vector<std::unique_ptr<MeshLib::IntegrationPointWriter>>
373
375
376protected:
381 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>
383
387 std::vector<BoundaryConditionCollection> _boundary_conditions;
388
389private:
393 std::vector<SourceTermCollection> _source_term_collections;
394
396};
397
398} // namespace ProcessLib
MathLib::SparsityPattern< GlobalIndexType > GlobalSparsityPattern
Global vector based on Eigen vector.
Definition: EigenVector.h:25
NumLib::LocalToGlobalIndexMap const & getDOFTable() const
NumLib::Extrapolator & getExtrapolator() const
std::string const name
Definition: Process.h:346
std::vector< NumLib::IndexValueVector< GlobalIndexType > > const * getKnownSolutions(double const t, GlobalVector const &x, int const process_id) const final
Definition: Process.h:141
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > const & getIntegrationPointWriters() const
Definition: Process.h:166
MeshLib::Mesh & getMesh() const
Definition: Process.h:153
virtual NumLib::IterationResult postIterationConcreteProcess(GlobalVector const &)
Definition: Process.h:286
virtual 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)=0
virtual Eigen::Vector3d getFlux(std::size_t, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &) const
Definition: Process.h:173
virtual void preTimestepConcreteProcess(std::vector< GlobalVector * > const &, const double, const double, const int)
Definition: Process.h:248
std::vector< BoundaryConditionCollection > _boundary_conditions
Definition: Process.h:387
std::vector< std::unique_ptr< MeshLib::IntegrationPointWriter > > _integration_point_writer
Definition: Process.h:372
SecondaryVariableCollection const & getSecondaryVariables() const
Definition: Process.h:160
virtual void preIterationConcreteProcess(const unsigned, GlobalVector const &)
Definition: Process.h:272
void constructDofTableOfSpecifiedProcessStaggeredScheme(const int specified_process_id)
Definition: Process.cpp:309
static PROCESSLIB_EXPORT const std::string constant_one_parameter_name
Definition: Process.h:46
void preAssemble(const double t, double const dt, GlobalVector const &x) final
Definition: Process.cpp:207
virtual void initializeBoundaryConditions(std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
Definition: Process.cpp:83
void initialize(std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
Definition: Process.cpp:97
void initializeProcessBoundaryConditionsAndSourceTerms(const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
Definition: Process.cpp:68
virtual void solveReactionEquation(std::vector< GlobalVector * > &, std::vector< GlobalVector * > const &, double const, double const, NumLib::EquationSystem &, int const)
Definition: Process.h:182
virtual void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order)=0
Process specific initialization called by initialize().
virtual bool isMonolithicSchemeUsed() const
Definition: Process.h:101
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
Definition: Process.cpp:237
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
Definition: Process.cpp:213
virtual void postTimestepConcreteProcess(std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, const double, int const)
Definition: Process.h:256
virtual void computeSecondaryVariableConcrete(double const, double const, std::vector< GlobalVector * > const &, GlobalVector const &, int const)
Definition: Process.h:277
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:350
void preOutput(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, int const process_id)
Definition: Process.cpp:459
MeshLib::Mesh & _mesh
Definition: Process.h:349
virtual 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)=0
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
Definition: Process.h:382
NumLib::LocalToGlobalIndexMap const & getSingleComponentDOFTable() const
Definition: Process.h:196
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const
Definition: Process.h:147
virtual void setInitialConditionsConcreteProcess(std::vector< GlobalVector * > &, double const, int const)
Definition: Process.h:224
NumLib::IterationResult postIteration(GlobalVector const &x) final
Definition: Process.cpp:453
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.
Definition: Process.cpp:398
SecondaryVariableCollection _secondary_variables
Definition: Process.h:354
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.
Definition: Process.cpp:388
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.
Definition: Process.cpp:434
void initializeExtrapolator()
Definition: Process.cpp:365
virtual void postNonLinearSolverConcreteProcess(std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, const double, double const, int const)
Definition: Process.h:265
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:155
virtual void constructDofTable()
Definition: Process.cpp:262
void constructMonolithicProcessDofTable()
Definition: Process.cpp:276
VectorMatrixAssembler _global_assembler
Definition: Process.h:359
virtual void preAssembleConcreteProcess(const double, double const, GlobalVector const &)
Definition: Process.h:231
void postNonLinearSolver(std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, const double t, double const dt, int const process_id)
Definition: Process.cpp:417
std::vector< GlobalIndexType > getIndicesOfResiduumWithoutInitialCompensation() const override
Definition: Process.cpp:471
unsigned const _integration_order
Definition: Process.h:366
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:352
std::vector< SourceTermCollection > _source_term_collections
Definition: Process.h:393
virtual void extrapolateIntegrationPointValuesToNodes(const double, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > &)
Definition: Process.h:106
std::unique_ptr< ProcessLib::AbstractJacobianAssembler > _jacobian_assembler
Definition: Process.h:358
void updateDeactivatedSubdomains(double const time, const int process_id)
Definition: Process.cpp:197
void computeSparsityPattern()
Definition: Process.cpp:382
MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const override
Definition: Process.cpp:189
void setInitialConditions(std::vector< GlobalVector * > &process_solutions, std::vector< GlobalVector * > const &process_solutions_prev, double const t, int const process_id)
Definition: Process.cpp:118
virtual void preOutputConcreteProcess(const double, double const, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > const &, int const)
Definition: Process.h:292
ExtrapolatorData _extrapolator_data
Definition: Process.h:395
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:191
GlobalSparsityPattern _sparsity_pattern
Definition: Process.h:374
void preIteration(const unsigned iter, GlobalVector const &x) final
Definition: Process.cpp:447
const bool _use_monolithic_scheme
Definition: Process.h:361
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const
Definition: Process.cpp:341
Handles configuration of several secondary variables from the project file.
IterationResult
Status flags telling the NonlinearSolver if an iteration succeeded.
#define PROCESSLIB_EXPORT