Fem

class Fem(*args)

GetFEM Fem object

This object represents a finite element method on a reference element.

General constructor for Fem objects

  • F = Fem('interpolated_fem', MeshFem mf_source, MeshIm mim_target, [ivec blocked_dofs[, bool caching]]) Build a special Fem which is interpolated from another MeshFem.

    Using this special finite element, it is possible to interpolate a given MeshFem mf_source on another mesh, given the integration method mim_target that will be used on this mesh.

    Note that this finite element may be quite slow or consume much memory depending if caching is used or not. By default caching is True

  • F = Fem('projected_fem', MeshFem mf_source, MeshIm mim_target, int rg_source, int rg_target[, ivec blocked_dofs[, bool caching]]) Build a special Fem which is interpolated from another MeshFem.

    Using this special finite element, it is possible to interpolate a given MeshFem mf_source on another mesh, given the integration method mim_target that will be used on this mesh.

    Note that this finite element may be quite slow or consume much memory depending if caching is used or not. By default caching is True

  • F = Fem(string fem_name) The fem_name should contain a description of the finite element method. Please refer to the GetFEM manual (especially the description of finite element and integration methods) for a complete reference. Here is a list of some of them:

    • FEM_PK(n,k) : classical Lagrange element Pk on a simplex of dimension n.

    • FEM_PK_DISCONTINUOUS(n,k[,alpha]) : discontinuous Lagrange element Pk on a simplex of dimension n.

    • FEM_QK(n,k) : classical Lagrange element Qk on quadrangles, hexahedrons etc.

    • FEM_QK_DISCONTINUOUS(n,k[,alpha]) : discontinuous Lagrange element Qk on quadrangles, hexahedrons etc.

    • FEM_Q2_INCOMPLETE(n) : incomplete Q2 elements with 8 and 20 dof (serendipity Quad 8 and Hexa 20 elements).

    • FEM_PK_PRISM(n,k) : classical Lagrange element Pk on a prism of dimension n.

    • FEM_PK_PRISM_DISCONTINUOUS(n,k[,alpha]) : classical discontinuous Lagrange element Pk on a prism.

    • FEM_PK_WITH_CUBIC_BUBBLE(n,k) : classical Lagrange element Pk on a simplex with an additional volumic bubble function.

    • FEM_P1_NONCONFORMING : non-conforming P1 method on a triangle.

    • FEM_P1_BUBBLE_FACE(n) : P1 method on a simplex with an additional bubble function on face 0.

    • FEM_P1_BUBBLE_FACE_LAG : P1 method on a simplex with an additional lagrange dof on face 0.

    • FEM_PK_HIERARCHICAL(n,k) : PK element with a hierarchical basis.

    • FEM_QK_HIERARCHICAL(n,k) : QK element with a hierarchical basis.

    • FEM_PK_PRISM_HIERARCHICAL(n,k) : PK element on a prism with a hierarchical basis.

    • FEM_STRUCTURED_COMPOSITE(Fem f,k) : Composite Fem f on a grid with k divisions.

    • FEM_PK_HIERARCHICAL_COMPOSITE(n,k,s) : Pk composite element on a grid with s subdivisions and with a hierarchical basis.

    • FEM_PK_FULL_HIERARCHICAL_COMPOSITE(n,k,s) : Pk composite element with s subdivisions and a hierarchical basis on both degree and subdivision.

    • FEM_PRODUCT(A,B) : tensorial product of two polynomial elements.

    • FEM_HERMITE(n) : Hermite element P3 on a simplex of dimension n = 1, 2, 3.

    • FEM_ARGYRIS : Argyris element P5 on the triangle.

    • FEM_HCT_TRIANGLE : Hsieh-Clough-Tocher element on the triangle (composite P3 element which is C1), should be used with IM_HCT_COMPOSITE() integration method.

    • FEM_QUADC1_COMPOSITE : Quadrilateral element, composite P3 element and C1 (16 dof).

    • FEM_REDUCED_QUADC1_COMPOSITE : Quadrilateral element, composite P3 element and C1 (12 dof).

    • FEM_RT0(n) : Raviart-Thomas element of order 0 on a simplex of dimension n.

    • FEM_NEDELEC(n) : Nedelec edge element of order 0 on a simplex of dimension n.

    Of course, you have to ensure that the selected fem is compatible with the geometric transformation: a Pk fem has no meaning on a quadrangle.

base_value(p)

Evaluate all basis functions of the FEM at point p.

p is supposed to be in the reference convex!

char()

Ouput a (unique) string representation of the Fem.

This can be used to perform comparisons between two different Fem objects.

dim()

Return the dimension (dimension of the reference convex) of the Fem.

display()

displays a short summary for a Fem object.

estimated_degree()

Return an estimation of the polynomial degree of the Fem.

This is an estimation for fem which are not polynomials.

grad_base_value(p)

Evaluate the gradient of all base functions of the Fem at point p.

p is supposed to be in the reference convex!

hess_base_value(p)

Evaluate the Hessian of all base functions of the Fem at point p.

p is supposed to be in the reference convex!

index_of_global_dof(cv)

Return the index of global dof for special fems such as interpolated fem.

is_equivalent()

Return 0 if the Fem is not equivalent.

Equivalent Fem are evaluated on the reference convex. This is the case of most classical Fem’s.

is_lagrange()

Return 0 if the Fem is not of Lagrange type.

is_polynomial()

Return 0 if the basis functions are not polynomials.

nbdof(cv=None)

Return the number of dof for the Fem.

Some specific Fem (for example ‘interpolated_fem’) may require a convex number cv to give their result. In most of the case, you can omit this convex number.

poly_str()

Return the polynomial expressions of its basis functions in the reference convex.

The result is expressed as a tuple of strings. Of course this will fail on non-polynomial Fem’s.

pts(cv=None)

Get the location of the dof on the reference element.

Some specific Fem may require a convex number cv to give their result (for example ‘interpolated_fem’). In most of the case, you can omit this convex number.

target_dim()

Return the dimension of the target space.

The target space dimension is usually 1, except for vector Fem.