Flow123d  master-f44eb46
element_values.hh
Go to the documentation of this file.
1 /*!
2  *
3  * Copyright (C) 2019 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 element_values.hh
15  * @brief Class ElementValues calculates data related to transformation
16  * of reference cell to actual cell (Jacobian, inverse Jacobian,
17  * determinant, point coordinates etc.).
18  * @author Jan Stebel
19  */
20 
21 #ifndef ELEMENT_VALUES_HH_
22 #define ELEMENT_VALUES_HH_
23 
24 // #include <string.h> // for memcpy
25 // #include <algorithm> // for swap
26 // #include <new> // for operator new[]
27 // #include <string> // for operator<<
28 #include <vector> // for vector
29 #include "mesh/ref_element.hh" // for RefElement
30 #include "mesh/accessors.hh" // for ElementAccessor
31 #include "fem/update_flags.hh" // for UpdateFlags
32 #include "fem/dh_cell_accessor.hh" // for DHCellAccessor, DHCellSide
33 
34 class Quadrature;
35 
36 
37 
38 
39 /**
40  * @brief Structure for storing the precomputed element data.
41  */
43 {
44 public:
45 
46  /// Resize vectors to size @p np.
47  RefElementData(unsigned int np);
48 
49  /// Barycentric coordinates of quadrature points.
51 
52  /// Quadrature weights.
54 
55  /// Number of quadrature points.
56  unsigned int n_points;
57 };
58 
59 
60 
61 
62 /**
63  * @brief Class ElementData holds the arrays of data computed by
64  * Mapping.
65  */
66 template<unsigned int spacedim = 3>
68 {
69 public:
70 
71  /**
72  * @brief Resize the data arrays.
73  * @param size Number of quadrature points.
74  * @param flags Update flags to be stores.
75  */
76  ElementData(unsigned int size, UpdateFlags flags, unsigned int dim);
77 
78  /// Print calculated data.
79  void print();
80 
81  /// Dimension of space of reference cell.
82  const unsigned int dim_;
83 
84  /**
85  * @brief Transformed quadrature weights.
86  *
87  * Values at quadrature points of the product of the Jacobian
88  * determinant of the mapping and the weight at the particular
89  * quadrature point.
90  */
92 
93  /// JxW values for sides.
95 
96  /// Jacobians (spacedim x dim) of the mapping at the quadrature points.
98 
99  /// Determinants of Jacobians at quadrature points.
101 
102  /// Inverse Jacobians (dim x spacedim) at the quadrature points.
104 
105  /// Coordinates (spacedim) of quadrature points in the actual cell coordinate system.
107 
108  /// Normal vectors (spacedim) to the element at the quadrature points lying on a side.
110 
111  /// Flags that indicate which finite element quantities are to be computed.
113 
114  /// Iterator to last updated cell.
116 
117  /// Iterator to last updated cell side.
119 
120 };
121 
122 
123 
124 /**
125  * @brief Class for computation of data on cell and side.
126  */
127 template<unsigned int spacedim>
129 {
130 public:
131 
132  /**
133  * @brief Constructor.
134  *
135  * Initializes structures and calculates
136  * cell-independent data. The quadrature can have dimension
137  * either dim or dim-1 (for values on side). In the later
138  * case also normal vectors and side jacobians are computed.
139  *
140  * @param _quadrature The quadrature rule.
141  * @param _flags The update flags.
142  * @param dim Dimension of space of reference cell.
143  */
144  ElementValues(Quadrature &_quadrature,
145  UpdateFlags _flags,
146  unsigned int dim);
147 
148  /// Correct deallocation of objects created by 'initialize' methods.
149  ~ElementValues();
150 
151 
152 
153  /**
154  * @brief Update cell-dependent data (gradients, Jacobians etc.)
155  *
156  * @param cell The actual cell.
157  */
159 
160  /**
161  * @brief Update side-dependent data (Jacobians etc.)
162  *
163  * @param cell_side The actual cell and side.
164  */
165  void reinit(const Side &cell_side);
166 
167  /**
168  * @brief Determine quantities to be recomputed on each cell.
169  *
170  * @param flags Flags that indicate what has to be recomputed.
171  */
173 
174  /**
175  * @brief Return the relative volume change of the cell (Jacobian determinant).
176  *
177  * If dim==spacedim then the sign may be negative, otherwise the
178  * result is a positive number.
179  *
180  * @param point_no Number of the quadrature point.
181  */
182  inline double determinant(const unsigned int point_no) const
183  {
184  ASSERT_LT(point_no, n_points_);
185  return data.determinants[point_no];
186  }
187 
188  /// Return Jacobian matrix at point @p point_no.
189  inline arma::mat jacobian(const unsigned int point_no) const
190  {
191  ASSERT_LT(point_no, n_points_);
192  return data.jacobians.arma_mat(point_no);
193  }
194 
195  /// Return inverse Jacobian matrix at point @p point_no.
196  inline arma::mat inverse_jacobian(const unsigned int point_no) const
197  {
198  ASSERT_LT(point_no, n_points_);
199  return data.inverse_jacobians.arma_mat(point_no);
200  }
201 
202  /**
203  * @brief Return the product of Jacobian determinant and the quadrature
204  * weight at given quadrature point.
205  *
206  * @param point_no Number of the quadrature point.
207  */
208  inline double JxW(const unsigned int point_no) const
209  {
210  ASSERT_LT(point_no, n_points_);
211  return data.JxW_values[point_no];
212  }
213 
214  /**
215  * @brief Return the product of side Jacobian determinant and the quadrature
216  * weight at given quadrature point.
217  *
218  * @param point_no Number of the quadrature point.
219  */
220  inline double side_JxW(const unsigned int point_no) const
221  {
222  ASSERT_LT(point_no, n_points_);
223  return data.side_JxW_values[point_no];
224  }
225 
226  /**
227  * @brief Return coordinates of the quadrature point in the actual cell system.
228  *
229  * @param point_no Number of the quadrature point.
230  */
231  inline arma::vec::fixed<spacedim> point(const unsigned int point_no) const
232  {
233  ASSERT_LT(point_no, n_points_);
234  return data.points.template vec<spacedim>(point_no);
235  }
236 
237  /// Return coordinates of all quadrature points in the actual cell system.
238  inline const Armor::array &point_list() const
239  {
240  return data.points;
241  }
242 
243 
244  /**
245  * @brief Returns the normal vector to a side at given quadrature point.
246  *
247  * @param point_no Number of the quadrature point.
248  */
249  inline arma::vec::fixed<spacedim> normal_vector(unsigned int point_no)
250  {
251  ASSERT_LT(point_no, n_points_);
252  return data.normal_vectors.template vec<spacedim>(point_no);
253  }
254 
255  /// Returns the number of quadrature points.
256  inline unsigned int n_points()
257  { return n_points_; }
258 
259  /// Return cell at which the values were reinited.
261  { return data.cell; }
262 
263  /// Return cell side where the values were reinited.
264  const Side &side() const
265  { return data.side; }
266 
267 
268 
269 protected:
270 
271  /// Precompute data on reference element.
273 
274  /// Compute data from reference cell and using MappingP1.
275  template<unsigned int dim>
276  void fill_data();
277 
278  /// Calculates the mapping data on a side of a cell.
279  template<unsigned int dim>
280  void fill_side_data();
281 
282 
283 
284  /// Dimension of space of reference cell.
285  const unsigned int dim_;
286 
287  /// Number of integration points.
288  const unsigned int n_points_;
289 
290  /// Number of sides in reference cell.
291  const unsigned int n_sides_;
292 
293  /// Data on reference element.
295 
296  /// Data on reference element (for each side ).
298 
299  /// Data computed by the mapping.
301 
302 };
303 
304 #endif /* ELEMENT_VALUES_HH_ */
#define ASSERT_LT(a, b)
Definition of comparative assert macro (Less Than) only for debug mode.
Definition: asserts.hh:301
arma::mat arma_mat(uint i) const
Definition: armor.hh:848
Class ElementData holds the arrays of data computed by Mapping.
UpdateFlags update_flags
Flags that indicate which finite element quantities are to be computed.
Armor::array inverse_jacobians
Inverse Jacobians (dim x spacedim) at the quadrature points.
std::vector< double > determinants
Determinants of Jacobians at quadrature points.
Armor::array jacobians
Jacobians (spacedim x dim) of the mapping at the quadrature points.
ElementAccessor< spacedim > cell
Iterator to last updated cell.
ElementData(unsigned int size, UpdateFlags flags, unsigned int dim)
Resize the data arrays.
void print()
Print calculated data.
std::vector< double > JxW_values
Transformed quadrature weights.
Armor::array normal_vectors
Normal vectors (spacedim) to the element at the quadrature points lying on a side.
Armor::array points
Coordinates (spacedim) of quadrature points in the actual cell coordinate system.
std::vector< double > side_JxW_values
JxW values for sides.
const unsigned int dim_
Dimension of space of reference cell.
Side side
Iterator to last updated cell side.
Class for computation of data on cell and side.
const ElementAccessor< spacedim > & cell() const
Return cell at which the values were reinited.
arma::vec::fixed< spacedim > normal_vector(unsigned int point_no)
Returns the normal vector to a side at given quadrature point.
arma::vec::fixed< spacedim > point(const unsigned int point_no) const
Return coordinates of the quadrature point in the actual cell system.
arma::mat jacobian(const unsigned int point_no) const
Return Jacobian matrix at point point_no.
double determinant(const unsigned int point_no) const
Return the relative volume change of the cell (Jacobian determinant).
double side_JxW(const unsigned int point_no) const
Return the product of side Jacobian determinant and the quadrature weight at given quadrature point.
const unsigned int n_sides_
Number of sides in reference cell.
RefElementData * ref_data
Data on reference element.
const unsigned int n_points_
Number of integration points.
ElementValues(Quadrature &_quadrature, UpdateFlags _flags, unsigned int dim)
Constructor.
const Armor::array & point_list() const
Return coordinates of all quadrature points in the actual cell system.
arma::mat inverse_jacobian(const unsigned int point_no) const
Return inverse Jacobian matrix at point point_no.
const Side & side() const
Return cell side where the values were reinited.
ElementData< spacedim > data
Data computed by the mapping.
RefElementData * init_ref_data(const Quadrature &q)
Precompute data on reference element.
void fill_side_data()
Calculates the mapping data on a side of a cell.
~ElementValues()
Correct deallocation of objects created by 'initialize' methods.
std::vector< RefElementData * > side_ref_data
Data on reference element (for each side ).
void reinit(const ElementAccessor< spacedim > &cell)
Update cell-dependent data (gradients, Jacobians etc.)
UpdateFlags update_each(UpdateFlags flags)
Determine quantities to be recomputed on each cell.
double JxW(const unsigned int point_no) const
Return the product of Jacobian determinant and the quadrature weight at given quadrature point.
void fill_data()
Compute data from reference cell and using MappingP1.
unsigned int n_points()
Returns the number of quadrature points.
const unsigned int dim_
Dimension of space of reference cell.
Base class for quadrature rules on simplices in arbitrary dimensions.
Definition: quadrature.hh:48
Structure for storing the precomputed element data.
RefElementData(unsigned int np)
Resize vectors to size np.
unsigned int n_points
Number of quadrature points.
std::vector< arma::vec > bar_coords
Barycentric coordinates of quadrature points.
std::vector< double > weights
Quadrature weights.
ArmaMat< double, N, M > mat
Definition: armor.hh:888
Class RefElement defines numbering of vertices, sides, calculation of normal vectors etc.
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