Journal of Applied Mathematics

Volume 2013 (2013), Article ID 239703, 10 pages

http://dx.doi.org/10.1155/2013/239703

Research Article

Piecewise Bivariate Hermite Interpolations for Large Sets of Scattered Data

Renzhong Feng1,2 and Yanan Zhang1,2

¹School of Mathematics and Systematic Science, Beijing University of Aeronautics and Astronautics, Beijing 100191, China

²Key Laboratory of Mathematics, Informatics and Behavioral Semantics, Ministry of Education, Beijing 100191, China

Received 12 January 2013; Accepted 13 March 2013

Academic Editor: Ray K. L. Su

Copyright © 2013 Renzhong Feng and Yanan Zhang. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

The requirements for interpolation of scattered data are high accuracy and high efficiency. In this paper, a piecewise bivariate Hermite interpolant satisfying these requirements is proposed. We firstly construct a triangulation mesh using the given scattered point set. Based on this mesh, the computational point () is divided into two types: interior point and exterior point. The value of Hermite interpolation polynomial on a triangle will be used as the approximate value if point () is an interior point, while the value of a Hermite interpolation function with the form of weighted combination will be used if it is an exterior point. Hermite interpolation needs the first-order derivatives of the interpolated function which is not directly given in scatted data, so this paper also gives the approximate derivative at every scatted point using local radial basis function interpolation. And numerical tests indicate that the proposed piecewise bivariate Hermite interpolations are economic and have good approximation capacity.

1. Introduction

The approximation to higher dimensional scattered data is one of the hot and difficult problems in approximation theory field. Due to the characteristics of scattered data, such as the large amount, irregularity, and high dimensionality, it is difficult to construct the approximation methods for them. However, the approximation method to scattered data has been widely applied in many fields, that is, estimating the parameter of input nonlinear system [1–3], solving partial differential equations, surface reconstruction in reverse engineering, data visualization, and so on. The most frequently-used approximation methods include interpolation by spline, interpolation by radial basis function, and the least square approximation. The interpolations by spline have two types, which are global interpolation and local interpolation, respectively. The global spline interpolation [4] is not able to deal with large scale of scattered data, while the local spline interpolation [5–7] needs smooth joining between piece and piece. The radial basis function interpolation [8, 9] requires solving a large scale of linear system, and the least square approximation [10, 11] also requires solving a certain scale of linear system. This is a problem that doesnot easily work out in the computation. It is well known that Hermite interpolation has higher approximation accuracy than interpolation that only interpolates function values. This is because it interpolates not only the function value, but also the derivative value. However, how to construct multivariate Hermite interpolation for large sets of scattered data is an important research consideration which is also a difficult task. At present, some results have been published [12–14]; however, these methods all require solving a certain scale of linear system. Therefore, they are not suitable for approximating to large scale of scattered data.

This paper proposes two piecewise bivariate Hermite interpolation methods for large sets of scattered data. One of them uses exact derivative in Hermite interpolation, and the other one does the approximate derivative. Both methods are economic, free us from solving any linear systems, and have better approximation capacity. Therefore, they are especially suitable for the approximation to large sets of scattered data. These interpolation methods firstly use Delaunay triangulation method [15] to make a triangulation mesh based on the given scattered point set and then divide the computational point into two types of interior point and exterior outer point. If point is an interior point, then the value of the Hermite interpolation function on the triangle which point lies on is regarded as the approximate value of the approximated function, while if point is an exterior point, then the value of a Hermite interpolation with a form of weighted combination is regarded as the approximate value. Hermite interpolation needs the first-order derivatives of the interpolated function; however, the derivative information is not directly given at scatted data, so this paper gives the approximate derivative at every scatted point by using local radial basis function interpolation.

The structure of this paper is organised as follows. The Hermite interpolation on triangle is presented in Section 2. The construction of piecewise bivariate Hermite interpolation is presented in Section 3. The estimation of partial derivative at every scattered point by local radial basis function interpolation is presented in Section 4. Some numerical tests and comparisons between constructed schemes and other interpolation schemes in accuracy and efficiency are carried out in Section 5. This paper ends with a section of brief conclusion.

2. Hermite Interpolation Based on Triangle

Given the values and partial derivatives of a bivariate function at three vertices of a triangle (their vertices 1, 2, 3 are arranged counter to clockwise), that is, ,   (see Figure 1). A bivariate Hermite interpolation polynomial is to be constructed satisfying the following interpolation conditions:

Figure 1: Interpolation triangle .

In order to infer the expression of Hermite interpolation from interpolation conditions (1), the barycenter coordinate of a triangle needs to be introduced first. It is assumed that is a point on triangle and its barycentric coordinates is denoted as . Then the relation between its cartesian coordinates and barycentric coordinates is

where is the area of the triangle and

Making use of transformation (2), the interpolation conditions (1) can be changed into the following form under barycentric coordinates:

A complete bivariate polynomial of three degrees has ten undetermined coefficients. Under barycentric coordinate its expression is

Since the Hermite interpolation conditions (1) or (4) only provide 9 conditions, we construct a Hermite interpolation polynomial with nine coefficients

which meets interpolation condition (4). Then the values of the coefficients in formula (6) can be obtained:

Combining (4), (6), and (7), the Hermite interpolation polynomial satisfying (4) can be written as

where

3. Bivariate Piecewise Hermite Interpolation for Large Sets of Scattered Data

Given a set of scattered data , which are assumed to be sampled from a function,. The partial derivatives of at the set of points are given as . Taking the scattered points as nodes, a triangulation mesh, , is constructed in domain using Delaunay triangulation method. If point lies on a triangle of , it is called as an interior point, otherwise, an exterior point. For example, points in Figure 2 are interior points, while point is an exterior point. Next we will give the approximate value of according to the category of point .

Figure 2: The exhibition of interior and exterior points.

3.1. Interior Point

Suppose that point lies on a triangle whose vertices are . We now construct a Hermite interpolation polynomial as (8) on triangle satisfying the following conditions:

Then let , so that .

3.2. Exterior Point

When point doesnot belong to any triangle within , we construct the following form of Hermite interpolation:

satisfying interpolation conditions

In this paper, we take as a Hermite interpolation polynomial with the same form of (8) on triangle with the least area among all triangles with vertex , satisfying interpolation conditions (10). is called as node basis function and function is called as weighted function. In order to satisfy interpolation conditions (12), these weighted functions are required for satisfying conditions

In the paper we take the weighted functions in (11) as those in the literature [16]:

where

The value of in (17) can be selected by user properly. It can be seen that when , , thereby ; when and , thereby . Thus satisfies condition (13). in formula (16) satisfies condition (14) obviously, and the proof of which can be found in literature [16]. If is determined by (17), then when . So in case of being appropriately selected ( is also called as influence radius of node basis function), there are only several nonzero terms in sum . Then it can be found that interpolation function (11) which takes functions (16) and (17) as its weight functions is a local Hermite interpolation. When applying (11), we can use different influence radius at every node, so we can use the same radius . The following computational formulation of uniform influence radius in (18) is from [16]:

is the number of points used inside a circle of radius .

3.3. Bivariate Piecewise Hermite Interpolation

Summing up Sections 3.1 and 3.2, a bivariate piecewise Hermite interpolation function is constructed on to approximate the function

where is a triangle within with node   as its vertices and is a Hermite interpolation polynomial on triangle .

We now summarize the description of the constructed interpolation function, written as Algorithm 1.

Algorithm 1. Consider the following.

Step  1. Generate a triangle mesh T in domain using Delaunay triangulation method based on given scattered point set .

Step  2. Judge the category of point : interior point or exterior point.

Step  3. If point is an interior point, then find out the triangle which includes point , compute the value of the interpolation function at point using formulas (2), (8), and (10), and take the value as the value of .

Step  4. If point is an exterior point, then compute and select to define . The default value of is set to 9, and it responds to an uniform radius. Find out the points from scattered point set , which belongs to the circle with radius and center , compute the values and , give the sum , and assign it to .

4. The Estimation of Partial Derivative

In the process of using bivariate piecewise Hermite interpolation , it demands the first-order derivatives of the interpolated function at every point, but the scattered data   doesnot provide the derivatives. Thus we use the first-order derivative of a local radial basis interpolation to approximate the derivative of the interpolated function and replace this exact derivative in (19) with the approximate derivative.

The Delaunay triangulation is a convex subdomain of . Every scattered point is either a boundary point or an interior node point of . We put the vertices of all triangles with point as into a set, written as , see Figure 3.

(a)

(b)

Figure 3: The local radial basis interpolation point set around point , (a) is a boundary point of and (b) is an interior node point of .

We take the point set as local radial basis function interpolation point set and multiquadric radial basis functions

as local interpolation basis functions.

By solving the linear system ,

we obtained a local radial basis interpolation function , and then compute the first order derivatives of at , that is,

The bivariate piecewise interpolation function using the previous approximate derivatives is written as

Algorithm 2 summarizes the computational process of interpolation .

Algorithm 2. Consider the following.

Step  1. Generate a triangle mesh in using Delaunay triangulation method based on the given scattered point set .

Step  2. Find out the local radial basis interpolation point set at every point , and compute the approximate derivatives at every point using interpolation formulation (21) and (22).

Step  3. Judge the category of point : interior point or exterior point.

Step  4. If point is an interior point, then find out the triangle which point falls in, compute the value of the interpolation function at point using formulas (2), (8), and (10), and assign the value to the function .

Step  5. If point is an exterior point, then compute and select to define . The default value of is set to 9, it responds to an uniform radius. Find out the points from scattered point set , which belongs to the circle with radius and center , compute the values and , give the sum , and assign it to .

5. Numerical Test

In this section, and are used to approximate Franke function

Firstly, a scattered data set is sampled from ; then, the interpolation functions and are generated using the sampled data. Aiming at different sampled data sets, mean square error (MSE) and maximum error (MME) of the two approximate functions are computed and their approximate capacities are also compared. Meanwhile, the approximate accuracy and computational efficiency of and are compared with Shepard interpolation [17] and Weighted Local RBF interpolant (abbr. WLRBF) [9]. Finally, the CPU time of these four methods is compared which are devised by MATLAB (7.12.0 R2011a) installed in a computer with Processor: Intel(R) Q8400 2.66 GHz, RAM: 4 GB. The computational results show that the methods introduced in this paper need less CPU time and have higher accuracy.

Figure 4 shows the 300 scattered points which are randomly selected from , and Figure 5 is Delaunay triangulation based on Figure 4. Figure 6 presents the graph of Franke function. Figures 7, 8, 9, and 10 describe the interpolation function graphs for the four methods of , , WLRBF, and Shepard which are also based on the 300 scattered points shown in Figure 4. It can be seen that , , and WLRBF can accurately approach Franke function and are better than Shepard method. Figures 11, 12, 13, and 14 display the error graphs of the four methods based on test points. Comparing Figure 11 with Figure 13, we can see that is much better than WLRBF method. Comparing Figures 12 and 13 for the interior points, is obviously better than WLRBF method. Besides, with the number of scattered data increasing, it is found that approximation capability of is closer to that of . In case of large set of scattered points, can be hence absolutely superior to WLRBF.

Figure 4: 300 scattered points randomly selected in .

Figure 5: Delaunay triangulation based on Figure 4.

Figure 6: The graph of Franke function.

Figure 7: The graph of based on Figure 4.

Figure 8: The graph of based on Figure 4.

Figure 9: The graph of WLRBF based on Figure 4.

Figure 10: The graph of Shepard interpolation based on Figure 4.

Figure 11: The error graph of .

Figure 12: The error graph of .

Figure 13: The error graph of WLRBF.

Figure 14: The error graph of Shepard interpolation.

It can be seen from Table 1 that the approximation accuracy of is higher than that of at the same scattered points. Meanwhile, MSE and MME of the two methods decrease with the increase of the point number. Besides, the approximating accuracies to of and get closer at the same time. These results indicate that the local radial function interpolation method can be used to estimate the partial derivative when the number of data points is large enough. Therefore, can be used to approximate when the partial derivative information of the scattered data is unknown.

Table 1: The comparison of MSE and MME of , and , for the same scattered points.

Number of scatter data

Influence radius

MSE
MME
MSE
MME

100
9
0.1818

0.0951

300
9
0.1179

0.0162

500
9
0.0918

0.0043

800
9
0.0705

0.0023

1000
9
0.0654

0.0022

It can be seen from Table 2 that the errors of approximation functions and increase with the radius of influence decreasing.

Table 2: The influence of influence radius on and .

Number of scatter data

Influence radius

MME
MME

300
4
0.0786
0.0489
0.0906

300
9
0.1179
0.0352
0.0747

300
12
0.1361
0.0325
0.0504

300
16
0.1572
0.0307
0.0481

300
25
0.1965
0.0293
0.0447

Tables 1 and 3 indicate that the approximating accuracy to using and is much more accurate than that using WLRBF method and Shepard method with the same number of scattered points.

Table 3: The comparison of MSE and MME of WLRBF method and Shepard method for the same scattered points.

Number of scatter data

Influence radius
WLRBF
Shepard’s Method

MSE
MME
MSE
MME

100
9
0.1818

300
9
0.1179

500
9
0.0918

800
9
0.0705

1000
9
0.0654

The data in Table 4 indicates that the computation efficiencies of and are much higher than those of WLRBF and Shepard methods.

Table 4: The computing time of the four methods on the grid points in the unit square.

Number of scatter data

WLRBF
Shepard’s method

100
5.6
15.46
185.64
52.54

300
28.59
30.35
399.38
113.99

500
42.31
13.77
586.63
180.94

800
64.21
68.04
779.32
225.97

1000
80.28
80.66
1102.64
271.22

6. Conclusion

In this paper, Delaunay triangulation based on planar point set is used to obtain the piecewise bivariate Hermite interpolation function in order to approximate three-dimensional scattered data sets. When point falls on a triangle of the triangulation, the Hermite interpolation on the triangle is used for approximate calculation. If point falls outside the triangulation area, the node basis function value weighted sum of the points that are closer to is used as the approximate value of . Since the constructed interpolation function needs the first derivative value of the approximated function, however, the given scattered data set does not provide such information, our interpolation scheme uses local radial basis interpolation function to estimate the first derivative of each scattered point. Numerical experiments show that our methods have strong approximation ability to scattered data and the estimation of derivatives by local radial basis interpolation has high accuracy. Owing to no demand for solving linear system and the weight functions with local support, our methods are easily implemented and have high computational efficiency, so it is better than B-Spline least square fitting which needs solving a big enough linear system. The use of local radial basis function interpolation to estimate the derivatives is still consuming time. Therefore, how to construct a simple and high-accurate numerical differential formula is one of our future works. Meanwhile, for nonuniform distributed scattered data, how to more reasonably approximate them is another task in further work. In addition, in the process of numerical experiments, we found that the triangulation based on scattered point set has a greater impact on the test results. Hence, how to construct a triangulation more suitable for approximate schemes also needs to be considered in the future work.

Acknowledgment

Supported by the National Natural Science Foundation of China (no. 11271041), National Natural Science Key Foundation of China (no. 10931004) and ISTCP of China grant (no. 2010DFR00700).

References

W. Xiong, W. Fan, and R. Ding, “Least-squares parameter estimation algorithm for a class of input nonlinear systems,” Journal of Applied Mathematics, vol. 2012, Article ID 684074, 14 pages, 2012.

F. Ding, H. Chen, and M. Li, “Multi-innovation least squares identification methods based on the auxiliary model for MISO systems,” Applied Mathematics and Computation, vol. 187, no. 2, pp. 658–668, 2007.

F. Ding, P. X. Liu, and G. Liu, “Multi-innovation least-squares identification for system modeling,” IEEE Transactions on Systems, Man, and Cybernetics B, vol. 40, no. 3, pp. 767–778, 2010.

R. Z. Feng and R. H. Wang, “Closed smooth surface defined from cubic triangular splines,” Journal of Computational Mathematics, vol. 23, no. 1, pp. 67–74, 2005.

R. Z. Feng and R. H. Wang, “Smooth spline surfaces over arbitrary topological triangular meshes,” Journal of Software, vol. 14, no. 4, pp. 830–837, 2003.

S. Lee, Ge. Wolberg, and S. Y. Shin, “Scattered data interpolation with multilevel B-splines,” IEEE Transactions on Visualization and Computer Graphics, vol. 3, no. 3, pp. 228–244, 1997.

W. Z. Xu, L. T. Guan, and Y. X. Xu, “Smoothing of space scattered data by polynomial natural splines,” Acta Scientiarum Naturalium Universitatis Sunyatseni, vol. 49, no. 6, pp. 20–30, 2010.

Z. M. Wu, “Radial basis functions in scattered data interpolation and the meshless method of numerical solution of PDEs,” Chinese Journal of Engineering Mathematics, vol. 19, no. 2, pp. 1–12, 2002.

D. Lazzaro and L. B. Montefusco, “Radial basis functions for the multivariate interpolation of large scattered data sets,” Journal of Computational and Applied Mathematics, vol. 140, no. 1-2, pp. 521–536, 2002.

R. Feng and L. Xu, “Large scattered data fitting based on radial basis functions,” Computer Aided Drafting, Design and Manufacturing, vol. 17, no. 1, pp. 66–72, 2007.

R. Franke and H. Hagen, “Least squares surface approximation using multiquadrics and parametric domain distortion,” Computer Aided Geometric Design, vol. 16, no. 3, pp. 177–196, 1999.

T. Sauer and Y. Xu, “On multivariate Hermite interpolation,” Advances in Computational Mathematics, vol. 4, no. 3, pp. 207–259, 1995.

Z. M. Wu, “Hermite-Birkhoff interpolation of scattered data by radial basis functions,” Approximation Theory and Its Applications, vol. 8, no. 2, pp. 1–10, 1992.

L. Zha and R. Feng, “A scattered hermite interpolation using radial basis Functions,” Journal of Information Computational Science, vol. 4, pp. 361–369, 2007.

B. Delaunay, “Sur la sphère vide,” Izvestia Akademii Nauk SSSR, Otdelenie Matematicheskikh i Estestvennykh Nauk, vol. 7, pp. 793–800, 1934.

R. Franke and G. Nielson, “Smooth interpolation of large sets of scattered data,” International Journal for Numerical Methods in Engineering, vol. 15, no. 11, pp. 1691–1704, 1980.

D. Shepard, “A two-dimensional interpolation function for irregularly spaced data,” in Proceedings of the 23rd ACM National Conference, pp. 517–524, ACM, 1968.