Interpolation Hermite Polynomials For Finite Element Method

We describe a new algorithm for analytic calculation of high-order Hermite interpolation polynomials of the simplex and give their classification. A typical example of triangle element, to be built in high accuracy finite element schemes, is given.


Introduction
For more than half a century, the finite element method (FEM) has won universal recognition as an efficient method for solving the most diverse problems of mathematical physics and engineering. In the multidimensional case the finite element grids of various shapes are used. The problem of constructing high-order interpolation polynomials for FEM has a simple solution only for simplex finite elements, such as the well known Lagrange interpolation polynomials (LIPs) [1]. Meanwhile, in the case of a d-dimensional simplex domain, the LIPs of the order p ≥ 1 are often sought by compiling and solving systems of (d + p )!/d!/p ! linear algebraic equations [2].
However, there are problems in which values of directional derivatives of the solutions are also necessary. They are of particular importance when high smoothness between the elements is required, or when highly accurate values of the gradient of the solution are necessary. The construction of such basis functions, referred to as Hermite interpolation polynomials (HIPs), is not possible on arbitrary meshes. This is one of the most important and difficult problems in the FEM and its applications in different fields, solved to date explicitly only for certain particular cases [2,3].
In this paper we report a new algorithm for the calculation of HIPs of the order p = κ max (p+1)−1 providing continuity of the approximating piecewise polynomial functions and of their directional derivatives up to an order κ along the normals to the boundaries of the simplex finite elements in the Euclidean space R d , which reduces to solving systems consisting of The conventional FEM implementation for a problem defined in terms of the set of coordinates z = (z 1 , . . . , z d ) ∈ R d performs all calculations in local (reference) coordinates z = (z 1 , . . . , z d ) ∈ R d , in which the d + 1 coordinates of the simplex vertices are [3]:ẑ j =(ẑ j1 , . . . ,ẑ jd ),ẑ jk =δ jk , j = 0, . . . , d, In the local coordinates of the d-dimensional simplex ∆, the LIP ϕ r (z ) of the order p = p which equates one at the node points ξ r = (n 1 /p, . . . , n d /p), n i ≥ 0, n 1 +· · ·+n d ≤ p, r = 1, . . . , (d + p)!/d!/p! and zero at the remaining node points ξ r , i.e., ϕ r (ξ r ) = δ rr , is determined by the formula: Step 1. To construct the HIPs in the local coordinates z , we define the set of auxiliary polynomials Here in contrast to the LIPs, the values of the functions themselves, and of their derivatives up to the order κ max −1 are specified at the node points ξ r . The explicit expressions of AP1 are given by where the coefficients a κ 1 ...κ d ,µ 1 ...µ d r are calculated from recurrence relations obtained by the substitution of Eq. (4) into the conditions (3).
Step 2. To enforce a uniquely defined polynomial basis, two types of auxiliary polynomials Q s (z) denoted respectively AP2 and AP3 are defined. AP2 and AP3 are linearly independent of AP1 from Eq. (4) and satisfy the following conditions at the node points ξ r of AP1: To provide the continuity of derivatives, the polynomials referred to as AP2 are asked to satisfy the condition Restriction of derivative order κ : 3pκ (κ + 1)/2 ≤ K where η s = (η s 1 , . . . , η s d ) are conveniently chosen points lying on the faces of various dimensionalities (from 1 to d − 1) of the d-dimensional simplex ∆ and do not coincide with the nodal points of HIP ξ r , where Eq. (3) is valid, ∂/∂n i(s) is the directional derivative along the vector n i , normal to the corresponding i-th face of the d-dimensional simplex ∆ q at the point η s in the physical frame, which is recalculated to the point η s of the face of the simplex ∆ in the local frame using the relations (1). Calculating the number T 1 (κ) of independent parameters required to provide the continuity of derivatives to the order κ, we determine its maximal value κ that can be obtained for the schemes with given p and κ max and, correspondingly, the additional conditions (6).
T 2 =K−T 1 (κ ) parameters remain independent and, correspondingly, T 2 additional conditions are added, necessary for the unique determination of the polynomials referred to as AP3, Q s (ζ s )=δ ss , s, s =T 1 (κ )+1, . . . , K, where ζ s = (ζ s 1 , . . . , ζ s d ) ∈ ∆ are the chosen points belonging to the simplex without the boundary, but not coincident with the node points of AP1 ξ r . The auxiliary polynomials AP2 and AP3 are given by the expression For AP2 k t = 1 if the point η s , at which the additional conditions (6) are specified, lies on the corresponding face of the simplex ∆, and k t = κ , otherwise, t = 0, . . . , d. For AP3 k t = κ , t = 0, . . . , d. The coefficients b j 1 ,..., j d ;s are determined from the uniquely solvable system of linear equations, obtained as a result of the substitution of Eq. (8) into the conditions (5)-(7).
For example, at d = 2 the derivatives ∂/∂n i along the direction n i , perpendicular to the appropriate face i = 0, 1, 2 in the physical frame are expressed in terms of the partial derivatives ∂/∂z j , j = 1, 2 in