Flow123d  DF_patch_fe_data_tables-419e950
generic_assembly.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 generic_assembly.hh
15  * @brief
16  */
17 
18 #ifndef GENERIC_ASSEMBLY_HH_
19 #define GENERIC_ASSEMBLY_HH_
20 
22 #include "fields/eval_subset.hh"
23 #include "fields/eval_points.hh"
25 #include "fem/fe_values.hh"
26 #include "fem/patch_fe_values.hh"
27 #include "tools/revertable_list.hh"
28 #include "system/sys_profiler.hh"
29 
30 
31 
32 /// Allow set mask of active integrals.
34  no_intg = 0,
35  bulk = 0x0001,
36  edge = 0x0002,
37  coupling = 0x0004,
38  boundary = 0x0008
39 };
40 
41 
42 /// Set of all used integral necessary in assemblation
44  std::array<std::shared_ptr<BulkIntegral>, 3> bulk_; ///< Bulk integrals of elements of dimensions 1, 2, 3
45  std::array<std::shared_ptr<EdgeIntegral>, 3> edge_; ///< Edge integrals between elements of dimensions 1, 2, 3
46  std::array<std::shared_ptr<CouplingIntegral>, 2> coupling_; ///< Coupling integrals between elements of dimensions 1-2, 2-3
47  std::array<std::shared_ptr<BoundaryIntegral>, 3> boundary_; ///< Boundary integrals betwwen elements of dimensions 1, 2, 3 and boundaries
48 };
49 
50 
51 /**
52  * Common interface class for all Assembly classes.
53  */
55 {
56 public:
57  /**
58  * Helper structzre holds data of cell (bulk) integral
59  *
60  * Data is specified by cell and subset index in EvalPoint object
61  */
63  /// Default constructor
65 
66  /// Constructor with data mebers initialization
67  BulkIntegralData(DHCellAccessor dhcell, unsigned int subset_idx)
68  : cell(dhcell), subset_index(subset_idx) {}
69 
70  /// Copy constructor
72  : cell(other.cell), subset_index(other.subset_index) {}
73 
74  DHCellAccessor cell; ///< Specified cell (element)
75  unsigned int subset_index; ///< Index (order) of subset in EvalPoints object
76  };
77 
78  /**
79  * Helper structzre holds data of edge integral
80  *
81  * Data is specified by side and subset index in EvalPoint object
82  */
84  /// Default constructor
87 
88  /// Copy constructor
91 
92  /// Constructor with data mebers initialization
94  : edge_side_range(range), subset_index(subset_idx) {}
95 
96  RangeConvert<DHEdgeSide, DHCellSide> edge_side_range; ///< Specified cell side (element)
97  unsigned int subset_index; ///< Index (order) of subset in EvalPoints object
98  };
99 
100  /**
101  * Helper structzre holds data of neighbour (coupling) integral
102  *
103  * Data is specified by cell, side and their subset indices in EvalPoint object
104  */
106  /// Default constructor
108 
109  /// Constructor with data mebers initialization
110  CouplingIntegralData(DHCellAccessor dhcell, unsigned int bulk_idx, DHCellSide dhside, unsigned int side_idx)
111  : cell(dhcell), bulk_subset_index(bulk_idx), side(dhside), side_subset_index(side_idx) {}
112 
113  /// Copy constructor
116 
118  unsigned int bulk_subset_index; ///< Index (order) of lower dim subset in EvalPoints object
119  DHCellSide side; ///< Specified cell side (higher dim element)
120  unsigned int side_subset_index; ///< Index (order) of higher dim subset in EvalPoints object
121  };
122 
123  /**
124  * Helper structzre holds data of boundary integral
125  *
126  * Data is specified by side and subset indices of side and appropriate boundary element in EvalPoint object
127  */
129  /// Default constructor
131 
132  /// Constructor with data mebers initialization
133  BoundaryIntegralData(unsigned int bdr_idx, DHCellSide dhside, unsigned int side_idx)
134  : bdr_subset_index(bdr_idx), side(dhside), side_subset_index(side_idx) {}
135 
136  /// Copy constructor
139 
140  // We don't need hold ElementAccessor of boundary element, side.cond().element_accessor() provides it.
141  unsigned int bdr_subset_index; ///< Index (order) of subset on boundary element in EvalPoints object
142  DHCellSide side; ///< Specified cell side (bulk element)
143  unsigned int side_subset_index; ///< Index (order) of subset on side of bulk element in EvalPoints object
144  };
145 
147 
149  virtual void assemble(std::shared_ptr<DOFHandlerMultiDim> dh) = 0;
150 
151  /// Getter to EvalPoints object
152  inline std::shared_ptr<EvalPoints> eval_points() const {
153  return eval_points_;
154  }
155 
156 protected:
157  AssemblyIntegrals integrals_; ///< Holds integral objects.
158  std::shared_ptr<EvalPoints> eval_points_; ///< EvalPoints object shared by all integrals
159  ElementCacheMap element_cache_map_; ///< ElementCacheMap according to EvalPoints
160 };
161 
162 
163 /**
164  * @brief Generic class of assemblation.
165  *
166  * Class
167  * - holds assemblation structures (EvalPoints, Integral objects, Integral data table).
168  * - associates assemblation objects specified by dimension
169  * - provides general assemble method
170  * - provides methods that allow construction of element patches
171  */
172 template < template<IntDim...> class DimAssembly>
174 {
175 public:
176  /// Constructor
177  GenericAssembly( typename DimAssembly<1>::EqFields *eq_fields, typename DimAssembly<1>::EqData *eq_data)
178  : use_patch_fe_values_(false),
179  multidim_assembly_(eq_fields, eq_data),
180  min_edge_sides_(2),
181  bulk_integral_data_(20, 10),
182  edge_integral_data_(12, 6),
185  {
186  initialize();
187  }
188 
189  /// Constructor
190  GenericAssembly( typename DimAssembly<1>::EqFields *eq_fields, typename DimAssembly<1>::EqData *eq_data, DOFHandlerMultiDim* dh)
191  : fe_values_(eq_data->quad_order(), dh->ds()->fe()),
192  use_patch_fe_values_(true),
193  multidim_assembly_(eq_fields, eq_data, &this->fe_values_),
194  min_edge_sides_(2),
195  bulk_integral_data_(20, 10),
196  edge_integral_data_(12, 6),
199  {
200  initialize();
201  }
202 
203  /// Getter to set of assembly objects
205  return multidim_assembly_;
206  }
207 
208  void set_min_edge_sides(unsigned int val) {
209  min_edge_sides_ = val;
210  }
211 
212  /**
213  * @brief General assemble methods.
214  *
215  * Loops through local cells and calls assemble methods of assembly
216  * object of each cells over space dimension.
217  *
218  * TODO:
219  * - make estimate of the cache fill for combination of (integral_type x element dimension)
220  * - add next cell to patch if current_patch_size + next_element_size <= fixed_cache_size
221  * - avoid reverting the integral data lists.
222  */
223  void assemble(std::shared_ptr<DOFHandlerMultiDim> dh) override {
224  START_TIMER( DimAssembly<1>::name() );
225  this->reallocate_cache();
226  multidim_assembly_[1_d]->begin();
227 
228  bool add_into_patch = false; // control variable
229  for(auto cell_it = dh->local_range().begin(); cell_it != dh->local_range().end(); )
230  {
231 
232  if (!add_into_patch) {
234  add_into_patch = true;
235  }
236 
237  START_TIMER("add_integrals_to_patch");
238  this->add_integrals_of_computing_step(*cell_it);
239  END_TIMER("add_integrals_to_patch");
240 
247  this->assemble_integrals();
248  add_into_patch = false;
249  } else {
256  this->assemble_integrals();
257  add_into_patch = false;
258  }
259  ++cell_it;
260  }
261  }
262  if (add_into_patch) {
263  this->assemble_integrals();
264  }
265 
266  multidim_assembly_[1_d]->end();
267  END_TIMER( DimAssembly<1>::name() );
268  }
269 
270  /// Return ElementCacheMap
271  inline const ElementCacheMap &cache_map() const {
272  return element_cache_map_;
273  }
274 
275 private:
276  /// Common part of GenericAssemblz constructors.
277  void initialize() {
278  eval_points_ = std::make_shared<EvalPoints>();
279  // first step - create integrals, then - initialize cache and initialize subobject of dimensions
280  multidim_assembly_[1_d]->create_integrals(eval_points_, integrals_);
281  multidim_assembly_[2_d]->create_integrals(eval_points_, integrals_);
282  multidim_assembly_[3_d]->create_integrals(eval_points_, integrals_);
284  multidim_assembly_[1_d]->initialize(&element_cache_map_);
285  multidim_assembly_[2_d]->initialize(&element_cache_map_);
286  multidim_assembly_[3_d]->initialize(&element_cache_map_);
287  active_integrals_ = multidim_assembly_[1_d]->n_active_integrals();
288  }
289 
290  /// Call assemblations when patch is filled
292  START_TIMER("create_patch");
294  END_TIMER("create_patch");
295  if (use_patch_fe_values_) {
296  START_TIMER("patch_reinit");
297  patch_reinit();
298  END_TIMER("patch_reinit");
299  }
300  START_TIMER("cache_update");
301  multidim_assembly_[1_d]->eq_fields_->cache_update(element_cache_map_); // TODO replace with sub FieldSet
302  END_TIMER("cache_update");
304 
305  {
306  START_TIMER("assemble_volume_integrals");
307  multidim_assembly_[1_d]->assemble_cell_integrals(bulk_integral_data_);
308  multidim_assembly_[2_d]->assemble_cell_integrals(bulk_integral_data_);
309  multidim_assembly_[3_d]->assemble_cell_integrals(bulk_integral_data_);
310  END_TIMER("assemble_volume_integrals");
311  }
312 
313  {
314  START_TIMER("assemble_fluxes_boundary");
315  multidim_assembly_[1_d]->assemble_boundary_side_integrals(boundary_integral_data_);
316  multidim_assembly_[2_d]->assemble_boundary_side_integrals(boundary_integral_data_);
317  multidim_assembly_[3_d]->assemble_boundary_side_integrals(boundary_integral_data_);
318  END_TIMER("assemble_fluxes_boundary");
319  }
320 
321  {
322  START_TIMER("assemble_fluxes_elem_elem");
323  multidim_assembly_[1_d]->assemble_edge_integrals(edge_integral_data_);
324  multidim_assembly_[2_d]->assemble_edge_integrals(edge_integral_data_);
325  multidim_assembly_[3_d]->assemble_edge_integrals(edge_integral_data_);
326  END_TIMER("assemble_fluxes_elem_elem");
327  }
328 
329  {
330  START_TIMER("assemble_fluxes_elem_side");
331  multidim_assembly_[2_d]->assemble_neighbour_integrals(coupling_integral_data_);
332  multidim_assembly_[3_d]->assemble_neighbour_integrals(coupling_integral_data_);
333  END_TIMER("assemble_fluxes_elem_side");
334  }
335  // clean integral data
341  if (use_patch_fe_values_) {
343  fe_values_.reset();
344  }
345  }
346 
347  void patch_reinit() {
350  multidim_assembly_[1_d]->add_patch_bulk_points(bulk_integral_data_);
351  multidim_assembly_[2_d]->add_patch_bulk_points(bulk_integral_data_);
352  multidim_assembly_[3_d]->add_patch_bulk_points(bulk_integral_data_);
353  }
355  multidim_assembly_[1_d]->add_patch_bdr_side_points(boundary_integral_data_);
356  multidim_assembly_[2_d]->add_patch_bdr_side_points(boundary_integral_data_);
357  multidim_assembly_[3_d]->add_patch_bdr_side_points(boundary_integral_data_);
358  }
360  multidim_assembly_[1_d]->add_patch_edge_points(edge_integral_data_);
361  multidim_assembly_[2_d]->add_patch_edge_points(edge_integral_data_);
362  multidim_assembly_[3_d]->add_patch_edge_points(edge_integral_data_);
363  }
365  multidim_assembly_[2_d]->add_patch_coupling_integrals(coupling_integral_data_);
366  multidim_assembly_[3_d]->add_patch_coupling_integrals(coupling_integral_data_);
367  }
368  this->fe_values_.reinit_patch();
369  }
370 
371  /**
372  * Add data of integrals to appropriate structure and register elements to ElementCacheMap.
373  *
374  * Types of used integrals must be set in data member \p active_integrals_.
375  */
378  if (cell.is_own()) { // Not ghost
379  this->add_volume_integral(cell);
380  }
381 
382  for( DHCellSide cell_side : cell.side_range() ) {
384  if (cell.is_own()) // Not ghost
385  if ( (cell_side.side().edge().n_sides() == 1) && (cell_side.side().is_boundary()) ) {
386  this->add_boundary_integral(cell_side);
387  continue;
388  }
390  if ( (cell_side.n_edge_sides() >= min_edge_sides_) && (cell_side.edge_sides().begin()->element().idx() == cell.elm_idx())) {
391  this->add_edge_integral(cell_side);
392  }
393  }
394 
396  bool add_low = true;
397  for( DHCellSide neighb_side : cell.neighb_sides() ) { // cell -> elm lower dim, neighb_side -> elm higher dim
398  if (cell.dim() != neighb_side.dim()-1) continue;
399  this->add_coupling_integral(cell, neighb_side, add_low);
400  add_low = false;
401  }
402  }
403  }
404 
405  /// Add data of volume integral to appropriate data structure.
406  inline void add_volume_integral(const DHCellAccessor &cell) {
407  uint subset_idx = integrals_.bulk_[cell.dim()-1]->get_subset_idx();
408  bulk_integral_data_.emplace_back(cell, subset_idx);
409 
410  unsigned int reg_idx = cell.elm().region_idx().idx();
411  table_sizes_.elem_sizes_[0][cell.dim()-1]++;
412  // Different access than in other integrals: We can't use range method CellIntegral::points
413  // because it passes element_patch_idx as argument that is not known during patch construction.
414  for (uint i=uint( eval_points_->subset_begin(cell.dim(), subset_idx) );
415  i<uint( eval_points_->subset_end(cell.dim(), subset_idx) ); ++i) {
416  element_cache_map_.add_eval_point(reg_idx, cell.elm_idx(), i, cell.local_idx());
417  table_sizes_.point_sizes_[0][cell.dim()-1]++;
418  }
419  }
420 
421  /// Add data of edge integral to appropriate data structure.
422  inline void add_edge_integral(const DHCellSide &cell_side) {
423  auto range = cell_side.edge_sides();
424  uint dim = range.begin()->dim();
425  edge_integral_data_.emplace_back(range, integrals_.edge_[dim-1]->get_subset_idx());
426  table_sizes_.elem_sizes_[1][dim-1]++;
427 
428  for( DHCellSide edge_side : range ) {
429  unsigned int reg_idx = edge_side.element().region_idx().idx();
430  for (auto p : integrals_.edge_[dim-1]->points(edge_side, &element_cache_map_) ) {
431  element_cache_map_.add_eval_point(reg_idx, edge_side.elem_idx(), p.eval_point_idx(), edge_side.cell().local_idx());
432  table_sizes_.point_sizes_[1][dim-1]++;
433  }
434  }
435  }
436 
437  /// Add data of coupling integral to appropriate data structure.
438  inline void add_coupling_integral(const DHCellAccessor &cell, const DHCellSide &ngh_side, bool add_low) {
439  coupling_integral_data_.emplace_back(cell, integrals_.coupling_[cell.dim()-1]->get_subset_low_idx(), ngh_side,
440  integrals_.coupling_[cell.dim()-1]->get_subset_high_idx());
441  table_sizes_.elem_sizes_[1][cell.dim()]++;
442  if (add_low) table_sizes_.elem_sizes_[0][cell.dim()-1]++;
443 
444  unsigned int reg_idx_low = cell.elm().region_idx().idx();
445  unsigned int reg_idx_high = ngh_side.element().region_idx().idx();
446  for (auto p : integrals_.coupling_[cell.dim()-1]->points(ngh_side, &element_cache_map_) ) {
447  element_cache_map_.add_eval_point(reg_idx_high, ngh_side.elem_idx(), p.eval_point_idx(), ngh_side.cell().local_idx());
448  table_sizes_.point_sizes_[1][cell.dim()]++;
449 
450  if (add_low) {
451  auto p_low = p.lower_dim(cell); // equivalent point on low dim cell
452  element_cache_map_.add_eval_point(reg_idx_low, cell.elm_idx(), p_low.eval_point_idx(), cell.local_idx());
453  table_sizes_.point_sizes_[0][cell.dim()-1]++;
454  }
455  }
456  }
457 
458  /// Add data of boundary integral to appropriate data structure.
459  inline void add_boundary_integral(const DHCellSide &bdr_side) {
460  boundary_integral_data_.emplace_back(integrals_.boundary_[bdr_side.dim()-1]->get_subset_low_idx(), bdr_side,
461  integrals_.boundary_[bdr_side.dim()-1]->get_subset_high_idx());
462 
463  unsigned int reg_idx = bdr_side.element().region_idx().idx();
464  table_sizes_.elem_sizes_[1][bdr_side.dim()-1]++;
465  for (auto p : integrals_.boundary_[bdr_side.dim()-1]->points(bdr_side, &element_cache_map_) ) {
466  element_cache_map_.add_eval_point(reg_idx, bdr_side.elem_idx(), p.eval_point_idx(), bdr_side.cell().local_idx());
467  table_sizes_.point_sizes_[1][bdr_side.dim()-1]++;
468 
469  BulkPoint p_bdr = p.point_bdr(bdr_side.cond().element_accessor()); // equivalent point on boundary element
470  unsigned int bdr_reg = bdr_side.cond().element_accessor().region_idx().idx();
471  // invalid local_idx value, DHCellAccessor of boundary element doesn't exist
472  element_cache_map_.add_eval_point(bdr_reg, bdr_side.cond().bc_ele_idx(), p_bdr.eval_point_idx(), -1);
473  }
474  }
475 
476  /// Calls cache_reallocate method on
477  inline void reallocate_cache() {
478  multidim_assembly_[1_d]->eq_fields_->cache_reallocate(this->element_cache_map_, multidim_assembly_[1_d]->used_fields_);
479  // DebugOut() << "Order of evaluated fields (" << DimAssembly<1>::name() << "):" << multidim_assembly_[1_d]->eq_fields_->print_dependency();
480  }
481 
482 
483  PatchFEValues<3> fe_values_; ///< Common FEValues object over all dimensions
484  bool use_patch_fe_values_; ///< Flag holds if common @p fe_values_ object is used in @p multidim_assembly_
486 
487  /// Holds mask of active integrals.
489 
490  /**
491  * Minimal number of sides on edge.
492  *
493  * Edge integral is created and calculated if number of sides is greater or equal than this value. Default value
494  * is 2 and can be changed
495  */
496  unsigned int min_edge_sides_;
497 
498  // Following variables hold data of all integrals depending of actual computed element.
499  // TODO sizes of arrays should be set dynamically, depend on number of elements in ElementCacheMap,
500  RevertableList<BulkIntegralData> bulk_integral_data_; ///< Holds data for computing bulk integrals.
501  RevertableList<EdgeIntegralData> edge_integral_data_; ///< Holds data for computing edge integrals.
502  RevertableList<CouplingIntegralData> coupling_integral_data_; ///< Holds data for computing couplings integrals.
503  RevertableList<BoundaryIntegralData> boundary_integral_data_; ///< Holds data for computing boundary integrals.
504 
505  /**
506  * Struct for pre-computing number of elements, sides, bulk points and side points on each dimension.
507  */
509 };
510 
511 
512 #endif /* GENERIC_ASSEMBLY_HH_ */
uint bc_ele_idx()
Definition: accessors.hh:374
ElementAccessor< 3 > element_accessor()
Base point accessor class.
Definition: eval_subset.hh:55
unsigned int eval_point_idx() const
Return index in EvalPoints object.
Definition: eval_subset.hh:84
static unsigned int get()
Return number of stored elements.
Cell accessor allow iterate over DOF handler cells.
RangeConvert< DHNeighbSide, DHCellSide > neighb_sides() const
Returns range of neighbour cell of lower dimension corresponding to cell of higher dimension.
bool is_own() const
Return true if accessor represents own element (false for ghost element)
unsigned int dim() const
Return dimension of element appropriate to cell.
Range< DHCellSide > side_range() const
Returns range of cell sides.
unsigned int elm_idx() const
Return serial idx to element of loc_ele_idx_.
ElementAccessor< 3 > elm() const
Return ElementAccessor to element of loc_ele_idx_.
unsigned int local_idx() const
Return local index to element (index of DOF handler).
Side accessor allows to iterate over sides of DOF handler cell.
unsigned int elem_idx() const
Boundary cond() const
const DHCellAccessor & cell() const
Return DHCellAccessor appropriate to the side.
RangeConvert< DHEdgeSide, DHCellSide > edge_sides() const
Returns range of all sides looped over common Edge.
unsigned int dim() const
Return dimension of element appropriate to the side.
ElementAccessor< 3 > element() const
Class allows to iterate over sides of edge.
Provides the numbering of the finite element degrees of freedom on the computational mesh.
Definition: dofhandler.hh:151
RegionIdx region_idx() const
Definition: accessors.hh:201
Directing class of FieldValueCache.
unsigned int get_simd_rounded_size()
Returns number of eval. points with addition of max simd duplicates due to regions.
void start_elements_update()
Start update of cache.
void finish_elements_update()
Finish update after reading data to cache.
void clear_element_eval_points_map()
Reset all items of elements_eval_points_map.
void create_patch()
Create patch of cached elements before reading data to cache.
RevertableList< EvalPointData > eval_point_data_
void add_eval_point(unsigned int i_reg, unsigned int i_ele, unsigned int i_eval_point, unsigned int dh_loc_idx)
void init(std::shared_ptr< EvalPoints > eval_points)
Init cache.
void make_paermanent_eval_points()
Mark eval_point_data_ as permanent.
ElementCacheMap element_cache_map_
ElementCacheMap according to EvalPoints.
std::shared_ptr< EvalPoints > eval_points_
EvalPoints object shared by all integrals.
std::shared_ptr< EvalPoints > eval_points() const
Getter to EvalPoints object.
AssemblyIntegrals integrals_
Holds integral objects.
virtual void assemble(std::shared_ptr< DOFHandlerMultiDim > dh)=0
Generic class of assemblation.
PatchFEValues< 3 > fe_values_
Common FEValues object over all dimensions.
MixedPtr< DimAssembly, 1 > multidim_assembly_
Assembly object.
RevertableList< EdgeIntegralData > edge_integral_data_
Holds data for computing edge integrals.
void assemble_integrals()
Call assemblations when patch is filled.
void assemble(std::shared_ptr< DOFHandlerMultiDim > dh) override
General assemble methods.
GenericAssembly(typename DimAssembly< 1 >::EqFields *eq_fields, typename DimAssembly< 1 >::EqData *eq_data, DOFHandlerMultiDim *dh)
Constructor.
unsigned int min_edge_sides_
void add_boundary_integral(const DHCellSide &bdr_side)
Add data of boundary integral to appropriate data structure.
const ElementCacheMap & cache_map() const
Return ElementCacheMap.
void add_volume_integral(const DHCellAccessor &cell)
Add data of volume integral to appropriate data structure.
void add_edge_integral(const DHCellSide &cell_side)
Add data of edge integral to appropriate data structure.
PatchFEValues< 3 >::TableSizes table_sizes_
RevertableList< BoundaryIntegralData > boundary_integral_data_
Holds data for computing boundary integrals.
int active_integrals_
Holds mask of active integrals.
void initialize()
Common part of GenericAssemblz constructors.
RevertableList< CouplingIntegralData > coupling_integral_data_
Holds data for computing couplings integrals.
void set_min_edge_sides(unsigned int val)
void reallocate_cache()
Calls cache_reallocate method on.
GenericAssembly(typename DimAssembly< 1 >::EqFields *eq_fields, typename DimAssembly< 1 >::EqData *eq_data)
Constructor.
bool use_patch_fe_values_
Flag holds if common fe_values_ object is used in multidim_assembly_.
void add_coupling_integral(const DHCellAccessor &cell, const DHCellSide &ngh_side, bool add_low)
Add data of coupling integral to appropriate data structure.
void add_integrals_of_computing_step(DHCellAccessor cell)
MixedPtr< DimAssembly, 1 > multidim_assembly() const
Getter to set of assembly objects.
RevertableList< BulkIntegralData > bulk_integral_data_
Holds data for computing bulk integrals.
void resize_tables(TableSizes table_sizes)
Resize tables of patch_point_vals_.
void reinit_patch()
Reinit data.
void reset()
Reset PatchpointValues structures.
IterConvert< ObjectIn, ObjectOut > begin()
Iterator to begin item of range.
unsigned int idx() const
Returns a global index of the region.
Definition: region.hh:81
Class FEValues calculates finite element data on the actual cells such as shape function values,...
Iter< Object > make_iter(Object obj)
ActiveIntegrals
Allow set mask of active integrals.
@ coupling
@ boundary
@ no_intg
@ bulk
@ edge
unsigned int uint
unsigned int IntDim
A dimension index type.
Definition: mixed.hh:19
Class FEValues calculates finite element data on the actual cells such as shape function values,...
Definitions of particular quadrature rules on simplices.
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.
unsigned int side_subset_index
Index (order) of subset on side of bulk element in EvalPoints object.
BoundaryIntegralData(unsigned int bdr_idx, DHCellSide dhside, unsigned int side_idx)
Constructor with data mebers initialization.
BoundaryIntegralData(const BoundaryIntegralData &other)
Copy constructor.
unsigned int bdr_subset_index
Index (order) of subset on boundary element in EvalPoints object.
DHCellSide side
Specified cell side (bulk element)
BulkIntegralData(DHCellAccessor dhcell, unsigned int subset_idx)
Constructor with data mebers initialization.
DHCellAccessor cell
Specified cell (element)
unsigned int subset_index
Index (order) of subset in EvalPoints object.
BulkIntegralData(const BulkIntegralData &other)
Copy constructor.
unsigned int bulk_subset_index
Index (order) of lower dim subset in EvalPoints object.
CouplingIntegralData(DHCellAccessor dhcell, unsigned int bulk_idx, DHCellSide dhside, unsigned int side_idx)
Constructor with data mebers initialization.
unsigned int side_subset_index
Index (order) of higher dim subset in EvalPoints object.
CouplingIntegralData(const CouplingIntegralData &other)
Copy constructor.
DHCellSide side
Specified cell side (higher dim element)
EdgeIntegralData(const EdgeIntegralData &other)
Copy constructor.
RangeConvert< DHEdgeSide, DHCellSide > edge_side_range
Specified cell side (element)
unsigned int subset_index
Index (order) of subset in EvalPoints object.
EdgeIntegralData(RangeConvert< DHEdgeSide, DHCellSide > range, unsigned int subset_idx)
Constructor with data mebers initialization.
std::size_t revert_temporary()
Erase temporary part of data.
std::size_t make_permanent()
Finalize temporary part of data.
void reset()
Clear the list.
std::size_t emplace_back(Args &&... args)
std::size_t permanent_size() const
Return permanent size of list.
#define END_TIMER(tag)
Ends a timer with specified tag.
#define START_TIMER(tag)
Starts a timer with specified tag.