GetFEM  5.4.2
getfem_fem.h
Go to the documentation of this file.
1 /* -*- c++ -*- (enables emacs c++ mode) */
2 /*===========================================================================
3 
4  Copyright (C) 1999-2020 Yves Renard
5 
6  This file is a part of GetFEM
7 
8  GetFEM is free software; you can redistribute it and/or modify it
9  under the terms of the GNU Lesser General Public License as published
10  by the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version along with the GCC Runtime Library
12  Exception either version 3.1 or (at your option) any later version.
13  This program is distributed in the hope that it will be useful, but
14  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  License and GCC Runtime Library Exception for more details.
17  You should have received a copy of the GNU Lesser General Public License
18  along with this program; if not, write to the Free Software Foundation,
19  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
20 
21  As a special exception, you may use this file as it is a part of a free
22  software library without restriction. Specifically, if other files
23  instantiate templates or use macros or inline functions from this file,
24  or you compile this file and link it with other files to produce an
25  executable, this file does not by itself cause the resulting executable
26  to be covered by the GNU Lesser General Public License. This exception
27  does not however invalidate any other reasons why the executable file
28  might be covered by the GNU Lesser General Public License.
29 
30 ===========================================================================*/
31 
32 /**@file getfem_fem.h
33  @author Yves Renard <Yves.Renard@insa-lyon.fr>
34  @date December 21, 1999.
35  @brief Definition of the finite element methods.
36 
37  This file defines the getfem::virtual_fem class, which is the common
38  base class of all FEM.
39 
40  @section fem_list List of FEM known by getfem::fem_descriptor :
41 
42  - "FEM_PK(N,K)" : classical Lagrange element PK on a simplex.
43 
44  - "FEM_PK_DISCONTINUOUS(N,K,alpha)" : discontinuous Lagrange
45  element PK on a simplex.
46 
47  - "FEM_QK(N,K)" : classical Lagrange element QK on a parellepiped.
48 
49  - "FEM_QK_DISCONTINUOUS(N,K,alpha)" : discontinuous Lagrange
50  element QK on a parallelepiped.
51 
52  - "FEM_Q2_INCOMPLETE(N)" : incomplete Q2 elements with 8 and 20 dof
53  (serendipity Quad 8 and Hexa 20 elements)
54 
55  - "FEM_Q2_INCOMPLETE_DISCONTINUOUS(N)" : discontinuous incomplete Q2
56  elements with 8 and 20 dof (serendipity Quad 8 and Hexa 20 elements)
57 
58  - "FEM_PRISM_PK(N,K)" : classical Lagrange element PK on a prism.
59 
60  - "FEM_PRISM_PK_DISCONTINUOUS(N,K,alpha)" : classical discontinuous
61  Lagrange element PK on a prism.
62 
63  - "FEM_PRISM_INCOMPLETE_P2" : Incomplete Lagrange element on a
64  quadratic 3D prism (serendipity, 15-node wedge element). Can be connected
65  toa standard P2 Lagrange on its triangular faces and a Q2_INCOMPLETE
66  Lagrange element on its quadrangular faces.
67 
68  - "FEM_PK_WITH_CUBIC_BUBBLE(N,K)" : classical Lagrange element PK
69  on a simplex with an additional volumic bubble function.
70 
71  - "FEM_PRODUCT(FEM1,FEM2)" : tensorial product of two polynomial
72  elements
73 
74  - "FEM_P1_NONCONFORMING" : Nonconforming P1 method on a
75  triangle.
76 
77  - "FEM_P1_BUBBLE_FACE(N)" : P1 method on a simplex with an
78  additional bubble function on face 0.
79 
80  - "FEM_P1_BUBBLE_FACE_LAG" : P1 method on a simplex with an
81  additional lagrange dof on face 0.
82 
83  - "FEM_HERMITE(N)" : Hermite element P3 on dimension N (1, 2 por 3).
84 
85  - "FEM_ARGYRIS" : Argyris element on the triangle.
86 
87  - "FEM_HCT_TRIANGLE" : Hsieh-Clough-Tocher element on the triangle
88  (composite P3 element which is C^1, 12 dof).
89 
90  - "FEM_REDUCED_HCT_TRIANGLE" : Hsieh-Clough-Tocher element on the triangle
91  (composite P3 element which is C^1, 9 dof).
92 
93  - "FEM_QUADC1_COMPOSITE" : quadrilateral element, composite P3 element
94  and C^1 (16 dof).
95 
96  - "FEM_REDUCED_QUADC1_COMPOSITE" : quadrilateral element, composite
97  P3 element and C^1 (12 dof).
98 
99  - "FEM_PK_HIERARCHICAL(N,K)" : PK element with a hierarchical basis.
100 
101  - "FEM_QK_HIERARCHICAL(N,K)" : QK element with a hierarchical basis.
102 
103  - "FEM_PRISM_PK_HIERARCHICAL(N,K)" : PK element on a prism with a
104  hierarchical basis.
105 
106  - "FEM_STRUCTURED_COMPOSITE(FEM, K)" : Composite fem on a grid with
107  K divisions.
108 
109  - "FEM_PK_HIERARCHICAL_COMPOSITE(N,K,S)" : PK composite element on
110  a grid with S subdivisions and with a hierarchical basis.
111 
112  - "FEM_PK_FULL_HIERARCHICAL_COMPOSITE(N,K,S)" : PK composite
113  element with S subdivisions and a hierarchical basis on both degree
114  and subdivision.
115 
116  - "FEM_PYRAMID_QK(K)" : Lagrange element on a 3D pyramid of degree
117  K=0, 1 or 2. Can be connected to a standard P1/P2 Lagrange element on its
118  triangular faces and a standard Q1/Q2 Lagrange element on its quadrangular
119  face.
120 
121  - "FEM_PYRAMID_QK_DISCONTINUOUS(K)" : Discontinuous Lagrange element
122  on a 3D pyramid of degree K = 0, 1 or 2.
123 
124  - "FEM_PYRAMID_Q2_INCOMPLETE" : Incomplete Lagrange element on a
125  quadratic 3D pyramid (serendipity, 13-node element). Can be connected to
126  a standard P2 Lagrange element on its triangular faces and a Q2_INCOMPLETE
127  Lagrange element on its quadrangular face.
128 
129  - "HHO(fem_interior, fem_face_1, ..., fem_face_n)" : Build a hybrid method
130  with "fem_interior" on the element itself and "fem_face_1", ...,
131  "fem_face_n" on each face. If only one method is given for the faces, it
132  is duplicated on each face.
133 
134 */
135 
136 #ifndef GETFEM_FEM_H__
137 #define GETFEM_FEM_H__
138 
140 #include "bgeot_geometric_trans.h"
141 #include "bgeot_poly_composite.h"
142 #include "getfem_integration.h"
143 #include "dal_naming_system.h"
144 #include <deque>
145 
146 namespace getfem {
147 
148  /** @defgroup dofdescr Dof description
149  * This set of functions gives a pointer to dof descriptions and
150  * tests the compatibility between two dof descriptions.
151  * A dof description describes a type of dof (Lagrange type,
152  * Hermite type ...) in order to be able to say wether two dofs are
153  * to be identified or not. The construction of dof type with
154  * the tensorial product is taken into
155  * account, making the dof_description structure a little bit complex
156  * (this structure will probably evoluate in the future).
157  * @{
158  */
159 
160  struct dof_description;
161 
162  /// Type representing a pointer on a dof_description
163  typedef dof_description *pdof_description;
164 
165  /** @brief Description of a unique dof of lagrange type (value at the node).
166  * @param d the dimension of the reference element (2 for triangles, 3 for tetrahedrons ...)
167  */
168  pdof_description lagrange_dof(dim_type d);
169 
170  /** Description of a unique dof of derivative type.
171  * (a derivative at the node).
172  * @param d the dimension of the reference element.
173  * @param r corresponds to the variable number for which the derivative is taken (0 <= r < d)
174  */
175  pdof_description derivative_dof(dim_type d, dim_type r);
176 
177  /** Description of a unique dof of second derivative type.
178  * @param d the dimension of the reference element.
179  * @param num_der1 corresponds to the variable number for which the first derivative is taken (0 <= r < d)
180  * @param num_der2 corresponds to the variable number for which the second derivative is taken (0 <= r < d)
181  */
183  dim_type num_der1,
184  dim_type num_der2);
185 
186  /** Description of a unique dof of normal derivative type
187  * (normal derivative at the node, regarding a face).
188  * @param d the dimension of the reference element.
189  */
191 
192  /** Description of a unique dof of mean value type.
193  * @param d the dimension of the reference element.
194  */
195  pdof_description mean_value_dof(dim_type d);
196 
197  pdof_description bubble1_dof(dim_type ct);
198 
199  /** Description of a global dof, i.e. a numbered dof having a global scope.
200  * @param d the dimension of the reference element.
201  */
202  pdof_description global_dof(dim_type d);
203 
204  /// Product description of the descriptions *pnd1 and *pnd2.
206 
207  pdof_description to_coord_dof(pdof_description p, dim_type ct);
208 
209  /// Description of a special dof for Xfem
211 
212  /// Returns the xfem_index of dof (0 for normal dof)
214 
215  size_type reserve_xfem_index();
216  dim_type coord_index_of_dof(pdof_description a);
217 
218  int dof_weak_compatibility(pdof_description a, pdof_description b);
219 
220  /** Gives a total order on the dof description compatible with the
221  * identification.
222  */
224 
225  /// Says if the dof is linkable.
227 
228  /// Says if the two dofs can be identified.
230 
231 
232  /* @} */
233 
234  /* ******************************************************************** */
235  /* Classes for description of a finite element. */
236  /* ******************************************************************** */
237 
238  /** @defgroup pfem Finite element description
239  *
240  * @{
241  */
242 
243 
244  class virtual_fem;
245 
246  /** type of pointer on a fem description @see getfem::virtual_fem */
247  typedef std::shared_ptr<const getfem::virtual_fem> pfem;
248 
249  class fem_precomp_;
250  typedef std::shared_ptr<const getfem::fem_precomp_> pfem_precomp;
251 
253 
254  /** @brief Base class for finite element description */
255  class virtual_fem : virtual public dal::static_stored_object,
256  public std::enable_shared_from_this<const virtual_fem> {
257  public :
258  enum vec_type { VECTORIAL_NOTRANSFORM_TYPE, VECTORIAL_PRIMAL_TYPE,
259  VECTORIAL_DUAL_TYPE };
260 
261  protected :
262 
263  mutable std::vector<pdof_description> dof_types_;
264  /* this convex structure is "owned" by the virtual_fem
265  (a deep copy is made when virtual_fems are copied)
266  But cvs_node has to be a pointer, as bgeot::convex_structure
267  inherits from dal::static_stored_object
268  */
269  std::shared_ptr<bgeot::convex_structure> cvs_node;
270  std::vector<std::vector<short_type>> face_tab; // face list for each dof
271  bgeot::convex<base_node> cv_node;
272  mutable bgeot::pstored_point_tab pspt;
273  mutable bool pspt_valid;
274  bgeot::pconvex_ref cvr; // reference element.
275  dim_type ntarget_dim; // dimension of the target space
276  mutable dim_type dim_; // dimension of the reference element
277  bool is_equiv; // true if the FEM is equivalent
278  bool is_lag; // true if the FEM is of Lagrange type
279  bool is_pol; // true if the FEM is polynomial
280  bool is_polycomp; // true if the FEM is polynomial composite
281  bool real_element_defined;
282  bool is_standard_fem; // is_equiv && !real_element_defined && scalar
283  short_type es_degree; // estimated polynomial degree of the FEM
284  short_type hier_raff; // hierarchical refinement of the FEM
285  vec_type vtype; // for vectorial elements, type of transformation
286  // from the reference element.
287  std::string debug_name_;
288 
289  public :
290  /** Number of degrees of freedom.
291 
292  @param cv the convex number for this FEM. This information is
293  rarely used, but is needed by some "special" FEMs, such as
294  getfem::interpolated_fem.
295  */
296  virtual size_type nb_dof(size_type /*cv*/) const
297  { return dof_types_.size(); }
298  /// Number of basis functions.
299  virtual size_type nb_base(size_type cv) const
300  { return nb_dof(cv); }
301  /// Number of components (nb_dof() * dimension of the target space).
303  { return nb_base(cv) * ntarget_dim; }
304  size_type nb_components(size_type cv) const
305  { return nb_dof(cv) * ntarget_dim; }
306  /// Get the array of pointer on dof description.
307  const std::vector<pdof_description> &dof_types() const
308  { return dof_types_; }
309  short_type hierarchical_raff() const { return hier_raff; }
310  /// dimension of the reference element.
311  dim_type dim() const { return dim_; }
312  dim_type &dim() { return dim_; }
313  /// dimension of the target space.
314  dim_type target_dim() const { return ntarget_dim; }
315  /// Type of vectorial element.
316  vec_type vectorial_type() const { return vtype; }
317  /// Return the convex of the reference element.
318  virtual bgeot::pconvex_ref ref_convex(size_type) const { return cvr; }
319  /// @internal
320  bgeot::pconvex_ref &mref_convex() { return cvr; }
321  /// Gives the convex of the reference element.
323  { return ref_convex(cv)->structure(); }
324  /// Gives the convex representing the nodes on the reference element.
326  { return cv_node; }
327  /// Gives the convex structure of the reference element nodes.
329  { return node_convex(cv).structure(); }
330  const std::string &debug_name() const { return debug_name_; }
331  std::string &debug_name() { return debug_name_; }
332  virtual bgeot::pstored_point_tab node_tab(size_type) const {
333  if (!pspt_valid) {
334  pspt = bgeot::store_point_tab(cv_node.points());
335  pspt_valid = true;
336  }
337  return pspt;
338  }
339  /** Gives the node corresponding to the dof i.
340  @param cv the convex number for this FEM. This information is
341  rarely used, by is needed by some "special" FEMs, such as
342  getfem::interpolated_fem.
343  @param i the local dof number (<tt>i < nb_dof(cv)</tt>)
344  */
345  const base_node &node_of_dof(size_type cv, size_type i) const
346  { return (*(node_tab(cv)))[i];}
347  virtual const std::vector<short_type> &
348  faces_of_dof(size_type /*cv*/, size_type i) const;
349  bool is_on_real_element() const { return real_element_defined; }
350  bool is_equivalent() const { return is_equiv; }
351  bool need_G() const
352  { return !(is_equivalent()) || real_element_defined; }
353  /// true if the base functions are such that @f$ \varphi_i(\textrm{node\_of\_dof(j)}) = \delta_{ij} @f$
354  bool is_lagrange() const { return is_lag; }
355  /// true if the base functions are polynomials
356  bool is_polynomial() const { return is_pol; }
357  bool is_polynomialcomp() const { return is_polycomp; }
358  bool is_standard() const { return is_standard_fem; }
359  bool &is_polynomialcomp() { return is_polycomp; }
360  bool &is_equivalent() { return is_equiv; }
361  bool &is_lagrange() { return is_lag; }
362  bool &is_polynomial() { return is_pol; }
363  bool &is_standard() { return is_standard_fem; }
364  short_type estimated_degree() const { return es_degree; }
365  short_type &estimated_degree() { return es_degree; }
366 
367  virtual void mat_trans(base_matrix &, const base_matrix &,
369  { GMM_ASSERT1(false, "This function should not be called."); }
370  /** Interpolate at an arbitrary point x given on the reference
371  element.
372 
373  @param c the fem_interpolation_context, should have been
374  suitably initialized for the point of evaluation.
375 
376  @param coeff is the vector of coefficient relatively to
377  the base functions, its length should be @c Qdim*this->nb_dof().
378 
379  @param val contains the interpolated value, on output (its
380  size should be @c Qdim*this->target_dim()).
381 
382  @param Qdim is the optional Q dimension, if the FEM is
383  considered as a "vectorized" one.
384  */
385  template<typename CVEC, typename VVEC>
386  void interpolation(const fem_interpolation_context& c,
387  const CVEC& coeff, VVEC &val, dim_type Qdim) const;
388 
389  /** Build the interpolation matrix for the interpolation at a
390  fixed point x, given on the reference element.
391 
392  The matrix @c M is filled, such that for a given @c coeff
393  vector, the interpolation is given by @c M*coeff.
394  */
395  template <typename MAT>
396  void interpolation(const fem_interpolation_context& c,
397  MAT &M, dim_type Qdim) const;
398 
399  /** Interpolation of the gradient. The output is stored in the @f$
400  Q\times N@f$ matrix @c val.
401  */
402  template<typename CVEC, typename VMAT>
403  void interpolation_grad(const fem_interpolation_context& c,
404  const CVEC& coeff, VMAT &val,
405  dim_type Qdim=1) const;
406 
407  /** Interpolation of the hessian. The output is stored in the @f$
408  Q\times (N^2)@f$ matrix @c val.
409  */
410  template<typename CVEC, typename VMAT>
411  void interpolation_hess(const fem_interpolation_context& c,
412  const CVEC& coeff, VMAT &val,
413  dim_type Qdim) const;
414 
415  /** Interpolation of the divergence. The output is stored in the
416  scalar @c val.
417  */
418  template<typename CVEC>
420  (const fem_interpolation_context& c, const CVEC& coeff,
421  typename gmm::linalg_traits<CVEC>::value_type &val) const;
422 
423  /** Give the value of all components of the base functions at the
424  * point x of the reference element. Basic function used essentially
425  * by fem_precomp.
426  */
427  virtual void base_value(const base_node &x, base_tensor &t) const = 0;
428 
429  /** Give the value of all gradients (on ref. element) of the components
430  * of the base functions at the point x of the reference element.
431  * Basic function used essentially by fem_precomp.
432  */
433  virtual void grad_base_value(const base_node &x, base_tensor &t) const = 0;
434 
435  /** Give the value of all hessians (on ref. element) of the components
436  * of the base functions at the point x of the reference element.
437  * Basic function used essentially by fem_precomp.
438  */
439  virtual void hess_base_value(const base_node &x, base_tensor &t) const = 0;
440 
441  /** Give the value of all components of the base functions at the
442  current point of the fem_interpolation_context. Used by
443  elementary computations. if withM is false the matrix M for
444  non tau-equivalent elements is not taken into account.
445  */
446  virtual void real_base_value(const fem_interpolation_context &c,
447  base_tensor &t, bool withM = true) const;
448 
449  /** Give the gradient of all components of the base functions at the
450  current point of the fem_interpolation_context. Used by
451  elementary computations. if withM is false the matrix M for
452  non tau-equivalent elements is not taken into account.
453  */
454  virtual void real_grad_base_value(const fem_interpolation_context &c,
455  base_tensor &t, bool withM = true) const;
456 
457  /** Give the hessian of all components of the base functions at the
458  current point of the fem_interpolation_context. Used by
459  elementary computations. if withM is false the matrix M for
460  non tau-equivalent elements is not taken into account.
461  */
462  virtual void real_hess_base_value(const fem_interpolation_context &c,
463  base_tensor &t, bool withM = true) const;
464 
465  virtual size_type index_of_global_dof(size_type, size_type) const
466  { GMM_ASSERT1(false, "internal error."); }
467 
468  /** internal function adding a node to an element for the creation
469  * of a finite element method. Important : the faces should be the faces
470  * on which the corresponding base function is non zero.
471  */
472  void add_node(const pdof_description &d, const base_node &pt,
473  const dal::bit_vector &faces);
474  void add_node(const pdof_description &d, const base_node &pt);
475  void init_cvs_node();
476  void unfreeze_cvs_node();
477 
478  virtual_fem &operator =(const virtual_fem &f) {
479  copy(f); return *this;
480  }
481 
482  virtual_fem() {
483  DAL_STORED_OBJECT_DEBUG_CREATED(this, "Fem");
484  ntarget_dim = 1; dim_ = 1;
485  is_equiv = is_pol = is_polycomp = is_lag = is_standard_fem = false;
486  pspt_valid = false; hier_raff = 0; real_element_defined = false;
487  es_degree = 5;
488  vtype = VECTORIAL_NOTRANSFORM_TYPE;
489  cvs_node = bgeot::new_convex_structure();
490  }
491  virtual_fem(const virtual_fem& f)
492  : dal::static_stored_object(),
493  std::enable_shared_from_this<const virtual_fem>()
494  { copy(f); DAL_STORED_OBJECT_DEBUG_CREATED(this, "Fem"); }
495  virtual ~virtual_fem() { DAL_STORED_OBJECT_DEBUG_DESTROYED(this, "Fem"); }
496  private:
497  void copy(const virtual_fem &f);
498  };
499 
500  /**
501  virtual_fem implementation as a vector of generic functions. The
502  class FUNC should provide "derivative" and "eval" member
503  functions (this is the case for bgeot::polynomial<T>).
504  */
505  template <class FUNC> class fem : public virtual_fem {
506  protected :
507  std::vector<FUNC> base_;
508  mutable std::vector<std::vector<FUNC>> grad_, hess_;
509  mutable bool grad_computed_ = false;
510  mutable bool hess_computed_ = false;
511 
512  void compute_grad_() const {
513  if (grad_computed_) return;
514  GLOBAL_OMP_GUARD
515  if (grad_computed_) return;
517  dim_type n = dim();
518  grad_.resize(R);
519  for (size_type i = 0; i < R; ++i) {
520  grad_[i].resize(n);
521  for (dim_type j = 0; j < n; ++j) {
522  grad_[i][j] = base_[i]; grad_[i][j].derivative(j);
523  }
524  }
525  grad_computed_ = true;
526  }
527 
528  void compute_hess_() const {
529  if (hess_computed_) return;
530  GLOBAL_OMP_GUARD
531  if (hess_computed_) return;
533  dim_type n = dim();
534  hess_.resize(R);
535  for (size_type i = 0; i < R; ++i) {
536  hess_[i].resize(n*n);
537  for (dim_type j = 0; j < n; ++j) {
538  for (dim_type k = 0; k < n; ++k) {
539  hess_[i][j+k*n] = base_[i];
540  hess_[i][j+k*n].derivative(j); hess_[i][j+k*n].derivative(k);
541  }
542  }
543  }
544  hess_computed_ = true;
545  }
546 
547  public :
548 
549  /// Gives the array of basic functions (components).
550  const std::vector<FUNC> &base() const { return base_; }
551  std::vector<FUNC> &base() { return base_; }
552  /** Evaluates at point x, all base functions and returns the result in
553  t(nb_base,target_dim) */
554  void base_value(const base_node &x, base_tensor &t) const {
555  bgeot::multi_index mi(2);
556  mi[1] = target_dim(); mi[0] = short_type(nb_base(0));
557  t.adjust_sizes(mi);
559  base_tensor::iterator it = t.begin();
560  for (size_type i = 0; i < R; ++i, ++it)
561  *it = bgeot::to_scalar(base_[i].eval(x.begin()));
562  }
563  /** Evaluates at point x, the gradient of all base functions w.r.t. the
564  reference element directions 0,..,dim-1 and returns the result in
565  t(nb_base,target_dim,dim) */
566  void grad_base_value(const base_node &x, base_tensor &t) const {
567  if (!grad_computed_) compute_grad_();
568  bgeot::multi_index mi(3);
569  dim_type n = dim();
570  mi[2] = n; mi[1] = target_dim(); mi[0] = short_type(nb_base(0));
571  t.adjust_sizes(mi);
573  base_tensor::iterator it = t.begin();
574  for (dim_type j = 0; j < n; ++j)
575  for (size_type i = 0; i < R; ++i, ++it)
576  *it = bgeot::to_scalar(grad_[i][j].eval(x.begin()));
577  }
578  /** Evaluates at point x, the hessian of all base functions w.r.t. the
579  reference element directions 0,..,dim-1 and returns the result in
580  t(nb_base,target_dim,dim,dim) */
581  void hess_base_value(const base_node &x, base_tensor &t) const {
582  if (!hess_computed_) compute_hess_();
583  bgeot::multi_index mi(4);
584  dim_type n = dim();
585  mi[3] = n; mi[2] = n; mi[1] = target_dim();
586  mi[0] = short_type(nb_base(0));
587  t.adjust_sizes(mi);
589  base_tensor::iterator it = t.begin();
590  for (dim_type k = 0; k < n; ++k)
591  for (dim_type j = 0; j < n; ++j)
592  for (size_type i = 0; i < R; ++i, ++it)
593  *it = bgeot::to_scalar(hess_[i][j+k*n].eval(x.begin()));
594  }
595 
596  };
597 
598  /** Classical polynomial FEM. */
600  /** Polynomial composite FEM */
602  /** Rational fration FEM */
604 
605  /** Give a pointer on the structures describing the classical
606  polynomial fem of degree k on a given convex type.
607 
608  @param pgt the geometric transformation (which defines the convex type).
609  @param k the degree of the fem.
610  @param complete a flag which requests complete Langrange polynomial
611  elements even if the provided pgt is an incomplete one (e.g. 8-node
612  quadrilateral or 20-node hexahedral).
613  @return a ppolyfem.
614  */
616  bool complete=false);
617 
618  /** Give a pointer on the structures describing the classical
619  polynomial discontinuous fem of degree k on a given convex type.
620 
621  @param pgt the geometric transformation (which defines the convex type).
622 
623  @param k the degree of the fem.
624 
625  @param alpha the "inset" factor for the dof nodes: with alpha =
626  0, the nodes are located as usual (i.e. with node on the convex border),
627  and for 0 < alpha < 1, they converge to the center of gravity of the
628  convex.
629 
630  @param complete a flag which requests complete Langrange polynomial
631  elements even if the provided pgt is an incomplete one (e.g. 8-node
632  quadrilateral or 20-node hexahedral).
633 
634  @return a ppolyfem.
635  */
637  scalar_type alpha=0, bool complete=false);
638 
639  /** get a fem descriptor from its string name. */
640  pfem fem_descriptor(const std::string &name);
641 
642  /** get the string name of a fem descriptor. */
643  std::string name_of_fem(pfem p);
644 
645  pfem PK_fem(size_type n, short_type k);
646  pfem QK_fem(size_type n, short_type k);
647  pfem PK_prism_fem(size_type n, short_type k);
648 
649  /**
650  Pre-computations on a fem (given a fixed set of points on the
651  reference convex, this object computes the value/gradient/hessian
652  of all base functions on this set of points and stores them.
653  */
654  class fem_precomp_ : virtual public dal::static_stored_object {
655  protected:
656  const pfem pf;
657  const bgeot::pstored_point_tab pspt;
658  mutable std::vector<base_tensor> c; // stored values of base functions
659  mutable std::vector<base_tensor> pc; // stored gradients of base functions
660  mutable std::vector<base_tensor> hpc; // stored hessians of base functions
661  public:
662  /// returns values of the base functions
663  inline const base_tensor &val(size_type i) const
664  { if (c.empty()) init_val(); return c[i]; }
665  /// returns gradients of the base functions
666  inline const base_tensor &grad(size_type i) const
667  { if (pc.empty()) init_grad(); return pc[i]; }
668  /// returns hessians of the base functions
669  inline const base_tensor &hess(size_type i) const
670  { if (hpc.empty()) init_hess(); return hpc[i]; }
671  inline pfem get_pfem() const { return pf; }
672  // inline const bgeot::stored_point_tab& get_point_tab() const
673  // { return *pspt; }
674  inline bgeot::pstored_point_tab get_ppoint_tab() const
675  { return pspt; }
676  fem_precomp_(const pfem, const bgeot::pstored_point_tab);
677  ~fem_precomp_() { DAL_STORED_OBJECT_DEBUG_DESTROYED(this, "Fem_precomp"); }
678  private:
679  void init_val() const;
680  void init_grad() const;
681  void init_hess() const;
682  };
683 
684 
685  /** @brief Handles precomputations for FEM. statically allocates a
686  fem-precomputation object, and returns a pointer to it. The
687  fem_precomp_ objects are "cached", i.e. they are stored in a
688  global pool and if this function is called two times with the
689  same arguments, a pointer to the same object will be returned.
690 
691  @param pf a pointer to the fem object. @param pspt a pointer to
692  a list of points in the reference convex.CAUTION: this array must
693  not be destroyed as long as the fem_precomp is used!!.
694 
695  Moreover pspt is supposed to identify uniquely the set of
696  points. This means that you should NOT alter its content at any
697  time after using this function.
698 
699  If you need a set of "temporary" getfem::fem_precomp_, create
700  them via a getfem::fem_precomp_pool structure. All memory will be
701  freed when this structure will be destroyed. */
702  pfem_precomp fem_precomp(pfem pf, bgeot::pstored_point_tab pspt,
703  dal::pstatic_stored_object dep);
704 
705  /** Request for the removal of a pfem_precomp */
706  inline void delete_fem_precomp(pfem_precomp pfp)
707  { dal::del_stored_object(pfp); }
708 
709 
710  /**
711  handle a pool (i.e. a set) of fem_precomp. The difference with
712  the global fem_precomp function is that these fem_precomp objects
713  are freed when the fem_precomp_pool is destroyed (they can eat
714  much memory). An example of use can be found in the
715  getfem::interpolation_solution functions of getfem_export.h
716  */
718  std::set<pfem_precomp> precomps;
719 
720  public :
721 
722  /** Request a pfem_precomp. If not already in the pool, the
723  pfem_precomp is computed, and added to the pool.
724 
725  @param pf a pointer to the fem object.
726  @param pspt a pointer to a list of points in the reference convex.
727 
728  CAUTION:
729  this array must not be destroyed as long as the fem_precomp is used!!
730 
731  Moreover pspt is supposed to identify uniquely the set of
732  points. This means that you should NOT alter its content until
733  the fem_precomp_pool is destroyed.
734  */
735  pfem_precomp operator()(pfem pf, bgeot::pstored_point_tab pspt) {
736  pfem_precomp p = fem_precomp(pf, pspt, 0);
737  precomps.insert(p);
738  return p;
739  }
740  void clear();
741  ~fem_precomp_pool() { clear(); }
742  };
743 
744 
745  /** structure passed as the argument of fem interpolation
746  functions. This structure can be partially filled (for example
747  the xreal will be computed if needed as long as pgp+ii is known).
748  */
751 
752  mutable base_matrix M_; // optional transformation matrix (for non tau-equivalent fems)
753  pfem pf_; // current fem
754  pfem_precomp pfp_; // optional fem_precomp_ (speed up the computations)
755  size_type convex_num_; // The element (convex) number
756  short_type face_num_; // Face number for boundary integration
757  int xfem_side_; // For the computation of a jump with fem_level_set only
758  public:
759  /// true if a fem_precomp_ has been supplied.
760  bool have_pfp() const { return pfp_ != 0; }
761  /// true if the pfem is available.
762  bool have_pf() const { return pf_ != 0; }
763  /// non tau-equivalent transformation matrix.
764  const base_matrix& M() const;
765  /** fill the tensor with the values of the base functions (taken
766  at point @c this->xref())
767  */
768  void base_value(base_tensor& t, bool withM = true) const;
769  // Optimized function for high level generic assembly
770  void pfp_base_value(base_tensor& t, const pfem_precomp &pfp__);
771  /** fill the tensor with the gradient of the base functions (taken
772  at point @c this->xref())
773  */
774  void grad_base_value(base_tensor& t, bool withM = true) const;
775  // Optimized function for high level generic assembly
776  void pfp_grad_base_value(base_tensor& t, const pfem_precomp &pfp__);
777  /** fill the tensor with the hessian of the base functions (taken
778  at point @c this->xref())
779  */
780  void hess_base_value(base_tensor& t, bool withM = true) const;
781  /** get the current FEM descriptor */
782  const pfem pf() const { return pf_; }
783  /** get the current convex number */
784  size_type convex_num() const;
785  bool is_convex_num_valid() const;
786  void invalid_convex_num() { convex_num_ = size_type(-1); }
787  /** set the current face number */
788  void set_face_num(short_type f) { face_num_ = f; }
789  /** get the current face number */
790  short_type face_num() const;
791  /** On a face ? */
792  bool is_on_face() const;
793  /** get the current fem_precomp_ */
794  pfem_precomp pfp() const { return pfp_; }
795  void set_pfp(pfem_precomp newpfp);
796  void set_pf(pfem newpf);
797  int xfem_side() const { return xfem_side_; }
798  void set_xfem_side(int side) { xfem_side_ = side; }
799  void change(bgeot::pgeotrans_precomp pgp__,
800  pfem_precomp pfp__, size_type ii__,
801  const base_matrix& G__, size_type convex_num__,
802  short_type face_num__ = short_type(-1)) {
803  bgeot::geotrans_interpolation_context::change(pgp__,ii__,G__);
804  convex_num_ = convex_num__; face_num_ = face_num__; xfem_side_ = 0;
805  set_pfp(pfp__);
806  }
807  void change(bgeot::pgeometric_trans pgt__,
808  pfem_precomp pfp__, size_type ii__,
809  const base_matrix& G__, size_type convex_num__,
810  short_type face_num__ = short_type(-1)) {
811  bgeot::geotrans_interpolation_context::change
812  (pgt__, pfp__->get_ppoint_tab(), ii__, G__);
813  convex_num_ = convex_num__; face_num_ = face_num__; xfem_side_ = 0;
814  set_pfp(pfp__);
815  }
816  void change(bgeot::pgeometric_trans pgt__,
817  pfem pf__, const base_node& xref__, const base_matrix& G__,
818  size_type convex_num__, short_type face_num__=short_type(-1)) {
819  bgeot::geotrans_interpolation_context::change(pgt__,xref__,G__);
820  pf_ = pf__; pfp_ = 0; convex_num_ = convex_num__; face_num_ = face_num__;
821  xfem_side_ = 0;
822  }
823  fem_interpolation_context()
824  : bgeot::geotrans_interpolation_context(),
825  convex_num_(size_type(-1)), face_num_(short_type(-1)), xfem_side_(0) {}
826  fem_interpolation_context(bgeot::pgeotrans_precomp pgp__,
827  pfem_precomp pfp__, size_type ii__,
828  const base_matrix& G__,
829  size_type convex_num__,
830  short_type face_num__ = short_type(-1))
831  : bgeot::geotrans_interpolation_context(pgp__,ii__,G__),
832  convex_num_(convex_num__), face_num_(face_num__), xfem_side_(0)
833  { set_pfp(pfp__); }
834  fem_interpolation_context(bgeot::pgeometric_trans pgt__,
835  pfem_precomp pfp__, size_type ii__,
836  const base_matrix& G__,
837  size_type convex_num__,
838  short_type face_num__ = short_type(-1))
839  : bgeot::geotrans_interpolation_context(pgt__,pfp__->get_ppoint_tab(),
840  ii__, G__),
841  convex_num_(convex_num__), face_num_(face_num__), xfem_side_(0)
842  { set_pfp(pfp__); }
843  fem_interpolation_context(bgeot::pgeometric_trans pgt__,
844  pfem pf__,
845  const base_node& xref__,
846  const base_matrix& G__,
847  size_type convex_num__,
848  short_type face_num__ = short_type(-1))
849  : bgeot::geotrans_interpolation_context(pgt__,xref__,G__),
850  pf_(pf__), pfp_(0), convex_num_(convex_num__), face_num_(face_num__),
851  xfem_side_(0) {}
852  };
853 
854  // IN : coeff(Qmult,nb_dof)
855  // OUT: val(Qdim), Qdim=target_dim*Qmult
856  // AUX: Z(nb_dof,target_dim)
857  template <typename CVEC, typename VVEC>
859  const CVEC& coeff, VVEC &val,
860  dim_type Qdim) const {
861  size_type Qmult = size_type(Qdim) / target_dim();
862  size_type nbdof = nb_dof(c.convex_num());
863  GMM_ASSERT1(gmm::vect_size(val) == Qdim, "dimensions mismatch");
864  GMM_ASSERT1(gmm::vect_size(coeff) == nbdof*Qmult,
865  "Wrong size for coeff vector");
866 
867  gmm::clear(val);
868  base_tensor Z; real_base_value(c, Z);
869 
870  for (size_type j = 0; j < nbdof; ++j) {
871  for (size_type q = 0; q < Qmult; ++q) {
872  typename gmm::linalg_traits<CVEC>::value_type co = coeff[j*Qmult+q];
873  for (size_type r = 0; r < target_dim(); ++r)
874  val[r + q*target_dim()] += co * Z[j + r*nbdof];
875  }
876  }
877  }
878 
879  template <typename MAT>
881  MAT &M, dim_type Qdim) const {
882  size_type Qmult = size_type(Qdim) / target_dim();
883  size_type nbdof = nb_dof(c.convex_num());
884  GMM_ASSERT1(gmm::mat_nrows(M) == Qdim && gmm::mat_ncols(M) == nbdof*Qmult,
885  "dimensions mismatch");
886 
887  gmm::clear(M);
888  base_tensor Z; real_base_value(c, Z);
889  for (size_type j = 0; j < nbdof; ++j) {
890  for (size_type q = 0; q < Qmult; ++q) {
891  for (size_type r = 0; r < target_dim(); ++r)
892  M(r+q*target_dim(), j*Qmult+q) = Z[j + r*nbdof];
893  }
894  }
895  }
896 
897 
898  // IN : coeff(Qmult,nb_dof)
899  // OUT: val(Qdim,N), Qdim=target_dim*Qmult
900  // AUX: t(nb_dof,target_dim,N)
901  template<typename CVEC, typename VMAT>
903  const CVEC& coeff, VMAT &val,
904  dim_type Qdim) const {
905  size_type N = c.N();
906  size_type nbdof = nb_dof(c.convex_num());
907  size_type Qmult = gmm::vect_size(coeff) / nbdof;
908  GMM_ASSERT1(gmm::mat_ncols(val) == N &&
909  gmm::mat_nrows(val) == target_dim()*Qmult &&
910  gmm::vect_size(coeff) == nbdof*Qmult,
911  "dimensions mismatch");
912  GMM_ASSERT1(Qdim == target_dim()*Qmult, // Qdim seems to be superfluous input, could be removed in the future
913  "dimensions mismatch");
914  base_tensor t;
915  real_grad_base_value(c, t); // t(nbdof,target_dim,N)
916 
917  gmm::clear(val);
918  for (size_type q = 0; q < Qmult; ++q) {
919  base_tensor::const_iterator it = t.begin();
920  for (size_type k = 0; k < N; ++k)
921  for (size_type r = 0; r < target_dim(); ++r)
922  for (size_type j = 0; j < nbdof; ++j, ++it)
923  val(r + q*target_dim(), k) += coeff[j*Qmult+q] * (*it);
924  }
925  }
926 
927 
928  template<typename CVEC, typename VMAT>
930  const CVEC& coeff, VMAT &val,
931  dim_type Qdim) const {
932  size_type Qmult = size_type(Qdim) / target_dim();
933  size_type N = c.N();
934  GMM_ASSERT1(gmm::mat_ncols(val) == N*N &&
935  gmm::mat_nrows(val) == Qdim, "dimensions mismatch");
936 
937  base_tensor t;
938  size_type nbdof = nb_dof(c.convex_num());
939 
940  gmm::clear(val);
941  real_hess_base_value(c, t);
942  for (size_type q = 0; q < Qmult; ++q) {
943  base_tensor::const_iterator it = t.begin();
944  for (size_type k = 0; k < N*N; ++k)
945  for (size_type r = 0; r < target_dim(); ++r)
946  for (size_type j = 0; j < nbdof; ++j, ++it)
947  val(r + q*target_dim(), k) += coeff[j*Qmult+q] * (*it);
948  }
949  }
950 
951 
952  // IN : coeff(Qmult,nb_dof)
953  // OUT: val
954  // AUX: t(nb_dof,target_dim,N), Qmult*target_dim == N
955  template<typename CVEC>
957  (const fem_interpolation_context& c, const CVEC& coeff,
958  typename gmm::linalg_traits<CVEC>::value_type &val) const {
959  size_type N = c.N();
960  size_type nbdof = nb_dof(c.convex_num());
961  size_type Qmult = gmm::vect_size(coeff) / nbdof;
962  GMM_ASSERT1(gmm::vect_size(coeff) == nbdof*Qmult , "dimensions mismatch");
963  GMM_ASSERT1(target_dim()*Qmult == N &&
964  (Qmult == 1 || target_dim() == 1),
965  "Dimensions mismatch. Divergence operator requires fem qdim equal to dim.");
966  base_tensor t;
967  real_grad_base_value(c, t); // t(nbdof,target_dim,N)
968  // for Qmult == 1 this is sub-optimal since it evaluates all (:,i,j)
969  // gradients instead of only the diagonal ones(:,i,i)
970 
971  val = scalar_type(0);
972  base_tensor::const_iterator it = t.begin();
973  if (Qmult == 1)
974  for (size_type k = 0; k < N; ++k) {
975  if (k) it += (N*nbdof + 1);
976  for (size_type j = 0; j < nbdof; ++j) {
977  if (j) ++it;
978  val += coeff[j] * (*it);
979  }
980  }
981  else // if (target_dim() == 1)
982  for (size_type k = 0; k < N; ++k) {
983  if (k) ++it;
984  for (size_type j = 0; j < nbdof; ++j) {
985  if (j) ++it;
986  val += coeff[j*N+k] * (*it);
987  }
988  }
989  }
990 
991  /**
992  Specific function for a HHO method to obtain the method in the interior.
993  If the method is not of composite type, return the argument.
994  */
996 
997 
998  /* Functions allowing the add of a finite element method outside
999  of getfem_fem.cc */
1000 
1001  typedef dal::naming_system<virtual_fem>::param_list fem_param_list;
1002 
1003  void inline read_poly(bgeot::base_poly &p, int d, const char *s)
1004  { p = bgeot::read_base_poly(short_type(d), s); }
1005 
1006  void add_fem_name(std::string name,
1008 
1009 
1010  /* @} */
1011 
1012 } /* end of namespace getfem. */
1013 
1014 
1015 #endif
getfem::virtual_fem::base_value
virtual void base_value(const base_node &x, base_tensor &t) const =0
Give the value of all components of the base functions at the point x of the reference element.
getfem::virtual_fem::interpolation_grad
void interpolation_grad(const fem_interpolation_context &c, const CVEC &coeff, VMAT &val, dim_type Qdim=1) const
Interpolation of the gradient.
Definition: getfem_fem.h:902
getfem::fem_interpolation_context::M
const base_matrix & M() const
non tau-equivalent transformation matrix.
Definition: getfem_fem.cc:43
getfem::virtual_fem::dim
dim_type dim() const
dimension of the reference element.
Definition: getfem_fem.h:311
getfem::interior_fem_of_hho_method
pfem interior_fem_of_hho_method(pfem hho_method)
Specific function for a HHO method to obtain the method in the interior.
Definition: getfem_fem_composite.cc:866
getfem::ppolycompfem
const typedef fem< bgeot::polynomial_composite > * ppolycompfem
Polynomial composite FEM.
Definition: getfem_fem.h:601
getfem::fem_precomp_pool
handle a pool (i.e.
Definition: getfem_fem.h:717
getfem::fem_precomp_::val
const base_tensor & val(size_type i) const
returns values of the base functions
Definition: getfem_fem.h:663
getfem::virtual_fem::nb_dof
virtual size_type nb_dof(size_type) const
Number of degrees of freedom.
Definition: getfem_fem.h:296
getfem::pdof_description
dof_description * pdof_description
Type representing a pointer on a dof_description.
Definition: getfem_fem.h:160
getfem::virtual_fem::dof_types
const std::vector< pdof_description > & dof_types() const
Get the array of pointer on dof description.
Definition: getfem_fem.h:307
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
getfem::virtual_fem::is_lagrange
bool is_lagrange() const
true if the base functions are such that
Definition: getfem_fem.h:354
getfem::fem::hess_base_value
void hess_base_value(const base_node &x, base_tensor &t) const
Evaluates at point x, the hessian of all base functions w.r.t.
Definition: getfem_fem.h:581
getfem::fem_interpolation_context::is_on_face
bool is_on_face() const
On a face ?
Definition: getfem_fem.cc:67
getfem::virtual_fem::ref_convex
virtual bgeot::pconvex_ref ref_convex(size_type) const
Return the convex of the reference element.
Definition: getfem_fem.h:318
getfem::virtual_fem::node_convex
virtual const bgeot::convex< base_node > & node_convex(size_type) const
Gives the convex representing the nodes on the reference element.
Definition: getfem_fem.h:325
getfem::virtual_fem::real_hess_base_value
virtual void real_hess_base_value(const fem_interpolation_context &c, base_tensor &t, bool withM=true) const
Give the hessian of all components of the base functions at the current point of the fem_interpolatio...
Definition: getfem_fem.cc:317
getfem::fem_precomp
pfem_precomp fem_precomp(pfem pf, bgeot::pstored_point_tab pspt, dal::pstatic_stored_object dep)
Handles precomputations for FEM.
Definition: getfem_fem.cc:4333
getfem::fem_interpolation_context::hess_base_value
void hess_base_value(base_tensor &t, bool withM=true) const
fill the tensor with the hessian of the base functions (taken at point this->xref())
Definition: getfem_fem.cc:255
dal
Dynamic Array Library.
Definition: dal_backtrace.cc:29
getfem_integration.h
Integration methods (exact and approximated) on convexes.
dal::naming_system
Associate a name to a method descriptor and store method descriptors.
Definition: dal_naming_system.h:56
getfem::virtual_fem::interpolation
void interpolation(const fem_interpolation_context &c, const CVEC &coeff, VVEC &val, dim_type Qdim) const
Interpolate at an arbitrary point x given on the reference element.
Definition: getfem_fem.h:858
getfem::virtual_fem::interpolation_hess
void interpolation_hess(const fem_interpolation_context &c, const CVEC &coeff, VMAT &val, dim_type Qdim) const
Interpolation of the hessian.
Definition: getfem_fem.h:929
bgeot::geotrans_interpolation_context
the geotrans_interpolation_context structure is passed as the argument of geometric transformation in...
Definition: bgeot_geometric_trans.h:411
bgeot::convex< base_node >
dal::del_stored_object
void del_stored_object(const pstatic_stored_object &o, bool ignore_unstored)
Delete an object and the object which depend on it.
Definition: dal_static_stored_objects.cc:369
getfem::classical_fem
pfem classical_fem(bgeot::pgeometric_trans pgt, short_type k, bool complete=false)
Give a pointer on the structures describing the classical polynomial fem of degree k on a given conve...
Definition: getfem_fem.cc:4141
getfem::fem_precomp_::grad
const base_tensor & grad(size_type i) const
returns gradients of the base functions
Definition: getfem_fem.h:666
getfem::virtual_fem
Base class for finite element description.
Definition: getfem_fem.h:255
getfem::fem_interpolation_context::convex_num
size_type convex_num() const
get the current convex number
Definition: getfem_fem.cc:52
dal_naming_system.h
Naming system.
bgeot::read_base_poly
base_poly read_base_poly(short_type n, std::istream &f)
read a base_poly on the stream ist.
Definition: bgeot_poly.cc:211
getfem::ppolyfem
const typedef fem< bgeot::base_poly > * ppolyfem
Classical polynomial FEM.
Definition: getfem_fem.h:599
getfem::fem_precomp_
Pre-computations on a fem (given a fixed set of points on the reference convex, this object computes ...
Definition: getfem_fem.h:654
bgeot::short_type
gmm::uint16_type short_type
used as the common short type integer in the library
Definition: bgeot_config.h:72
getfem::virtual_fem::structure
bgeot::pconvex_structure structure(size_type cv) const
Gives the convex structure of the reference element nodes.
Definition: getfem_fem.h:328
getfem::dof_compatibility
bool dof_compatibility(pdof_description, pdof_description)
Says if the two dofs can be identified.
Definition: getfem_fem.cc:618
getfem
GEneric Tool for Finite Element Methods.
Definition: getfem_accumulated_distro.h:46
getfem::product_dof
pdof_description product_dof(pdof_description pnd1, pdof_description pnd2)
Product description of the descriptions *pnd1 and *pnd2.
Definition: getfem_fem.cc:552
getfem::virtual_fem::basic_structure
bgeot::pconvex_structure basic_structure(size_type cv) const
Gives the convex of the reference element.
Definition: getfem_fem.h:322
getfem::fem_interpolation_context
structure passed as the argument of fem interpolation functions.
Definition: getfem_fem.h:749
getfem::virtual_fem::add_node
void add_node(const pdof_description &d, const base_node &pt, const dal::bit_vector &faces)
internal function adding a node to an element for the creation of a finite element method.
Definition: getfem_fem.cc:648
getfem::fem_interpolation_context::face_num
short_type face_num() const
get the current face number
Definition: getfem_fem.cc:61
getfem::fem
virtual_fem implementation as a vector of generic functions.
Definition: getfem_fem.h:505
getfem::second_derivative_dof
pdof_description second_derivative_dof(dim_type d, dim_type num_der1, dim_type num_der2)
Description of a unique dof of second derivative type.
Definition: getfem_fem.cc:496
getfem::virtual_fem::target_dim
dim_type target_dim() const
dimension of the target space.
Definition: getfem_fem.h:314
bgeot::polynomial
This class deals with plain polynomials with several variables.
Definition: bgeot_poly.h:180
getfem::virtual_fem::vectorial_type
vec_type vectorial_type() const
Type of vectorial element.
Definition: getfem_fem.h:316
getfem::derivative_dof
pdof_description derivative_dof(dim_type d, dim_type r)
Description of a unique dof of derivative type.
Definition: getfem_fem.cc:487
getfem::virtual_fem::hess_base_value
virtual void hess_base_value(const base_node &x, base_tensor &t) const =0
Give the value of all hessians (on ref.
getfem::delete_fem_precomp
void delete_fem_precomp(pfem_precomp pfp)
Request for the removal of a pfem_precomp.
Definition: getfem_fem.h:706
getfem::virtual_fem::nb_base
virtual size_type nb_base(size_type cv) const
Number of basis functions.
Definition: getfem_fem.h:299
dal_static_stored_objects.h
Stores interdependent getfem objects.
getfem::dof_linkable
bool dof_linkable(pdof_description)
Says if the dof is linkable.
Definition: getfem_fem.cc:615
getfem::pfem
std::shared_ptr< const getfem::virtual_fem > pfem
type of pointer on a fem description
Definition: getfem_fem.h:244
getfem::fem_interpolation_context::have_pf
bool have_pf() const
true if the pfem is available.
Definition: getfem_fem.h:762
getfem::virtual_fem::is_polynomial
bool is_polynomial() const
true if the base functions are polynomials
Definition: getfem_fem.h:356
getfem::mean_value_dof
pdof_description mean_value_dof(dim_type d)
Description of a unique dof of mean value type.
Definition: getfem_fem.cc:534
bgeot_geometric_trans.h
Geometric transformations on convexes.
getfem::virtual_fem::node_of_dof
const base_node & node_of_dof(size_type cv, size_type i) const
Gives the node corresponding to the dof i.
Definition: getfem_fem.h:345
getfem::normal_derivative_dof
pdof_description normal_derivative_dof(dim_type d)
Description of a unique dof of normal derivative type (normal derivative at the node,...
Definition: getfem_fem.cc:507
getfem::fem_precomp_pool::operator()
pfem_precomp operator()(pfem pf, bgeot::pstored_point_tab pspt)
Request a pfem_precomp.
Definition: getfem_fem.h:735
getfem::virtual_fem::nb_base_components
size_type nb_base_components(size_type cv) const
Number of components (nb_dof() * dimension of the target space).
Definition: getfem_fem.h:302
getfem::fem_interpolation_context::have_pfp
bool have_pfp() const
true if a fem_precomp_ has been supplied.
Definition: getfem_fem.h:760
getfem::lagrange_dof
pdof_description lagrange_dof(dim_type d)
Description of a unique dof of lagrange type (value at the node).
Definition: getfem_fem.cc:390
getfem::prationalfracfem
const typedef fem< bgeot::base_rational_fraction > * prationalfracfem
Rational fration FEM.
Definition: getfem_fem.h:603
getfem::virtual_fem::grad_base_value
virtual void grad_base_value(const base_node &x, base_tensor &t) const =0
Give the value of all gradients (on ref.
bgeot
Basic Geometric Tools.
Definition: bgeot_convex_ref.cc:27
bgeot::pconvex_structure
std::shared_ptr< const convex_structure > pconvex_structure
Pointer on a convex structure description.
Definition: bgeot_convex_structure.h:54
getfem::fem_interpolation_context::pfp
pfem_precomp pfp() const
get the current fem_precomp_
Definition: getfem_fem.h:794
dal::static_stored_object
base class for static stored objects
Definition: dal_static_stored_objects.h:206
getfem::fem::grad_base_value
void grad_base_value(const base_node &x, base_tensor &t) const
Evaluates at point x, the gradient of all base functions w.r.t.
Definition: getfem_fem.h:566
getfem::dof_xfem_index
size_type dof_xfem_index(pdof_description)
Returns the xfem_index of dof (0 for normal dof)
Definition: getfem_fem.cc:621
bgeot::pgeometric_trans
std::shared_ptr< const bgeot::geometric_trans > pgeometric_trans
pointer type for a geometric transformation
Definition: bgeot_geometric_trans.h:186
getfem::name_of_fem
std::string name_of_fem(pfem p)
get the string name of a fem descriptor.
Definition: getfem_fem.cc:4241
getfem::fem_interpolation_context::set_face_num
void set_face_num(short_type f)
set the current face number
Definition: getfem_fem.h:788
getfem::fem_interpolation_context::pf
const pfem pf() const
get the current FEM descriptor
Definition: getfem_fem.h:782
getfem::fem::base_value
void base_value(const base_node &x, base_tensor &t) const
Evaluates at point x, all base functions and returns the result in t(nb_base,target_dim)
Definition: getfem_fem.h:554
getfem::fem_interpolation_context::grad_base_value
void grad_base_value(base_tensor &t, bool withM=true) const
fill the tensor with the gradient of the base functions (taken at point this->xref())
Definition: getfem_fem.cc:202
getfem::xfem_dof
pdof_description xfem_dof(pdof_description p, size_type ind)
Description of a special dof for Xfem.
Definition: getfem_fem.cc:433
getfem::virtual_fem::interpolation_diverg
void interpolation_diverg(const fem_interpolation_context &c, const CVEC &coeff, typename gmm::linalg_traits< CVEC >::value_type &val) const
Interpolation of the divergence.
Definition: getfem_fem.h:957
getfem::fem::base
const std::vector< FUNC > & base() const
Gives the array of basic functions (components).
Definition: getfem_fem.h:550
getfem::fem_interpolation_context::base_value
void base_value(base_tensor &t, bool withM=true) const
fill the tensor with the values of the base functions (taken at point this->xref())
Definition: getfem_fem.cc:120
bgeot_poly_composite.h
Handle composite polynomials.
getfem::classical_discontinuous_fem
pfem classical_discontinuous_fem(bgeot::pgeometric_trans pg, short_type k, scalar_type alpha=0, bool complete=false)
Give a pointer on the structures describing the classical polynomial discontinuous fem of degree k on...
Definition: getfem_fem.cc:4146
getfem::virtual_fem::real_grad_base_value
virtual void real_grad_base_value(const fem_interpolation_context &c, base_tensor &t, bool withM=true) const
Give the gradient of all components of the base functions at the current point of the fem_interpolati...
Definition: getfem_fem.cc:313
getfem::fem_precomp_::hess
const base_tensor & hess(size_type i) const
returns hessians of the base functions
Definition: getfem_fem.h:669
getfem::fem_descriptor
pfem fem_descriptor(const std::string &name)
get a fem descriptor from its string name.
Definition: getfem_fem.cc:4232
getfem::dof_description_compare
int dof_description_compare(pdof_description a, pdof_description b)
Gives a total order on the dof description compatible with the identification.
Definition: getfem_fem.cc:606
getfem::global_dof
pdof_description global_dof(dim_type d)
Description of a global dof, i.e.
Definition: getfem_fem.cc:542
getfem::virtual_fem::real_base_value
virtual void real_base_value(const fem_interpolation_context &c, base_tensor &t, bool withM=true) const
Give the value of all components of the base functions at the current point of the fem_interpolation_...
Definition: getfem_fem.cc:309