Hermes2D  3.0
discrete_problem_dg_assembler.h
1 
16 #ifndef __H2D_DISCRETE_PROBLEM_DG_ASSEMBLER_H
17 #define __H2D_DISCRETE_PROBLEM_DG_ASSEMBLER_H
18 
19 #include "hermes_common.h"
20 #include "forms.h"
21 #include "weakform/weakform.h"
22 #include "function/function.h"
23 #include "neighbor_search.h"
24 #include "refinement_selectors/selector.h"
25 #include "exceptions.h"
26 #include "mixins2d.h"
27 #include "multimesh_dg_neighbor_tree.h"
28 #include "discrete_problem/discrete_problem_selective_assembler.h"
29 
30 namespace Hermes
31 {
32  namespace Hermes2D
33  {
38  template<typename Scalar>
39  class HERMES_API DiscreteProblemDGAssembler
40  {
41  public:
43  DiscreteProblemDGAssembler(DiscreteProblemThreadAssembler<Scalar>* threadAssembler, const std::vector<SpaceSharedPtr<Scalar> > spaces, std::vector<MeshSharedPtr>& meshes);
44 
47 
49  void init_assembling_one_state(Traverse::State* current_state_);
51  void assemble_one_state();
53  void deinit_assembling_one_state();
54 
55  static unsigned int dg_order;
56  private:
58  bool DG_matrix_forms_present;
59 
61  bool DG_vector_forms_present;
62 
64  void init_assembling_one_neighbor();
66  void assemble_one_neighbor(bool edge_processed, unsigned int neighbor_i, NeighborSearch<Scalar>** neighbor_searches);
68  void deinit_assembling_one_neighbor();
69 
71  DiscontinuousFunc<Scalar>** init_ext_fns(std::vector<MeshFunctionSharedPtr<Scalar> > ext,
72  NeighborSearch<Scalar>** neighbor_searches, int order);
73 
75  bool init_neighbors(NeighborSearch<Scalar>** neighbor_searches, Traverse::State* current_state);
77  void deinit_neighbors(NeighborSearch<Scalar>** neighbor_searches, Traverse::State* current_state);
78 
79  NeighborSearch<Scalar>*** neighbor_searches;
80  unsigned int* num_neighbors;
81  bool** processed;
82 
83  // Neighbor psss, refmaps.
85  RefMap ** nrefmaps;
86 
88  RefMap** refmaps;
89  Solution<Scalar>** u_ext;
90  AsmList<Scalar>* als;
91  std::vector<Transformable *> fns;
93  int spaces_size;
94  bool nonlinear;
96 
97  SparseMatrix<Scalar>* current_mat;
98  Vector<Scalar>* current_rhs;
99 
100  Traverse::State* current_state;
102  Scalar local_stiffness_matrix[H2D_MAX_LOCAL_BASIS_SIZE * H2D_MAX_LOCAL_BASIS_SIZE * 4];
103 
104  const std::vector<SpaceSharedPtr<Scalar> > spaces;
105  const std::vector<MeshSharedPtr>& meshes;
106 
107  template<typename T> friend class DiscreteProblem;
108  template<typename T> friend class DiscreteProblemIntegrationOrderCalculator;
109 
111  static NeighborSearch<Scalar>* get_neighbor_search_ext(WeakFormSharedPtr<Scalar> wf, NeighborSearch<Scalar>** neighbor_searches, int index);
112 
113 #ifdef DEBUG_DG_ASSEMBLING
114  void debug();
115 #endif
116  };
117  }
118 }
119 #endif
PrecalcShapeset variant for fast assembling.
Definition: precalc.h:109
Definition: adapt.h:24
This class represents a function with jump discontinuity on an interface of two elements.
Definition: forms.h:335
Provides capabilities to (re-)assemble a matrix / vector only where necessary. See also Solver::keep_...
Used to pass the instances of Space around.
Definition: space.h:34
Used to pass the instances of WeakForm around.
Definition: weakform.h:55
This class is a one-thread (non-DG) assembly worker.
Represents the reference mapping.
Definition: refmap.h:40
This class characterizes a neighborhood of a given edge in terms of adjacent elements and provides me...
Represents the solution of a PDE.
Definition: api2d.h:35