Hermes2D
2.0
|
#include <discrete_problem.h>
Classes | |
class | CacheRecordPerElement |
Caching. More... | |
class | CacheRecordPerSubIdx |
Public Member Functions | |
DiscreteProblem (const WeakForm< Scalar > *wf, Hermes::vector< const Space< Scalar > * > spaces) | |
Constructor for multiple components / equations. | |
DiscreteProblem (const WeakForm< Scalar > *wf, const Space< Scalar > *space) | |
Constructor for one equation. | |
virtual bool | isOkay () const |
State querying helpers. | |
virtual std::string | getClassName () const |
Get class name, for the purpose of messaging. | |
void | set_fvm () |
Set this problem to Finite Volume (no integration order calculation). | |
virtual void | set_spaces (Hermes::vector< const Space< Scalar > * > spaces) |
Sets new spaces for the instance. More... | |
virtual void | set_space (const Space< Scalar > *space) |
DiscreteProblem () | |
Non-parameterized constructor. | |
virtual | ~DiscreteProblem () |
Destuctor. | |
void | set_do_not_use_cache () |
If the cache should not be used for any reason. | |
const WeakForm< Scalar > * | get_weak_formulation () const |
Get the weak forms. | |
void | set_weak_formulation (const WeakForm< Scalar > *wf) |
Set the weak forms. | |
virtual Hermes::vector< const Space< Scalar > * > | get_spaces () const |
Get all spaces as a Hermes::vector. | |
int | get_num_dofs () const |
Get the number of unknowns. | |
bool | is_matrix_free () const |
Get info about presence of a matrix. | |
virtual void | set_time (double time) |
set time information for time-dependent problems. | |
virtual void | set_time_step (double time_step) |
void | delete_cache () |
void | assemble (Scalar *coeff_vec, SparseMatrix< Scalar > *mat, Vector< Scalar > *rhs=NULL, bool force_diagonal_blocks=false, Table *block_weights=NULL) |
void | assemble (Scalar *coeff_vec, Vector< Scalar > *rhs=NULL, bool force_diagonal_blocks=false, Table *block_weights=NULL) |
virtual void | assemble (SparseMatrix< Scalar > *mat, Vector< Scalar > *rhs=NULL, bool force_diagonal_blocks=false, Table *block_weights=NULL) |
void | assemble (Vector< Scalar > *rhs=NULL, bool force_diagonal_blocks=false, Table *block_weights=NULL) |
![]() | |
virtual const Space< Scalar > * | get_space (int n) const |
![]() | |
void | check () const |
Method to handle the state. | |
Static Public Member Functions | |
static int | init_geometry_points (RefMap *reference_mapping, int order, Geom< double > *&geometry, double *&jacobian_x_weights) |
static int | init_surface_geometry_points (RefMap *reference_mapping, int &order, Traverse::State *current_state, Geom< double > *&geometry, double *&jacobian_x_weights) |
Protected Member Functions | |
void | init_assembling (Scalar *coeff_vec, PrecalcShapeset ***pss, PrecalcShapeset ***spss, RefMap ***refmaps, Solution< Scalar > ***u_ext, AsmList< Scalar > ***als, WeakForm< Scalar > **weakforms) |
void | deinit_assembling (PrecalcShapeset ***pss, PrecalcShapeset ***spss, RefMap ***refmaps, Solution< Scalar > ***u_ext, AsmList< Scalar > ***als, WeakForm< Scalar > **weakforms) |
bool | form_to_be_assembled (MatrixForm< Scalar > *form, Traverse::State *current_state) |
The form will be assembled. | |
bool | form_to_be_assembled (MatrixFormVol< Scalar > *form, Traverse::State *current_state) |
bool | form_to_be_assembled (MatrixFormSurf< Scalar > *form, Traverse::State *current_state) |
bool | form_to_be_assembled (VectorForm< Scalar > *form, Traverse::State *current_state) |
bool | form_to_be_assembled (VectorFormVol< Scalar > *form, Traverse::State *current_state) |
bool | form_to_be_assembled (VectorFormSurf< Scalar > *form, Traverse::State *current_state) |
double | block_scaling_coeff (MatrixForm< Scalar > *form) const |
void | create_sparse_structure () |
void | create_sparse_structure (SparseMatrix< Scalar > *mat, Vector< Scalar > *rhs=NULL) |
void | set_RK (int original_spaces_count) |
Set the special handling of external functions of Runge-Kutta methods, including information how many spaces were there in the original problem. | |
bool | state_needs_recalculation (AsmList< Scalar > **current_als, Traverse::State *current_state) |
void | calculate_cache_records (PrecalcShapeset **current_pss, PrecalcShapeset **current_spss, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, AsmList< Scalar > **current_als, Traverse::State *current_state, AsmList< Scalar > **current_alsSurface, WeakForm< Scalar > *current_wf) |
Calculate cache records for this set of parameters. | |
void | assemble_one_state (PrecalcShapeset **current_pss, PrecalcShapeset **current_spss, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, AsmList< Scalar > **current_als, Traverse::State *current_state, WeakForm< Scalar > *current_wf) |
Assemble one state. | |
void | adjust_order_to_refmaps (Form< Scalar > *form, int &order, Hermes::Ord *o, RefMap **current_refmaps) |
Adjusts order to refmaps. | |
int | calc_order_matrix_form (MatrixForm< Scalar > *mfv, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, Traverse::State *current_state) |
Matrix volumetric forms - calculate the integration order. | |
virtual void | assemble_matrix_form (MatrixForm< Scalar > *form, int order, Func< double > **base_fns, Func< double > **test_fns, Func< Scalar > **ext, Func< Scalar > **u_ext, AsmList< Scalar > *current_als_i, AsmList< Scalar > *current_als_j, Traverse::State *current_state, int n_quadrature_points, Geom< double > *geometry, double *jacobian_x_weights) |
Matrix volumetric forms - assemble the form. | |
int | calc_order_vector_form (VectorForm< Scalar > *mfv, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, Traverse::State *current_state) |
Vector volumetric forms - calculate the integration order. | |
void | assemble_vector_form (VectorForm< Scalar > *form, int order, Func< double > **test_fns, Func< Scalar > **ext, Func< Scalar > **u_ext, AsmList< Scalar > *current_als, Traverse::State *current_state, int n_quadrature_points, Geom< double > *geometry, double *jacobian_x_weights) |
Vector volumetric forms - assemble the form. | |
void | init_ext_orders (Form< Scalar > *form, Func< Hermes::Ord > **oi, Func< Hermes::Ord > **oext, Solution< Scalar > **current_u_ext, Traverse::State *current_state) |
void | deinit_ext_orders (Form< Scalar > *form, Func< Hermes::Ord > **oi, Func< Hermes::Ord > **oext) |
void | init () |
Init function. Common code for the constructors. | |
bool | is_up_to_date () const |
Matrix structure as well as spaces and weak formulation is up-to-date. | |
void | assemble_one_DG_state (PrecalcShapeset **current_pss, PrecalcShapeset **current_spss, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, AsmList< Scalar > **current_als, Traverse::State *current_state, Hermes::vector< MatrixFormDG< Scalar > * > current_mfDG, Hermes::vector< VectorFormDG< Scalar > * > current_vfDG, Transformable **fn) |
Assemble DG forms. More... | |
void | assemble_DG_one_neighbor (bool edge_processed, unsigned int neighbor_i, PrecalcShapeset **current_pss, PrecalcShapeset **current_spss, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, AsmList< Scalar > **current_als, Traverse::State *current_state, Hermes::vector< MatrixFormDG< Scalar > * > current_mfDG, Hermes::vector< VectorFormDG< Scalar > * > current_vfDG, Transformable **fn, std::map< unsigned int, PrecalcShapeset * > npss, std::map< unsigned int, PrecalcShapeset * > nspss, std::map< unsigned int, RefMap * > nrefmap, LightArray< NeighborSearch< Scalar > * > &neighbor_searches, unsigned int min_dg_mesh_seq) |
Assemble one DG neighbor. | |
void | assemble_DG_matrix_forms (PrecalcShapeset **current_pss, PrecalcShapeset **current_spss, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, AsmList< Scalar > **current_als, Traverse::State *current_state, MatrixFormDG< Scalar > **current_mfDG, std::map< unsigned int, PrecalcShapeset * > npss, std::map< unsigned int, PrecalcShapeset * > nspss, std::map< unsigned int, RefMap * > nrefmap, LightArray< NeighborSearch< Scalar > * > &neighbor_searches) |
Assemble DG matrix forms. | |
void | assemble_DG_vector_forms (PrecalcShapeset **current_spss, RefMap **current_refmaps, Solution< Scalar > **current_u_ext, AsmList< Scalar > **current_als, Traverse::State *current_state, VectorFormDG< Scalar > **current_vfDG, std::map< unsigned int, PrecalcShapeset * > nspss, std::map< unsigned int, RefMap * > nrefmap, LightArray< NeighborSearch< Scalar > * > &neighbor_searches) |
Assemble DG vector forms. | |
DiscontinuousFunc< Hermes::Ord > * | init_ext_fn_ord (NeighborSearch< Scalar > *ns, MeshFunction< Scalar > *fu) |
int | calc_order_dg_matrix_form (MatrixFormDG< Scalar > *mfDG, Hermes::vector< Solution< Scalar > * > u_ext, PrecalcShapeset *fu, PrecalcShapeset *fv, RefMap *ru, SurfPos *surf_pos, bool neighbor_supp_u, bool neighbor_supp_v, LightArray< NeighborSearch< Scalar > * > &neighbor_searches, int neighbor_index_u) |
Calculates integration order for DG matrix forms. | |
Scalar | eval_dg_form (MatrixFormDG< Scalar > *mfDG, Hermes::vector< Solution< Scalar > * > u_ext, PrecalcShapeset *fu, PrecalcShapeset *fv, RefMap *ru_central, RefMap *ru_actual, RefMap *rv, bool neighbor_supp_u, bool neighbor_supp_v, SurfPos *surf_pos, LightArray< NeighborSearch< Scalar > * > &neighbor_searches, int neighbor_index_u, int neighbor_index_v) |
Evaluates DG matrix forms on an edge between elements identified by ru_actual, rv. | |
int | calc_order_dg_vector_form (VectorFormDG< Scalar > *vfDG, Hermes::vector< Solution< Scalar > * > u_ext, PrecalcShapeset *fv, RefMap *ru, SurfPos *surf_pos, LightArray< NeighborSearch< Scalar > * > &neighbor_searches, int neighbor_index_v) |
Calculates integration order for DG vector forms. | |
Scalar | eval_dg_form (VectorFormDG< Scalar > *vfDG, Hermes::vector< Solution< Scalar > * > u_ext, PrecalcShapeset *fv, RefMap *rv, SurfPos *surf_pos, LightArray< NeighborSearch< Scalar > * > &neighbor_searches, int neighbor_index_v) |
Evaluates DG vector forms on an edge between elements identified by ru_actual, rv. | |
Func< Hermes::Ord > ** | init_ext_fns_ord (Hermes::vector< MeshFunction< Scalar > * > &ext, LightArray< NeighborSearch< Scalar > * > &neighbor_searches) |
Initialize orders of external functions for DG forms. | |
Func< Scalar > ** | init_ext_fns (Hermes::vector< MeshFunction< Scalar > * > &ext, LightArray< NeighborSearch< Scalar > * > &neighbor_searches, int order, unsigned int min_dg_mesh_seq) |
Initialize external functions for DG forms. | |
bool | init_neighbors (LightArray< NeighborSearch< Scalar > * > &neighbor_searches, Traverse::State *current_state, unsigned int min_dg_mesh_seq) |
Initialize neighbors. | |
void | build_multimesh_tree (NeighborNode *root, LightArray< NeighborSearch< Scalar > * > &neighbor_searches) |
Initialize the tree for traversing multimesh neighbors. | |
void | insert_into_multimesh_tree (NeighborNode *node, unsigned int *transformations, unsigned int transformation_count) |
Recursive insertion function into the tree. | |
Hermes::vector< Hermes::vector < unsigned int > * > | get_multimesh_neighbors_transformations (NeighborNode *multimesh_tree) |
Return a global (unified list of central element transformations representing the neighbors on the union mesh. | |
void | traverse_multimesh_tree (NeighborNode *node, Hermes::vector< Hermes::vector< unsigned int > * > &running_transformations) |
Traverse the multimesh tree. Used in the function get_multimesh_neighbors_transformations(). | |
void | update_neighbor_search (NeighborSearch< Scalar > *ns, NeighborNode *multimesh_tree) |
Update the NeighborSearch according to the multimesh tree. | |
NeighborNode * | find_node (unsigned int *transformations, unsigned int transformation_count, NeighborNode *node) |
int | update_ns_subtree (NeighborSearch< Scalar > *ns, NeighborNode *node, unsigned int ith_neighbor) |
void | traverse_multimesh_subtree (NeighborNode *node, Hermes::vector< Hermes::vector< unsigned int > * > &running_central_transformations, Hermes::vector< Hermes::vector< unsigned int > * > &running_neighbor_transformations, const typename NeighborSearch< Scalar >::NeighborEdgeInfo &edge_info, const int &active_edge, const int &mode) |
Traverse the multimesh subtree. Used in the function update_ns_subtree(). | |
Protected Attributes | |
const WeakForm< Scalar > * | wf |
Weak formulation. | |
Hermes::vector< const Space < Scalar > * > | spaces |
Space instances for all equations in the system. | |
Hermes::vector< unsigned int > | spaces_first_dofs |
int * | sp_seq |
Seq numbers of Space instances in spaces. | |
int | ndof |
Number of DOFs of all Space instances in spaces. | |
Geom< Hermes::Ord > | geom_ord |
Instance of the class Geom used in the calculation of integration order. | |
bool | is_fvm |
bool | is_linear |
Internal. | |
bool | have_matrix |
bool | DG_matrix_forms_present |
There is a matrix form set on DG_INNER_EDGE area or not. | |
bool | DG_vector_forms_present |
There is a vector form set on DG_INNER_EDGE area or not. | |
bool | RungeKutta |
Turn on Runge-Kutta specific handling of external functions. | |
int | RK_original_spaces_count |
Number of spaces in the original problem in a Runge-Kutta method. | |
SparseMatrix< Scalar > * | current_mat |
Storing assembling info. | |
Vector< Scalar > * | current_rhs |
bool | current_force_diagonal_blocks |
Table * | current_block_weights |
std::map< uint64_t, CacheRecordPerSubIdx * > *** | cache_records_sub_idx |
CacheRecordPerElement *** | cache_records_element |
bool ** | cache_element_stored |
int | cache_size |
bool | do_not_use_cache |
Hermes::Exceptions::Exception * | caughtException |
Exception caught in a parallel region. | |
Static Protected Attributes | |
static double | fake_wt = 1.0 |
Fake weight used in the calculation of integration order. | |
Discrete problem class.
This class does assembling into external matrix / vector structures.
Definition at line 63 of file discrete_problem.h.
void Hermes::Hermes2D::DiscreteProblem< Scalar >::assemble | ( | Scalar * | coeff_vec, |
SparseMatrix< Scalar > * | mat, | ||
Vector< Scalar > * | rhs = NULL , |
||
bool | force_diagonal_blocks = false , |
||
Table * | block_weights = NULL |
||
) |
Assembling. General assembling procedure for nonlinear problems. coeff_vec is the previous Newton vector. If force_diagonal_block == true, then (zero) matrix antries are created in diagonal blocks even if corresponding matrix weak forms do not exist. This is useful if the matrix is later to be merged with a matrix that has nonzeros in these blocks. The Table serves for optional weighting of matrix blocks in systems. The parameter add_dir_lift decides whether Dirichlet lift will be added while coeff_vec is converted into Solutions.
Definition at line 1031 of file discrete_problem.cpp.
void Hermes::Hermes2D::DiscreteProblem< Scalar >::assemble | ( | Scalar * | coeff_vec, |
Vector< Scalar > * | rhs = NULL , |
||
bool | force_diagonal_blocks = false , |
||
Table * | block_weights = NULL |
||
) |
Assembling. Without the matrix.
Definition at line 1204 of file discrete_problem.cpp.
|
virtual |
Light version passing NULL for the coefficient vector. External solutions are initialized with zeros.
Reimplemented in Hermes::Hermes2D::DiscreteProblemLinear< Scalar >.
Definition at line 866 of file discrete_problem.cpp.
void Hermes::Hermes2D::DiscreteProblem< Scalar >::assemble | ( | Vector< Scalar > * | rhs = NULL , |
bool | force_diagonal_blocks = false , |
||
Table * | block_weights = NULL |
||
) |
Light version passing NULL for the coefficient vector. External solutions are initialized with zeros. Without the matrix.
Definition at line 874 of file discrete_problem.cpp.
|
protected |
|
protected |
Preassembling. Precalculate matrix sparse structure. If force_diagonal_block == true, then (zero) matrix antries are created in diagonal blocks even if corresponding matrix weak forms do not exist. This is useful if the matrix is later to be merged with a matrix that has nonzeros in these blocks. The Table serves for optional weighting of matrix blocks in systems.
Definition at line 671 of file discrete_problem.cpp.
|
protected |
{*, *} Cleans up after init_ext_orders.
Definition at line 2156 of file discrete_problem.cpp.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::calc_order_matrix_form(), and Hermes::Hermes2D::DiscreteProblem< Scalar >::calc_order_vector_form().
|
protected |
Finds a node in the multimesh tree that corresponds to the array transformations, with the length of transformation_count, starting to look for it in the NeighborNode node.
Definition at line 2903 of file discrete_problem.cpp.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::update_neighbor_search().
|
protected |
{*, *} Calculates orders for external functions.
Definition at line 2118 of file discrete_problem.cpp.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::calc_order_matrix_form(), and Hermes::Hermes2D::DiscreteProblem< Scalar >::calc_order_vector_form().
|
static |
{*, *} Init geometry, jacobian * weights, return the number of integration points.
Definition at line 2079 of file discrete_problem.cpp.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::calculate_cache_records().
|
virtual |
Sets new spaces for the instance.
Implements Hermes::Hermes2D::Mixins::SettableSpaces< Scalar >.
Definition at line 333 of file discrete_problem.cpp.
|
protected |
Assemble one state - needs recalculation?
Definition at line 1266 of file discrete_problem.cpp.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::assemble_one_state().
|
protected |
Updates the NeighborSearch ns according to the subtree of NeighborNode node. Returns 0 if no neighbor was deleted, -1 otherwise.
Definition at line 2930 of file discrete_problem.cpp.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::update_neighbor_search().
|
protected |
Matrix structure can be reused. If other conditions apply.
Definition at line 246 of file discrete_problem.h.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::DiscreteProblem().
|
protected |
If the problem has only constant test functions, there is no need for order calculation, which saves time.
Definition at line 239 of file discrete_problem.h.
Referenced by Hermes::Hermes2D::DiscreteProblem< Scalar >::calc_order_matrix_form(), Hermes::Hermes2D::DiscreteProblem< Scalar >::calc_order_vector_form(), and Hermes::Hermes2D::DiscreteProblem< Scalar >::DiscreteProblem().