.TH "basis" 2rheolef "Sat Mar 13 2021" "Version 7.1" "rheolef" \" -*- nroff -*- .ad l .nh .SH NAME basis \- finite element method (rheolef-7\&.1) .PP .SH "DESCRIPTION" .PP The class constructor takes a string that characterizes the finite element method, e\&.g\&. \fC'P0'\fP, \fC'P1'\fP or \fC'P2'\fP\&. The letters characterize the finite element family and the number is the polynomial degree in this family\&. The finite element also depends upon the reference element, e\&.g\&. triangle, square, tetrahedron, etc\&. See \fBreference_element(6)\fP for more details\&. For instance, on the square reference element, the \fC'P1'\fP string designates the common \fCQ1\fP four-nodes finite element\&. .SH "AVAILABLE FINITE ELEMENTS" .PP \fCPk\fP .PP .RS 4 The classical Lagrange finite element family, where \fCk\fP is any polynomial degree greater or equal to 1\&. .RE .PP \fCPkd\fP .PP .RS 4 The discontinuous Galerkin finite element family, where \fCk\fP is any polynomial degree greater or equal to 0\&. For convenience, \fCP0d\fP could also be written as \fCP0\fP\&. .RE .PP \fCbubble\fP .PP .RS 4 The product of baycentric coordinates\&. Only simplicial elements are supported: \fBedge(6)\fP, \fBtriangle(6)\fP and \fBtetrahedron(6)\fP\&. .RE .PP \fCRTk\fP .br \fCRTkd\fP .PP .RS 4 The vector-valued Raviart-Thomas family, where \fCk\fP is any polynomial degree greater or equal to 0\&. The \fCRTkd\fP piecewise discontinuous version is fully implemented for the \fBtriangle(6)\fP\&. Its \fCRTk\fP continuous counterpart is still under development\&. .RE .PP \fCBk\fP .PP .RS 4 The Bernstein finite element family, where \fCk\fP is any polynomial degree greater or equal to 1\&. This \fCbasis\fP was initially introduced by Bernstein (Comm\&. Soc\&. Math\&. Kharkov, 2th series, 1912) and more recently used in the context of finite elements\&. This feature is still under development\&. .RE .PP \fCSk\fP .PP .RS 4 The spectral finite element family, as proposed by Sherwin and Karniadakis (Oxford University Press, 2nd ed\&., 2005)\&. Here, \fCk\fP is any polynomial degree greater or equal to 1\&. This feature is still under development\&. .RE .PP .SH "CONTINUOUS FEATURE" .PP Some finite element families could support either continuous or discontinuous junction at inter-element boundaries\&. For instance, the \fCPk\fP Lagrange finite element family, with k >= 1, has a continuous interelements junction: it defines a piecewise polynomial and globally continuous finite element functional \fBspace(2)\fP\&. Conversely, the \fCPkd\fP Lagrange finite element family, with k >= 0, has a discontinuous interelements junction: it defines a piecewise polynomial and globally discontinuous finite element functional \fBspace(2)\fP\&. .SH "OPTIONS" .PP The \fCbasis\fP class supports some options associated to each finite element method\&. These options are appended to the string bewteen square braces, and are separated by a coma, e\&.g\&. \fC'P5[feteke,bernstein]'\fP\&. See \fBbasis(1)\fP for some examples\&. .SH "LAGRANGE NODES OPTION" .PP \fCequispaced\fP .PP .RS 4 Nodes are equispaced\&. This is the default\&. .RE .PP \fCfekete\fP .PP .RS 4 Node are non-equispaced: for high-order polynomial degree, e\&.g\&. greater or equal to 3, their placement is fully optimized, as proposed by Taylor, Wingate and Vincent, 2000, SIAM J\&. Numer\&. Anal\&. With this choice, the interpolation error is dramatically decreased for high order polynomials\&. This feature is still restricted to the triangle reference element and to polynomial degree lower or equal to 19\&. Otherwise, it fall back to the following \fCwarburton\fP node set\&. .RE .PP \fCwarburton\fP .PP .RS 4 Node are non-equispaced: for high-order polynomial degree, e\&.g\&. greater or equal to 3, their placement is optimized thought an heuristic solution, as proposed by Warburton, 2006, J\&. Eng\&. Math\&. With this choice, the interpolation error is dramatically decreased for high order polynomials\&. This feature applies to any reference element and polynomial degree\&. .RE .PP .SH "RAW POLYNOMIAL BASIS" .PP The raw (or initial) polynomial basis is used for building the dual basis, associated to degrees of freedom, via the Vandermonde matrix and its inverse\&. Changing the raw basis do not change the definition of the FEM basis, but only the way it is constructed\&. There are three possible raw basis: .PP \fCmonomial\fP .PP .RS 4 The monomial basis is the simplest choice\&. While it is suitable for low polynomial degree (less than five), for higher polynomial degree, the Vandermonde matrix becomes ill-conditioned and its inverse leads to errors in double precision\&. .RE .PP \fCdubiner\fP .PP .RS 4 The Dubiner basis (see Dubiner, 1991 J\&. Sci\&. Comput\&.) leads to much better condition number\&. This is the default\&. .RE .PP \fCbernstein\fP .PP .RS 4 The Bernstein basis could also be an alternative raw basis\&. .RE .PP .SH "INTERNALS: EVALUATE AND INTERPOLATE" .PP The \fCbasis\fP class defines member functions that evaluates all the polynomial basis functions of a finite element and their derivatives at a point or a set of point\&. .PP The basis polynomial functions are evaluated by the \fCevaluate\fP member function\&. This member function is called during the assembly process of the \fBintegrate(3)\fP function\&. .PP The interpolation nodes used by the \fBinterpolate(3)\fP function are available by the member function \fChat_node\fP\&. Conversely, the member function \fCcompute_dofs\fP compute the degrees of freedom\&. .SH "IMPLEMENTATION" .PP This documentation has been generated from file fem/lib/basis\&.h .PP The \fCbasis\fP class is simply an alias to the \fC\fBbasis_basic\fP\fP class .PP .PP .nf typedef basis_basic basis; .fi .PP \fB\fP .RS 4 .RE .PP The \fC\fBbasis_basic\fP\fP class provides an interface to a data container: .PP .PP .nf template class basis_basic : public smart_pointer_nocopy >, public persistent_table > { public: // typedefs: typedef basis_rep rep; typedef smart_pointer_nocopy base; typedef typename rep::size_type size_type; typedef typename rep::value_type value_type; typedef typename rep::valued_type valued_type; // allocators: basis_basic (std::string name = ""); void reset (std::string& name); void reset_family_index (size_type k); // accessors: bool is_initialized() const { return base::operator->() != 0; } size_type degree() const; size_type family_index() const; std::string family_name() const; std::string name() const; size_type ndof (reference_element hat_K) const; size_type nnod (reference_element hat_K) const; bool is_continuous() const; bool is_discontinuous() const; bool is_nodal() const; bool have_continuous_feature() const; bool have_compact_support_inside_element() const; bool is_hierarchical() const; size_type size() const; const basis_basic& operator[] (size_type i_comp) const; bool have_index_parameter() const; const basis_option& option() const; valued_type valued_tag() const; const std::string& valued() const; const piola_fem& get_piola_fem() const; .fi .PP .PP .nf }; .fi .PP .SH AUTHOR Pierre Saramito .SH COPYRIGHT Copyright (C) 2000-2018 Pierre Saramito GPLv3+: GNU GPL version 3 or later . This is free software: you are free to change and redistribute it. There is NO WARRANTY, to the extent permitted by law.