GetFEM  5.4.2
Finite element description

Classes

class  getfem::virtual_fem
 Base class for finite element description. More...
 
class  getfem::fem< FUNC >
 virtual_fem implementation as a vector of generic functions. More...
 
class  getfem::fem_precomp_
 Pre-computations on a fem (given a fixed set of points on the reference convex, this object computes the value/gradient/hessian of all base functions on this set of points and stores them. More...
 
class  getfem::fem_precomp_pool
 handle a pool (i.e. More...
 
class  getfem::fem_interpolation_context
 structure passed as the argument of fem interpolation functions. More...
 

Typedefs

typedef std::shared_ptr< const getfem::virtual_femgetfem::pfem
 type of pointer on a fem description More...
 

Functions

virtual size_type getfem::virtual_fem::nb_dof (size_type) const
 Number of degrees of freedom. More...
 
virtual size_type getfem::virtual_fem::nb_base (size_type cv) const
 Number of basis functions.
 
size_type getfem::virtual_fem::nb_base_components (size_type cv) const
 Number of components (nb_dof() * dimension of the target space).
 
const std::vector< pdof_description > & getfem::virtual_fem::dof_types () const
 Get the array of pointer on dof description.
 
dim_type getfem::virtual_fem::dim () const
 dimension of the reference element.
 
dim_type getfem::virtual_fem::target_dim () const
 dimension of the target space.
 
vec_type getfem::virtual_fem::vectorial_type () const
 Type of vectorial element.
 
virtual bgeot::pconvex_ref getfem::virtual_fem::ref_convex (size_type) const
 Return the convex of the reference element.
 
bgeot::pconvex_structure getfem::virtual_fem::basic_structure (size_type cv) const
 Gives the convex of the reference element.
 
virtual const bgeot::convex< base_node > & getfem::virtual_fem::node_convex (size_type) const
 Gives the convex representing the nodes on the reference element.
 
bgeot::pconvex_structure getfem::virtual_fem::structure (size_type cv) const
 Gives the convex structure of the reference element nodes.
 
const base_node & getfem::virtual_fem::node_of_dof (size_type cv, size_type i) const
 Gives the node corresponding to the dof i. More...
 
bool getfem::virtual_fem::is_lagrange () const
 true if the base functions are such that $ \varphi_i(\textrm{node\_of\_dof(j)}) = \delta_{ij} $
 
bool getfem::virtual_fem::is_polynomial () const
 true if the base functions are polynomials
 
template<typename CVEC , typename VVEC >
void getfem::virtual_fem::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. More...
 
template<typename MAT >
void getfem::virtual_fem::interpolation (const fem_interpolation_context &c, MAT &M, dim_type Qdim) const
 Build the interpolation matrix for the interpolation at a fixed point x, given on the reference element. More...
 
template<typename CVEC , typename VMAT >
void getfem::virtual_fem::interpolation_grad (const fem_interpolation_context &c, const CVEC &coeff, VMAT &val, dim_type Qdim=1) const
 Interpolation of the gradient. More...
 
template<typename CVEC , typename VMAT >
void getfem::virtual_fem::interpolation_hess (const fem_interpolation_context &c, const CVEC &coeff, VMAT &val, dim_type Qdim) const
 Interpolation of the hessian. More...
 
template<typename CVEC >
void getfem::virtual_fem::interpolation_diverg (const fem_interpolation_context &c, const CVEC &coeff, typename gmm::linalg_traits< CVEC >::value_type &val) const
 Interpolation of the divergence. More...
 
virtual void getfem::virtual_fem::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. More...
 
virtual void getfem::virtual_fem::grad_base_value (const base_node &x, base_tensor &t) const =0
 Give the value of all gradients (on ref. More...
 
virtual void getfem::virtual_fem::hess_base_value (const base_node &x, base_tensor &t) const =0
 Give the value of all hessians (on ref. More...
 
virtual void getfem::virtual_fem::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_context. More...
 
virtual void getfem::virtual_fem::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_interpolation_context. More...
 
virtual void getfem::virtual_fem::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_interpolation_context. More...
 
void getfem::virtual_fem::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. More...
 
const std::vector< FUNC > & getfem::fem< FUNC >::base () const
 Gives the array of basic functions (components).
 
void getfem::fem< FUNC >::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)
 
void getfem::fem< FUNC >::grad_base_value (const base_node &x, base_tensor &t) const
 Evaluates at point x, the gradient of all base functions w.r.t. More...
 
void getfem::fem< FUNC >::hess_base_value (const base_node &x, base_tensor &t) const
 Evaluates at point x, the hessian of all base functions w.r.t. More...
 
pfem getfem::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 convex type. More...
 
pfem getfem::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 a given convex type. More...
 
pfem getfem::fem_descriptor (const std::string &name)
 get a fem descriptor from its string name.
 
std::string getfem::name_of_fem (pfem p)
 get the string name of a fem descriptor.
 
const base_tensor & getfem::fem_precomp_::val (size_type i) const
 returns values of the base functions
 
const base_tensor & getfem::fem_precomp_::grad (size_type i) const
 returns gradients of the base functions
 
const base_tensor & getfem::fem_precomp_::hess (size_type i) const
 returns hessians of the base functions
 
pfem_precomp getfem::fem_precomp (pfem pf, bgeot::pstored_point_tab pspt, dal::pstatic_stored_object dep)
 Handles precomputations for FEM. More...
 
void getfem::delete_fem_precomp (pfem_precomp pfp)
 Request for the removal of a pfem_precomp.
 
pfem_precomp getfem::fem_precomp_pool::operator() (pfem pf, bgeot::pstored_point_tab pspt)
 Request a pfem_precomp. More...
 
bool getfem::fem_interpolation_context::have_pfp () const
 true if a fem_precomp_ has been supplied.
 
bool getfem::fem_interpolation_context::have_pf () const
 true if the pfem is available.
 
const base_matrix & getfem::fem_interpolation_context::M () const
 non tau-equivalent transformation matrix.
 
void getfem::fem_interpolation_context::base_value (base_tensor &t, bool withM=true) const
 fill the tensor with the values of the base functions (taken at point this->xref())
 
void getfem::fem_interpolation_context::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())
 
void getfem::fem_interpolation_context::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())
 
const pfem getfem::fem_interpolation_context::pf () const
 get the current FEM descriptor
 
size_type getfem::fem_interpolation_context::convex_num () const
 get the current convex number
 
void getfem::fem_interpolation_context::set_face_num (short_type f)
 set the current face number
 
short_type getfem::fem_interpolation_context::face_num () const
 get the current face number
 
bool getfem::fem_interpolation_context::is_on_face () const
 On a face ?
 
pfem_precomp getfem::fem_interpolation_context::pfp () const
 get the current fem_precomp_
 
pfem getfem::interior_fem_of_hho_method (pfem hho_method)
 Specific function for a HHO method to obtain the method in the interior. More...
 

Variables

const typedef fem< bgeot::base_poly > * getfem::ppolyfem
 Classical polynomial FEM.
 
const typedef fem< bgeot::polynomial_composite > * getfem::ppolycompfem
 Polynomial composite FEM.
 
const typedef fem< bgeot::base_rational_fraction > * getfem::prationalfracfem
 Rational fration FEM.
 

Detailed Description

Typedef Documentation

◆ pfem

typedef std::shared_ptr<const getfem::virtual_fem> getfem::pfem

type of pointer on a fem description

See also
getfem::virtual_fem

Definition at line 244 of file getfem_fem.h.

Function Documentation

◆ nb_dof()

virtual size_type getfem::virtual_fem::nb_dof ( size_type  ) const
inlinevirtual

Number of degrees of freedom.

Parameters
cvthe convex number for this FEM. This information is rarely used, but is needed by some "special" FEMs, such as getfem::interpolated_fem.

Reimplemented in getfem::interpolated_fem, getfem::projected_fem, and getfem::fem_global_function.

Definition at line 296 of file getfem_fem.h.

◆ node_of_dof()

const base_node& getfem::virtual_fem::node_of_dof ( size_type  cv,
size_type  i 
) const
inline

Gives the node corresponding to the dof i.

Parameters
cvthe convex number for this FEM. This information is rarely used, by is needed by some "special" FEMs, such as getfem::interpolated_fem.
ithe local dof number (i < nb_dof(cv))

Definition at line 345 of file getfem_fem.h.

◆ interpolation() [1/2]

template<typename CVEC , typename VVEC >
void getfem::virtual_fem::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.

Parameters
cthe fem_interpolation_context, should have been suitably initialized for the point of evaluation.
coeffis the vector of coefficient relatively to the base functions, its length should be Qdim*this->nb_dof().
valcontains the interpolated value, on output (its size should be Qdim*this->target_dim()).
Qdimis the optional Q dimension, if the FEM is considered as a "vectorized" one.

Definition at line 858 of file getfem_fem.h.

◆ interpolation() [2/2]

template<typename MAT >
void getfem::virtual_fem::interpolation ( const fem_interpolation_context c,
MAT &  M,
dim_type  Qdim 
) const

Build the interpolation matrix for the interpolation at a fixed point x, given on the reference element.

The matrix M is filled, such that for a given coeff vector, the interpolation is given by M*coeff.

Definition at line 880 of file getfem_fem.h.

◆ interpolation_grad()

template<typename CVEC , typename VMAT >
void getfem::virtual_fem::interpolation_grad ( const fem_interpolation_context c,
const CVEC &  coeff,
VMAT &  val,
dim_type  Qdim = 1 
) const

Interpolation of the gradient.

The output is stored in the $ Q\times N$ matrix val.

Definition at line 902 of file getfem_fem.h.

◆ interpolation_hess()

template<typename CVEC , typename VMAT >
void getfem::virtual_fem::interpolation_hess ( const fem_interpolation_context c,
const CVEC &  coeff,
VMAT &  val,
dim_type  Qdim 
) const

Interpolation of the hessian.

The output is stored in the $ Q\times (N^2)$ matrix val.

Definition at line 929 of file getfem_fem.h.

◆ interpolation_diverg()

template<typename CVEC >
void getfem::virtual_fem::interpolation_diverg ( const fem_interpolation_context c,
const CVEC &  coeff,
typename gmm::linalg_traits< CVEC >::value_type &  val 
) const

Interpolation of the divergence.

The output is stored in the scalar val.

Definition at line 957 of file getfem_fem.h.

◆ base_value()

virtual void getfem::virtual_fem::base_value ( const base_node &  x,
base_tensor &  t 
) const
pure virtual

Give the value of all components of the base functions at the point x of the reference element.

Basic function used essentially by fem_precomp.

Implemented in getfem::fem< FUNC >, getfem::fem< base_poly >, getfem::fem< bgeot::polynomial_composite >, getfem::fem_level_set, getfem::interpolated_fem, getfem::projected_fem, getfem::fem_global_function, and getfem::torus_fem.

◆ grad_base_value() [1/2]

virtual void getfem::virtual_fem::grad_base_value ( const base_node &  x,
base_tensor &  t 
) const
pure virtual

Give the value of all gradients (on ref.

element) of the components of the base functions at the point x of the reference element. Basic function used essentially by fem_precomp.

Implemented in getfem::fem< FUNC >, getfem::fem< base_poly >, getfem::fem< bgeot::polynomial_composite >, getfem::fem_level_set, getfem::interpolated_fem, getfem::projected_fem, getfem::fem_global_function, and getfem::torus_fem.

◆ hess_base_value() [1/2]

virtual void getfem::virtual_fem::hess_base_value ( const base_node &  x,
base_tensor &  t 
) const
pure virtual

Give the value of all hessians (on ref.

element) of the components of the base functions at the point x of the reference element. Basic function used essentially by fem_precomp.

Implemented in getfem::fem< FUNC >, getfem::fem< base_poly >, getfem::fem< bgeot::polynomial_composite >, getfem::fem_level_set, getfem::interpolated_fem, getfem::projected_fem, getfem::fem_global_function, and getfem::torus_fem.

◆ real_base_value()

void getfem::virtual_fem::real_base_value ( const fem_interpolation_context c,
base_tensor &  t,
bool  withM = true 
) const
virtual

Give the value of all components of the base functions at the current point of the fem_interpolation_context.

Used by elementary computations. if withM is false the matrix M for non tau-equivalent elements is not taken into account.

Reimplemented in getfem::interpolated_fem, getfem::projected_fem, getfem::fem_level_set, getfem::fem_global_function, and getfem::torus_fem.

Definition at line 309 of file getfem_fem.cc.

◆ real_grad_base_value()

void getfem::virtual_fem::real_grad_base_value ( const fem_interpolation_context c,
base_tensor &  t,
bool  withM = true 
) const
virtual

Give the gradient of all components of the base functions at the current point of the fem_interpolation_context.

Used by elementary computations. if withM is false the matrix M for non tau-equivalent elements is not taken into account.

Reimplemented in getfem::interpolated_fem, getfem::projected_fem, getfem::fem_level_set, getfem::fem_global_function, and getfem::torus_fem.

Definition at line 313 of file getfem_fem.cc.

◆ real_hess_base_value()

void getfem::virtual_fem::real_hess_base_value ( const fem_interpolation_context c,
base_tensor &  t,
bool  withM = true 
) const
virtual

Give the hessian of all components of the base functions at the current point of the fem_interpolation_context.

Used by elementary computations. if withM is false the matrix M for non tau-equivalent elements is not taken into account.

Reimplemented in getfem::fem_level_set, getfem::interpolated_fem, getfem::projected_fem, getfem::fem_global_function, and getfem::torus_fem.

Definition at line 317 of file getfem_fem.cc.

◆ add_node()

void getfem::virtual_fem::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.

Important : the faces should be the faces on which the corresponding base function is non zero.

Definition at line 648 of file getfem_fem.cc.

◆ grad_base_value() [2/2]

template<class FUNC >
void getfem::fem< FUNC >::grad_base_value ( const base_node &  x,
base_tensor &  t 
) const
inlinevirtual

Evaluates at point x, the gradient of all base functions w.r.t.

the reference element directions 0,..,dim-1 and returns the result in t(nb_base,target_dim,dim)

Implements getfem::virtual_fem.

Definition at line 566 of file getfem_fem.h.

◆ hess_base_value() [2/2]

template<class FUNC >
void getfem::fem< FUNC >::hess_base_value ( const base_node &  x,
base_tensor &  t 
) const
inlinevirtual

Evaluates at point x, the hessian of all base functions w.r.t.

the reference element directions 0,..,dim-1 and returns the result in t(nb_base,target_dim,dim,dim)

Implements getfem::virtual_fem.

Definition at line 581 of file getfem_fem.h.

◆ classical_fem()

pfem getfem::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 convex type.

Parameters
pgtthe geometric transformation (which defines the convex type).
kthe degree of the fem.
completea flag which requests complete Langrange polynomial elements even if the provided pgt is an incomplete one (e.g. 8-node quadrilateral or 20-node hexahedral).
Returns
a ppolyfem.

Definition at line 4141 of file getfem_fem.cc.

◆ classical_discontinuous_fem()

pfem getfem::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 a given convex type.

Parameters
pgtthe geometric transformation (which defines the convex type).
kthe degree of the fem.
alphathe "inset" factor for the dof nodes: with alpha = 0, the nodes are located as usual (i.e. with node on the convex border), and for 0 < alpha < 1, they converge to the center of gravity of the convex.
completea flag which requests complete Langrange polynomial elements even if the provided pgt is an incomplete one (e.g. 8-node quadrilateral or 20-node hexahedral).
Returns
a ppolyfem.

Definition at line 4146 of file getfem_fem.cc.

◆ fem_precomp()

pfem_precomp getfem::fem_precomp ( pfem  pf,
bgeot::pstored_point_tab  pspt,
dal::pstatic_stored_object  dep 
)

Handles precomputations for FEM.

statically allocates a fem-precomputation object, and returns a pointer to it. The fem_precomp_ objects are "cached", i.e. they are stored in a global pool and if this function is called two times with the same arguments, a pointer to the same object will be returned.

Parameters
pfa pointer to the fem object.
pspta pointer to a list of points in the reference convex.CAUTION: this array must not be destroyed as long as the fem_precomp is used!!.

Moreover pspt is supposed to identify uniquely the set of points. This means that you should NOT alter its content at any time after using this function.

If you need a set of "temporary" getfem::fem_precomp_, create them via a getfem::fem_precomp_pool structure. All memory will be freed when this structure will be destroyed.

Definition at line 4333 of file getfem_fem.cc.

◆ operator()()

pfem_precomp getfem::fem_precomp_pool::operator() ( pfem  pf,
bgeot::pstored_point_tab  pspt 
)
inline

Request a pfem_precomp.

If not already in the pool, the pfem_precomp is computed, and added to the pool.

Parameters
pfa pointer to the fem object.
pspta pointer to a list of points in the reference convex.

CAUTION: this array must not be destroyed as long as the fem_precomp is used!!

Moreover pspt is supposed to identify uniquely the set of points. This means that you should NOT alter its content until the fem_precomp_pool is destroyed.

Definition at line 735 of file getfem_fem.h.

◆ interior_fem_of_hho_method()

pfem getfem::interior_fem_of_hho_method ( pfem  hho_method)

Specific function for a HHO method to obtain the method in the interior.

If the method is not of composite type, return the argument.

Definition at line 866 of file getfem_fem_composite.cc.