OGS
Process.h
Go to the documentation of this file.
1
11#pragma once
12
13#include <string>
14#include <tuple>
15
27#include "ProcessVariable.h"
29
30namespace MeshLib
31{
32class Mesh;
33}
34
35namespace ProcessLib
36{
37struct CoupledSolutionsForStaggeredScheme;
38
40 : public NumLib::ODESystem< // TODO: later on use a simpler ODE system
41 NumLib::ODESystemTag::FirstOrderImplicitQuasilinear,
42 NumLib::NonlinearSolverTag::Newton>
43{
44public:
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, const double t,
65 const double delta_t, int const process_id);
66
69 void postNonLinearSolver(GlobalVector const& x, GlobalVector const& xdot,
70 const double t, double const dt,
71 int const process_id);
72
73 void preIteration(const unsigned iter, GlobalVector const& x) final;
74
76 void computeSecondaryVariable(double const t,
77 double const dt,
78 std::vector<GlobalVector*> const& x,
79 GlobalVector const& x_dot,
80 int const process_id);
81
83
84 void initialize();
85
87 std::vector<GlobalVector*>& process_solutions,
88 std::vector<GlobalVector*> const& process_solutions_prev,
89 double const t,
90 int const process_id);
91
93 const int process_id) const override;
94
96 CoupledSolutionsForStaggeredScheme* const coupled_solutions)
97 {
98 _coupled_solutions = coupled_solutions;
99 }
100
101 void updateDeactivatedSubdomains(double const time, const int process_id);
102
105 int const /*process_id*/)
106 {
107 }
108
110 const double /*t*/,
111 std::vector<GlobalVector*> const& /*integration_point_values_vectors*/,
112 std::vector<GlobalVector*>& /*nodal_values_vectors*/)
113 {
114 }
115
116 void preAssemble(const double t, double const dt,
117 GlobalVector const& x) final;
118
119 void assemble(const double t, double const dt,
120 std::vector<GlobalVector*> const& x,
121 std::vector<GlobalVector*> const& xdot, int const process_id,
122 GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b) final;
123
124 void assembleWithJacobian(const double t, double const dt,
125 std::vector<GlobalVector*> const& x,
126 std::vector<GlobalVector*> const& xdot,
127 int const process_id, GlobalMatrix& M,
129 GlobalMatrix& Jac) final;
130
131 std::vector<NumLib::IndexValueVector<GlobalIndexType>> const*
132 getKnownSolutions(double const t, GlobalVector const& x,
133 int const process_id) const final
134 {
135 return _boundary_conditions[process_id].getKnownSolutions(t, x);
136 }
137
139 const int /*process_id*/) const
140 {
142 }
143
144 MeshLib::Mesh& getMesh() const { return _mesh; }
145 std::vector<std::reference_wrapper<ProcessVariable>> const&
146 getProcessVariables(const int process_id) const
147 {
148 return _process_variables[process_id];
149 }
150
152 {
154 }
155
156 std::vector<std::unique_ptr<IntegrationPointWriter>> const&
158 {
160 }
161
162 // Used as a call back for CalculateSurfaceFlux process.
163
164 virtual Eigen::Vector3d getFlux(
165 std::size_t /*element_id*/,
166 MathLib::Point3d const& /*p*/,
167 double const /*t*/,
168 std::vector<GlobalVector*> const& /*x*/) const
169 {
170 return Eigen::Vector3d{};
171 }
172
174 std::vector<GlobalVector*>& /*x*/,
175 std::vector<GlobalVector*> const& /*x_prev*/, double const /*t*/,
176 double const /*dt*/, NumLib::EquationSystem& /*ode_sys*/,
177 int const /*process_id*/)
178 {
179 }
180
181protected:
183 {
185 }
186
188 {
190 }
191
198 const NumLib::LocalToGlobalIndexMap& dof_table, const int process_id);
199
200private:
203 NumLib::LocalToGlobalIndexMap const& dof_table,
204 MeshLib::Mesh const& mesh,
205 unsigned const integration_order) = 0;
206
209 virtual void initializeBoundaryConditions();
210
212 std::vector<GlobalVector*>& /*x*/,
213 double const /*t*/,
214 int const /*process_id*/)
215 {
216 }
217
218 virtual void preAssembleConcreteProcess(const double /*t*/,
219 double const /*dt*/,
220 GlobalVector const& /*x*/)
221 {
222 }
223
224 virtual void assembleConcreteProcess(const double t, double const dt,
225 std::vector<GlobalVector*> const& x,
226 std::vector<GlobalVector*> const& xdot,
227 int const process_id, GlobalMatrix& M,
228 GlobalMatrix& K, GlobalVector& b) = 0;
229
231 const double t, double const dt, std::vector<GlobalVector*> const& x,
232 std::vector<GlobalVector*> const& xdot, int const process_id,
234 GlobalMatrix& Jac) = 0;
235
237 std::vector<GlobalVector*> const& /*x*/,
238 const double /*t*/,
239 const double /*dt*/,
240 const int /*process_id*/)
241 {
242 }
243
245 std::vector<GlobalVector*> const& /*x*/,
246 const double /*t*/,
247 const double /*dt*/,
248 int const /*process_id*/)
249 {
250 }
251
253 GlobalVector const& /*x*/, GlobalVector const& /*xdot*/,
254 const double /*t*/, double const /*dt*/, int const /*process_id*/)
255 {
256 }
257
258 virtual void preIterationConcreteProcess(const unsigned /*iter*/,
259 GlobalVector const& /*x*/)
260 {
261 }
262
264 double const /*t*/,
265 double const /*dt*/,
266 std::vector<GlobalVector*> const& /*x*/,
267 GlobalVector const& /*x_dot*/,
268 int const /*process_id*/)
269 {
270 }
271
273 GlobalVector const& /*x*/)
274 {
276 }
277
278protected:
285 virtual void constructDofTable();
286
292
299 const int specified_prosess_id);
300
309 virtual std::tuple<NumLib::LocalToGlobalIndexMap*, bool>
311
314 std::vector<GlobalIndexType>
316
317private:
319
323
324public:
325 std::string const name;
326
327protected:
329 std::unique_ptr<MeshLib::MeshSubset const> _mesh_subset_all_nodes;
330
331 std::unique_ptr<NumLib::LocalToGlobalIndexMap> _local_to_global_index_map;
332
334
336
338
342
346 unsigned const _integration_order;
347
351 std::vector<std::unique_ptr<IntegrationPointWriter>>
353
355
356protected:
361 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>
363
367 std::vector<BoundaryConditionCollection> _boundary_conditions;
368
369private:
373 std::vector<SourceTermCollection> _source_term_collections;
374
376};
377
378} // namespace ProcessLib
MathLib::SparsityPattern< GlobalIndexType > GlobalSparsityPattern
Global vector based on Eigen vector.
Definition: EigenVector.h:28
NumLib::LocalToGlobalIndexMap const & getDOFTable() const
NumLib::Extrapolator & getExtrapolator() const
std::string const name
Definition: Process.h:325
std::vector< NumLib::IndexValueVector< GlobalIndexType > > const * getKnownSolutions(double const t, GlobalVector const &x, int const process_id) const final
Definition: Process.h:132
virtual void initializeBoundaryConditions()
Definition: Process.cpp:81
MeshLib::Mesh & getMesh() const
Definition: Process.h:144
virtual NumLib::IterationResult postIterationConcreteProcess(GlobalVector const &)
Definition: Process.h:272
virtual Eigen::Vector3d getFlux(std::size_t, MathLib::Point3d const &, double const, std::vector< GlobalVector * > const &) const
Definition: Process.h:164
virtual void preTimestepConcreteProcess(std::vector< GlobalVector * > const &, const double, const double, const int)
Definition: Process.h:236
std::vector< BoundaryConditionCollection > _boundary_conditions
Definition: Process.h:367
void constructDofTableOfSpecifiedProcessStaggeredScheme(const int specified_prosess_id)
Definition: Process.cpp:305
virtual void postTimestepConcreteProcess(std::vector< GlobalVector * > const &, const double, const double, int const)
Definition: Process.h:244
void postNonLinearSolver(GlobalVector const &x, GlobalVector const &xdot, const double t, double const dt, int const process_id)
Definition: Process.cpp:404
SecondaryVariableCollection const & getSecondaryVariables() const
Definition: Process.h:151
virtual void preIterationConcreteProcess(const unsigned, GlobalVector const &)
Definition: Process.h:258
virtual 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)=0
void assembleWithJacobian(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) final
Definition: Process.cpp:233
void preAssemble(const double t, double const dt, GlobalVector const &x) final
Definition: Process.cpp:203
virtual void solveReactionEquation(std::vector< GlobalVector * > &, std::vector< GlobalVector * > const &, double const, double const, NumLib::EquationSystem &, int const)
Definition: Process.h:173
virtual void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order)=0
Process specific initialization called by initialize().
void initializeProcessBoundaryConditionsAndSourceTerms(const NumLib::LocalToGlobalIndexMap &dof_table, const int process_id)
Definition: Process.cpp:67
void computeSecondaryVariable(double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, int const process_id)
compute secondary variables for the coupled equations or for output.
Definition: Process.cpp:413
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:22
virtual void computeSecondaryVariableConcrete(double const, double const, std::vector< GlobalVector * > const &, GlobalVector const &, int const)
Definition: Process.h:263
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:329
MeshLib::Mesh & _mesh
Definition: Process.h:328
std::vector< std::vector< std::reference_wrapper< ProcessVariable > > > _process_variables
Definition: Process.h:362
NumLib::LocalToGlobalIndexMap const & getSingleComponentDOFTable() const
Definition: Process.h:187
virtual void postNonLinearSolverConcreteProcess(GlobalVector const &, GlobalVector const &, const double, double const, int const)
Definition: Process.h:252
virtual NumLib::LocalToGlobalIndexMap const & getDOFTable(const int) const
Definition: Process.h:138
virtual 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)=0
virtual void setInitialConditionsConcreteProcess(std::vector< GlobalVector * > &, double const, int const)
Definition: Process.h:211
NumLib::IterationResult postIteration(GlobalVector const &x) final
Definition: Process.cpp:432
std::vector< std::unique_ptr< IntegrationPointWriter > > _integration_point_writer
Definition: Process.h:352
SecondaryVariableCollection _secondary_variables
Definition: Process.h:333
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:341
void setCoupledSolutionsForStaggeredScheme(CoupledSolutionsForStaggeredScheme *const coupled_solutions)
Definition: Process.h:95
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:384
void initializeExtrapolator()
Definition: Process.cpp:361
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:146
virtual void constructDofTable()
Definition: Process.cpp:258
void constructMonolithicProcessDofTable()
Definition: Process.cpp:272
VectorMatrixAssembler _global_assembler
Definition: Process.h:335
virtual void preAssembleConcreteProcess(const double, double const, GlobalVector const &)
Definition: Process.h:218
std::vector< GlobalIndexType > getIndicesOfResiduumWithoutInitialCompensation() const override
Definition: Process.cpp:439
unsigned const _integration_order
Definition: Process.h:346
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:331
std::vector< SourceTermCollection > _source_term_collections
Definition: Process.h:373
virtual void extrapolateIntegrationPointValuesToNodes(const double, std::vector< GlobalVector * > const &, std::vector< GlobalVector * > &)
Definition: Process.h:109
void updateDeactivatedSubdomains(double const time, const int process_id)
Definition: Process.cpp:193
void computeSparsityPattern()
Definition: Process.cpp:378
MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const override
Definition: Process.cpp:185
void setInitialConditions(std::vector< GlobalVector * > &process_solutions, std::vector< GlobalVector * > const &process_solutions_prev, double const t, int const process_id)
Definition: Process.cpp:114
bool isMonolithicSchemeUsed() const
Definition: Process.h:103
void assemble(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) final
Definition: Process.cpp:209
void postTimestep(std::vector< GlobalVector * > const &x, const double t, const double delta_t, int const process_id)
Postprocessing after a complete timestep.
Definition: Process.cpp:394
ExtrapolatorData _extrapolator_data
Definition: Process.h:375
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:182
GlobalSparsityPattern _sparsity_pattern
Definition: Process.h:354
virtual void setCoupledTermForTheStaggeredSchemeToLocalAssemblers(int const)
Definition: Process.h:104
std::vector< std::unique_ptr< IntegrationPointWriter > > const & getIntegrationPointWriters() const
Definition: Process.h:157
void preIteration(const unsigned iter, GlobalVector const &x) final
Definition: Process.cpp:426
const bool _use_monolithic_scheme
Definition: Process.h:337
virtual std::tuple< NumLib::LocalToGlobalIndexMap *, bool > getDOFTableForExtrapolatorData() const
Definition: Process.cpp:337
Handles configuration of several secondary variables from the project file.
IterationResult
Status flags telling the NonlinearSolver if an iteration succeeded.
static const double t