Flow123d  master-f44eb46
assembly_base.hh
Go to the documentation of this file.
1 /*!
2  *
3  * Copyright (C) 2015 Technical University of Liberec. All rights reserved.
4  *
5  * This program is free software; you can redistribute it and/or modify it under
6  * the terms of the GNU General Public License version 3 as published by the
7  * Free Software Foundation. (http://www.gnu.org/licenses/gpl-3.0.en.html)
8  *
9  * This program is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12  *
13  *
14  * @file assembly_base.hh
15  * @brief
16  */
17 
18 #ifndef ASSEMBLY_BASE_HH_
19 #define ASSEMBLY_BASE_HH_
20 
21 
24 #include "fields/eval_points.hh"
26 #include "fem/update_flags.hh"
27 
28 
29 
30 /**
31  * Base class define empty methods, these methods can be overwite in descendants.
32  */
33 template <unsigned int dim>
35 {
36 public:
41 
42  /// Constructor
43  AssemblyBase(unsigned int quad_order) {
44  quad_ = new QGauss(dim, 2*quad_order);
45  quad_low_ = new QGauss(dim-1, 2*quad_order);
46  }
47 
48  /// Destructor
49  virtual ~AssemblyBase() {
50  delete quad_;
51  delete quad_low_;
52  }
53 
54  /// Assembles the volume integrals on cell.
55  virtual inline void cell_integral(FMT_UNUSED DHCellAccessor cell, FMT_UNUSED unsigned int element_patch_idx) {}
56 
57  /// Assembles the fluxes on the boundary.
58  virtual inline void boundary_side_integral(FMT_UNUSED DHCellSide cell_side) {}
59 
60  /// Assembles the fluxes between sides on the edge.
61  virtual inline void edge_integral(FMT_UNUSED RangeConvert<DHEdgeSide, DHCellSide> edge_side_range) {}
62 
63  /// Assembles the fluxes between elements of different dimensions.
64  virtual inline void dimjoin_intergral(FMT_UNUSED DHCellAccessor cell_lower_dim, FMT_UNUSED DHCellSide neighb_side) {}
65 
66  /// Method prepares object before assemblation (e.g. balance, ...).
67  virtual void begin() {}
68 
69  /// Method finishes object after assemblation (e.g. balance, ...).
70  virtual void end() {}
71 
72  /// Getter of active_integrals.
73  inline int n_active_integrals() const {
74  return active_integrals_;
75  }
76 
77  /// Create integrals according to dim of assembly object
78  void create_integrals(std::shared_ptr<EvalPoints> eval_points, AssemblyIntegrals &integrals) {
80  ASSERT_PERMANENT_PTR(quad_).error("Data member 'quad_' must be initialized if you use bulk integral!\n");
81  integrals_.bulk_ = eval_points->add_bulk<dim>(*quad_);
82  integrals.bulk_[dim-1] = integrals_.bulk_;
83  }
85  ASSERT_PERMANENT_PTR(quad_low_).error("Data member 'quad_low_' must be initialized if you use edge integral!\n");
86  integrals_.edge_ = eval_points->add_edge<dim>(*quad_low_);
87  integrals.edge_[dim-1] = integrals_.edge_;
88  }
89  if ((dim>1) && (active_integrals_ & ActiveIntegrals::coupling)) {
90  ASSERT_PERMANENT_PTR(quad_).error("Data member 'quad_' must be initialized if you use coupling integral!\n");
91  ASSERT_PERMANENT_PTR(quad_low_).error("Data member 'quad_low_' must be initialized if you use coupling integral!\n");
92  integrals_.coupling_ = eval_points->add_coupling<dim>(*quad_low_);
93  integrals.coupling_[dim-2] = integrals_.coupling_;
94  }
96  ASSERT_PERMANENT_PTR(quad_).error("Data member 'quad_' must be initialized if you use boundary integral!\n");
97  ASSERT_PERMANENT_PTR(quad_low_).error("Data member 'quad_low_' must be initialized if you use boundary integral!\n");
98  integrals_.boundary_ = eval_points->add_boundary<dim>(*quad_low_);
99  integrals.boundary_[dim-1] = integrals_.boundary_;
100  }
101  }
102 
103  /// Return BulkPoint range of appropriate dimension
104  inline Range< BulkPoint > bulk_points(unsigned int element_patch_idx) const {
105  return integrals_.bulk_->points(element_patch_idx, element_cache_map_);
106  }
107 
108  /// Return EdgePoint range of appropriate dimension
109  inline Range< EdgePoint > edge_points(const DHCellSide &cell_side) const {
110  ASSERT( cell_side.dim() > 0 ).error("Invalid cell dimension, must be 1, 2 or 3!\n");
111  return integrals_.edge_->points(cell_side, element_cache_map_);
112  }
113 
114  /// Return CouplingPoint range of appropriate dimension
115  inline Range< CouplingPoint > coupling_points(const DHCellSide &cell_side) const {
116  ASSERT( cell_side.dim() > 1 ).error("Invalid cell dimension, must be 2 or 3!\n");
117  return integrals_.coupling_->points(cell_side, element_cache_map_);
118  }
119 
120  /// Return BoundaryPoint range of appropriate dimension
121  inline Range< BoundaryPoint > boundary_points(const DHCellSide &cell_side) const {
122  ASSERT( cell_side.dim() > 0 ).error("Invalid cell dimension, must be 1, 2 or 3!\n");
123  return integrals_.boundary_->points(cell_side, element_cache_map_);
124  }
125 
126  /// Assembles the cell integrals for the given dimension.
127  virtual inline void assemble_cell_integrals(const RevertableList<BulkIntegralData> &bulk_integral_data) {
128  for (unsigned int i=0; i<bulk_integral_data.permanent_size(); ++i) {
129  if (bulk_integral_data[i].cell.dim() != dim) continue;
130  this->cell_integral(bulk_integral_data[i].cell, element_cache_map_->position_in_cache(bulk_integral_data[i].cell.elm_idx()));
131  }
132  // Possibly optimization but not so fast as we would assume (needs change interface of cell_integral)
133  /*for (unsigned int i=0; i<element_cache_map_->n_elements(); ++i) {
134  unsigned int elm_start = element_cache_map_->element_chunk_begin(i);
135  if (element_cache_map_->eval_point_data(elm_start).i_eval_point_ != 0) continue;
136  this->cell_integral(i, element_cache_map_->eval_point_data(elm_start).dh_loc_idx_);
137  }*/
138  }
139 
140  /// Assembles the boundary side integrals for the given dimension.
141  inline void assemble_boundary_side_integrals(const RevertableList<BoundaryIntegralData> &boundary_integral_data) {
142  for (unsigned int i=0; i<boundary_integral_data.permanent_size(); ++i) {
143  if (boundary_integral_data[i].side.dim() != dim) continue;
144  this->boundary_side_integral(boundary_integral_data[i].side);
145  }
146  }
147 
148  /// Assembles the edge integrals for the given dimension.
149  inline void assemble_edge_integrals(const RevertableList<EdgeIntegralData> &edge_integral_data) {
150  for (unsigned int i=0; i<edge_integral_data.permanent_size(); ++i) {
151  auto range = edge_integral_data[i].edge_side_range;
152  if (range.begin()->dim() != dim) continue;
153  this->edge_integral(edge_integral_data[i].edge_side_range);
154  }
155  }
156 
157  /// Assembles the neighbours integrals for the given dimension.
158  inline void assemble_neighbour_integrals(const RevertableList<CouplingIntegralData> &coupling_integral_data) {
159  for (unsigned int i=0; i<coupling_integral_data.permanent_size(); ++i) {
160  if (coupling_integral_data[i].side.dim() != dim) continue;
161  this->dimjoin_intergral(coupling_integral_data[i].cell, coupling_integral_data[i].side);
162  }
163  }
164 
165 protected:
166  /// Set of integral of given dimension necessary in assemblation
167  struct DimIntegrals {
168  std::shared_ptr<BulkIntegral> bulk_; ///< Bulk integrals of elements
169  std::shared_ptr<EdgeIntegral> edge_; ///< Edge integrals between elements of same dimensions
170  std::shared_ptr<CouplingIntegral> coupling_; ///< Coupling integrals between elements of dimensions dim and dim-1
171  std::shared_ptr<BoundaryIntegral> boundary_; ///< Boundary integrals betwwen side and boundary element of dim-1
172  };
173 
174  /**
175  * Default constructor.
176  *
177  * Be aware if you use this constructor. Quadrature objects must be initialized manually in descendant.
178  */
180  : quad_(nullptr), quad_low_(nullptr) {}
181 
182  /// Print update flags to string format.
183  std::string print_update_flags(UpdateFlags u) const {
184  std::stringstream s;
185  s << u;
186  return s.str();
187  }
188 
189  Quadrature *quad_; ///< Quadrature used in assembling methods.
190  Quadrature *quad_low_; ///< Quadrature used in assembling methods (dim-1).
191  int active_integrals_; ///< Holds mask of active integrals.
192  DimIntegrals integrals_; ///< Set of used integrals.
193  ElementCacheMap *element_cache_map_; ///< ElementCacheMap shared with GenericAssembly object.
194 };
195 
196 
197 #endif /* ASSEMBLY_BASE_HH_ */
#define ASSERT(expr)
Definition: asserts.hh:351
#define ASSERT_PERMANENT_PTR(ptr)
Definition of assert macro checking non-null pointer (PTR)
Definition: asserts.hh:337
virtual void assemble_cell_integrals(const RevertableList< BulkIntegralData > &bulk_integral_data)
Assembles the cell integrals for the given dimension.
Range< BulkPoint > bulk_points(unsigned int element_patch_idx) const
Return BulkPoint range of appropriate dimension.
virtual void boundary_side_integral(FMT_UNUSED DHCellSide cell_side)
Assembles the fluxes on the boundary.
virtual void end()
Method finishes object after assemblation (e.g. balance, ...).
virtual ~AssemblyBase()
Destructor.
virtual void cell_integral(FMT_UNUSED DHCellAccessor cell, FMT_UNUSED unsigned int element_patch_idx)
Assembles the volume integrals on cell.
Range< EdgePoint > edge_points(const DHCellSide &cell_side) const
Return EdgePoint range of appropriate dimension.
GenericAssemblyBase::CouplingIntegralData CouplingIntegralData
Quadrature * quad_
Quadrature used in assembling methods.
void create_integrals(std::shared_ptr< EvalPoints > eval_points, AssemblyIntegrals &integrals)
Create integrals according to dim of assembly object.
void assemble_edge_integrals(const RevertableList< EdgeIntegralData > &edge_integral_data)
Assembles the edge integrals for the given dimension.
DimIntegrals integrals_
Set of used integrals.
virtual void begin()
Method prepares object before assemblation (e.g. balance, ...).
std::string print_update_flags(UpdateFlags u) const
Print update flags to string format.
Quadrature * quad_low_
Quadrature used in assembling methods (dim-1).
GenericAssemblyBase::BulkIntegralData BulkIntegralData
void assemble_neighbour_integrals(const RevertableList< CouplingIntegralData > &coupling_integral_data)
Assembles the neighbours integrals for the given dimension.
Range< CouplingPoint > coupling_points(const DHCellSide &cell_side) const
Return CouplingPoint range of appropriate dimension.
virtual void dimjoin_intergral(FMT_UNUSED DHCellAccessor cell_lower_dim, FMT_UNUSED DHCellSide neighb_side)
Assembles the fluxes between elements of different dimensions.
virtual void edge_integral(FMT_UNUSED RangeConvert< DHEdgeSide, DHCellSide > edge_side_range)
Assembles the fluxes between sides on the edge.
int active_integrals_
Holds mask of active integrals.
Range< BoundaryPoint > boundary_points(const DHCellSide &cell_side) const
Return BoundaryPoint range of appropriate dimension.
AssemblyBase(unsigned int quad_order)
Constructor.
void assemble_boundary_side_integrals(const RevertableList< BoundaryIntegralData > &boundary_integral_data)
Assembles the boundary side integrals for the given dimension.
GenericAssemblyBase::EdgeIntegralData EdgeIntegralData
int n_active_integrals() const
Getter of active_integrals.
ElementCacheMap * element_cache_map_
ElementCacheMap shared with GenericAssembly object.
GenericAssemblyBase::BoundaryIntegralData BoundaryIntegralData
Cell accessor allow iterate over DOF handler cells.
Side accessor allows to iterate over sides of DOF handler cell.
unsigned int dim() const
Return dimension of element appropriate to the side.
Directing class of FieldValueCache.
unsigned int position_in_cache(unsigned mesh_elm_idx, bool bdr=false) const
Return position of element stored in ElementCacheMap.
Symmetric Gauss-Legendre quadrature formulae on simplices.
Base class for quadrature rules on simplices in arbitrary dimensions.
Definition: quadrature.hh:48
Range helper class.
@ coupling
@ boundary
@ bulk
@ edge
#define FMT_UNUSED
Definition: posix.h:75
Definitions of particular quadrature rules on simplices.
Set of integral of given dimension necessary in assemblation.
std::shared_ptr< EdgeIntegral > edge_
Edge integrals between elements of same dimensions.
std::shared_ptr< CouplingIntegral > coupling_
Coupling integrals between elements of dimensions dim and dim-1.
std::shared_ptr< BulkIntegral > bulk_
Bulk integrals of elements.
std::shared_ptr< BoundaryIntegral > boundary_
Boundary integrals betwwen side and boundary element of dim-1.
Set of all used integral necessary in assemblation.
std::array< std::shared_ptr< BulkIntegral >, 3 > bulk_
Bulk integrals of elements of dimensions 1, 2, 3.
std::array< std::shared_ptr< CouplingIntegral >, 2 > coupling_
Coupling integrals between elements of dimensions 1-2, 2-3.
std::array< std::shared_ptr< BoundaryIntegral >, 3 > boundary_
Boundary integrals betwwen elements of dimensions 1, 2, 3 and boundaries.
std::array< std::shared_ptr< EdgeIntegral >, 3 > edge_
Edge integrals between elements of dimensions 1, 2, 3.
std::size_t permanent_size() const
Return permanent size of list.
Enum type UpdateFlags indicates which quantities are to be recomputed on each finite element cell.
UpdateFlags
Enum type UpdateFlags indicates which quantities are to be recomputed on each finite element cell.
Definition: update_flags.hh:68