Flow123d  master-f44eb46
fe_values.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 fe_values.hh
15  * @brief Class FEValues calculates finite element data on the actual
16  * cells such as shape function values, gradients, Jacobian of
17  * the mapping from the reference cell etc.
18  * @author Jan Stebel
19  */
20 
21 #ifndef FE_VALUES_HH_
22 #define FE_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 "fem/element_values.hh" // for ElementValues
30 #include "fem/fe_values_views.hh" // for FEValuesViews
31 #include "mesh/ref_element.hh" // for RefElement
32 #include "mesh/accessors.hh"
33 #include "fem/update_flags.hh" // for UpdateFlags
34 #include "tools/mixed.hh"
36 
37 class Quadrature;
38 template<unsigned int dim> class FiniteElement;
39 
40 template<unsigned int dim> class MapScalar;
41 template<unsigned int dim> class MapPiola;
42 template<unsigned int dim> class MapContravariant;
43 template<unsigned int dim> class MapVector;
44 template<unsigned int dim> class MapTensor;
45 template<unsigned int dim> class MapSystem;
46 
47 
48 
49 
50 
51 
52 
53 
54 
55 /**
56  * @brief Calculates finite element data on the actual cell.
57  *
58  * FEValues takes care of the calculation of finite element data on
59  * the actual cell or on its side, (values of shape functions at quadrature
60  * points, gradients of shape functions, Jacobians of the mapping
61  * from the reference cell etc.).
62  * @param spacedim Dimension of the Euclidean space where the actual
63  * cell lives.
64  */
65 template<unsigned int spacedim = 3>
66 class FEValues
67 {
68 private:
69 
70  // internal structure that stores all possible views
71  // for scalar and vector-valued components of the FE
72  struct ViewsCache {
76 
77  template<unsigned int DIM>
78  void initialize(const FEValues &fv, const FiniteElement<DIM> &fe);
79  };
80 
81 public:
82 
83  /// Default constructor with postponed initialization.
84  FEValues();
85 
86 
87  /// Constructor with initialization of data structures
88  /// (see initialize() for description of parameters).
89  template<unsigned int DIM>
90  FEValues(Quadrature &_quadrature,
91  FiniteElement<DIM> &_fe,
92  UpdateFlags _flags)
93  : FEValues()
94  {
95  initialize(_quadrature, _fe, _flags);
96  }
97 
98 
99  /// Correct deallocation of objects created by 'initialize' methods.
100  ~FEValues();
101 
102 
103  /**
104  * @brief Allocates space for computed data.
105  *
106  * @param n_points Number of quadrature points.
107  * @param _fe The finite element.
108  * @param flags The update flags.
109  */
110  template<unsigned int DIM>
111  void allocate(unsigned int n_points,
112  FiniteElement<DIM> &_fe,
113  UpdateFlags flags);
114 
115  /**
116  * @brief Initialize structures and calculates cell-independent data.
117  *
118  * @param _quadrature The quadrature rule for the cell associated
119  * to given finite element or for the cell side.
120  * @param _fe The finite element.
121  * @param _flags The update flags.
122  */
123  template<unsigned int DIM>
124  void initialize(Quadrature &_quadrature,
125  FiniteElement<DIM> &_fe,
126  UpdateFlags _flags);
127 
128  /**
129  * @brief Update cell-dependent data (gradients, Jacobians etc.)
130  *
131  * @param cell The actual cell.
132  */
133  void reinit(const ElementAccessor<spacedim> &cell);
134 
135  /**
136  * @brief Update cell side-dependent FE data (values, gradients).
137  *
138  * @param cell_side Accessor to cell side.
139  */
140  void reinit(const Side &cell_side);
141 
142  /**
143  * @brief Return the value of the @p function_no-th shape function at
144  * the @p point_no-th quadrature point.
145  *
146  * @param function_no Number of the shape function.
147  * @param point_no Number of the quadrature point.
148  */
149  inline double shape_value(const unsigned int function_no, const unsigned int point_no) const
150  {
151  ASSERT_LT(function_no, n_dofs_);
152  ASSERT_LT(point_no, n_points_);
153  return shape_values[point_no][function_no];
154  }
155 
156 
157  /**
158  * @brief Return the gradient of the @p function_no-th shape function at
159  * the @p point_no-th quadrature point.
160  *
161  * @param function_no Number of the shape function.
162  * @param point_no Number of the quadrature point.
163  */
164  inline arma::vec::fixed<spacedim> shape_grad(const unsigned int function_no, const unsigned int point_no) const
165  {
166  ASSERT_LT(function_no, n_dofs_);
167  ASSERT_LT(point_no, n_points_);
168  return shape_gradients[point_no][function_no];
169  }
170 
171  /**
172  * @brief Return the value of the @p function_no-th shape function at
173  * the @p point_no-th quadrature point.
174  *
175  * For vectorial finite elements.
176  *
177  * @param function_no Number of the shape function.
178  * @param point_no Number of the quadrature point.
179  */
180  inline double shape_value_component(const unsigned int function_no,
181  const unsigned int point_no,
182  const unsigned int comp) const
183  {
184  ASSERT_LT(function_no, n_dofs_);
185  ASSERT_LT(point_no, n_points_);
186  ASSERT_LT(comp, n_components_);
187  return shape_values[point_no][function_no*n_components_+comp];
188  }
189 
190  /**
191  * @brief Return the gradient of the @p function_no-th shape function at
192  * the @p point_no-th quadrature point.
193  *
194  * For vectorial finite elements.
195  *
196  * @param function_no Number of the shape function.
197  * @param point_no Number of the quadrature point.
198  */
199  arma::vec::fixed<spacedim> shape_grad_component(const unsigned int function_no,
200  const unsigned int point_no,
201  const unsigned int comp) const;
202 
203  /**
204  * @brief Return the relative volume change of the cell (Jacobian determinant).
205  *
206  * If dim_==spacedim then the sign may be negative, otherwise the
207  * result is a positive number.
208  *
209  * @param point_no Number of the quadrature point.
210  */
211  inline double determinant(const unsigned int point_no)
212  {
213  ASSERT_LT(point_no, n_points_);
214  return elm_values->determinant(point_no);
215  }
216 
217  /**
218  * @brief Return the product of Jacobian determinant and the quadrature
219  * weight at given quadrature point.
220  *
221  * @param point_no Number of the quadrature point.
222  */
223  inline double JxW(const unsigned int point_no)
224  {
225  ASSERT_LT(point_no, n_points_);
226  // TODO: This is temporary solution to distinguish JxW on element and side_JxW on side.
227  // In future we should call the appropriate method in elm_values.
228  return (elm_values->cell().is_valid()) ? elm_values->JxW(point_no) : elm_values->side_JxW(point_no);
229  }
230 
231  /**
232  * @brief Return coordinates of the quadrature point in the actual cell system.
233  *
234  * @param point_no Number of the quadrature point.
235  */
236  inline arma::vec::fixed<spacedim> point(const unsigned int point_no)
237  {
238  ASSERT_LT(point_no, n_points_);
239  return elm_values->point(point_no);
240  }
241 
242  /**
243  * @brief Return coordinates of all quadrature points in the actual cell system.
244  *
245  */
246  inline const Armor::array &point_list() const
247  {
248  return elm_values->point_list();
249  }
250 
251 
252  /**
253  * @brief Returns the normal vector to a side at given quadrature point.
254  *
255  * @param point_no Number of the quadrature point.
256  */
257  inline arma::vec::fixed<spacedim> normal_vector(unsigned int point_no)
258  {
259  ASSERT_LT(point_no, n_points_);
260  return elm_values->normal_vector(point_no);
261  }
262 
263  /**
264  * @brief Accessor to scalar values of multicomponent FE.
265  * @param i Index of scalar component.
266  */
267  const FEValuesViews::Scalar<spacedim> &scalar_view(unsigned int i) const
268  {
269  ASSERT_LT(i, views_cache_.scalars.size());
270  return views_cache_.scalars[i];
271  }
272 
273  /**
274  * @brief Accessor to vector values of multicomponent FE.
275  * @param i Index of first vector component.
276  */
277  const FEValuesViews::Vector<spacedim> &vector_view(unsigned int i) const
278  {
279  ASSERT_LT(i, views_cache_.vectors.size());
280  return views_cache_.vectors[i];
281  }
282 
283  /**
284  * @brief Accessor to tensor values of multicomponent FE.
285  * @param i Index of first tensor component.
286  */
287  const FEValuesViews::Tensor<spacedim> &tensor_view(unsigned int i) const
288  {
289  ASSERT_LT(i, views_cache_.tensors.size());
290  return views_cache_.tensors[i];
291  }
292 
293  /**
294  * @brief Returns the number of quadrature points.
295  */
296  inline unsigned int n_points() const
297  { return n_points_; }
298 
299  /**
300  * @brief Returns the number of shape functions.
301  */
302  inline unsigned int n_dofs() const
303  {
304  return n_dofs_;
305  }
306 
307  /// Return dimension of reference space.
308  inline unsigned int dim() const
309  { return dim_; }
310 
311 
312 protected:
313 
314  /// Structure for storing the precomputed finite element data.
316  {
317  public:
318 
319  FEInternalData(unsigned int np, unsigned int nd);
320 
321  /// Create a new instance of FEInternalData for a FESystem component or subvector.
322  FEInternalData(const FEInternalData &fe_system_data,
323  const std::vector<unsigned int> &dof_indices,
324  unsigned int first_component_idx,
325  unsigned int ncomponents = 1);
326 
327  /**
328  * @brief Precomputed values of basis functions at the quadrature points.
329  *
330  * Dimensions: (no. of quadrature points)
331  * x (no. of dofs)
332  * x (no. of components in ref. cell)
333  */
335 
336  /**
337  * @brief Precomputed gradients of basis functions at the quadrature points.
338  *
339  * Dimensions: (no. of quadrature points)
340  * x (no. of dofs)
341  * x ((dim_ of. ref. cell)x(no. of components in ref. cell))
342  */
344 
345  /// Number of quadrature points.
346  unsigned int n_points;
347 
348  /// Number of dofs (shape functions).
349  unsigned int n_dofs;
350  };
351 
352 
353 
354 
355 
356  /// Precompute finite element data on reference element.
357  template<unsigned int DIM>
358  std::shared_ptr<FEInternalData> init_fe_data(const FiniteElement<DIM> &fe, const Quadrature &q);
359 
360  /**
361  * @brief Computes the shape function values and gradients on the actual cell
362  * and fills the FEValues structure.
363  *
364  * @param fe_data Precomputed finite element data.
365  */
367 
368  /**
369  * @brief Computes the shape function values and gradients on the actual cell
370  * and fills the FEValues structure. Specialized variant of previous method for
371  * different FETypes given by template parameter.
372  */
373  template<class MapType>
375 
376 
377  /// Dimension of reference space.
378  unsigned int dim_;
379 
380  /// Number of integration points.
381  unsigned int n_points_;
382 
383  /// Number of finite element dofs.
384  unsigned int n_dofs_;
385 
386  /// Type of finite element (scalar, vector, tensor).
388 
389  /// Dof indices of FESystem sub-elements.
391 
392  /// Numbers of components of FESystem sub-elements in reference space.
394 
395  /// Numbers of components of FESystem sub-elements in real space.
397 
398  /// Shape functions evaluated at the quadrature points.
400 
401  /// Gradients of shape functions evaluated at the quadrature points.
402  /// Each row of the matrix contains the gradient of one shape function.
404 
405  /// Flags that indicate which finite element quantities are to be computed.
407 
408  /// Auxiliary object for calculation of element-dependent data.
409  std::shared_ptr<ElementValues<spacedim> > elm_values;
410 
411  /// Vector of FEValues for sub-elements of FESystem.
413 
414  /// Number of components of the FE.
415  unsigned int n_components_;
416 
417  /// Auxiliary storage of FEValuesViews accessors.
419 
420  /// Precomputed finite element data.
421  std::shared_ptr<FEInternalData> fe_data;
422 
423  /// Precomputed FE data (shape functions on reference element) for all side quadrature points.
425 
426  friend class MapScalar<spacedim>;
427  friend class MapPiola<spacedim>;
428  friend class MapContravariant<spacedim>;
429  friend class MapVector<spacedim>;
430  friend class MapTensor<spacedim>;
431  friend class MapSystem<spacedim>;
432 };
433 
434 
435 
436 
437 
438 
440  QGauss::array &quadrature,
442  UpdateFlags flags);
443 
444 
445 
446 
447 
448 
449 
450 #endif /* FE_VALUES_HH_ */
#define ASSERT_LT(a, b)
Definition of comparative assert macro (Less Than) only for debug mode.
Definition: asserts.hh:301
Class for computation of data on cell and side.
Structure for storing the precomputed finite element data.
Definition: fe_values.hh:316
std::vector< std::vector< arma::vec > > ref_shape_values
Precomputed values of basis functions at the quadrature points.
Definition: fe_values.hh:334
unsigned int n_dofs
Number of dofs (shape functions).
Definition: fe_values.hh:349
FEInternalData(unsigned int np, unsigned int nd)
Definition: fe_values.cc:41
std::vector< std::vector< arma::mat > > ref_shape_grads
Precomputed gradients of basis functions at the quadrature points.
Definition: fe_values.hh:343
unsigned int n_points
Number of quadrature points.
Definition: fe_values.hh:346
Calculates finite element data on the actual cell.
Definition: fe_values.hh:67
double determinant(const unsigned int point_no)
Return the relative volume change of the cell (Jacobian determinant).
Definition: fe_values.hh:211
~FEValues()
Correct deallocation of objects created by 'initialize' methods.
Definition: fe_values.cc:131
std::shared_ptr< ElementValues< spacedim > > elm_values
Auxiliary object for calculation of element-dependent data.
Definition: fe_values.hh:409
std::vector< unsigned int > fe_sys_n_space_components_
Numbers of components of FESystem sub-elements in real space.
Definition: fe_values.hh:396
std::shared_ptr< FEInternalData > init_fe_data(const FiniteElement< DIM > &fe, const Quadrature &q)
Precompute finite element data on reference element.
double JxW(const unsigned int point_no)
Return the product of Jacobian determinant and the quadrature weight at given quadrature point.
Definition: fe_values.hh:223
std::vector< std::vector< double > > shape_values
Shape functions evaluated at the quadrature points.
Definition: fe_values.hh:399
FEValues()
Default constructor with postponed initialization.
Definition: fe_values.cc:123
const FEValuesViews::Vector< spacedim > & vector_view(unsigned int i) const
Accessor to vector values of multicomponent FE.
Definition: fe_values.hh:277
arma::vec::fixed< spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int comp) const
Return the gradient of the function_no-th shape function at the point_no-th quadrature point.
Definition: fe_values.cc:277
FEType fe_type_
Type of finite element (scalar, vector, tensor).
Definition: fe_values.hh:387
std::vector< unsigned int > fe_sys_n_components_
Numbers of components of FESystem sub-elements in reference space.
Definition: fe_values.hh:393
const FEValuesViews::Tensor< spacedim > & tensor_view(unsigned int i) const
Accessor to tensor values of multicomponent FE.
Definition: fe_values.hh:287
std::vector< shared_ptr< FEInternalData > > side_fe_data
Precomputed FE data (shape functions on reference element) for all side quadrature points.
Definition: fe_values.hh:424
ViewsCache views_cache_
Auxiliary storage of FEValuesViews accessors.
Definition: fe_values.hh:418
unsigned int n_dofs_
Number of finite element dofs.
Definition: fe_values.hh:384
double shape_value(const unsigned int function_no, const unsigned int point_no) const
Return the value of the function_no-th shape function at the point_no-th quadrature point.
Definition: fe_values.hh:149
std::vector< std::vector< unsigned int > > fe_sys_dofs_
Dof indices of FESystem sub-elements.
Definition: fe_values.hh:390
unsigned int dim_
Dimension of reference space.
Definition: fe_values.hh:378
void initialize(Quadrature &_quadrature, FiniteElement< DIM > &_fe, UpdateFlags _flags)
Initialize structures and calculates cell-independent data.
Definition: fe_values.cc:137
const Armor::array & point_list() const
Return coordinates of all quadrature points in the actual cell system.
Definition: fe_values.hh:246
std::vector< std::vector< arma::vec::fixed< spacedim > > > shape_gradients
Definition: fe_values.hh:403
FEValues(Quadrature &_quadrature, FiniteElement< DIM > &_fe, UpdateFlags _flags)
Definition: fe_values.hh:90
void reinit(const ElementAccessor< spacedim > &cell)
Update cell-dependent data (gradients, Jacobians etc.)
Definition: fe_values.cc:537
unsigned int n_points() const
Returns the number of quadrature points.
Definition: fe_values.hh:296
UpdateFlags update_flags
Flags that indicate which finite element quantities are to be computed.
Definition: fe_values.hh:406
arma::vec::fixed< spacedim > normal_vector(unsigned int point_no)
Returns the normal vector to a side at given quadrature point.
Definition: fe_values.hh:257
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int comp) const
Return the value of the function_no-th shape function at the point_no-th quadrature point.
Definition: fe_values.hh:180
void fill_data(const ElementValues< spacedim > &elm_values, const FEInternalData &fe_data)
Computes the shape function values and gradients on the actual cell and fills the FEValues structure.
Definition: fe_values.cc:494
unsigned int n_components_
Number of components of the FE.
Definition: fe_values.hh:415
const FEValuesViews::Scalar< spacedim > & scalar_view(unsigned int i) const
Accessor to scalar values of multicomponent FE.
Definition: fe_values.hh:267
unsigned int n_points_
Number of integration points.
Definition: fe_values.hh:381
arma::vec::fixed< spacedim > point(const unsigned int point_no)
Return coordinates of the quadrature point in the actual cell system.
Definition: fe_values.hh:236
void fill_data_specialized(const ElementValues< spacedim > &elm_values, const FEInternalData &fe_data)
Computes the shape function values and gradients on the actual cell and fills the FEValues structure....
unsigned int dim() const
Return dimension of reference space.
Definition: fe_values.hh:308
std::shared_ptr< FEInternalData > fe_data
Precomputed finite element data.
Definition: fe_values.hh:421
void allocate(unsigned int n_points, FiniteElement< DIM > &_fe, UpdateFlags flags)
Allocates space for computed data.
Definition: fe_values.cc:180
std::vector< FEValues< spacedim > > fe_values_vec
Vector of FEValues for sub-elements of FESystem.
Definition: fe_values.hh:412
arma::vec::fixed< spacedim > shape_grad(const unsigned int function_no, const unsigned int point_no) const
Return the gradient of the function_no-th shape function at the point_no-th quadrature point.
Definition: fe_values.hh:164
unsigned int n_dofs() const
Returns the number of shape functions.
Definition: fe_values.hh:302
Abstract class for the description of a general finite element on a reference simplex in dim dimensio...
Helper class allows update values and gradients of FEValues of FEVectorContravariant type.
Helper class allows update values and gradients of FEValues of FEVectorPiola type.
Helper class allows update values and gradients of FEValues of FEScalar type.
Helper class allows update values and gradients of FEValues of FEMixedSystem type.
Helper class allows update values and gradients of FEValues of FETensor type.
Helper class allows update values and gradients of FEValues of FEVector type.
std::array< QGauss, 4 > array
Base class for quadrature rules on simplices in arbitrary dimensions.
Definition: quadrature.hh:48
Class ElementValues calculates data related to transformation of reference cell to actual cell (Jacob...
std::vector< FEValues< 3 > > mixed_fe_values(QGauss::array &quadrature, MixedPtr< FiniteElement > fe, UpdateFlags flags)
Definition: fe_values.cc:570
FEType
Definitions of particular quadrature rules on simplices.
Class RefElement defines numbering of vertices, sides, calculation of normal vectors etc.
vector< FEValuesViews::Tensor< spacedim > > tensors
Definition: fe_values.hh:75
vector< FEValuesViews::Vector< spacedim > > vectors
Definition: fe_values.hh:74
vector< FEValuesViews::Scalar< spacedim > > scalars
Definition: fe_values.hh:73
void initialize(const FEValues &fv, const FiniteElement< DIM > &fe)
Definition: fe_values.cc:69
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