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} // NumLib
22
23namespace ProcessLib
24{
25struct CoupledSolutionsForStaggeredScheme;
26struct LocalCoupledSolutions;
27
34{
35public:
36 virtual ~LocalAssemblerInterface() = default;
37
38 void setInitialConditions(std::size_t const mesh_item_id,
39 NumLib::LocalToGlobalIndexMap const& dof_table,
40 GlobalVector const& x, double const t,
41 bool const use_monolithic_scheme,
42 int const process_id);
43
44 virtual void initialize(std::size_t const mesh_item_id,
45 NumLib::LocalToGlobalIndexMap const& dof_table);
46
47 virtual void preAssemble(double const /*t*/, double const /*dt*/,
48 std::vector<double> const& /*local_x*/){};
49
50 virtual void assemble(double const t, double const dt,
51 std::vector<double> const& local_x,
52 std::vector<double> const& local_xdot,
53 std::vector<double>& local_M_data,
54 std::vector<double>& local_K_data,
55 std::vector<double>& local_b_data);
56
57 virtual void assembleForStaggeredScheme(double const t, double const dt,
58 Eigen::VectorXd const& local_x,
59 Eigen::VectorXd const& local_xdot,
60 int const process_id,
61 std::vector<double>& local_M_data,
62 std::vector<double>& local_K_data,
63 std::vector<double>& local_b_data);
64
65 virtual void assembleWithJacobian(double const t, double const dt,
66 std::vector<double> const& local_x,
67 std::vector<double> const& local_xdot,
68 std::vector<double>& local_M_data,
69 std::vector<double>& local_K_data,
70 std::vector<double>& local_b_data,
71 std::vector<double>& local_Jac_data);
72
74 double const t, double const dt, Eigen::VectorXd const& local_x,
75 Eigen::VectorXd const& local_xdot, int const process_id,
76 std::vector<double>& local_M_data, std::vector<double>& local_K_data,
77 std::vector<double>& local_b_data, std::vector<double>& local_Jac_data);
78
79 virtual void computeSecondaryVariable(
80 std::size_t const mesh_item_id,
81 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
82 double const t, double const dt, std::vector<GlobalVector*> const& x,
83 GlobalVector const& x_dot, int const process_id);
84
85 virtual void preTimestep(std::size_t const mesh_item_id,
86 NumLib::LocalToGlobalIndexMap const& dof_table,
87 GlobalVector const& x, double const t,
88 double const delta_t);
89
90 virtual void postTimestep(
91 std::size_t const mesh_item_id,
92 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
93 std::vector<GlobalVector*> const& x, double const t, double const dt);
94
95 void postNonLinearSolver(std::size_t const mesh_item_id,
96 NumLib::LocalToGlobalIndexMap const& dof_table,
97 GlobalVector const& x, GlobalVector const& xdot,
98 double const t, double const dt,
99 bool const use_monolithic_scheme,
100 int const process_id);
101
105 virtual Eigen::Vector3d getFlux(
106 MathLib::Point3d const& /*p_local_coords*/,
107 double const /*t*/,
108 std::vector<double> const& /*local_x*/) const
109 {
110 return Eigen::Vector3d{};
111 }
112
114 virtual Eigen::Vector3d getFlux(
115 MathLib::Point3d const& /*p_local_coords*/,
116 double const /*t*/,
117 std::vector<std::vector<double>> const& /*local_xs*/) const
118 {
119 return Eigen::Vector3d{};
120 }
121
122private:
124 std::vector<double> const& /*local_x*/, double const /*t*/,
125 bool const /*use_monolithic_scheme*/, int const /*process_id*/)
126 {
127 }
128
129 virtual void initializeConcrete() {}
130
131 virtual void preTimestepConcrete(std::vector<double> const& /*local_x*/,
132 double const /*t*/, double const /*dt*/)
133 {
134 }
135
136 virtual void postTimestepConcrete(Eigen::VectorXd const& /*local_x*/,
137 double const /*t*/, double const /*dt*/)
138 {
139 }
140
142 std::vector<double> const& /*local_x*/,
143 std::vector<double> const& /*local_xdot*/, double const /*t*/,
144 double const /*dt*/, bool const /*use_monolithic_scheme*/,
145 int const /*process_id*/)
146 {
147 }
148
150 double const /*t*/,
151 double const /*dt*/,
152 Eigen::VectorXd const& /*local_x*/,
153 Eigen::VectorXd const& /*local_x_dot*/)
154 {
155 }
156
158 double const /*t*/, std::vector<std::vector<double>> const&
159 /*coupled_local_solutions*/)
160 {
161 }
162};
163
164} // namespace ProcessLib
Definition of the Point3d class.
Global vector based on Eigen vector.
Definition: EigenVector.h:25
virtual void preTimestepConcrete(std::vector< double > const &, double const, double const)
virtual ~LocalAssemblerInterface()=default
virtual void assembleForStaggeredScheme(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
virtual void setInitialConditionsConcrete(std::vector< double > const &, double const, bool const, int const)
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 preAssemble(double const, double const, std::vector< double > const &)
virtual void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, 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 postNonLinearSolverConcrete(std::vector< double > const &, std::vector< double > const &, 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, double const t, double const dt)
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_xdot, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_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_dot, int const process_id)
virtual void assembleWithJacobianForStaggeredScheme(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, 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)
void postNonLinearSolver(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, GlobalVector const &xdot, double const t, double const dt, bool const use_monolithic_scheme, int const process_id)
virtual void postTimestepConcrete(Eigen::VectorXd const &, double const, double const)