OGS
LocalAssemblerInterface.cpp
Go to the documentation of this file.
1
12
13#include <cassert>
14
18
19namespace ProcessLib
20{
22 double const /*t*/,
23 double const /*dt*/,
24 std::vector<double> const& /*local_x*/,
25 std::vector<double> const& /*local_xdot*/,
26 std::vector<double>& /*local_M_data*/,
27 std::vector<double>& /*local_K_data*/,
28 std::vector<double>& /*local_b_data*/)
29{
31 "The assemble() function is not implemented in the local assembler.");
32}
33
35 double const /*t*/, double const /*dt*/, Eigen::VectorXd const& /*local_x*/,
36 Eigen::VectorXd const& /*local_xdot*/, int const /*process_id*/,
37 std::vector<double>& /*local_M_data*/,
38 std::vector<double>& /*local_K_data*/,
39 std::vector<double>& /*local_b_data*/)
40{
42 "The assembleForStaggeredScheme() function is not implemented in the "
43 "local assembler.");
44}
45
47 double const /*t*/, double const /*dt*/,
48 std::vector<double> const& /*local_x*/,
49 std::vector<double> const& /*local_xdot*/,
50 std::vector<double>& /*local_M_data*/,
51 std::vector<double>& /*local_K_data*/,
52 std::vector<double>& /*local_b_data*/,
53 std::vector<double>& /*local_Jac_data*/)
54{
56 "The assembleWithJacobian() function is not implemented in the local "
57 "assembler.");
58}
59
61 double const /*t*/, double const /*dt*/, Eigen::VectorXd const& /*local_x*/,
62 Eigen::VectorXd const& /*local_xdot*/, int const /*process_id*/,
63 std::vector<double>& /*local_M_data*/,
64 std::vector<double>& /*local_K_data*/,
65 std::vector<double>& /*local_b_data*/,
66 std::vector<double>& /*local_Jac_data*/)
67{
69 "The assembleWithJacobianForStaggeredScheme() function is not "
70 "implemented in the local assembler.");
71}
72
74 std::size_t const mesh_item_id,
75 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
76 double const t, double const dt, std::vector<GlobalVector*> const& x,
77 GlobalVector const& x_dot, int const process_id)
78{
79 std::vector<double> local_x_vec;
80
81 auto const n_processes = x.size();
82 for (std::size_t process_id = 0; process_id < n_processes; ++process_id)
83 {
84 auto const indices =
85 NumLib::getIndices(mesh_item_id, *dof_tables[process_id]);
86 assert(!indices.empty());
87 auto const local_solution = x[process_id]->get(indices);
88 local_x_vec.insert(std::end(local_x_vec), std::begin(local_solution),
89 std::end(local_solution));
90 }
91 auto const local_x = MathLib::toVector(local_x_vec);
92
93 // Todo: A more decent way is to directly pass x_dots as done for x
94 auto const indices =
95 NumLib::getIndices(mesh_item_id, *dof_tables[process_id]);
96 auto const local_x_dot_vec = x_dot.get(indices);
97 auto const local_x_dot = MathLib::toVector(local_x_dot_vec);
98
99 computeSecondaryVariableConcrete(t, dt, local_x, local_x_dot);
100}
101
103 std::size_t const mesh_item_id,
104 NumLib::LocalToGlobalIndexMap const& dof_table, GlobalVector const& x,
105 double const t, bool const use_monolithic_scheme, int const process_id)
106{
107 auto const indices = NumLib::getIndices(mesh_item_id, dof_table);
108 auto const local_x = x.get(indices);
109
110 setInitialConditionsConcrete(local_x, t, use_monolithic_scheme, process_id);
111}
112
114 std::size_t const /*mesh_item_id*/,
115 NumLib::LocalToGlobalIndexMap const& /*dof_table*/)
116{
118}
119
121 std::size_t const mesh_item_id,
122 NumLib::LocalToGlobalIndexMap const& dof_table, GlobalVector const& x,
123 double const t, double const delta_t)
124{
125 auto const indices = NumLib::getIndices(mesh_item_id, dof_table);
126 auto const local_x = x.get(indices);
127
128 preTimestepConcrete(local_x, t, delta_t);
129}
130
132 std::size_t const mesh_item_id,
133 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
134 std::vector<GlobalVector*> const& x, double const t, double const dt)
135{
136 std::vector<double> local_x_vec;
137
138 auto const n_processes = x.size();
139 for (std::size_t process_id = 0; process_id < n_processes; ++process_id)
140 {
141 auto const indices =
142 NumLib::getIndices(mesh_item_id, *dof_tables[process_id]);
143 assert(!indices.empty());
144 auto const local_solution = x[process_id]->get(indices);
145 local_x_vec.insert(std::end(local_x_vec), std::begin(local_solution),
146 std::end(local_solution));
147 }
148 auto const local_x = MathLib::toVector(local_x_vec);
149
150 postTimestepConcrete(local_x, t, dt);
151}
152
154 std::size_t const mesh_item_id,
155 NumLib::LocalToGlobalIndexMap const& dof_table, GlobalVector const& x,
156 GlobalVector const& xdot, double const t, double const dt,
157 bool const use_monolithic_scheme, int const process_id)
158{
159 auto const indices = NumLib::getIndices(mesh_item_id, dof_table);
160 auto const local_x = x.get(indices);
161 auto const local_xdot = xdot.get(indices);
162
163 postNonLinearSolverConcrete(local_x, local_xdot, t, dt,
164 use_monolithic_scheme, process_id);
165}
166
167} // namespace ProcessLib
#define OGS_FATAL(...)
Definition: Error.h:26
Global vector based on Eigen vector.
Definition: EigenVector.h:28
double get(IndexType rowId) const
get entry
Definition: EigenVector.h:61
virtual void preTimestepConcrete(std::vector< double > const &, double const, double const)
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 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 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)
Eigen::Map< const Vector > toVector(std::vector< double > const &data, Eigen::VectorXd::Index size)
Creates an Eigen mapped vector from the given data vector.
static const double t
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)