Macroelement Analysis in T-Patches Using Lagrange Polynomials

1. Introduction

Lagrange polynomials have been extensively used from the very beginning of the Finite Element Method (FEM) [1-3]. First, linear and then quadratic elements were used in conjunction with the isoparametric concept [4]. In the 1980s, a widely accepted assumption was that the complexity in the implementation of elements of higher degree increases considerably, and therefore less attention was devoted to high-order discretization techniques (e.g., Ref. [5], p. 309). However, in the advent of powerful computers, the idea of using polynomials of higher degrees gradually gained acceptance and popularity [6,7]. Two alternative methodologies dealing with high polynomial degrees are: (i) Higher-order Finite Element Method (p-FEM) [8-10] and (ii) Spectral Finite Element Methods (SEM) [11]. The p-method (polynomial order refinement) in FEM increases the degree (p) of the polynomial basis functions used in each element while keeping the mesh fixed. Spectral methods also use high-order polynomial approximations but differ in that they use global basis functions rather than piecewise-defined basis functions. In the paper at hand, we concentrate on a variation of SEM as the topic is elements based on Lagrangian interpolation polynomials.

The opposition in the use of Lagrange polynomials is based on the well-known Runge phenomenon, where a specific non-polynomial function is gradually badly approximated by Lagrange polynomials as far as their degree increases [11-13]. Since the Lagrange polynomials obtain large values near the edges of the interval, it was found that the use of nonuniform node distributions (usually based on Gauß-Lobatto-Chebyshev (GLC) or Gauß-Lobatto-Legendre (GLL) points, which are dense near the ends), lead to excellent results [7,14,15].

Aiming at the easy manipulation and coupling of dissimilar geometric entities, since early 1970s transfinite elements were developed and used in automotive industry [16-18]. In Ref. [18] these elements have been called ‘macroelements’. In the same spirit, within the context of global approximation, which is highly related to SEM, FEM macroelements which occupy the entire or a large portion of the domain have been also successfully used outside the strict context of ideal tensor products of Lagrange polynomials, for example using transfinite interpolation (including Coons interpolation [19]) in conjunction with a few sets of trial functions such as piecewise-linear, B-splines, etc. [20,21]. As it is overviewed in a monograph [22], these elements perform well for the transition between unequally meshed opposite edges of a patch, for the coupling between dissimilar patches, as well as for the treatment of a structured set of internal nodes which do not exactly correspond to the boundary ones [15,23,24]. Among others, typical works in which an entire circle or quadrilateral patch is successfully modeled using a single transfinite element (based on trial functions in the form of Lagrange polynomials) are [25,26].

Unlike structured grids, in 1974, a pilot study demonstrated that it is possible for a rectangular polygon to be decomposed into rectangular cells, where locally blended cubic splines or Hermite interpolation may be applied [27]. Almost forty years later, in 2003, nets of control points which form T-junctions between perpendicular isolines ($\xi$=const. and $\eta$=const.), were efficiently treated by Sederberg et al. [28] in the context of computer-aided geometric design (CAGD)-surfaces using B-splines and NURBS and thus were called ‘T-splines’. Since 2010, the application of isogeometric analysis (IGA) using T-splines has been a matter of intensive research to date [29-33].

Based on the above discussion, in the paper at hand, the subdivision of a rectangular domain into non-overlapping quadrilateral elements such that (i) Elements are quadrilaterals (axis-aligned), (ii) T-junctions are allowed (meaning that not all edges need to extend completely across the mesh), and (iii) it is typically hierarchical (allowing local refinements by inserting new elements without disturbing the entire structure), is called ‘T-mesh’.

A recent report has shown that simple T-meshes may be treated using the transfinite interpolation concept, because the auxiliary nodes are automatically eliminated [34]. However, the general T-spline has not been investigated yet. Within this context, the paper at hand investigates the capability of transfinite elements to treat complicated T-meshed patches and eventually proposes a slight extension or an alternative technique using constraints to systematically achieve it. The numerical performance of this approach is verified through a boundary value problem in which the analytical solution of the Laplace equation is of non-polynomial form.

2. Tensor Product Elements

2.1. General Expressions

We consider a reference patch (unit square) in which the horizontal $x\left(\xi \right)$-axis ($0\le \xi \le 1$) is uniformly subdivided into $m$ segments, and thus the $\left(m+1\right)$ associated Lagrange polynomials (of degree $m$) form the column vector ${\mathbf{L}}_{m,x}$:

$${\mathbf{L}}_{m,x}={\left[\begin{array}{ccc}{L}_{1,m}\left(\xi \right)& \cdots & L_{m+1,m}\left(\xi \right)\end{array}\right]}^{\mathit{{\rm T}}}.$$

(1)

Similarly, the vertical $y\left(\eta \right)$-axis ($0\le \eta \le 1$) is uniformly subdivided into $n$ segments, and thus the $\left(n+1\right)$ associated Lagrange polynomials (of degree $n$) form the column vector ${\mathbf{L}}_{n,y}$:

$${\mathbf{L}}_{n,y}={\left[\begin{array}{ccc}{L}_{1,n}\left(\eta \right)& \cdots & L_{n+1,n}\left(\eta \right)\end{array}\right]}^{\mathit{{\rm T}}}.$$

(2)

In a tensor product of degree ($m,n$), using $\left(m+1\right)$ and $\left(n+1\right)$ nodal points along the axes $x$ and $y$, respectively (i.e., totally $\left(m+1\right)\times \left(n+1\right)$ nodal points) in the interval $0\le \xi ,\eta \le 1$, we can write (see [12]):

$$U\left(\xi ,\eta \right)={\mathbf{L}}_{m,x}^T\mathbf{U}{\mathbf{L}}_{n,y},$$

(3)

where the matrix $\mathbf{U}$ includes all the nodal values in a structured way (see [34]).

2.2. Quadratic Interpolation

As an example, for a quadratic tensor product element of Lagrange type, the univariate set ${\mathbf{L}}_{2,x}$ is defined by the nodes shown in Figure 1a, the element in which $m=n=2$ is shown in Figure 1b, and thus we have:

$$L_{2,x}=\left[\begin{array}{c}{L}_{1,2}(\xi )\\ L_{2,2}(\xi )\\ L_{3,2}(\xi )\end{array}\right]=\left[\begin{array}{c}2(\xi -1)(\xi -1/2)\\ -4\xi (\xi -1)\\ 2\xi (\xi -1/2)\end{array}\right]=\left[\begin{array}{c}2{\xi }^2-3\xi +1\\ -4{\xi }^2+4\xi \\ 2{\xi }^2-\xi \end{array}\right],$$

(4)

$$L_{2,y}=\left[\begin{array}{c}{L}_{1,2}(\eta )\\ L_{2,2}(\eta )\\ L_{3,2}(\eta )\end{array}\right]=\left[\begin{array}{c}2(\eta -1)(\eta -1/2)\\ -4\eta (\eta -1)\\ 2\eta (\eta -1/2)\end{array}\right]=\left[\begin{array}{c}2{\eta }^2-3\eta +1\\ -4{\eta }^2+4\eta \\ 2{\eta }^2-\eta \end{array}\right],$$

(5)

and

$$U=\left[\begin{array}{ccc}{U}_1& U_4& U_7\\ U_2& U_5& U_8\\ U_3& U_6& U_9\end{array}\right].$$

(6)

Substituting Equations (4) to (6) into Equation (3), results in the well-known bi-quadratic shape functions:

$$\begin{array}{l}{\varphi }_1(\xi ,\eta )=L_{1,2}(\xi )L_{1,2}(\eta ),\\ {\varphi }_2(\xi ,\eta )=L_{2,2}(\xi )L_{1,2}(\eta ),\\ {\varphi }_3(\xi ,\eta )=L_{3,2}(\xi )L_{1,2}(\eta ),\\ {\varphi }_4(\xi ,\eta )=L_{1,2}(\xi )L_{2,2}(\eta ),\\ {\varphi }_5(\xi ,\eta )=L_{2,2}(\xi )L_{2,2}(\eta ),\\ {\varphi }_6(\xi ,\eta )=L_{3,2}(\xi )L_{2,2}(\eta ),\\ {\varphi }_7(\xi ,\eta )=L_{1,2}(\xi )L_{3,2}(\eta ),\\ {\varphi }_8(\xi ,\eta )=L_{2,2}(\xi )L_{3,2}(\eta ),\\ {\varphi }_9(\xi ,\eta )=L_{3,2}(\xi )L_{3,2}(\eta ).\end{array}$$

(7)

3. Traditional Transfinite Elements

There are many higher order macroelements that cannot be described through a tensor product (discussed in Section 2), and some of them have been discussed in Ref. [34]. To prepare the reader for the next section, which is the main topic of the paper at hand, we first focus on traditional structured transfinite elements, which are characterized by some horizontal and vertical stations that are subdivided by several nodes, like those in the 21-node element shown in Figure 2. In Ref. [35], it has been discussed how the individual projectors ($P_{\xi },P_{\eta },P_{\xi \eta }$) of the Boolean sum formulation can be constructed for quadratic blending functions (Lagrange polynomials) and quartic trial functions (Lagrange polynomials. To shorten their presentation, here we prefer an alternative form exploiting matrix notation as follows:

$$P_{\xi }=\underset{\left(L_{4,y}^T\right)}{\underbrace{\left[\begin{array}{ccccc}{L}_{1,4}(\eta )& L_{2,4}(\eta )& L_{3,4}(\eta )& L_{4,4}(\eta )& L_{5,4}(\eta )\end{array}\right]}}\cdot \underset{\left[U_{\xi }\right]}{\underbrace{\left[\begin{array}{ccc}{U}_1& U_3& U_5\\ U_6& U_7& U_8\\ U_9& U_{11}& U_{13}\\ U_{14}& U_{15}& U_{16}\\ U_{17}& U_{19}& U_{21}\end{array}\right]}}\cdot \underset{\left(L_{2,x}\right)}{\underbrace{\left[\begin{array}{c}{L}_{1,2}(\xi )\\ L_{2,2}(\xi )\\ L_{3,2}(\xi )\end{array}\right]}},$$

(8)

$$P_{\eta }=\underset{\left(L_{4,x}^T\right)}{\underbrace{\left[\begin{array}{ccccc}{L}_{1,4}(\xi )& L_{2,4}(\xi )& L_{3,4}(\xi )& L_{4,4}(\xi )& L_{5,4}(\xi )\end{array}\right]}}\cdot \underset{\left[U_{\eta }\right]}{\underbrace{\left[\begin{array}{ccc}{U}_1& U_9& U_{17}\\ U_2& U_{10}& U_{18}\\ U_3& U_{11}& U_{19}\\ U_4& U_{12}& U_{20}\\ U_5& U_{13}& U_{21}\end{array}\right]}}\cdot \underset{\left(L_{2,y}\right)}{\underbrace{\left[\begin{array}{c}{L}_{1,2}(\eta )\\ L_{2,2}(\eta )\\ L_{3,2}(\eta )\end{array}\right]}},$$

(9)

and

$$P_{\xi \eta }=\underset{\left(L_{2,x}^T\right)}{\underbrace{\left[\begin{array}{ccc}{L}_{1,2}(\xi )& L_{2,2}(\xi )& L_{3,2}(\xi )\end{array}\right]}}\cdot \underset{\left[U_{\xi \eta }\right]}{\underbrace{\left[\begin{array}{ccc}{U}_1& U_9& U_{17}\\ U_3& U_{11}& U_{19}\\ U_5& U_{13}& U_{21}\end{array}\right]}}\cdot \underset{\left(L_{2,y}\right)}{\underbrace{\left[\begin{array}{c}{L}_{1,2}(\eta )\\ L_{2,2}(\eta )\\ L_{3,2}(\eta )\end{array}\right]}}.$$

(10)

According to the standard theory, the transfinite interpolation of bivariate functions expressed be $U(\xi ,\eta )$ is given by the following Boolean sum:

$$U(\xi ,\eta )=P_{\xi }\oplus P_{\eta }=P_{\xi }+P_{\eta }-P_{\xi \eta }.$$

(11)

Substituting the three projectors, i.e., Equations (8) to (10) into Equation (11), extraction of coefficients gives the expression:

$$U(\xi ,\eta )={\displaystyle \sum _{j=1}^{21}{\varphi }_j}(\xi ,\eta )\cdot U_j,$$

(12)

where the global bivariate shape functions ${\varphi }_j(\xi ,\eta )$ will be given in terms of univariate Lagrange polynomials by:

$$\begin{array}{l}{\varphi }_1(\xi ,\eta )=L_{1,2}(\xi )L_{1,4}(\eta )+L_{1,4}(\xi )L_{1,2}(\eta )-L_{1,2}(\xi )L_{1,2}(\eta ),\\ {\varphi }_2(\xi ,\eta )=L_{2,4}(\xi )L_{1,2}(\eta ),\\ {\varphi }_3(\xi ,\eta )=L_{2,2}(\xi )L_{1,4}(\eta )+L_{3,4}(\xi )L_{1,2}(\eta )-L_{2,2}(\xi )L_{1,2}(\eta ),\\ {\varphi }_4(\xi ,\eta )=L_{4,4}(\xi )L_{1,2}(\eta ),\\ {\varphi }_5(\xi ,\eta )=L_{3,2}(\xi )L_{1,4}(\eta )+L_{5,4}(\xi )L_{1,2}(\eta )-L_{3,2}(\xi )L_{1,2}(\eta ),\\ {\varphi }_6(\xi ,\eta )=L_{1,2}(\xi )L_{2,4}(\eta ),\\ {\varphi }_7(\xi ,\eta )=L_{2,2}(\xi )L_{2,4}(\eta ),\\ {\varphi }_8(\xi ,\eta )=L_{3,2}(\xi )L_{2,4}(\eta ),\\ {\varphi }_9(\xi ,\eta )=L_{1,2}(\xi )L_{3,4}(\eta )+L_{1,4}(\xi )L_{2,2}(\eta )-L_{1,2}(\xi )L_{2,2}(\eta ),\\ {\varphi }_{10}(\xi ,\eta )=L_{2,4}(\xi )L_{2,2}(\eta ),\\ {\varphi }_{11}(\xi ,\eta )=L_{2,2}(\xi )L_{3,4}(\eta )+L_{3,4}(\xi )L_{2,2}(\eta )-L_{2,2}(\xi )L_{2,2}(\eta ),\\ {\varphi }_{12}(\xi ,\eta )=L_{4,4}(\xi )L_{2,2}(\eta ),\\ {\varphi }_{13}(\xi ,\eta )=L_{3,2}(\xi )L_{3,4}(\eta )+L_{5,4}(\xi )L_{2,2}(\eta )-L_{3,2}(\xi )L_{2,2}(\eta ),\\ {\varphi }_{14}(\xi ,\eta )=L_{1,2}(\xi )L_{4,4}(\eta ),\\ {\varphi }_{15}(\xi ,\eta )=L_{2,2}(\xi )L_{4,4}(\eta ),\\ {\varphi }_{16}(\xi ,\eta )=L_{3,2}(\xi )L_{4,4}(\eta ),\\ {\varphi }_{17}(\xi ,\eta )=L_{1,2}(\xi )L_{5,4}(\eta )+L_{1,4}(\xi )L_{3,2}(\eta )-L_{1,2}(\xi )L_{3,2}(\eta ),\\ {\varphi }_{18}(\xi ,\eta )=L_{2,4}(\xi )L_{3,2}(\eta ),\\ {\varphi }_{19}(\xi ,\eta )=L_{2,2}(\xi )L_{5,4}(\eta )+L_{3,4}(\xi )L_{3,2}(\eta )-L_{2,2}(\xi )L_{3,2}(\eta ),\\ {\varphi }_{20}(\xi ,\eta )=L_{4,4}(\xi )L_{3,2}(\eta ),\\ {\varphi }_{21}(\xi ,\eta )=L_{3,2}(\xi )L_{5,4}(\eta )+L_{5,4}(\xi )L_{3,2}(\eta )-L_{3,2}(\xi )L_{3,2}(\eta ).\end{array}$$

(13)

As has been previously discussed in [35] and elsewhere, the form of the bivariate shape functions in Equation (13) is generic and can be easily computerized. In brief, nodal points at the intersections of stations (i.e., nodes 1,3,5,9,11,13,17,19,21 in Figure 2) are associated with shape functions (${\varphi }_1,{\varphi }_3,{\varphi }_5,\dots ,{\varphi }_{19},{\varphi }_{21}$) which are influenced by all the three projectors, whereas the rest are influenced by only one projector (the one perpendicular to the isoline the node belongs to). It is worth mentioning that each projector is well-defined by the actual nodes of the mesh and does not require any auxiliary points to be constructed. The set of the 21 shape functions is illustrated in Figure 3.

4. T-Meshed Patch Elements: Constraints and Elimination

In this section, we deal with higher order macroelements that are neither described by a tensor product nor have the form of the traditional structured elements discussed in Section 3. The key characteristic of these elements is that they come from a tensor product in which some nodes are missing. As shown in Figure 4, in general we could categorize these elements in the following three classes:

• Elements where the missing nodes belong to the boundary of the patch (e.g., nodes G and F as illustrated in Figure 4a).
• Elements with hanging internal nodes (denoted by a red cross $\times$) that belong to the extension of isolines in only one direction (e.g., in Figure 4b these isolines are directed toward the $\eta$-direction).
• Elements with hanging internal nodes that belong to the intersection of two isolines (e.g., node H as illustrated in Figure 4c).

In all formulations of the paper at hand, the missing nodes to form a complete tensor product are initially filled by auxiliary points which are also called ‘artificial’.

Regarding the first category, where the missing nodes belong to the boundary (shown by letter G in Figure 4), obviously they do not influence the interpolation along the corresponding edge, but only influence the interpolation in the vertical direction of the isoline to which they belong. For example, the interpolation along the bottom edges will be of degree $p=2, 3$ and 2 for Figures 4a, 4b, and 4c, respectively, although the artificial node G was introduced. As we shall see below, the bivariate shape functions of these elements associated with boundary nodes on incomplete edges, may be directly obtained applying transfinite interpolation. This happens because each artificial boundary node belongs to the projector having as subscript the direction being vertical to the corresponding edge. For example, nodes G and F (Figure 4a) are involved in the projectors $P_{\xi }$ and $P_{\eta }$, respectively, but eventually are automatically eliminated because they appear in the subtracted projector $P_{\xi \eta }$ as well.

Concerning the second category, the artificial nodes eventually appear in the projector whose index indicates the direction in which to extend the isoline. This is necessary to connect the internal hanging nodes to the boundary. For example, in Figure 4b we need to extend the two vertical isolines from inside to the bottom edge, and thus the artificial nodes belong to the projector $P_{\xi }$, which is further cancelled by the corresponding terms from the subtracted projector $P_{\xi \eta }$. Note that this is possible, because both projectors are based on the blending functions. Again, with respect to Figure 4b, it should become clear that, if the two nodes that belong to the horizontal isoline 5-6-7-8 had been used to define the trial functions along 5-8 (as a uniform polynomial of six nodes and thus degree 5, instead of as a nonuniform polynomial consisting of nodes 5-6-7-8 and thus of degree 3), then they would have been included in the projector $P_{\eta }$ as well, and thus could not be eliminated.

Regarding the third category, the artificial (auxiliary) node H (Figure 4c) belongs to the projector $P_{\xi \eta }$, because it lies at the intersection of the horizontal (4-5) and vertical (8-12-16) stations. For reasons of equal treatment, we cannot prefer a specific direction. Therefore, H will belong to $P_{\xi }$ as well as to $P_{\eta }$, resulting in the inability to delete it within the algebraic sum $P_{\xi }+P_{\eta }-P_{\xi \eta }$. In general, there are two possible techniques:

• Technique 1: Eliminate H by considering it within the projector $P_{\xi }$ but ignoring it within the projector $P_{\eta }$, i.e., considering the trial function along the horizontal isoline $\eta ={\eta }_4={\eta }_5$ to be defined by the nodes (4, 5, and F). However, this trick cannot be applied to complicated T-meshes.
• Technique 2: Interpolate the nodal value at H once along a horizontal and another time along the vertical isoline passing through H and then consider the mean average value of these interpolations. Therefore, it is possible to eliminate the auxiliary node H, as later shown in Subsection 5.5.

5. Constraints

Regarding the first category where the missing nodes belong to the boundary, an indicative collection of elements in which a single boundary node (designated by G) is missing on the bottom edge, is shown in Figure 5. The polynomial interpolation along the bottom edge (to which G belongs) is of degree $p_a=1$, $p_b=2$ and $p_c=3$ for Figures 5a, 5b, and 5c, respectively. As we will discuss below in Subsections 5.1 to 5.3, the same interpolation is produced when starting from an initial interpolation of one unit higher (${p\text{'}}_a=2$, ${p\text{'}}_b=3$ and ${p\text{'}}_c=4$) in which the proper constraint is imposed along the nodal points of the bottom edge. The next subsections (5.1 to 5.3) are consistent with the discretization of the bottom edge in the three cases shown in Figure 5. Note that the entire 8-node transfinite element shown in Figure 5a is discussed in detail within Section 5.4.

5.1. Quadratic Polynomials and Linear Constraints

Let us consider a polynomial of second degree, which can be written either in the power form: $f_1\left(\xi \right)=c_2{\xi }^2+c_1\xi +c_0$, or as a sum of three quadratic Lagrange polynomials:

$$f_2\left(\xi \right)=\left(2\xi -1\right)\left(\xi -1\right)U_1+4\xi \left(1-\xi \right)U_2+\xi \left(2\xi -1\right)U_3,$$

(14)

where $\left(U_1,U_2,U_3\right)$ are the nodal values associated with the nodes (1,2,3) shown in Figure 1a (in local numbering), as well as for the horizontal isolines of Figure 6a, of which the bottom edge is of interest.

If we assume that node 2 is located in the middle of the interval [0,1] (Figure 6a), having the average value of the two ends:

$$U_2=(U_1+U_3)/2,$$

(15)

after performing trivial manual operations (replace $f_1\left(\xi \right)$ and $f_2\left(\xi \right)$ with $U\left(\xi \right)$), we find that the second-order polynomial becomes linear:

$$\begin{array}{l}U(\xi )=\left[(2\xi -1)(\xi -1)\right]U_1+\left[4\xi (1-\xi )\right]U_2+\left[\xi (2\xi -1)\right]U_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left[(2\xi -1)(\xi -1)\right]U_1+\left[4\xi (1-\xi )\right]\left({\displaystyle \frac{U_1+U_3}{2}}\right)+\left[\xi (2\xi -1)\right]U_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=(1-\xi )U_1+\xi U_3.\end{array}$$

(16)

In other words, the constraint $U_2=(U_1+U_3)/2$ forces the full quadratic polynomial to degenerate to a linear one.

5.2. Cubic Polynomials and Quadratic Constraints

Now we test the case in which the interpolation is quadratic, however the intermediate (internal) node ‘2’ is not in the middle as usual, but at one-third measured from the left end (i.e., is located at ${\xi }_2=1/3$), and thus the node ‘3’ at position ${\xi }_3=2/3$ is missing (shown in red color in the right side of Figure 6b). Clearly, we have a unit length with three non-uniform nodes: 1 (at $\xi =0$), 2 (at $\xi =1/3$), and 4 (at $\xi =1$), and thus non-uniform Lagrange polynomials of degree two are created. In the previously mentioned sequence, node ‘3’ was omitted because it is eliminated from the uniform sequence of four nodes: 1-2-3-4, and thus only the nodes 1-2-4 remain (right part of Figure 6b).

Since node ‘3’ is missing, this in turn means that if we consider the value $U_3$ produced by the non-uniform polynomial of degree 2 at node 3, as a constraint in terms of ($U_1,U_2,U_4$), if the latter is introduced into the series expansion of four (cubic) Lagrange polynomials, we shall identically derive the expression of the non-uniform quadratic polynomial.

Actually, considering the three non-uniform Lagrange polynomials (based on the nodal values, ${\xi }_1=0,{\xi }_2=1/3,{\xi }_4=1$), we obtain:

$$U(\xi )=\left[3(\xi -1)(\xi -1/3)\right]U_1+\left[-{\displaystyle \frac{9}{2}}\xi (\xi -1)\right]U_2+\left[{\displaystyle \frac{3}{2}}\xi (\xi -{\displaystyle \frac{1}{3}})\right]U_4.$$

(17)

Setting in Equation (17) the value (${\xi =\xi }_3=2/3$), we derive $U_3=U\left({\xi }_3\right)$, and thus the constraint becomes:

$$U_3=-{\displaystyle \frac{1}{3}}{U}_1+U_2+{\displaystyle \frac{1}{3}}{U}_4.$$

(18)

Using the full cubic Lagrange polynomials based on the uniform interpolation points at (${\xi }_1=0,{\xi }_2={\displaystyle \frac{1}{3}},{\xi }_3={\displaystyle \frac{2}{3}},{\xi }_4=1$), the univariate function $U\left(\xi \right)$ is written in terms of the associated nodal values ($U_1,U_2,U_3$, and $U_4$) as:

$$\begin{array}{l}U(\xi )=\left[-{\displaystyle \frac{9}{2}}(\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\left(\xi -{\displaystyle \frac{2}{3}}\right)\right]U_1+\left[{\displaystyle \frac{27}{2}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{2}{3}}\right)\right]U_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[-{\displaystyle \frac{27}{2}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\right]U_3+\left[{\displaystyle \frac{9}{2}}\xi \left(\xi -{\displaystyle \frac{1}{3}}\right)\left(\xi -{\displaystyle \frac{2}{3}}\right)\right]U_4.\end{array}$$

(19)

Substituting the constraint described by Equation (18) into the general Equation (19), we receive:

$$\begin{array}{l}U(\xi )=\left[-{\displaystyle \frac{9}{2}}(\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\left(\xi -{\displaystyle \frac{2}{3}}\right)\right]U_1+\left[{\displaystyle \frac{27}{2}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{2}{3}}\right)\right]U_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[-{\displaystyle \frac{27}{2}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\right]\left(-{\displaystyle \frac{1}{3}}{U}_1+U_2+{\displaystyle \frac{1}{3}}{U}_4\right)+\left[{\displaystyle \frac{9}{2}}\xi \left(\xi -{\displaystyle \frac{1}{3}}\right)\left(\xi -{\displaystyle \frac{2}{3}}\right)\right]U_4\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left[-{\displaystyle \frac{9}{2}}(\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\left(\xi -{\displaystyle \frac{2}{3}}\right)+{\displaystyle \frac{27}{6}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\right]U_1\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[{\displaystyle \frac{27}{2}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{2}{3}}\right)-{\displaystyle \frac{27}{2}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\right]U_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[{\displaystyle \frac{9}{2}}\xi \left(\xi -{\displaystyle \frac{1}{3}}\right)\left(\xi -{\displaystyle \frac{2}{3}}\right)-{\displaystyle \frac{27}{6}}\xi (\xi -1)\left(\xi -{\displaystyle \frac{1}{3}}\right)\right]U_4\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left[3(\xi -1)(\xi -1/3)\right]U_1+\left[-{\displaystyle \frac{9}{2}}\xi (\xi -1)\right]U_2+\left[{\displaystyle \frac{3}{2}}\xi (\xi -{\displaystyle \frac{1}{3}})\right]U_4.\end{array}$$

(20)

One may observe that the last equality in Equation (20) is identical with the expression of Equation (17). In other words, the imposition of the linear constraint by Equation (18) (which refers to a non-uniform set of Lagrange polynomials of degree $p=2$) into the full expression by Equation (19) (which refers to a uniform set of Lagrange polynomials of degree $p\mathrm{\text{'}}=3$), results in the same approximation as that of the non-uniform set of Lagrange polynomials.

5.3. Quartic Polynomials and Cubic Constraints

In this section, we consider the nodes with global and local numbering (1-2-3-4-5) along the bottom edge (Figure 6c), where node 2 (at $\xi =1/4$) is to be eliminated, so as the initial quartic polynomial degenerates to a cubic one (a similar elimination is valid for the bottom edge in Figure 5c).

Since the local node ‘2’ is missing, this in turn means that if we consider the value $U_2$ being produced by the non-uniform cubic polynomial (based on nodes 1-3-4-5) at the point 2, as a constraint in terms of ($U_1,U_3,U_4,U_5$) and the latter is introduced into the series expansion of five (quartic) Lagrange polynomials (based on nodes 1-2-3-4-5), we shall identically derive the expression of the non-uniform cubic polynomial.

Therefore, using the explicit expressions of the four non-uniform Lagrange polynomials of degree $p=3$ (based on the nodal values, ${\xi }_1=0,{\xi }_3=2/4,{{\xi }_4=3/4,\xi }_5=1$):

$$\begin{array}{l}U(\xi )=\left[-{\displaystyle \frac{8}{3}}(\xi -1)(\xi -1/2)(\xi -3/4)\right]U_1+\left[16\xi (\xi -1)(\xi -3/4)\right]U_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[-{\displaystyle \frac{64}{3}}\xi (\xi -1)(\xi -{\displaystyle \frac{1}{2}})\right]U_4+\left[8\xi (\xi -1/2)(\xi -3/4)\right]U_5,\end{array}$$

(21)

when setting the value${\xi =\xi }_2=1/4$ in Equation (21), we obtain the constraint:

$$U_2={\displaystyle \frac{1}{4}}{U}_1+{\displaystyle \frac{3}{2}}{U}_3-U_4+{\displaystyle \frac{1}{4}}{U}_5.$$

(22)

Now we consider the set of uniform Lagrange polynomials of degree $p\text{'}=4$, in which the univariate function $U\left(\xi \right)$ is approximated as follows:

$$\begin{array}{l}U(\xi )=\left[{\displaystyle \frac{32}{3}}(\xi -1)(\xi -1/2)(\xi -1/4)(\xi -3/4)\right]U_1+\left[\xi (\xi -1)(\xi -1/2)(\xi -3/4)\right]U_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[64\xi (\xi -1)(\xi -{\displaystyle \frac{1}{4}})(\xi -{\displaystyle \frac{3}{4}})\right]U_3+\left[-{\displaystyle \frac{128}{3}}\xi (\xi -1)(\xi -1/2)(\xi -1/4)\right]U_4\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[{\displaystyle \frac{32}{3}}\xi (\xi -1/2)(\xi -{\displaystyle \frac{1}{4}})(\xi -{\displaystyle \frac{3}{4}})\right]U_5.\end{array}$$

(23)

Substituting the current constraint (i.e., Equation (22)) into Equation (23), results in:

$$\begin{array}{l}U(\xi )=\left[{\displaystyle \frac{32}{3}}(\xi -1)(\xi -1/2)(\xi -1/4)(\xi -3/4)\right]U_1\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[\xi (\xi -1)(\xi -1/2)(\xi -3/4)\right]\left({\displaystyle \frac{1}{4}}{U}_1+{\displaystyle \frac{3}{2}}{U}_3-U_4+{\displaystyle \frac{1}{4}}{U}_5\right)\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[64\xi (\xi -1)(\xi -{\displaystyle \frac{1}{4}})(\xi -{\displaystyle \frac{3}{4}})\right]U_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[-{\displaystyle \frac{128}{3}}\xi (\xi -1)(\xi -1/2)(\xi -1/4)\right]U_4\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[{\displaystyle \frac{32}{3}}\xi (\xi -1/2)(\xi -{\displaystyle \frac{1}{4}})(\xi -{\displaystyle \frac{3}{4}})\right]U_5\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left\{\left[{\displaystyle \frac{32}{3}}(\xi -1)(\xi -1/2)(\xi -1/4)(\xi -3/4)\right]+{\displaystyle \frac{1}{4}}\left[\xi (\xi -1)(\xi -1/2)(\xi -3/4)\right]\right\}{U}_1\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left\{\left[64\xi (\xi -1)(\xi -{\displaystyle \frac{1}{4}})(\xi -{\displaystyle \frac{3}{4}})\right]+{\displaystyle \frac{3}{2}}\left[\xi (\xi -1)(\xi -1/2)(\xi -3/4)\right]\right\}{U}_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left\{\left[-{\displaystyle \frac{128}{3}}\xi (\xi -1)(\xi -1/2)(\xi -1/4)\right]-\left[\xi (\xi -1)(\xi -1/2)(\xi -3/4)\right]\right\}{U}_4\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left\{\left[{\displaystyle \frac{32}{3}}\xi (\xi -1/2)(\xi -{\displaystyle \frac{1}{4}})(\xi -{\displaystyle \frac{3}{4}})\right]+{\displaystyle \frac{1}{4}}\left[\xi (\xi -1)(\xi -1/2)(\xi -3/4)\right]\right\}{U}_5\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left[-{\displaystyle \frac{8}{3}}(\xi -1)(\xi -1/2)(\xi -3/4)\right]U_1+\left[16\xi (\xi -1)(\xi -3/4)\right]U_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[-{\displaystyle \frac{64}{3}}\xi (\xi -1)(\xi -{\displaystyle \frac{1}{2}})\right]U_4+\left[8\xi (\xi -1/2)(\xi -3/4)\right]U_5.\end{array}$$

(24)

One may observe that the last equality of Equation (24) is identical to Equation (21). In other words, as previously happened, the implementation of the linear constraint forces the full set of uniform Lagrange polynomials of degree $p$’ to degenerate to the same set as of non-uniform ones with degree $p=p^{\mathrm{\text{'}}}-1$.

5.4. Implementation of an Eight-Node Transfinite Element

This subsection discusses the construction of a transfinite element, which is produced by a tensor product with missing boundary nodes. Although the procedure is general, for the sake of brevity we focus on one missing boundary node and a tensor product of quadratic interpolation ($p=2$).

Therefore, let us consider the 8-node element shown in Figure 7 by adding the inactive auxiliary node G. There are at least two approaches to derive the shape functions of this element, as follows:

5.4.1. Approach 1: Transfinite Interpolation

We construct three projectors, ($P_{\xi },P_{\eta },P_{\xi \eta }$). To produce $P_{\xi }$, along the middle vertical station at $\xi =1/2$, we need to introduce the artificial (auxiliary) node G (Figure 7), which is everywhere denoted in red color:

$$\begin{array}{l}{P}_{\xi }=\left[L_{1,2}(\eta )U_1+L_{2,2}(\eta )U_3+L_{3,2}(\eta )U_6\right]L_{1,2}(\xi )\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{1,2}(\eta )U_G+L_{2,2}(\eta )U_4+L_{3,2}(\eta )U_7\right]L_{2,2}(\xi )\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{1,2}(\eta )U_2+L_{2,2}(\eta )U_5+L_{3,2}(\eta )U_8\right]L_{3,2}(\xi ).\end{array}$$

(25)

In contrast, the artificial node G is not included in the projector $P_{\eta }$, because the representation of $U(\xi ,\eta )$ along the bottom boundary edge is a-priori linear, and thus we have ($L_{1,1}(\xi )=1-\xi ,L_{2,1}(\xi )=\xi$):

$$\begin{array}{l}{P}_{\eta }=\left[L_{1,1}(\xi )U_1+\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}0\cdot U_G\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+L_{2,1}(\xi )U_2\right]L_{1,2}(\eta )\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{1,2}(\xi )U_3+L_{2,2}(\xi )U_4+L_{3,2}(\xi )U_5\right]L_{2,2}(\eta )\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{1,2}(\xi )U_6+L_{2,2}(\xi )U_7+L_{3,2}(\xi )U_8\right]L_{3,2}(\eta ).\end{array}$$

(26)

The third projector is:

$$\begin{array}{l}{P}_{\xi \eta }=L_{1,2}(\xi )L_{1,2}(\eta )U_1+L_{2,2}(\xi )L_{1,2}(\eta )U_G+L_{3,2}(\xi )L_{1,2}(\eta )U_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+L_{1,2}(\xi )L_{2,2}(\eta )U_3+L_{2,2}(\xi )L_{2,2}(\eta )U_4+L_{3,2}(\xi )L_{2,2}(\eta )U_5\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+L_{1,2}(\xi )L_{3,2}(\eta )U_6+L_{2,2}(\xi )L_{3,2}(\eta )U_7+L_{3,2}(\xi )L_{3,2}(\eta )U_8.\end{array}$$

(27)

One may observe that the artificial node ‘G’ exists in only two out of the three projectors (for which we have $P_{\xi }=P_{\xi \eta }$), and thus when the Boolean sum of Equation (11) is implemented, the variable $U_G$ is eventually eliminated. This in turn means that the dependent variable $U(\xi ,\eta )$ is written in terms of the following eight shape functions:

$$\begin{array}{l}U(\xi ,\eta )=\bcancel{P_{\xi }}+P_{\eta }-\bcancel{P_{\xi \eta }}=P_{\eta }\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left[L_{1,1}(\xi )L_{1,2}(\eta )\right]U_1\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{2,1}(\xi )L_{1,2}(\eta )\right]U_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{1,2}(\xi )L_{2,2}(\eta )\right]U_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{2,2}(\xi )L_{2,2}(\eta )\right]U_4\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{3,2}(\xi )L_{2,2}(\eta )\right]U_5\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{1,2}(\xi )L_{3,2}(\eta )\right]U_6\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{2,2}(\xi )L_{3,2}(\eta )\right]U_7\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[L_{3,2}(\xi )L_{3,2}(\eta )\right]U_8.\end{array}$$

(28)

Regarding the bottom edge 1-2 which is included in the projector $P_{\eta }$, one may observe in the first square brackets of Equation (26) that the univariate function $U(\xi ,0)$ is directly approximated by two polynomials of first degree $\left(U(\xi ,0)=L_{1,1}(\xi )U_1+L_{2,1}(\xi )U_2\right)$, and thus no constraint is considered.

5.4.2. Approach 2: Successive Node Elimination from an Initial Tensor Product

In this approach we start with the full tensor product, of which node ‘G’ (Figure 7) is vital part:

$$\begin{array}{l}U(\xi ,\eta )=\left[\begin{array}{ccc}{L}_{1,2}(\xi )& L_{2,2}(\xi )& L_{3,2}(\xi )\end{array}\right]\cdot \left[\begin{array}{ccc}{U}_1& U_3& U_6\\ U_G& U_4& U_7\\ U_2& U_5& U_8\end{array}\right]\cdot \left[\begin{array}{c}{L}_{1,2}(\eta )\\ L_{2,2}(\eta )\\ L_{3,2}(\eta )\end{array}\right]\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=L_{1,2}(\xi )L_{1,2}(\eta )U_1+L_{2,2}(\xi )L_{1,2}(\eta )U_G+L_{3,2}(\xi )L_{1,2}(\eta )U_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+L_{1,2}(\xi )L_{2,2}(\eta )U_3+L_{2,2}(\xi )L_{2,2}(\eta )U_4+L_{3,2}(\xi )L_{2,2}(\eta )U_5\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+L_{1,2}(\xi )L_{3,2}(\eta )U_6+L_{2,2}(\xi )L_{3,2}(\eta )U_7+L_{3,2}(\xi )L_{3,2}(\eta )U_8.\end{array}$$

(29)

Then, in consistency to Section 5.1, we impose the constraint:

$$U_G=\left(U_1+U_2\right)/2.$$

(30)

Substituting Equation (30) into Equation (29), the action of the basis function $L_{2,2}(\xi )L_{1,2}(\eta )$, which is the factor of the eliminated variable $U_G$, is related to two DOFs (i.e., $U_1$ and $U_2$), and thus Equation (29) becomes:

$$\begin{array}{l}\widehat{U}(\xi ,\eta )=\left\{\stackrel{(1-\xi )}{\overbrace{\left[L_{1,2}(\xi )+{\displaystyle \frac{1}{2}}{L}_{2,2}(\xi )\right]}}{L}_{1,2}(\eta )\right\}{U}_1+0\cdot U_G\hspace{0.17em}+\left\{\stackrel{\xi }{\overbrace{\left[L_{3,2}(\xi )+{\displaystyle \frac{1}{2}}{L}_{2,2}(\xi )\right]}}{L}_{1,2}(\eta )\right\}{U}_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+L_{1,2}(\xi )L_{2,2}(\eta )U_3+L_{2,2}(\xi )L_{2,2}(\eta )U_4+L_{3,2}(\xi )L_{2,2}(\eta )U_5\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+L_{1,2}(\xi )L_{3,2}(\eta )U_6\hspace{0.17em}+L_{2,2}(\xi )L_{3,2}(\eta )U_7+L_{3,2}(\xi )L_{3,2}(\eta )U_8.\end{array}$$

(31)

Using the explicit form of the quadratic Lagrange polynomials, in $\xi$:

$$L_{1,2}(\xi )=\left[(2\xi -1)(\xi -1)\right],\hspace{0.33em}{L}_{2,2}(\xi )=\left[4\xi (1-\xi )\right],\hspace{0.33em}{L}_{3,2}(\xi )=\left[\xi (2\xi -1)\right],$$

(32)

it can be verified that the two sums involved into Equation (31) simplify to:

$$L_{1,2}(\xi )+{\displaystyle \frac{1}{2}}{L}_{2,2}(\xi )=1-\xi =L_{1,1}(\xi ),\hspace{1em}\mathrm{and}\hspace{1em}{L}_{3,2}(\xi )+{\displaystyle \frac{1}{2}}{L}_{2,2}(\xi )=\xi =L_{2,1}(\xi ).$$

(33)

Substituting Equation (33) into Equation (31), one may observe that the basis functions, which constitute Approach 2, are identical with those into Equation (28), which constitutes Approach 1. In other words, whether transfinite interpolation is directly applied, or the constraint is implemented into the initial background tensor product, the results are identical.

Interestingly, Equations (28) and (31) show that the shape function associated with each nodal point equals the local tensor product provided linear interpolation along the bottom edge is assumed.

5.5. Elimination of Internal Nodes

In general, the elimination of an internal nodal point is a task that is more difficult than eliminating a boundary node. For example, let us consider the element in Figure 8a in which node 18 (at point H) is to be eliminated. Whatever approach is followed (i.e., transfinite interpolation or background tensor product), the latter point may be considered to belong to the vertical isoline 19-16 as well as to the horizontal isoline 4-20. Obviously, since the value $U_H$ appears in all three projectors ($P_{\xi },P_{\eta },P_{\xi \eta }$) of the transfinite interpolation, it cannot be eliminated by itself in the Boolean sum (Equation (11)) and thus Approach 1 seems to fail as a first single step; therefore, a different treatment is required by extending Approach 1 through a second step. Alternatively, we can deal with Approach 2.

In both approaches, we begin the interpolation through point 18 and seek a straightforward way to eliminate $U_{18}$. Since node 18 belongs to the horizontal isoline 4-5-20, it is reasonable to interpolate it along this isoline as a nonuniform quadratic polynomial which comes from the elimination of the third node, like the local nodal value $U_3$ in Equation (18), and thus we have:

$${\left(U_{18}\right)}_{\xi }=-{\displaystyle \frac{1}{3}}{U}_4+U_5+{\displaystyle \frac{1}{3}}{U}_{20}.$$

(34)

On the other hand, node 18 also belongs to the vertical isoline 19-18-8-12-16, and thus plays the role of $U_2$ in Equation (22), whence:

$${\left(U_{18}\right)}_{\eta }={\displaystyle \frac{1}{4}}{U}_{19}+{\displaystyle \frac{3}{2}}{U}_8-U_{12}+{\displaystyle \frac{1}{4}}{U}_{16}.$$

(35)

The obvious handling is to take the average of these two constraints, and thus the actual constraint becomes:

$$U_{18}={\displaystyle \frac{{\left(U_{18}\right)}_{\xi }+{\left(U_{18}\right)}_{\eta }}{2}}={\displaystyle \frac{1}{2}}\left[\left(-{\displaystyle \frac{1}{3}}{U}_4+U_5+{\displaystyle \frac{1}{3}}{U}_{20}\right)+\left({\displaystyle \frac{1}{4}}{U}_{19}+{\displaystyle \frac{3}{2}}{U}_8-U_{12}+{\displaystyle \frac{1}{4}}{U}_{16}\right)\right].$$

(36)

It is reminded that Equation (36) is useful for Approach 1 and Approach 2 as well.

6. Construction of the 17-Node and 18-Node T-Mesh Elements

In this section, we will demonstrate the construction of 17-node element (illustrated in Figure 8a) and the 18-node element which would be generated by including node 18 (H) in the formulation.

We recall that having in hand the interpolation of the internal node given by Equation (36), we can follow two alternative but equivalent approaches, i.e., transfinite interpolation or background tensor product formulation. The lowest level of these two approaches were demonstrated in Section 5.4.1 and 5.4.2, respectively.

6.1. Approach 1 to Derive T-Mesh Elements

6.1.1. General Remarks

In Ref. [34], it was shown that in some special cases it is possible to start with the transfinite interpolation formula and eventually lead to a T-spline expression. This is the case when the inserted artificial (auxiliary) nodes, which complete the background (tensor-product) mesh, appear only in two out of the three projectors, ${\{P}_{\xi },P_{\eta }$ and $P_{\xi \eta }\}$, i.e., either in ${\{P}_{\xi },P_{\xi \eta }\}$ or in ${\{P}_{\eta },P_{\xi \eta }\}$, and thus can be eventually eliminated. For example, if the projector $P_{\xi }$ includes all the nodal points of the tensor product background mesh in such a way that $P_{\xi }=P_{\xi \eta }$ (e.g., Figure 4b), then the bivariate variable will be approximated by $U\left(\xi ,\eta \right)=P_{\eta }$, which will not include any artificial DOF. Moreover, in the paper at hand, it was made clear that missing nodes at the boundary are easily treated through a simple linear constraint of the nodes along the edge to which the missing node belongs.

Nevertheless, given a quadrilateral patch in which several ordered nodes $\left({\xi }_i,{\eta }_j\right)$ are interconnected in the form of a T-mesh, we shall see that transfinite interpolation is not always capable of immediately determining the associated shape functions. To make this point clear, below we refer to the 17-node and 18-node elements (Figure 8a) and then in Section 7 we continue with a larger T-mesh.

6.1.2 17- and 18-Node Elements

We consider the T-mesh of seventeen nodal points (1 to 17) shown in Figure 8a, where a rather large gap exists due to three missing points (F, G, H) in the background tensor-product representation (Figure 8b). Nodes 1 to 17 are called primary nodes, whereas the auxiliary points (F, G, H) are called artificial or secondary nodes.

We begin implementing the transfinite interpolation using Lagrange polynomials throughout, and thus all the initial seventeen nodal points, as well as the three artificial ones (F, G, and H created by extension of existing lines), contribute to the formation of the three projectors $(P_{\xi },P_{\eta },P_{\xi \eta })$ involved.

Within this patch, the variable of the problem $U(\xi ,\eta )$ is interpolated by the Boolean sum according to Equation (11). First, we start with the 18-node element and then continue with the 17-node element.

Transfinite Elements Using 18 Nodal Points

Assuming that the blending and trial functions are Lagrange polynomials, the three projectors are given as

$$\begin{array}{l}{P}_{\xi }\left\{U\right\}=E_1(\xi )\left[E_1(\eta )U_1+E_2(\eta )U_4+E_3(\eta )U_6+E_4(\eta )U_{10}+E_5(\eta )U_{14}\right]\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+E_2(\xi )\left[E_1(\eta )U_2+E_2(\eta )U_5+E_3(\eta )U_7+E_4(\eta )U_{11}+E_5(\eta )U_{15}\right]\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+E_3(\xi )\left[E_1(\eta )U_G+E_2(\eta )U_H+E_3(\eta )U_8+E_4(\eta )U_{12}+E_5(\eta )U_{16}\right]\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+E_4(\xi ){\left[L_1^{(3)}(\eta )U_3+L_2^{(3)}(\eta )U_9+L_3^{(3)}(\eta )U_{13}+L_4^{(3)}(\eta )U_{17}\right]}_{BC}\end{array},$$

(37)

$$\begin{array}{l}{P}_{\eta }\left\{U\right\}=E_1(\eta ){\left[L_1^{(2)}(\xi )U_1+L_2^{(2)}(\xi )U_2+L_3^{(2)}(\xi )U_3\right]}_{AB}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+E_2(\eta )\left[E_1(\xi )U_4+E_2(\xi )U_5+E_3(\xi )U_H+E_4(\xi )U_F\right]\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+E_3(\eta )\left[E_1(\xi )U_6+E_2(\xi )U_7+E_3(\xi )U_8+E_4(\xi )U_9\right]\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+E_4(\eta )\left[E_1(\xi )U_{10}+E_2(\xi )U_{11}+E_3(\xi )U_{12}+E_4(\xi )U_{13}\right]\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+E_5(\eta )\left[E_1(\xi )U_{14}+E_2(\xi )U_{15}+E_3(\xi )U_{16}+E_4(\xi )U_{17}\right]\end{array},$$

(38)

and

$$\begin{array}{l}{P}_{\xi \eta }\left\{U\right\}=\left[E_1(\xi )E_1(\eta )\right]U_1+\left[E_2(\xi )E_1(\eta )\right]U_2+\left[E_3(\xi )E_1(\eta )\right]U_G+\left[E_4(\xi )E_1(\eta )\right]U_3\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[E_1(\xi )E_2(\eta )\right]U_4+\left[E_2(\xi )E_2(\eta )\right]U_5+\left[E_3(\xi )E_2(\eta )\right]U_H+\left[E_4(\xi )E_2(\eta )\right]U_F\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[E_1(\xi )E_3(\eta )\right]U_6+\left[E_2(\xi )E_3(\eta )\right]U_7+\left[E_3(\xi )E_3(\eta )\right]U_8+\left[E_4(\xi )E_3(\eta )\right]U_9\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[E_1(\xi )E_4(\eta )\right]U_{10}+\left[E_2(\xi )E_4(\eta )\right]U_{11}+\left[E_3(\xi )E_4(\eta )\right]U_{12}+\left[E_4(\xi )E_4(\eta )\right]U_{13}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left[E_1(\xi )E_5(\eta )\right]U_{14}+\left[E_2(\xi )E_5(\eta )\right]U_{15}+\left[E_3(\xi )E_5(\eta )\right]U_{16}+\left[E_4(\xi )E_5(\eta )\right]U_{17}.\end{array}$$

(39)

In more detail, the interpolation along the third vertical station at $(\xi =2/3)$ utilizes the artificial nodes H and G, whereas along the fourth vertical station measured from left (i.e., the edge BC at $\xi =1$) uses only the four primary nodes (3, 9, 13, 17) and thus nonuniform Lagrange polynomials $L_j^{(3)}(\eta ),j=1,2,3,4$ of degree 3 (i.e., the artificial point F(20) is ignored).

Moreover, the interpolation along the bottom horizontal station utilizes the three primary nodes (1, 2, and 3) and nonuniform Lagrange polynomials $L_i^{(2)}(\xi ),i=1,2,3$ of degree 2 (i.e., the artificial point G(19) is ignored). Furthermore, regarding the second horizontal station measured from the bottom (at $\eta =1/4$), the interpolation is performed in terms of the blending functions on the primary (4, 5) continued by the artificial (H, F) nodes.

Also, the tensor product $P_{\xi \eta }\left\{U\right\}$ includes all the points involved, i.e., the 17 real nodes and the 3 artificial ones.

Substituting Equations (37) to (39) into the Boolean sum described by Equation (11), one may observe that:

• the two artificial values on the boundary ($U_G,U_F$) are eliminated,
• the one in the interior ($U_H$) remains.

Therefore, eighteen shape functions appear in the approximation of the variable

$$U(\xi ,\eta )={\displaystyle \sum _{i=1}^{18}{N}_j(\xi ,\eta )}\hspace{0.17em}{U}_j, \mathrm{with} U_{18}=U_H$$

(40)

with:

$$\begin{array}{l}{N}_1(\xi ,\eta )=L_1^{(2)}(\xi )\cdot E_1(\eta ),\\ N_2(\xi ,\eta )=L_2^{(2)}(\xi )\cdot E_1(\eta ),\\ N_3(\xi ,\eta )=E_4(\xi )L_1^{(3)}(\eta )+L_3^{(2)}(\xi )E_1(\eta )-E_4(\xi )E_1(\eta ),\\ N_4(\xi ,\eta )=E_1(\xi )E_2(\eta ),\\ N_5(\xi ,\eta )=E_2(\xi )E_2(\eta ),\\ N_6(\xi ,\eta )=E_1(\xi )E_3(\eta ),\\ N_7(\xi ,\eta )=E_2(\xi )E_3(\eta ),\\ N_8(\xi ,\eta )=E_3(\xi )E_3(\eta ),\\ N_9(\xi ,\eta )=E_4(\xi )L_2^{(3)}(\eta ),\\ N_{10}(\xi ,\eta )=E_1(\xi )E_4(\eta ),\\ N_{11}(\xi ,\eta )=E_2(\xi )E_4(\eta ),\\ N_{12}(\xi ,\eta )=E_3(\xi )E_4(\eta ),\\ N_{13}(\xi ,\eta )=E_4(\xi )L_3^{(3)}(\eta ),\\ N_{14}(\xi ,\eta )=E_1(\xi )E_5(\eta ),\\ N_{15}(\xi ,\eta )=E_2(\xi )E_5(\eta ),\\ N_{16}(\xi ,\eta )=E_3(\xi )E_5(\eta ),\\ N_{17}(\xi ,\eta )=E_4(\xi )L_4^{(3)}(\eta ),\\ N_{18}(\xi ,\eta )=N_H(\xi ,\eta )=E_3(\xi )E_2(\eta ).\end{array}$$

(41)

It can be verified that the eighteen global shape functions in Equation (41) fulfil the Partition of Unity Property:

$${\displaystyle \sum _{i=1}^{18}{N}_j(\xi ,\eta )}\hspace{0.17em}=1,\hspace{1em}\forall (\xi ,\eta )\in {[0,1]}^2.$$

(42)

The shape functions $(N_j(\xi ,\eta ),j=1,\dots ,18)$ of Equation (41) are illustrated in Figure 9b.

One may also observe that all the shape functions of Equation (41) associated with the first seventeen DOFs vanish at the position $({\xi }_H,{\eta }_H)=(2/3,1/4)$ of point H, i.e., we have $N_i(\frac{2}{3},\frac{1}{4})=0,\hspace{0.33em}i=1,\dots ,17$. The latter means that a nodal point is necessary to be put at the position of H (called node ‘18’), otherwise the function $U(\frac{2}{3},\frac{1}{4})$ would vanish at it and the Partition of Unity Property does not hold. In other words, the direct application of the transfinite interpolation to the element of Figure 8a gives the 18-node element.

Transfinite Elements Using 17 Nodal Points

Despite the limitation of transfinite interpolation to automatically derive the shape functions of the 17-node element, Equation (41) concerning the 18-node element (in which the secondary node H(18) exists) may become a useful start for further processing. In more detail, we interpolate $U_{18}$ using the rest of the nodal values along the horizontal and vertical isolines passing through it, take the mean average and receive again Equation (36). However, the latter is still not applicable because it contains the unknown artificial DOFs $U_{19}$ and$U_{20}$.

Since node 19 is constrained by the node sequence 1-2-3, by virtue of Equation (18), in global numbering we have:

$$U_{19}=-{\displaystyle \frac{1}{3}}{U}_1+U_2+{\displaystyle \frac{1}{3}}{U}_3.$$

(43)

Similarly, since node 20 is constrained by the node sequence 3-9-13-17, by virtue of Equation (22), in global numbering we have:

$$U_{20}={\displaystyle \frac{1}{4}}{U}_3+{\displaystyle \frac{3}{2}}{U}_9-U_{13}+{\displaystyle \frac{1}{4}}{U}_{17}.$$

(44)

Substituting Equations (43) and (44) into Equation (36), we eventually receive the following constraint for internal node 18:

$$U_{18}=-{\displaystyle \frac{U_1}{24}}+{\displaystyle \frac{U_2}{8}}+{\displaystyle \frac{U_3}{12}}-{\displaystyle \frac{U_4}{6}}+{\displaystyle \frac{U_5}{2}}+{\displaystyle \frac{3U_8}{4}}+{\displaystyle \frac{U_9}{4}}-{\displaystyle \frac{U_{12}}{2}}-{\displaystyle \frac{U_{13}}{6}}+{\displaystyle \frac{U_{16}}{8}}+{\displaystyle \frac{U_{17}}{24}}.$$

(45)

One may observe that nodal value $U_{18}$ is expressed in terms of 11 out of the 17 actual nodal values. Therefore, the substitution of Equation (45) into the expansion of the bivariate function $U(\xi ,\eta )$ results in:

$$\begin{array}{l}U(\xi ,\eta )=N_1(\xi ,\eta )U_1+\cdots +N_{17}(\xi ,\eta )U_{17}+N_{18}(\xi ,\eta )U_{18}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=N_1(\xi ,\eta )U_1+\cdots +N_{17}(\xi ,\eta )U_{17}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+N_{18}(\xi ,\eta )\left(-{\displaystyle \frac{U_1}{24}}+{\displaystyle \frac{U_2}{8}}+{\displaystyle \frac{U_3}{12}}-{\displaystyle \frac{U_4}{6}}+{\displaystyle \frac{U_5}{2}}+{\displaystyle \frac{3U_8}{4}}+{\displaystyle \frac{U_9}{4}}-{\displaystyle \frac{U_{12}}{2}}-{\displaystyle \frac{U_{13}}{6}}+{\displaystyle \frac{U_{16}}{8}}+{\displaystyle \frac{U_{17}}{24}}\right).\end{array}$$

(46)

Equation (46) shows that all the 11 DOFs $(U_1,U_2,U_3,U_4,U_5,U_8,U_9,U_{12},U_{13},U_{16},U_{17})$ are affected by $N_{18}(\xi ,\eta )\equiv N_H(\xi ,\eta )=E_3(\xi )E_2(\eta )$ (i.e., the shape function of the artificial node H(18), which is described by the lowest equality in Equation (41)). Therefore, after factorization in Equation (46), the set of shape functions of the 17-node element will be given by:

$$\begin{array}{l}{N}_1(\xi ,\eta )=L_1^{(2)}(\xi )\cdot E_1(\eta )-1/24N_H(\xi ,\eta ),\\ N_2(\xi ,\eta )=L_2^{(2)}(\xi )\cdot E_1(\eta )+1/8N_H(\xi ,\eta ),\\ N_3(\xi ,\eta )=E_4(\xi )L_1^{(3)}(\eta )+L_3^{(2)}(\xi )E_1(\eta )-E_4(\xi )E_1(\eta )+1/12N_H(\xi ,\eta ),\\ N_4(\xi ,\eta )=E_1(\xi )E_2(\eta )-1/6N_H(\xi ,\eta ),\\ N_5(\xi ,\eta )=E_2(\xi )E_2(\eta )+1/2N_H(\xi ,\eta ),\\ N_6(\xi ,\eta )=E_1(\xi )E_3(\eta ),\\ N_7(\xi ,\eta )=E_2(\xi )E_3(\eta ),\\ N_8(\xi ,\eta )=E_3(\xi )E_3(\eta )+3/4N_H(\xi ,\eta ),\\ N_9(\xi ,\eta )=E_4(\xi )L_2^{(3)}(\eta )+1/4N_H(\xi ,\eta ),\\ N_{10}(\xi ,\eta )=E_1(\xi )E_4(\eta ),\\ N_{11}(\xi ,\eta )=E_2(\xi )E_4(\eta ),\\ N_{12}(\xi ,\eta )=E_3(\xi )E_4(\eta )-1/2N_H(\xi ,\eta ),\\ N_{13}(\xi ,\eta )=E_4(\xi )L_3^{(3)}(\eta )-1/6N_H(\xi ,\eta ),\\ N_{14}(\xi ,\eta )=E_1(\xi )E_5(\eta ),\\ N_{15}(\xi ,\eta )=E_2(\xi )E_5(\eta ),\\ N_{16}(\xi ,\eta )=E_3(\xi )E_5(\eta )+1/8N_H(\xi ,\eta ),\\ N_{17}(\xi ,\eta )=E_4(\xi )L_4^{(3)}(\eta )+1/24N_H(\xi ,\eta ),\\ \mathrm{with}\hspace{1em}{N}_H(\xi ,\eta )=E_3(\xi )E_2(\eta ).\end{array}$$

(47)

The shape functions $(N_j(\xi ,\eta ),j=1,\dots ,17)$ of Equation (47), which refer to the 17-node element shown in Figure 8a, are illustrated in Figure 9a.

6.2. Approach 2 to Derive the 17-Node and 18-Node T-Mesh Elements

6.2.1. The 18-Node Element

We start with a background tensor product mesh element, of size 4 × 5 (i.e., 20-node element), shown in Figure 8b, in which the bivariate function is approximated by:

$$U(\xi ,\eta )={\displaystyle \sum _{k=1}^{20}{\varphi }_k(\xi ,\eta )U_k}={\displaystyle \sum _{i=1}^4{\displaystyle \sum _{j=1}^5{L}_{i,3}(\xi )L_{j,4}(\eta ){\tilde{U}}_{ij}}},$$

(48)

where $U_k$ is the $k$-th nodal value whereas ${\tilde{U}}_{ij}$ denotes the same in terms of the associated univariate Lagrange polynomials. It is noted that the set of Lagrange polynomials $L_{i,3}(\xi )$, which are involved in Equation (48), coincides with the set of blending functions $E_i(\xi ),i=1,2,3,4$, whereas the set of Lagrange polynomials $L_{j,4}(\eta )$ coincides with the set of blending functions $E_j(\eta ),j=1,2,3,4,5$, by definition. Thus, we can write $L_{i,3}(\xi )\equiv E_i(\xi )$ and $L_{j,4}(\eta )\equiv E_j(\eta )$, which means that the symbols $L$ and $E$can be used interchangeably.

Regarding the 18-node element, it may be produced from the tensor product mesh, in which the two boundary nodes (G, F) are secondary (not independent) and thus must be eliminated. These nodes simply lie inside the intervals made by the nonuniform nodes 1-2-3 (at $\xi =2/3$) and 3-9-13-17 (at $\eta =1/4$), respectively. Therefore, based on these nonuniform Lagrange polynomials, and applying Equations (18) and (22), we can put the following constraints on the points G and F along the edges AB and BC, respectively:

$$\mathrm{For} \xi ={\displaystyle \frac{2}{3}}:\hspace{1em}{U}_G=-{\displaystyle \frac{1}{3}}{U}_1+U_2+{\displaystyle \frac{1}{3}}{U}_3,$$

(49)

and

$$\mathrm{For} \eta ={\displaystyle \frac{1}{4}}:\hspace{1em}{U}_F={\displaystyle \frac{1}{4}}{U}_3+{\displaystyle \frac{3}{2}}{U}_9-U_{13}+{\displaystyle \frac{1}{4}}{U}_{17}.$$

(50)

Substituting Equations (49) and (50) into Equation (48), and then extracting the common factors, we receive:

$$\begin{array}{l}U(\xi ,\eta )=\stackrel{{\tilde{N}}_1(\xi ,\eta )}{\overbrace{\left[E_1(\xi )E_1(\eta )-{\displaystyle \frac{1}{3}}{E}_3(\xi )E_1(\eta )\right]}}{U}_1+\stackrel{{\tilde{N}}_2(\xi ,\eta )}{\overbrace{\left[E_2(\xi )E_1(\eta )+E_3(\xi )E_1(\eta )\right]}}{U}_2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\stackrel{{\tilde{N}}_3(\xi ,\eta )}{\overbrace{\left[E_4(\xi )E_1(\eta )+{\displaystyle \frac{1}{3}}{E}_3(\xi )E_1(\eta )+{\displaystyle \frac{1}{4}}{E}_4(\xi )E_2(\eta )\right]}}{U}_3+\stackrel{{\tilde{N}}_4(\xi ,\eta )}{\overbrace{\left[E_1(\xi )E_2(\eta )\right]}}{U}_4\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\stackrel{{\tilde{N}}_5(\xi ,\eta )}{\overbrace{\left[E_2(\xi )E_2(\eta )\right]}}{U}_5+\stackrel{{\tilde{N}}_6(\xi ,\eta )}{\overbrace{\left[E_1(\xi )E_3(\eta )\right]}}{U}_6+\stackrel{{\tilde{N}}_7(\xi ,\eta )}{\overbrace{\left[E_2(\xi )E_3(\eta )\right]}}{U}_7+\stackrel{{\tilde{N}}_8(\xi ,\eta )}{\overbrace{\left[E_3(\xi )E_3(\eta )\right]}}{U}_8\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\stackrel{{\tilde{N}}_9(\xi ,\eta )}{\overbrace{\left[E_4(\xi )E_3(\eta )+{\displaystyle \frac{3}{2}}{E}_4(\xi )E_2(\eta )\right]}}{U}_9+\stackrel{{\tilde{N}}_{10}(\xi ,\eta )}{\overbrace{\left[E_1(\xi )E_4(\eta )\right]}}{U}_{10}+\stackrel{{\tilde{N}}_{11}(\xi ,\eta )}{\overbrace{\left[E_2(\xi )E_4(\eta )\right]}}{U}_{11}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\stackrel{{\tilde{N}}_{12}(\xi ,\eta )}{\overbrace{\left[E_3(\xi )E_4(\eta )\right]}}{U}_{12}+\stackrel{{\tilde{N}}_{13}(\xi ,\eta )}{\overbrace{\left[E_4(\xi )E_4(\eta )-E_4(\xi )E_2(\eta )\right]}}{U}_{13}+\stackrel{{\tilde{N}}_{14}(\xi ,\eta )}{\overbrace{\left[E_1(\xi )E_5(\eta )\right]}}{U}_{14}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\stackrel{{\tilde{N}}_{15}(\xi ,\eta )}{\overbrace{\left[E_2(\xi )E_5(\eta )\right]}}{U}_{15}+\stackrel{{\tilde{N}}_{16}(\xi ,\eta )}{\overbrace{\left[E_3(\xi )E_5(\eta )\right]}}{U}_{16}+\stackrel{{\tilde{N}}_{17}(\xi ,\eta )}{\overbrace{\left[E_4(\xi )E_5(\eta )+{\displaystyle \frac{1}{4}}{E}_4(\xi )E_2(\eta )\right]}}{U}_{17}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\stackrel{{\tilde{N}}_{18}(\xi ,\eta )}{\overbrace{\left[E_3(\xi )E_2(\eta )\right]}}{U}_H\hspace{1em}\hspace{1em}\square \end{array}$$

(51)

One may observe in Equation (51) that six out of the total eighteen DOFs, i.e. the nodal values $\{U_1,U_2,U_3,U_9,U_{13},U_{17}\}$ (located along the edges AB and BC), are not associated with a simple tensor product of the blending functions. This fact may be also observed in Equation (41), in which the associated shape functions include the non-uniform Lagrange polynomials $L_i^{(2)}(\xi ),i=1,2,3$ (along the edge AB) and $L_j^{(3)}(\eta ),j=1,2,3,4$ (along the edge BC).

Of course, it is trivial to prove that each of the above non-tensor product shape functions are identical (i.e., $N_i\equiv {\tilde{N}}_i,i=1,\dots .18$), in both Equation (41) and Equation (51). For example, the modified shape function of the nodal point ‘1’ is written as:

$$\begin{array}{l}{\tilde{N}}_1(\xi ,\eta )=E_1(\xi )E_1(\eta )-1/3E_3(\xi )E_1(\eta )=\left[E_1(\xi )-1/3E_3(\xi )\right]\cdot E_1(\eta )\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left[{\displaystyle \frac{\left(\xi -1/3\right)\left(\xi -2/3\right)\left(\xi -1\right)}{\left(0-1/3\right)\left(0-2/3\right)\left(0-1\right)}}-{\displaystyle \frac{1}{3}}\cdot {\displaystyle \frac{\left(\xi -0\right)\left(\xi -1/3\right)\left(\xi -1\right)}{\left(2/3-0\right)\left(2/3-1/3\right)\left(2/3-1\right)}}\right]\cdot E_1(\eta )\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\left[{\displaystyle \frac{\left(\xi -1/3\right)\left(\xi -1\right)}{\left(0-1/3\right)\left(0-1\right)}}\right]\cdot E_1(\eta )\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=L_1^{AB}(\xi )\cdot E_1(\eta )\triangleq L_1^{(2)}(\xi )\cdot E_1(\eta ),\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\mathrm{Q}.\mathrm{E}.\mathrm{D}.\end{array}$$

(52)

and, eventually, becomes identical to the nonuniform Lagrange polynomial $L_1^{(2)}(\xi )$ involved in the first equality of Equation (41).

6.2.2. The 17-Node Element

Substituting Equation (45) into Equation (51), we obtain:

$$\begin{array}{l}{\overline{N}}_1(\xi ,\eta )=\left[E_1(\xi )-{\displaystyle \frac{1}{3}}{E}_3(\xi )\right]E_1(\eta )-{\displaystyle \frac{1}{24}}{E}_3(\xi )E_2(\eta ),\\ {\overline{N}}_2(\xi ,\eta )=\left[E_2(\xi )+E_3(\xi )\right]E_1(\eta )+{\displaystyle \frac{1}{8}}{E}_3(\xi )E_2(\eta ),\\ {\overline{N}}_3(\xi ,\eta )=E_3(\xi )\left[{\displaystyle \frac{1}{3}}{E}_1(\eta )+{\displaystyle \frac{1}{12}}{E}_2(\eta )\right]+E_4(\xi )\left[E_1(\eta )+{\displaystyle \frac{1}{4}}{E}_2(\eta )\right],\\ {\overline{N}}_4(\xi ,\eta )=E_1(\xi )E_2(\eta )-{\displaystyle \frac{1}{6}}{E}_3(\xi )E_2(\eta ),\\ {\overline{N}}_5(\xi ,\eta )=E_2(\xi )E_2(\eta )+{\displaystyle \frac{1}{2}}{E}_3(\xi )E_2(\eta ),\\ {\overline{N}}_6(\xi ,\eta )=E_1(\xi )E_3(\eta ),\\ {\overline{N}}_7(\xi ,\eta )=E_2(\xi )E_3(\eta ),\\ {\overline{N}}_8(\xi ,\eta )={\displaystyle \frac{3}{4}}{E}_3(\xi )E_2(\eta )+E_3(\xi )E_3(\eta ), \\ {\overline{N}}_9(\xi ,\eta )=\left[{\displaystyle \frac{1}{4}}{E}_3(\xi )+{\displaystyle \frac{3}{2}}{E}_4(\xi )\right]E_2(\eta )+E_4(\xi )E_3(\eta ), \\ {\overline{N}}_{10}(\xi ,\eta )=E_1(\xi )E_4(\eta ),\\ {\overline{N}}_{11}(\xi ,\eta )=E_2(\xi )E_4(\eta ),\\ {\overline{N}}_{12}(\xi ,\eta )=E_3(\xi )E_4(\eta )-{\displaystyle \frac{1}{2}}{E}_3(\xi )E_2(\eta ),\\ {\overline{N}}_{13}(\xi ,\eta )=E_4(\xi )\left[E_4(\eta )-E_2(\eta )\right]-{\displaystyle \frac{1}{6}}{E}_3(\xi )E_2(\eta ),\\ {\overline{N}}_{14}(\xi ,\eta )=E_1(\xi )E_5(\eta ),\\ {\overline{N}}_{15}(\xi ,\eta )=E_2(\xi )E_5(\eta ),\\ {\overline{N}}_{16}(\xi ,\eta )={\displaystyle \frac{1}{8}}{E}_3(\xi )E_2(\eta )+E_3(\xi )E_5(\eta ),\\ {\overline{N}}_{17}(\xi ,\eta )=\left[{\displaystyle \frac{1}{24}}{E}_3(\xi )+{\displaystyle \frac{1}{4}}{E}_4(\xi )\right]E_2(\eta )+E_4(\xi )E_5(\eta ).\end{array}$$

(53)

One may observe that all the seventeen shape functions in Equation (53) have been expressed in terms of the blending functions only, and thus the accurate representation of the boundary conditions is not obvious yet. Nevertheless, it can be easily verified that the shape functions of Equation (53), which are associated with boundary nodes, give the same result as that produced by the nonuniform Lagrange polynomials in each edge. For example, the shape function ${\overline{N}}_1(\xi ,\eta )$, which is associated with the corner node A(1) of Figure 8a, was drawn along the mutually vertical edges AD and AB (shown in Figure 8b), and was found to coincide with the polynomials (determined by the uniform node sequence {0,1/4,2/4,3/4,1} and the nonuniform one {0,1/3,1}, respectively). Alternatively, by substituting the analytical expressions of the equal-spaced univariate Lagrange polynomials $E_1(\xi ),E_3(\xi ),E_1(\eta ),$ and $E_2(\eta )$ in the first equality of Equation (53), the reader may prove the same thing.

Moreover, despite their slightly different analytical forms, one may easily show (theoretically and numerically) that each of the 17 equalities in Equation (53) is identical to the corresponding equality in Equation (47).

In conclusion, the desired 17-node element can be derived either using Approach 1 (transfinite interpolation) accompanied by the imposition of a constraint for only the internal node H(18) or following Approach 2 (background tensor product) in which the constraint for internal node 18 is followed by the constraints at the boundary nodes G(19) and F(20).

A MATLAB® code which deals with all the three T-mesh elements (19-, 18-, and 17-node), including the background tensor product, is cited in Appendix A.

7.1. Description of the Stations

Let us consider the index space of the T-mesh shown in Figure 10 (with parameters $0\le \xi ,\eta \le 1$), which is inspired by Ref. [29](p. 244), but was purposely extended for a Lagrange based formulation adding one extra layer along the entire boundary.

One reason for the above choice is that it is clearly written in Ref. [29](p. 243) that this index space does not ensure a-priori the Partition of Unit Property (PUP), whereas another reason is its complexity and thus its treatment is instructive for judgement of the previously presented methodology in practical applications.

One may observe that the $n_p=46$ nodes, numbered from 1 to 46, belong to 9 horizontal and 9 vertical stations, which overall determine a uniform background mesh (tensor product required to form the projector $P_{\xi \eta }$) of 9 × 9 = 81 points in total. Again, all the 81 background nodes need to be considered in the beginning, because all of them contribute to the projector $P_{\xi \eta }$. Therefore, since only 46 out of the 81 total points are primary nodes associated with 46 degrees of freedom, there will be $n_s=35$ secondary nodes as well. Moreover, since 11 out of the 35 secondary nodes belong to the boundary, there will be only 24 secondary nodes in the interior of the element. This in turn means that if Approach 1 is adopted (i.e., transfinite interpolation), we need only 24 constraints, one for each internal secondary node. Alternatively, if Approach 2 is adopted (background tensor product), we need 35 constraints, because the secondary boundary nodes will be equally treated.

In this specific T-element, the complexity is high, because among the total 18 stations only one of them is uniformly meshed (i.e., the third vertical one, measured from left) whereas in all the other 17 stations at least one out of the local 9 background nodes is missing.

In a tensor product with missing nodes, all empty node positions must first be filled in. Regarding the boundary, the artificial nodes (51,61) on bottom edge AB belong to only two projectors, i.e. $P_{\xi }$ (non-uniform polynomials with respect to parameter $\eta$ along isolines 51-43 and 61-44, respectively) and $P_{\xi \eta }$. Nevertheless, the associated shape functions $N_{51}(\xi ,\eta )$ and $N_{61}(\xi ,\eta )$ are not eliminated, because $P_{\xi \eta }$ consists of uniform blending functions whereas $P_{\xi }$ consists of non-uniform polynomials (elimination restricts to the boundary only). Similarly, the artificial nodes of BC (72,75,78) belong to $P_{\eta }$ (non-uniform polynomial with respect to parameter $\xi$ along isolines 69-72, 28-75, and 79-78, respectively) and are not eliminated as well (as previously, elimination restricts to the boundary only). For each completed node, it is a matter of choosing in which isoline direction it will be considered. Here, a question arises:

Is it possible to number the auxiliary secondary nodes in a clever way so as they appear only twice in the Boolean sum?

For example, secondary nodes 47-64-48-49-50 may be chosen to belong to the vertical isoline $\xi =3/8$, where they create a uniform partition of it, and thus this part of the projector $P_{\xi }$ is cancelled by the corresponding part of the projector ${{\rm P}}_{\xi \eta }$. On the other hand, the secondary nodes 52-to-54 present the following difficulty. To have uniformity along the isoline $\xi =4/8$, node 81 would have to be included. However, node 81 is also necessary (along with 64) for the formation of ${{\rm P}}_{\eta }$ along the uniformly partitioned isoline $\eta =5/8$.

From the above discussion it becomes clear that for a general T-meshed patch (T-element) like that of Figure 10, it is not possible to find a sequence of nodal points along all the isolines, which will ensure that secondary nodes belong to only two out of the three projectors. Therefore, no tricks will be adopted, and each station will be described mainly by its corresponding primary nodes. Below, the terms secondary, auxiliary, or artificial, will be equally used.

Let us again consider interpolation using Lagrange polynomials for the entire element. Then, the blending functions of the transfinite interpolation will be uniform Lagrange polynomials of degree 8, in either of the two directions $\xi$ and $\eta$. Regarding the formation of the three projectors, each station must be taken with its individual set of trial functions (based on the primary nodes and the secondary nodes, when there is no support). To facilitate understanding the procedure, the horizontal isolines and the associated auxiliary (secondary) points have been drawn in blue color, whereas the vertical ones in red color. In addition, auxiliary points at intersection have been drawn in green color.

In more detail, regarding the horizontal stations we have the following trial functions ${\tilde{B}}_i(\xi ):$

• Station H1 at $\eta =0$, which is defined by the seven nodes $(U_1,\dots ,U_7)$, is described by a set of 7 nonuniform Lagrange polynomials. All the involved nodes are primary.
• Station H2 at $\eta =1/8$, which is defined by the seven nodes ($U_8,\dots ,U_{14}$), is again described by a set of 7 nonuniform Lagrange polynomials. All the involved nodes are primary.
• Station H3 at $\eta =2/8$, which is defined by the seven nodes $(U_{65},U_{66},U_{15},U_{16},U_{17},U_{18},U_{19})$, is described by a set of 7 nonuniform Lagrange polynomials. Note that the first two nodes $(U_{65},U_{66})$ are secondary and are required to complete the support on the left part of station H3, whereas the auxiliary (background) nodes 53 and 63 are not required because they are between existing supports.
• Station H4 at $\eta =3/8$, which is defined by the seven nodes$(U_{67},U_{68},U_{20},U_{21},U_{22},U_{23},U_{24})$, is described by a set of 7 nonuniform Lagrange polynomials. Note that the first two nodes $(U_{67},U_{68})$ are secondary.
• Station H5 at $\eta =4/8$, which is defined by the seven nodes $(U_{69},U_{70},U_{25},U_{26},U_{27},U_{71},U_{72})$, is described by a set of 7 nonuniform Lagrange polynomials. Note that here, in addition to the first secondary nodes $(U_{69},U_{70})$, additional auxiliary nodes $(U_{71},U_{72})$ are used to complete the support on the right side of the station.
• Station H6 at $\eta =5/8$, which is defined by the seven nodes $(U_{28},U_{29},U_{30},U_{81},U_{73},U_{74},U_{75})$ of which the last four are secondary, is described by a set of 7 nonuniform Lagrange polynomials.
• Station H7 at $\eta =6/8$, which is defined by the nodes $(U_{79},U_{80},U_{31},U_{32},U_{76},U_{77},U_{78})$, of which five are secondary, is described by a set of 7 nonuniform Lagrange polynomials.
• Station H8 at $\eta =7/8$, which is defined by the nodes $(U_{33},U_{34},U_{35},U_{36},U_{37},U_{38},U_{39})$, is described by a set of 7 nonuniform Lagrange polynomials. All the involved nodes are primary.
• Station H9 at $\eta =1$, which is defined by the seven nodes $(U_{40},U_{41},U_{42},U_{43},U_{44},U_{45},U_{46})$, is described by a set of 7 nonuniform Lagrange polynomials. All the involved nodes are primary.
• One may observe that all the above trial functions ${\tilde{B}}_i(\xi ),i=1,\dots ,7$ are nonuniform polynomials of the same degree 7, but this is an accidental event (which does not occur for the vertical stations).
• Regarding the vertical stations we have the following trial functions ${\tilde{B}}_j(\eta ):$
• Station V1 at $\xi =0$, which is defined by the five nodes $(U_1,U_8,U_{28},U_{33},U_{40})$, is described by a set of 5 nonuniform Lagrange polynomials. All the involved nodes are primary.
• Station V2 at $\xi =1/8$, which is defined by the five nodes $(U_2,U_9,U_{29},U_{34},U_{41})$, is described by a set of 5 nonuniform Lagrange polynomials. All the involved nodes are primary.
• Station V3 at $\xi =2/8$, which is defined by the nine nodes $(U_3,U_{10},U_{15},U_{20},U_{25},U_{30},U_{31},U_{35},U_{42})$, is described by a set of 9 uniform Lagrange polynomials. All the involved nodes are primary, and it is the only station which involves the maximum allowable number of nine primary nodes.
• Station V4 at $\xi =3/8$, which is defined by the eight nodes $(U_4,U_{11},U_{16},U_{21},U_{47},U_{48},U_{49},U_{50})$ of which the last four are secondary, is described by a set of 8 nonuniform Lagrange polynomials.
• Station V5 at $\xi =4/8$, which is defined by the eight nodes $(U_{51},U_{52},U_{53},U_{54},U_{26},U_{32},U_{36},U_{43})$ of which the first four are secondary, is described by a set of 8 nonuniform Lagrange polynomials.
• Station V6 at $\xi =5/8$, which is defined by the seven nodes $(U_5,U_{12},U_{17},U_{55},U_{56},U_{59},U_{60})$ of which the last four are secondary, is described by a set of 7 nonuniform Lagrange polynomials. Note the big gap between the secondary nodes $U_{56}$ and $U_{59}$.
• Station V7 at $\xi =6/8$, which is defined by the seven nodes $(U_{61},U_{62},U_{63},U_{22},U_{27},U_{37},U_{44})$ of which the first three are secondary, is described by a set of 7 nonuniform Lagrange polynomials. Note the big gap between $U_{27}$ and $U_{37}$.
• Station V8 at $\xi =7/8$, which is defined by the six nodes $(U_6,U_{13},U_{18},U_{23},U_{38},U_{45})$, is described by a set of 6 nonuniform Lagrange polynomials. Note that all the involved nodes are primary, while there is a big gap between $U_{23}$ and $U_{38}$.
• Station V9 at $\xi =1$, which is defined by the six nodes $(U_7,U_{14},U_{19},U_{24},U_{39},U_{46})$, is described by a set of 6 nonuniform Lagrange polynomials. Note that all the involved nodes are primary, while there is a big gap between $U_{24}$ and $U_{39}$.

Based on the above blending and trial functions, the Boolean sum of Equation (11) was found to include shape functions associated with all the primary nodes 1 to 46 as well as all the artificial nodes (81 nonzero shape functions in total). In other words, for the T-mesh shown in Figure 10, transfinite interpolation is not capable by itself of producing a set of 46 global shape functions possessing the Partition of Unity Property. Therefore, a different approach is required to resolve this problem.

7.2. Imposition of Linear Constraints at Secondary Nodes

As was discussed earlier, we may follow either Approach 1 or Approach 2.

7.2.1. Approach 1:

We apply transfinite interpolation, initially considering all the 81 nodes (primary and secondary). The 11 secondary nodes on the boundary do not contribute because the boundary is interpolated through the primary nodal values only, whereas each secondary node on the boundary is involved in only two projectors and thus are eventually eliminated. Therefore, in this approach we deal with only the remaining 24 secondary nodes in the interior, of which 21 belong to isolines and thus are easily eliminated applying simple linear constraints like those of Equations (15), (18) and (22). Regarding the secondary internal nodes 57, 58 and 64, we apply the mean average procedure. Thus, overall, we must deal with a total of 24 linear constraints in the general form:

$$U_s=S{U}_p,$$

(54)

where $U_s$ is the column vector (of size 35 × 1) which includes the secondary DOFs, $U_p$ is the column vector (of size 46 × 1) which includes the primary DOFs, and $S$ is the transformation matrix (of size 35 × 46) which relates the secondary DOFs in terms of the primary ones. When Equation (54) is considered, the tensor products associated with the secondary DOFs $U_s$ are embodied into the primary DOFs $U_p$, and thus a set of 46 bivariate shape functions is eventually produced.

7.2.2. Approach 2

We blindly start with the complete tensor product of the 81 nodes in a background mesh, without being interested whether a node lies on the boundary or not:

$$U(\xi ,\eta )={\displaystyle \sum _{j=1}^{81}{N}_j(\xi ,\eta )U_j}.$$

(55)

In this initial state, if the nodal value $U_j$ is associated with the intersection of the $k-\mathrm{th}$ row and $l-\mathrm{th}$ column in the T-mesh, the global shape function $N_j(\xi ,\eta )$ is given by the tensor product as

$$N_j(\xi ,\eta )=E_k(\xi )E_l(\eta );\hspace{0.33em}k,l=1,\dots ,9;\hspace{0.33em}j=1,\dots ,81.$$

(56)

As we discussed in the previous section, the univariate trial function at each station was determined in terms of the primary DOFs ($U_i, i=1,\dots ,46$). The next step is to interpolate the auxiliary DOFs ($U_{47}$ to $U_{81}$) in terms of the primary DOFs ($U_1$ to $U_{46}$). To this purpose we distinguish two cases:

• If the auxiliary DOF $U_s$ clearly belongs to a single station ($H_i$ or $V_j$), then we interpolate $U_s$ in terms of the nodal values associated with this specific station. All the auxiliary nodes except of three ones belong to this category (see below).
• If the auxiliary DOF $U_s$ does not clearly belong to a single station but to the intersection of two sections, then we interpolate $U_s$, in both directions, using the mean average of the two values (one for the horizontal and the other for the vertical section). For the configuration of Figure 10, the relevant nodes (for averaging like Equation (36)) are 57, 58 and 64.
• Therefore, for the 46node T-element of this section, we impose 32 simple interpolations along either horizontal or vertical isolines and 3 averaging procedures using both mutually perpendicular isolines passing through the intersected nodes (57, 58 and 64).

Practically, we create a function which calculates the values of all the Lagrange polynomials associated with any given node sequence. Using this function, we identify the section to which a secondary point belongs, and thus we calculate the 35 nodal values $U_j,j=47,\dots ,81$ in terms of the primary ones ($U_i,i=1,\dots ,46$) and may be in terms of some secondary values as well. Using symbol manipulation software, we express the secondary DOFs in terms of the primary ones, and in this way (elimination process) we eventually determine closed-form analytical expressions for all the global shape functions, $N_i(\xi ,\eta ),i=1,\dots ,46$.

• For example, considering the nonuniform Lagrange polynomials based on the nodal points (1,2,3,4,5,6,7) on the bottom edge AB of the quadrilateral patch, the nodal values at the points 51 and 61 are eliminated in terms of the primary nodal values of the same edge by the linear relationships:

  $$\left[\begin{array}{c}{U}_{51}\\ U_{61}\end{array}\right]=\left[\begin{array}{ccccccc}-{\displaystyle \frac{3}{70}}& {\displaystyle \frac{2}{7}}& -{\displaystyle \frac{4}{5}}& {\displaystyle \frac{6}{5}}& {\displaystyle \frac{2}{5}}& -{\displaystyle \frac{2}{35}}& {\displaystyle \frac{1}{70}}\\ {\displaystyle \frac{1}{14}}& -{\displaystyle \frac{3}{7}}& 1& -1& 1& {\displaystyle \frac{3}{7}}& -{\displaystyle \frac{1}{14}}\end{array}\right]\left[\begin{array}{c}{U}_1\\ U_2\\ U_3\\ U_4\\ U_5\\ U_6\\ U_7\end{array}\right]$$

  (57)
• Similar constraints such as those of Equation (57) are obtained for the 17 (i.e., not all the 18) iso-lines and for all the secondary nodes except for (57, 58, and 64) which belong to two isolines simultaneously. As already mentioned, for the latter three secondary nodes, we take the mean average of the two constraints as previously shown in Equation (36).
• Since the secondary nodes are a linear combination of the primary ones (a portion given by Equation (57)), we can eventually find a numerical matrix $\mathbf{S}$ (of size 35 × 46), according to Equation (54). Obviously, when Equation (54) is considered in the tensor product of the 81 terms, this will be eventually expressed in terms of the 46 primary variables only (i.e., the vector ${\mathbf{U}}_p$). The procedure is very similar to that for a smaller mesh (see Subsection 6.2.2), which was described in detail through a MATLAB® computer program in Appendix A.

Implementing the above approach, the 46 shape functions are illustrated in Figure 11. It was validated that they fulfil the Partition of Unity Property. The analytical expressions are given in Appendix B.

8. Numerical Verification

In addition to the above theoretical presentation regarding the construction of the shape functions in T-elements, below we demonstrate the numerical performance of the above elements in the following boundary value problem (BVP) that is characterized by different polynomial degrees in each direction and a non-polynomial exact solution.

Example : Heat-flow in rectangular sheet

Let us consider a square sheet of dimensions ($a=b=1$) in which a Laplace equation governs physics. The boundary conditions are partially Dirichlet and Neumann type, as shown in Figure 12a. The temperature along the top edge is given as

$$U(x,y=b)=U_m\cos \left({\displaystyle \frac{\pi x}{2a}}\right),\hspace{1em}0\le x\le a.$$

(58)

The exact solution is given as:

$$U(x,y)=U_m{\displaystyle \frac{\sinh \left(\frac{\pi y}{2a}\right)}{\sinh \left(\frac{\pi b}{2a}\right)}}\cos \left({\displaystyle \frac{\pi x}{2a}}\right),\hspace{1em}0\le x\le a,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}0\le y\le b.$$

(59)

whereas the error norm $L_2$, is defined as:

$$L_2={\left[{\displaystyle \frac{\int {\left(u_{\mathrm{c}\mathrm{a}\mathrm{l}\mathrm{c}\mathrm{u}\mathrm{l}\mathrm{a}\mathrm{t}\mathrm{e}\mathrm{d}}-u_{\mathrm{e}\mathrm{x}\mathrm{a}\mathrm{c}\mathrm{t}}\right)}^2\mathrm{d}\Omega }{\int {\left(u_{\mathrm{e}\mathrm{x}\mathrm{a}\mathrm{c}\mathrm{t}}\right)}^2\mathrm{d}\Omega }}\right]}^{1/2}\times 100(\%)$$

(60)

Solution: The entire patch ABCD was discretized as a single tensor product element of 20 degrees of freedom (DOF) shown in Figure 12b (with polynomial degrees $p_x=3$ and $p_y=4$), of which only nine DOF are unrestrained. Using Lagrange polynomials, the error norm was found to be $L_2=0.2025\%$.

It is worth mentioning that next to the above tensor product, we further used non-uniform Lagrange polynomials toward the y-direction, with nodal points at $y=0.0, 0.3, 0.6, 0.8 (\mathrm{o}\mathrm{r} 0.85), 1.0$ (i.e., slightly refined close to the top edge DC). It was found that, in the formulation using nonuniform Lagrange polynomials, the L₂-norm was practically insensitive, since it changed with respect to the above uniform solution only after the 13^th decimal point. Despite the non-uniform location of nodal points in the y-direction, the linear mapping was preserved (i.e., $y(\xi ,\eta )=\eta$ and $x(\xi ,\eta )=\xi$) and (closely related) the determinant of the Jacobian was always constant, equal to unity.

For all the elements developed in the paper at hand, which model the entire unit square patch as a single macroelement, the calculated L₂ error norm is shown in Table 1. One may observe that:

• The accuracy of the 8-node transfinite element is acceptable.
• The accuracy of the traditional 21-node transfinite element is excellent.
• The accuracy in the 17-, 18-, and 20-node transfinite elements is acceptable and very similar.
• The additional constraints to come from a 20 to an 18 and 17 node element decrease the accuracy
• The constraints on internal secondary nodes are worse for the accuracy.
• The accuracy of the complicated 46-node T-element is the best of all.

9. Discussion

In computational mechanics, one disadvantage of the tensor product is that it requires the same number and the same relative location of nodes on opposite edges of a large element or patch. Therefore, these nodes must be arranged along discrete isolines. This in turn leads to connectivity problems when two dissimilar patches must be joined together. As a result, the construction of special ‘transition’ elements is imperative. Such elements can be automatically generated using the transfinite interpolation.

Moreover, the existence of high gradients requires local mesh refinement leading to ‘hanging nodes’, but this issue is easily treated through the transfinite interpolation, using a larger number of stations close to the singularity. The hanging nodes are treated by implicitly extending the isolines to which they belong until they reach the boundary. However, the produced nodes are ‘artificial’, in the sense that they are not eventually involved in the final expressions. This happened in all test cases investigated by the authors so far (before the paper at hand), just because the isolines were in no case intersected.

Beyond the above state-of-the-art, in cases where some isolines are intersected at points where no primary node exists, the present paper proposes the substitution of this nodal value at the intersection (associated with a secondary DOF) through a linear constraint which is based on the primary nodal points along these two isolines. Although this procedure would be sufficient and complete this research, the paper continues and shows that, in principle, it would be also possible to start with a background tensor product and then impose linear constraints between secondary and primary DOFs along the boundary and the interior. In this way, since secondary DOFs are related with the primary DOFs through a (generally non-square) matrix (i.e., ${\mathbf{U}}_s=\mathbf{S}{\mathbf{U}}_p$), it is a trivial task to reshape the initial tensor product into a set of global shape functions associated with only the primary DOFs ${\mathbf{U}}_p$. For 17-, 18-, and 19-node elements, the procedure was fully explained in Appendix A.

10. Conclusions

It was shown that in T-elements with mutual vertical isolines which may be intersected at points where no primary node exists, the transfinite interpolation may not directly eliminate the nodal values of primary degrees of freedom. A remedy was proposed, and it consists of the elimination of the secondary degrees of freedom at the intersection of these isolines using the primary nodal values lying along these isolines. Alternatively, it was shown that an equivalent approach is to begin with a tensor product and then successively impose the linear constraints which constitute the elimination of the secondary in terms of the primary degrees of freedom.