Generating bicubic B-spline surfaces by a sixth order PDE

1.   Introduction

PDE based modelling approach is especially suitable for satisfying surface boundary constraints (^[1]). Using PDE to design surfaces can be traced back to the method of Bloor and Wilson (^[2]), which is based upon the notion that a surface can be created as the solution to a PDE. In ^[3,4,5], Bloor and Wilson indicate the ways to create surfaces satisfying a wide range of functional requirements by a suitable choice of PDE and boundary conditions. Furthermore, they describe a method whereby PDE surfaces may be obtained in closed form, even for the case of general boundary conditions (^[6]). In ^[7], Ugail et al. show how surfaces of practical significance can be constructed interactively in real time.

The minimal surfaces can be related to the harmonic equation. In ^[8,9,10], Cosin and Monterde obtain polynomial approximations to minimal Bézier surfaces. Xu and Wang present the problems related to sextic and quintic parametric polynomial minimal surfaces in ^[11] and ^[12]. Xu et al. present the quasi-harmonic Bézier surfaces to approximate the minimal surfaces in ^[13].

The biharmonic equation can be used to generate surfaces with minimal blending energy and design the tangent-continuous surfaces. In ^[14], Monterde and Ugail present a method to generate Bézier surfaces based on a general 4th-order PDE. On the existence of biharmonic tensor-product Bézier surface patches are presented in ^[15].

B-spline curves and surfaces, especially bicubic B-spline surfaces are widespread throughout computer aided design (CAD), CAGD and numerical solutions of PDE. In ^[16], a method of approximating the surface in terms of B-spline is proposed. In ^[17], the bicubic B-spline surfaces with fixed boundary conditions by minimizing functions are constructed. In ^[18,19,20], cubic B-splines are widely applied to solve PDE problems. The optimization method of the quadratic functional in a bicubic splines functional space is presented in ^[21].

The method of mask is used to construct a surface that satisfies a certain PDE, which makes the formation of the surface simple and stable. In ^[22], Farin and Hansford formulate discrete Coons patches. Representation of piecewise biharmonic surfaces using biquadratic B-splines is proposed by Han and Han in ^[23]. A main result is that they construct biquadratic B-spline surfaces by applying biharmonic PDE. Afterwards, they present the non-uniform bicubic B-spline surfaces subjected to the biharmonic PDE in ^[24].

The sixth order PDE is generally formulated as

where $a, b, c, d, e, f, g \in R$ are constants, whose corresponding quadratic functional is

where $\Omega$ is the parametric domain of $\mathit{\boldsymbol{x}}(u, v)$. Particularly, the triharmonic functional $L_{1, 0, 3, 0, 3, 0, 1}$ related to the Laplacian gradient energy of surface $\mathit{\boldsymbol{x}}(u, v)$ have found their way into various areas of application, such as surface design and geometric mesh smoothing (^[25]). In fluid mechanics, the triharmonic equation is used to describe two-dimensional slowly rotating highly viscous fluid flow in small cavities (^[26]). In ^[27], Kapl and Vitrih present a framework for solving the triharmonic equation over bilinearly parameterized planar multi-patch domains by means of isogeometric analysis (IGA). Wu and Zhu construct the tensor product Bézier surfaces which satisfies the triharmonic equation in ^[28]. The result is that the internal control points of the resulting triharmonic Bézier surface can be obtained uniquely by the given boundary condition.

In ^[29,30], the general sixth order PDEs are proposed to solve the problems of $C^{2}$ continuity in fair surface modeling and surface blending. A framework for real-time freeform modeling based on the sixth order PDE is proposed in ^[31]. In ^[32], Zhang and Xu apply a general sixth order geometric PDE from minimizing a curvature integral functional to solve several surface modeling problems. In ^[33], a general formulation of sixth order geometric flows is proposed.

The sixth order Eq (1.1) provides enough degrees of freedom not only to accommodate tangent, but also curvature boundary conditions, and offers more shape control parameters to serve as user controls for the manipulation of surface shapes. Furthermore, the main control of the surface is to gain through the manipulation of the boundary conditions, which makes it have a good application in interactive design. The use of sixth order surfaces would enable surfaces to be generated and controlled using a set of control points, where the interior of the resulting surface can be controlled by the boundary condition.

It is worth noting that if the surface $\mathit{\boldsymbol{x}}(u, v)$ is bicubic B-spline and $d\neq 0$, then the general sixth order PDE (1.1) is reduced to the particular PDE

The aim of this paper is to analyze the theory and present effective algorithms for generating surfaces from the sixth order Eq (1.3). As we know, high order continuity surfaces can be used to efficiently solve problems of surface processing and modeling. And bicubic B-spline surfaces are a very useful tool in surface modeling and scientific computing. If we apply the operator (1.3) to a bicubic B-spline surface, only one of terms is nonzero, and the others are always zeros. This makes the rest of the parameters no impact and then leads the problem easy to solve. If we apply (1.3) to B-spline surfaces of higher degrees, then the corresponding linear systems will be very complex. In this paper, inspired by Han's methods (^[23,24]) and the idea from ^[28], we apply the sixth order operator to generate bicubic B-spline surfaces with the uniform and the quasi-uniform knots, and the algorithms are also provided. We provide the mask of $4\times4$, which is a generalization of Farin's mask (^[22]), that is, each equation in the linear system is related to 16 control points.

The paper is organised as follows. In the Section 2, we recall some notations related to bicubic B-spline surfaces, and then we study the general sixth order PDE operator acting on every patch of bicubic B-spline to form a mask. In the Section 3, bicubic B-spline surfaces satisfying sixth order PDE are constructed with the uniform knots and the quasi-uniform knots respectively. In the Section 4, we provide some examples for construction of bicubic B-spline surfaces based on sixth order PDE to verify our methods.

2.   Formulations of bicubic B-spline surfaces and sixth order PDEs

2.1.   Notations related to bicubic B-spline surfaces

A bicubic B-spline surface $\mathit{\boldsymbol{x}}(u, v)$ with the knot vectors $U = (u_{0}, u_{1}, \cdots, u_{n+4})$ and $V = (v_{0}, v_{1}, \cdots, v_{n+4})$ is defined as (^[1])

where $N_{i}(u)$ and $N_{j}(v)$ are cubic B-spline basis functions defined over $U$ and $V$ respectively, $\{\mathit{\boldsymbol{P}}_{i, j}\}_{i, j = 0}^{n}\in \mathbb{R}^{3}$ denote the control points. The basis functions $N_{i}(u)$, $i = 0, 1, 2, \cdots, n$, defined over $U$ can be represented as follow (^[1])

And the basis functions $N_{j}(v)$, $j = 1, 2, \cdots, n$, defined over $V$ can be represented similarly.

For the convenience, we make some notations for knot vector $U$, which is similar to the one used by Han in ^[23]. When the knots do not coincide with each other, we introduce the symbols,

Thereupon, from (2.2), we get

In the same way, we set $\tilde{\alpha}_{j}, \tilde{\beta}_{j}, \tilde{\gamma}_{j}$ and $\tilde{\lambda}_{j}$ for the knot vector $V$ and obtain $N_{j}^{(3)}(v)$.

2.2.   Bicubic B-spline solutions to the general sixth order PDEs

We consider a sixth order PDE operator (1.3) to act on the bicubic B-spline surface (2.1). For a patch of $\mathit{\boldsymbol{x}}(u, v)$ on $[u_{i}, u_{i+1}]\times [v_{j}, v_{j+1}]$, we can apply the condition Eq (1.3) and get

Therefore, according to (2.4), we get

In order to simplify the Eq (2.5), let us remark that

where

Hence, according to these notations, Eq (2.5) holds if and only if

where $t_{i, k+3-i}\tilde{t}_{j, l+3-j}$ ($k = i-3, \cdots, i$, $l = j-3, \cdots, j$) are the coefficients associated with the knot vectors $U$ and $V$.

Considering the scheme with the control points as follow,

According to Eq (2.9), the coefficients of mask in (2.10) are obtained,

Specially, if the knot vectors for $U$ and $V$ are uniform ($\lambda = u_{i+1}-u_{i}$, $\tilde{\lambda} = v_{j+1}-v_{j}$), then, the coefficients in (2.11) can be rewritten in the below form,

3.   Construction of bicubic B-spline surfaces satisfying sixth order PDEs

The goal is to give the boundary control points whereby the choice of the remaining control points will be done by a linear system. In this paper, the boundary condition is

Let's consider applying the mask (2.11) to the boundary condition (3.1), and then the linear system related to the unknown inner control points as $'\ast'$ shown in (3.1) can be obtained by

where

and

where $H = \{h_{i, j}\}$ and $\tilde{E} = \{e_{i, j}\}$ are diagonal matrices.

The solution to the linear system (3.2) which subjects to boundary condition (3.1) and the knot vectors for $U$ and $V$. Moreover, the coefficient matrix $T$ of the corresponding linear system can be divided into the form of the Kronecker product of two matrices. we substitute $T_{(n-2)^{2}}(n\geq 3)$ for $T$ as

where

and

where the Kronecker product of matrices is

for any $K = (k_{i, j})_{m\times n}$, $L = (l_{i, j})_{r\times s}$, which satisfy the following conditions,

Note that it is difficult to prove the invertibility of the above matrix $M_{n-2}$ or $N_{n-2}$. Since the B-spline functions and B-spline surfaces with the uniform or the quasi-uniform knots offer more advantages, such as simple structure, interpolation property of boundary and corners and easy implementation of geometric algorithms, so that they are the most commonly used tools in CAGD and other application fields. In this paper, we will prove the invertibility of the above matrices for the uniform and the quasi-uniform knot vectors, and present the algorithms of surface construction respectively.

3.1.   Uniform knots

Suppose bicubic B-spline surface defined over knot vectors are uniform, then we make

where $\lambda = u_{i+1}-u_{i}, \tilde{\lambda} = v_{j+1}-v_{j}, i, j = 0, 1, \cdots, n+3$.

Theorem 1. If

is a bicubic B-spline surface satisfying the sixth order PDE (1.3) with control net $\{\mathit{\boldsymbol{P}}_{i, j}\}_{i, j = 0}^{n}$ for the uniform knot vectors $U$ and $V$, then all the inner control points $\{\mathit{\boldsymbol{P}}_{i, j}\}_{i, j = 2}^{n-1}$ are determined by the remaining control points.

Proof. Since $U$ and $V$ are the uniform knot vectors, the mask (2.12) can be moved by the boundary condition (3.1), then we can get a linear system related to the internal control points. The coefficient matrix (3.6) is simplified as

where $\bigotimes$ is the Kronecker product, and $S_{n-2}$ is a quad-diagonal matrix

Now, we define $L_{n}$ as the determinant of the matrix $S_{n}$. Then, recurrence relation is

and the initial values are

Therefore, we get

Obviously, matrix $S_{n-2}$ is invertible. Furthermore, the matrix $T_{(n-2)^{2}}$ is non-singular. So that all the inner control points $\{\mathit{\boldsymbol{P}}_{i, j}\}_{i, j = 2}^{n-1}$ are obtained from the linear system (3.2).

3.2.   Quasi-uniform knots

The knot vector $U$ of $p$ degree B-spline surfaces is called quasi-uniform if we add $p+1$ multiple knots at the endpoints $u_{0} = u_{1} = \cdots = u_{p}, u_{n+1} = u_{n+2} = \cdots = u_{n+p+1}$ and $u_{i+1}-u_{i} = \lambda, i = p, \cdots, n.$ This type of B-spline surfaces is also called clamped B-spline surfaces. Since the clamped B-spline surface interpolates four corner control points and has better boundary property, it is widely used in practical applications.

Suppose bicubic B-spline surface defined over knot vectors are quasi-uniform, for $p = 3$, we make

where

and

Thus, (2.3) can be rewritten as follows

Theorem 2. If

is a bicubic B-spline surface satisfying the sixth order PDE (1.3) with control net $\{\mathit{\boldsymbol{P}}_{i, j}\}_{i, j = 0}^{n}$ for the quasi-uniform knot vectors $U$ and $V$, then all the inner control points $\{\mathit{\boldsymbol{P}}_{i, j}\}_{i, j = 2}^{n-1}$ are determined by the remaining control points.

The proof is similar to the proof of Theorem 1, and then we will omit it here. Next, we will show how to generate the bicubic B-spline surface satisfying the sixth order PDE (1.3) for different $n$ with the quasi-uniform knots from given boundary condition (3.1).

Let us define a set of vectors for $U$,

where $\{\alpha_{i}\}_{i = 3}^{n}$, $\{\beta_{i}\}_{i = 2}^{n}$ and $\{\gamma_{i}\}_{i = 1}^{n}$ can be taken from (3.10), thus,

Now the model can be obtained by Eq (3.11) as

Similarly, $\tilde{\alpha}, \tilde{\beta}, \tilde{\gamma}$ and $\tilde{\lambda}$ for the knot vector $V$ can be obtained in the same way.

Now, let's consider that the surface (2.1) satisfying the sixth order PDE over the quasi-uniform knots, then we get the following results.

(1) For $n = 3$, $(u, v)\in [u_{3}, u_{4}]\times [v_{3}, v_{4}]$. According to (3.10) and (3.11), we get

there is the unique mask (2.11) satisfying the condition of sixth order PDE (1.3) exactly, and the coefficient matrix is

This mask is suitable for boundary condition (3.1) with $n = 3$, so the unknown point $\mathit{\boldsymbol{P}}_{2, 2}$ can be obtained.

(2) For $n = 4$, $(u, v)\in [u_{i}, u_{i+1}]\times [v_{j}, v_{j+1}]$ $(i, j = 3, 4)$, there are four different masks moving on the control points. According to and (3.10) and (3.11), we get

and then the linear system related to the internal unknown control points is obtained. By simple calculation, the coefficient matrix of the linear system is given as

where

Obviously, by the boundary condition (3.1) with $n = 4$, the all the unknown points can be obtained.

(3) For $n\geq 5$, we will give some notations in order to write $a_{i}, b_{i}, c_{i}, d_{i}$ in matrix forms,

where $a_{i}, b_{i}, c_{i}, d_{i}$ can be taken by (2.8) and (3.10). Then, we get below matrices

where $A = \{a_{i, j}\}_{i, j = 3}^{n}$, $B = \{b_{i, j}\}_{i, j = 3}^{n}$, $C = \{c_{i, j}\}_{i, j = 3}^{n}$ and $D = \{d_{i, j}\}_{i, j = 3}^{n}$.

Given boundary condition (3.1) with $n\geq3$, the bicubic B-spline surface (2.1) satisfying the sixth order PDE with the quasi-uniform knots $U$ and $V$ can be constructed by above analysis, and the algorithm is listed as below.

4.   Examples

In this section, we will illustrate some examples for constructing bicubic B-spline surfaces satisfying the sixth order PDEs from different boundary conditions. By applying the general sixth order PDE operator, the result is that the surfaces are curvature continuous. Furthermore, the asymmetry boundary conditions can be an advantage for solving blending problems (see Example 4 and Example 5).

Example 1. From the boundary condition (3.1), a bicubic B-spline surface satisfying the sixth order PDE (1.3) of $n = 5$ with the uniform knots is constructed (see Figure 1). First, the boundary condition is given (see Figure 1(a)), and then all the control points are calculated according to Algorithm 1 (see the red color dots in Figure 1(b)). Finally, the resulting surface is constructed in Figure 1(c). Obviously, we can adjust the shape of the surface by the given control points.

Example 2. The known boundary control points (see Figure 2(a)) are taken from points on the catenoid

whose parameter intervals are $[0, \pi]$ and $[0, acosh(2)]$. By Algorithm 1, the all the control points and the resulting surface with the uniform knots has been taken as Figure 2(b) and Figure 2(c) shows.

Example 3. The known boundary control points (see Figure 3(a)) are taken from points on the surface

where $(u, v)\in[0, \pi]\times[0, 2\pi]$. The control points and the resulting surface with the quasi-uniform knots are obtained by Algorithm 2 (see Figure 3(b), Figure 3(c) and Figure 3(d)).

Surface blending is a useful technique in CAGD and mesh reconstruction. According to condition (3.1), by the boundary interpolation property with the quasi-uniform knots, we can blend two bicubic B-spline surfaces along common boundary or four bicubic B-spline surfaces shared one point, if the boundary condition of them are selected properly. Since it is easy to implement, we will not show this type of examples here. As we know, for given four boundary curves only, we can't generate a single bicubic B-spline surface satisfying the sixth order PDE (1.3) to interpolate these boundary curves uniquely. Based on our proposed method (Algorithm 2), we can construct four piecewise bicubic B-spline surfaces to satisfy the sixth order PDEs with $C^{1}/G^{1}$ continuity to interpolate the given boundary curves. Furthermore, for given closed boundary curves, we can generate blending bicubic B-spline surfaces to satisfy the sixth order PDEs with $C^{1}/G^{1}$ continuity to interpolate whose boundary curves in the same way.

Example 4. In this example, we construct blending surface to satisfy the sixth order PDE (1.3) with $C^{1}$ continuity interpolating given boundary curves. The idea of the method is: We split every boundary B-spline curve into two parts by knot insertion firstly (see Figure 4(a)). Then we assume two curves with one corner is the boundary of one patch of the resulting surfaces and we choose the other boundary condition according (3.1) for every of these surfaces, such that they satisfy $C^{1}$ continuity (see Figure 4(b)). Next, we can obtain all the control points of every part by Algorithm 2 (see Figure 4(c)). Finally, four piecewise bicubic B-spline surfaces satisfying a sixth order PDE are generated, which satisfy $C^{1}$ continuity and interpolate the initial given boundary curves (see Figure 4(e)).

Next, we will explain all the figure in detail. Figure 4(a) shows the initial given boundary curves, in which every boundary has been split into to two parts. We get boundary condition (in blue) and chosen boundary condition (in black) for one part for as Figure 4(b) shown, and then all the unknown control points for this part with the quasi-uniform knots are obtained in Figure 4(c) (in red) by Algorithm 2. And then, we get the remaining three sets of control points (see Figure 4(d)) similarly. The blending bicubic B-spline surfaces are constructed in Figure 4(e). This means the resulting four piecewise bicubic B-spline surfaces interpolate the four given boundary curves in Figure 4(a), and satisfy the sixth order PDE and $C^{1}$ continuity.

Figure 4(f) shows the relations of all the control points in a grid of different color dots on the plane. The blue, cyan, yellow and green dots represents the control points of the initial given boundary curves. The black dots are selected by the $C^{1}$ conditions of every two connected parts. The red dots are unknown control points, which can be obtained uniquely by Algorithm 2.

Example 5. Similar to the above example, we construct a closed blending surface from the initial given boundary closed curves in Figure 5(a). We split every boundary B-spline curve into two parts by knot insertion (seeFigure 5(a)), and then choose boundary conditions for these two parts by $C^{1}$ continuity. Figure 5(b) shows the initial boundary condition and chosen boundary condition for one part. All the control points for this part with the quasi-uniform knots are obtained in Figure 5(c) by Algorithm 2. We get the other set of control points (see Figure 5(c)) similarly. Figure 5(d) and Figure 5(e) show the resulting piecewise bicubic B-spline surface satisfying the sixth order PDE with $C^{1}$ continuity and interpolating initial boundary curves. From this example, we find that the asymmetry of boundary conditions brings more changes and flexibility to geometric modelings.

5.   Conclusion

In this paper, we propose a method to generate bicubic B-spline surfaces by a class of sixth order PDE (1.3) with boundary conditions (3.1). The proposed method not only generalizes the triharmonic PDE surfaces, but also deepen the theory of the polynomial solutions for PDEs. Furthermore, if we choose some control points properly, then the $C^{1}/G^{1}$ piecewise bicubic B-spline surfaces satisfying the sixth order PDE (1.3) and the given boundary curves can be constructed by our method (see Figure 4(e) and Figure 5(e)). Compared with Bézier surface case, the surface based on bicubic B-spline is more stable and more suitable interactive design.

Note that the selection of some control points in Example 4 (see black dots in Figure 4(f)) and Example 5 can be fulfilled by energy, shape, or other constrains, together with $C^{1}/G^{1}$ continuity condition. This leads to solve a constrained optimization problem and we will consider this problem in the future.

Acknowledgments

The authors appreciate the valuable comments and suggestions from the anonymous reviewers, which improve the clarity of the paper. This work is partly supported by National Natural Science Foundation of China (Nos. 12071057, 11671068).

Conflict of interest

There are no conflicts of interest.