OGS
LocalAssemblerInterface.h
Go to the documentation of this file.
1
11#pragma once
12
13#include <Eigen/Core>
14
15#include "MathLib/Point3d.h"
17
18namespace NumLib
19{
20class LocalToGlobalIndexMap;
21} // namespace NumLib
22
23namespace ProcessLib
24{
25
32{
33public:
34 virtual ~LocalAssemblerInterface() = default;
35
36 void setInitialConditions(std::size_t const mesh_item_id,
37 NumLib::LocalToGlobalIndexMap const& dof_table,
38 GlobalVector const& x, double const t,
39 bool const use_monolithic_scheme,
40 int const process_id);
41
42 virtual void initialize(std::size_t const mesh_item_id,
43 NumLib::LocalToGlobalIndexMap const& dof_table);
44
45 virtual void preAssemble(double const /*t*/, double const /*dt*/,
46 std::vector<double> const& /*local_x*/){};
47
48 virtual void assemble(double const t, double const dt,
49 std::vector<double> const& local_x,
50 std::vector<double> const& local_x_prev,
51 std::vector<double>& local_M_data,
52 std::vector<double>& local_K_data,
53 std::vector<double>& local_b_data);
54
55 virtual void assembleForStaggeredScheme(double const t, double const dt,
56 Eigen::VectorXd const& local_x,
57 Eigen::VectorXd const& local_x_prev,
58 int const process_id,
59 std::vector<double>& local_M_data,
60 std::vector<double>& local_K_data,
61 std::vector<double>& local_b_data);
62
63 virtual void assembleWithJacobian(double const t, double const dt,
64 std::vector<double> const& local_x,
65 std::vector<double> const& local_x_prev,
66 std::vector<double>& local_M_data,
67 std::vector<double>& local_K_data,
68 std::vector<double>& local_b_data,
69 std::vector<double>& local_Jac_data);
70
72 double const t, double const dt, Eigen::VectorXd const& local_x,
73 Eigen::VectorXd const& local_x_prev, int const process_id,
74 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
75 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
76
77 virtual void computeSecondaryVariable(
78 std::size_t const mesh_item_id,
79 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
80 double const t, double const dt, std::vector<GlobalVector*> const& x,
81 GlobalVector const& x_prev, int const process_id);
82
83 virtual void preTimestep(std::size_t const mesh_item_id,
84 NumLib::LocalToGlobalIndexMap const& dof_table,
85 GlobalVector const& x, double const t,
86 double const delta_t);
87
88 virtual void postTimestep(
89 std::size_t const mesh_item_id,
90 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
91 std::vector<GlobalVector*> const& x,
92 std::vector<GlobalVector*> const& x_prev, double const t,
93 double const dt, bool const use_monolithic_scheme,
94 int const process_id);
95
96 void postNonLinearSolver(std::size_t const mesh_item_id,
97 NumLib::LocalToGlobalIndexMap const& dof_table,
98 std::vector<GlobalVector*> const& x,
99 std::vector<GlobalVector*> const& x_prev,
100 double const t, double const dt,
101 bool const use_monolithic_scheme,
102 int const process_id);
103
107 virtual Eigen::Vector3d getFlux(
108 MathLib::Point3d const& /*p_local_coords*/,
109 double const /*t*/,
110 std::vector<double> const& /*local_x*/) const
111 {
112 return Eigen::Vector3d{};
113 }
114
116 virtual Eigen::Vector3d getFlux(
117 MathLib::Point3d const& /*p_local_coords*/,
118 double const /*t*/,
119 std::vector<std::vector<double>> const& /*local_xs*/) const
120 {
121 return Eigen::Vector3d{};
122 }
123
124private:
126 std::vector<double> const& /*local_x*/, double const /*t*/,
127 bool const /*use_monolithic_scheme*/, int const /*process_id*/)
128 {
129 }
130
131 virtual void initializeConcrete() {}
132
133 virtual void preTimestepConcrete(std::vector<double> const& /*local_x*/,
134 double const /*t*/, double const /*dt*/)
135 {
136 }
137
138 virtual void postTimestepConcrete(Eigen::VectorXd const& /*local_x*/,
139 Eigen::VectorXd const& /*local_x_prev*/,
140 double const /*t*/, double const /*dt*/,
141 bool const /*use_monolithic_scheme*/,
142 int const /*process_id*/)
143 {
144 }
145
147 std::vector<double> const& /*local_x*/,
148 std::vector<double> const& /*local_x_prev*/, double const /*t*/,
149 double const /*dt*/, bool const /*use_monolithic_scheme*/,
150 int const /*process_id*/)
151 {
152 }
153
155 double const /*t*/,
156 double const /*dt*/,
157 Eigen::VectorXd const& /*local_x*/,
158 Eigen::VectorXd const& /*local_x_prev*/)
159 {
160 }
161
163 double const /*t*/, std::vector<std::vector<double>> const&
164 /*coupled_local_solutions*/)
165 {
166 }
167};
168
169} // namespace ProcessLib
Definition of the Point3d class.
Global vector based on Eigen vector.
Definition: EigenVector.h:25
virtual void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
virtual void preTimestepConcrete(std::vector< double > const &, double const, double const)
virtual ~LocalAssemblerInterface()=default
virtual void postTimestepConcrete(Eigen::VectorXd const &, Eigen::VectorXd const &, double const, double const, bool const, int const)
virtual void setInitialConditionsConcrete(std::vector< double > const &, double const, bool const, int const)
virtual void postTimestep(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, bool const use_monolithic_scheme, int const process_id)
virtual void computeSecondaryVariableConcrete(double const, double const, Eigen::VectorXd const &, Eigen::VectorXd const &)
virtual Eigen::Vector3d getFlux(MathLib::Point3d const &, double const, std::vector< double > const &) const
virtual Eigen::Vector3d getFlux(MathLib::Point3d const &, double const, std::vector< std::vector< double > > const &) const
Fits to staggered scheme.
virtual void assembleWithJacobianForStaggeredScheme(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
virtual void computeSecondaryVariable(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
virtual void preAssemble(double const, double const, std::vector< double > const &)
virtual void assembleForStaggeredScheme(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
virtual void postNonLinearSolverConcrete(std::vector< double > const &, std::vector< double > const &, double const, double const, bool const, int const)
virtual void preTimestep(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, double const delta_t)
virtual void initialize(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
void setInitialConditions(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, bool const use_monolithic_scheme, int const process_id)
virtual void computeSecondaryVariableWithCoupledProcessConcrete(double const, std::vector< std::vector< double > > const &)
virtual void assemble(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
void postNonLinearSolver(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, bool const use_monolithic_scheme, int const process_id)