16 #ifndef __H2D_DISCRETE_PROBLEM_H
17 #define __H2D_DISCRETE_PROBLEM_H
19 #include "hermes_common.h"
20 #include "weakform/weakform.h"
21 #include "function/function.h"
22 #include "exceptions.h"
24 #include "discrete_problem_helpers.h"
25 #include "discrete_problem_thread_assembler.h"
31 class PrecalcShapeset;
36 template<
typename Scalar>
38 public Hermes::Mixins::Loggable,
39 public Hermes::Mixins::TimeMeasurable,
41 public Hermes::Mixins::StateQueryable,
44 public Hermes::Mixins::IntegrableWithGlobalOrder,
69 DiscreteProblem(
bool linear =
false,
bool dirichlet_lift_accordingly =
true);
74 bool assemble(Scalar*& coeff_vec, SparseMatrix<Scalar>* mat, Vector<Scalar>* rhs =
nullptr);
77 bool assemble(Scalar*& coeff_vec, Vector<Scalar>* rhs =
nullptr);
80 bool assemble(SparseMatrix<Scalar>* mat, Vector<Scalar>* rhs =
nullptr);
84 bool assemble(Vector<Scalar>* rhs);
87 void set_time(
double time);
88 void set_time_step(
double time_step);
98 std::vector<SpaceSharedPtr<Scalar> > get_spaces();
101 typedef void(*reassembled_states_reuse_linear_system_fn)(
Traverse::State**& states,
unsigned int& num_states, SparseMatrix<Scalar>* mat, Vector<Scalar>* rhs, Vector<Scalar>* dirichlet_lift_rhs, Scalar*& coeff_vec);
102 void set_reassembled_states_reuse_linear_system_fn(reassembled_states_reuse_linear_system_fn fn) {
103 this->reassembled_states_reuse_linear_system = fn;
105 reassembled_states_reuse_linear_system_fn reassembled_states_reuse_linear_system;
106 void set_reusable_DOFs(
bool **reusable_DOFs,
bool **reusable_Dirichlet)
108 for (
int i = 0; i < this->num_threads_used; i++)
110 this->threadAssembler[i]->reusable_DOFs = reusable_DOFs;
111 this->threadAssembler[i]->reusable_Dirichlet = reusable_Dirichlet;
116 virtual void set_verbose_output(
bool to_set);
120 void init_assembling(Traverse::State**& states,
unsigned int& num_states, std::vector<MeshSharedPtr>& meshes);
121 void deinit_assembling(Traverse::State** states,
unsigned int num_states);
124 void set_RK(
int original_spaces_count,
bool force_diagonal_blocks =
nullptr, Table* block_weights =
nullptr);
128 inline std::string getClassName()
const {
return "DiscreteProblem"; }
131 void init(
bool linear,
bool dirichlet_lift_accordingly);
134 std::vector<SpaceSharedPtr<Scalar> >
spaces;
144 bool set_matrix(SparseMatrix<Scalar>* mat);
145 bool set_rhs(Vector<Scalar>* rhs);
146 void invalidate_matrix();
154 template<
typename T>
friend class Solver;
156 template<
typename T,
typename S>
friend class AdaptSolver;
157 template<
typename T>
friend class NonlinearSolver;
Class utilizes parallel calculation.
std::vector< SpaceSharedPtr< Scalar > > spaces
Space instances for all equations in the system.
Provides capabilities to (re-)assemble a matrix / vector only where necessary. See also Solver::keep_...
Used to pass the instances of Space around.
::xsd::cxx::tree::time< char, simple_type > time
C++ type corresponding to the time XML Schema built-in type.
DiscreteProblemSelectiveAssembler< Scalar > selectiveAssembler
Select the right things to assemble.
DiscreteProblemThreadAssembler< Scalar > ** threadAssembler
Assembly data.
::xsd::cxx::tree::string< char, simple_type > string
C++ type corresponding to the string XML Schema built-in type.
This class is a one-thread (non-DG) assembly worker.
Vector< Scalar > * dirichlet_lift_rhs
Dirichlet lift rhs part.