Flow123d  master-f44eb46
finite_element.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 finite_element.hh
15  * @brief Abstract class for description of finite elements.
16  * @author Jan Stebel
17  */
18 
19 #ifndef FINITE_ELEMENT_HH_
20 #define FINITE_ELEMENT_HH_
21 
22 #include <armadillo>
23 #include <stdio.h> // for sprintf
24 #include <string.h> // for memcpy
25 #include <algorithm> // for max, min
26 
27 #include <new> // for operator new[]
28 #include <ostream> // for operator<<
29 #include <string> // for operator<<
30 #include <vector> // for vector
31 #include "fem/update_flags.hh" // for operator&, operator|=
32 #include "system/exceptions.hh" // for ExcAssertMsg::~ExcAsse...
33 
34 template<unsigned int dim> class FESystem;
35 template<unsigned int spacedim> class FEValues;
36 template<unsigned int dim> class FE_P_disc;
37 
38 
39 
40 
41 
42 // Possible types are: value, gradient, cell integral, ...
43 enum DofType { Value = 1 };
44 
45 /**
46  * Class Dof is a general description of functionals (degrees of freedom)
47  * determining the finite element. We assume that the Dof is defined as
48  * a linear combination of components of function value at a given point:
49  *
50  * Dof_value = a_1.f_1(x) + ... + a_n.f_n(x),
51  *
52  * where (a_1, ... , a_n) are given by @p coefs, x is the support point
53  * given by barycentric @p coords and (f_1, ..., f_n) is a generally
54  * vector-valued function. For the simplest Dof, i.e. value of a scalar
55  * function at point p, we set
56  *
57  * @p coords = p, @p coefs = { 1 }.
58  * The member @p dim denotes the affiliation of Dof to n-face:
59  * Nodal dofs have dim = 0,
60  * Dofs on lines: dim = 1,
61  * Dofs on triangles: dim = 2,
62  * Dofs in tetrahedron: dim = 3.
63  * It means that when a node, line or triangle is shared by adjacent cells,
64  * also the Dofs on this n-face are shared by the cells. Therefore
65  * for DG finite elements we set for all dofs the highest possible dimension.
66  *
67  * The class implements the method evaluate() which computes the Dof value
68  * for a basis function from given FunctionSpace.
69  */
70 class Dof {
71 public:
72 
73  Dof(unsigned int dim_,
74  unsigned int n_face_idx_,
75  arma::vec coords_,
76  arma::vec coefs_,
77  DofType type_)
78 
79  : dim(dim_),
80  n_face_idx(n_face_idx_),
81  coords(coords_),
82  coefs(coefs_),
83  type(type_)
84  {}
85 
86  /// Evaulate dof for basis function of given function space.
87  template<class FS>
88  double evaluate(const FS &function_space,
89  unsigned int basis_idx) const;
90 
91  /// Association to n-face of given dimension (point, line, triangle, tetrahedron.
92  unsigned int dim;
93 
94  /// Index of n-face to which the dof is associated.
95  unsigned int n_face_idx;
96 
97  /// Barycentric coordinates.
99 
100  /// Coefficients of linear combination of function value components.
102 
103  /**
104  * Currently we consider only type=Value, possibly we could have Gradient,
105  * CellIntegral, FaceIntegral or other types.
106  */
108 };
109 
110 
111 /**
112  * FunctionSpace is an abstract class that is used to describe finite elements.
113  * It is determined by the dimension of the field of definition (@p space_dim_),
114  * by the dimension of the range (@p n_components_) and by the values and
115  * gradients of a basis functions.
116  *
117  * Combining FunctionSpace with Dof(s), the FiniteElement class constructs the
118  * shape functions, i.e. basis of FunctionSpace for which the Dof(s) attain
119  * the value 0 or 1.
120  */
122 public:
123 
124  FunctionSpace(unsigned int space_dim,
125  unsigned int n_components)
126 
129  {}
130 
131 
132  /**
133  * @brief Value of the @p i th basis function at point @p point.
134  * @param basis_index Index of the basis function.
135  * @param point Point coordinates.
136  * @param comp_index Index of component (>0 for vector-valued functions).
137  */
138  virtual double basis_value(unsigned int basis_index,
139  const arma::vec &point,
140  unsigned int comp_index = 0
141  ) const = 0;
142 
143  /**
144  * @brief Gradient of the @p i th basis function at point @p point.
145  * @param basis_index Index of the basis function.
146  * @param point Point coordinates.
147  * @param comp_index Index of component (>0 for vector-valued functions).
148  */
149  virtual const arma::vec basis_grad(unsigned int basis_index,
150  const arma::vec &point,
151  unsigned int comp_index = 0
152  ) const = 0;
153 
154  /// Dimension of function space (number of basis functions).
155  virtual unsigned int dim() const = 0;
156 
157  /// Getter for space dimension.
158  unsigned int space_dim() const { return space_dim_; }
159 
160  /// Getter for number of components.
161  unsigned int n_components() const { return n_components_; }
162 
163  virtual ~FunctionSpace() {}
164 
165 protected:
166 
167  /// Space dimension of function arguments (i.e. 1, 2 or 3).
168  unsigned int space_dim_;
169 
170  /// Number of components of function values.
171  unsigned int n_components_;
172 };
173 
174 
175 /**
176  * Types of FiniteElement: scalar, vector-valued, tensor-valued or mixed system.
177  *
178  * The type also indicates how the shape functions and their gradients are transformed
179  * from reference element to arbitrary element. In particular:
180  *
181  * TYPE OBJECT EXPRESSION
182  * -----------------------------------------------------------
183  * FEScalar, FEVector, value ref_value
184  * FETensor grad J^{-T} * ref_grad
185  * -----------------------------------------------------------
186  * FEVectorContravariant value J * ref_value
187  * grad J^{-T} * ref_grad * J^T
188  * -----------------------------------------------------------
189  * FEVectorPiola value J * ref_value / |J|
190  * grad J^{-T} * ref_grad * J^T / |J|
191  * -----------------------------------------------------------
192  * FEMixedSystem value depends on sub-elements
193  * grad depends on sub-elements
194  *
195  * The transformation itself is done in FEValuesBase::fill_..._data() methods.
196  *
197  * Note that we use columnwise gradients, i.e. gradient of each component is a column vector.
198  */
199 enum FEType {
200  FEScalar = 0,
201  FEVector = 1,
204  FETensor = 4,
205  FEMixedSystem = 5
206 };
207 
208 
209 
210 /**
211  * @brief Abstract class for the description of a general finite element on
212  * a reference simplex in @p dim dimensions.
213  *
214  * The finite element is determined by a @p function_space_ and a set
215  * of @p dofs_. Further it must be set whether the finite element
216  * @p is_primitive_, which means that even if the functions in
217  * @p function_space_ are vector-valued, the basis functions have
218  * only one nonzero component (this is typical for tensor-product FE,
219  * e.g. vector-valued polynomial FE, but not for Raviart-Thomas FE).
220  *
221  * Description of dofs:
222  *
223  * The reference cell consists of lower dimensional entities (points,
224  * lines, triangles). Each dof is associated to one of these
225  * entities. If the entity is shared by 2 or more neighbouring cells
226  * in the mesh then this dof is shared by the finite elements on all
227  * of these cells. If a dof is associated to the cell itself then it
228  * is not shared with neighbouring cells.
229  *
230  *
231  * Shape functions:
232  *
233  * Sometimes it is convenient to describe the function space using
234  * a basis (called the raw basis) that is different from the set of
235  * shape functions for the finite element (the actual basis). For
236  * this reason we define the support points which play the role of
237  * nodal functionals associated to the particular dofs. To convert
238  * between the two bases one can use the @p node_matrix, which is
239  * constructed by the method compute_node_matrix(). In the case of
240  * non-Lagrangean finite elements the dofs are not associated to
241  * nodal functionals but e.g. to derivatives or averages. For that
242  * reason we distinguish between the unit support points which are
243  * uniquely associated to the dofs and the generalized support
244  * points that are auxiliary for the calculation of the dof
245  * functionals.
246  *
247  *
248  */
249 template<unsigned int dim>
251 public:
252 
253  /**
254  * @brief Constructor.
255  */
256  FiniteElement();
257 
258  /**
259  * @brief Returns the number of degrees of freedom needed by the finite
260  * element.
261  */
262  inline unsigned int n_dofs() const
263  { return dofs_.size(); }
264 
265  /**
266  * @brief Calculates the value of the @p comp-th component of
267  * the @p i-th shape function at the
268  * point @p p on the reference element.
269  *
270  * @param i Number of the shape function.
271  * @param p Point of evaluation.
272  * @param comp Number of vector component.
273  */
274  double shape_value(const unsigned int i,
275  const arma::vec::fixed<dim> &p,
276  const unsigned int comp = 0) const;
277 
278  /**
279  * @brief Calculates the @p comp-th component of the gradient
280  * of the @p i-th shape function at the point @p p on the
281  * reference element.
282  *
283  * @param i Number of the shape function.
284  * @param p Point of evaluation.
285  * @param comp Number of vector component.
286  */
287  arma::vec::fixed<dim> shape_grad(const unsigned int i,
288  const arma::vec::fixed<dim> &p,
289  const unsigned int comp = 0) const;
290 
291  /// Returns numer of components of the basis function.
292  inline unsigned int n_components() const
293  { return function_space_->n_components(); }
294 
295  /// Returns @p i -th degree of freedom.
296  inline const Dof &dof(unsigned int i) const
297  { return dofs_[i]; }
298 
299  /// Number of components of FE in a mapped space with dimension @p spacedim.
300  unsigned int n_space_components(unsigned int spacedim);
301 
302  /// Get barycentric coordinates of the points on the reference element associated with the dofs.
303  /// Used in BDDC for unknown reason.
305 
306  /**
307  * @brief Destructor.
308  */
309  virtual ~FiniteElement() {};
310 
311 protected:
312 
313  /**
314  * @brief Clears all internal structures.
315  */
316  void init(bool primitive = true,
317  FEType type = FEScalar);
318 
319  /**
320  * @brief Initialize vectors with information about components of basis functions.
321  */
322  void setup_components();
323 
324  /**
325  * @brief Decides which additional quantities have to be computed
326  * for each cell.
327  *
328  * @param flags Computed update flags.
329  */
330  virtual UpdateFlags update_each(UpdateFlags flags);
331 
332  /**
333  * @brief Initializes the @p node_matrix for computing the coefficients
334  * of the shape functions in the raw basis of @p functions_space_.
335  * This is done by evaluating the @p dofs_ for the basis function
336  * and by inverting the resulting matrix.
337  */
338  virtual void compute_node_matrix();
339 
340  /**
341  * @brief Indicates whether the basis functions have one or more
342  * nonzero components (scalar FE spaces are always primitive).
343  */
344  inline bool is_primitive() const
345  { return is_primitive_; }
346 
347  /**
348  * @brief Returns the component index for vector valued finite elements.
349  * @param sys_idx Index of shape function.
350  */
351  inline unsigned int system_to_component_index(unsigned sys_idx) const
352  { return component_indices_[sys_idx]; }
353 
354  /**
355  * @brief Returns the mask of nonzero components for given basis function.
356  * @param sys_idx Index of basis function.
357  */
358  inline const std::vector<bool> &get_nonzero_components(unsigned int sys_idx) const
359  { return nonzero_components_[sys_idx]; }
360 
361 
362 
363  /// Type of FiniteElement.
365 
366  /**
367  * @brief Primitive FE is using componentwise shape functions,
368  * i.e. only one component is nonzero for each shape function.
369  */
371 
372  /// Indices of nonzero components of shape functions (for primitive FE).
374 
375  /// Footprints of nonzero components of shape functions.
377 
378  /**
379  * @brief Matrix that determines the coefficients of the raw basis
380  * functions from the values at the support points.
381  */
383 
384  /// Function space defining the FE.
385  std::shared_ptr<FunctionSpace> function_space_;
386 
387  /// Set of degrees of freedom (functionals) defining the FE.
389 
390 
391  friend class FESystem<dim>;
392  friend class FEValues<3>;
393  friend class FE_P_disc<dim>;
394  friend class SubDOFHandlerMultiDim;
395 };
396 
397 
398 
399 
400 #endif /* FINITE_ELEMENT_HH_ */
double evaluate(const FS &function_space, unsigned int basis_idx) const
Evaulate dof for basis function of given function space.
arma::vec coefs
Coefficients of linear combination of function value components.
Dof(unsigned int dim_, unsigned int n_face_idx_, arma::vec coords_, arma::vec coefs_, DofType type_)
unsigned int n_face_idx
Index of n-face to which the dof is associated.
DofType type
arma::vec coords
Barycentric coordinates.
unsigned int dim
Association to n-face of given dimension (point, line, triangle, tetrahedron.
Compound finite element on dim dimensional simplex.
Definition: fe_system.hh:102
Calculates finite element data on the actual cell.
Definition: fe_values.hh:67
Discontinuous Lagrangean finite element on dim dimensional simplex.
Definition: fe_p.hh:108
Abstract class for the description of a general finite element on a reference simplex in dim dimensio...
unsigned int n_space_components(unsigned int spacedim)
Number of components of FE in a mapped space with dimension spacedim.
virtual UpdateFlags update_each(UpdateFlags flags)
Decides which additional quantities have to be computed for each cell.
std::shared_ptr< FunctionSpace > function_space_
Function space defining the FE.
const Dof & dof(unsigned int i) const
Returns i -th degree of freedom.
unsigned int n_dofs() const
Returns the number of degrees of freedom needed by the finite element.
bool is_primitive_
Primitive FE is using componentwise shape functions, i.e. only one component is nonzero for each shap...
arma::mat node_matrix
Matrix that determines the coefficients of the raw basis functions from the values at the support poi...
FEType type_
Type of FiniteElement.
const std::vector< bool > & get_nonzero_components(unsigned int sys_idx) const
Returns the mask of nonzero components for given basis function.
virtual std::vector< arma::vec::fixed< dim+1 > > dof_points() const
void init(bool primitive=true, FEType type=FEScalar)
Clears all internal structures.
unsigned int system_to_component_index(unsigned sys_idx) const
Returns the component index for vector valued finite elements.
double shape_value(const unsigned int i, const arma::vec::fixed< dim > &p, const unsigned int comp=0) const
Calculates the value of the comp-th component of the i-th shape function at the point p on the refere...
std::vector< std::vector< bool > > nonzero_components_
Footprints of nonzero components of shape functions.
virtual ~FiniteElement()
Destructor.
bool is_primitive() const
Indicates whether the basis functions have one or more nonzero components (scalar FE spaces are alway...
void setup_components()
Initialize vectors with information about components of basis functions.
virtual void compute_node_matrix()
Initializes the node_matrix for computing the coefficients of the shape functions in the raw basis of...
arma::vec::fixed< dim > shape_grad(const unsigned int i, const arma::vec::fixed< dim > &p, const unsigned int comp=0) const
Calculates the comp-th component of the gradient of the i-th shape function at the point p on the ref...
unsigned int n_components() const
Returns numer of components of the basis function.
std::vector< Dof > dofs_
Set of degrees of freedom (functionals) defining the FE.
std::vector< unsigned int > component_indices_
Indices of nonzero components of shape functions (for primitive FE).
FiniteElement()
Constructor.
unsigned int n_components() const
Getter for number of components.
unsigned int space_dim() const
Getter for space dimension.
virtual unsigned int dim() const =0
Dimension of function space (number of basis functions).
virtual ~FunctionSpace()
virtual double basis_value(unsigned int basis_index, const arma::vec &point, unsigned int comp_index=0) const =0
Value of the i th basis function at point point.
virtual const arma::vec basis_grad(unsigned int basis_index, const arma::vec &point, unsigned int comp_index=0) const =0
Gradient of the i th basis function at point point.
unsigned int n_components_
Number of components of function values.
unsigned int space_dim_
Space dimension of function arguments (i.e. 1, 2 or 3).
FunctionSpace(unsigned int space_dim, unsigned int n_components)
FEType
@ FEScalar
@ FEMixedSystem
@ FEVectorPiola
@ FETensor
@ FEVectorContravariant
@ FEVector
DofType
@ Value
ArmaMat< double, N, M > mat
Definition: armor.hh:888
ArmaVec< double, N > vec
Definition: armor.hh:885
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