Applied Mathematics and Computation 217 (2011) 4644–4653
Contents lists available at ScienceDirect
Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc
A weighted bivariate blending rational interpolation based on function values Wei-Xian Huang, Guo-Jin Wang ⇑ Department of Mathematics, Zhejiang University, Hangzhou 310027, China State Key Laboratory of CAD and CG, Zhejiang University, Hangzhou 310027, China
a r t i c l e
i n f o
Keywords: Computer aided geometric design Rational spline interpolation C1-continuous Bivariate interpolation Cubic spline
a b s t r a c t This paper presents a new weighted bivariate blending rational spline interpolation based on function values. This spline interpolation has the following advantages: firstly, it can modify the shape of the interpolating surface by changing the parameters under the condition that the values of the interpolating nodes are fixed; secondly, the interpolating function is C1-continuous for any positive parameters; thirdly, the interpolating function has a simple and explicit mathematical representation; fourthly, the interpolating function only depends on the values of the function being interpolated, so the computation is simple. In addition, this paper discusses some properties of the interpolating function, such as the bases of the interpolating function, the matrix representation, the bounded property, the error between the interpolating function and the function being interpolated. Ó 2010 Elsevier Inc. All rights reserved.
1. Introduction Spline interpolation is a useful tool for designing curves and surfaces in computer aided geometric design, hence in the past decade, there has been considerable interest in spline interpolation. For instance, literature [1] gives a method that interpolating the spatial points by a uniform cubic B-spline, and literatures [2,3] realize the smooth interpolation for surfaces. But this kind of spline interpolation has one disadvantage that given the interpolating data, the interpolating function is unique. In order to overcome this disadvantage, in recent years the rational spline, especially the rational cubic spline with linear and quadratic denominators, and its application to shape control, have received attention in literatures [4–10]. This type of interpolation has the advantage that it can modify the shape of curves and surfaces by changing the parameters for the given interpolating data. Literature [4] presents a bivariate rational interpolation based on function values with linear denominators. Literature [5] presents a univariate rational cubic interpolation based on function values with linear denominators. Literature [6] gives a bivariate rational interpolation based on function values and derivatives with linear denominators. Literature [7] discusses the bounded property and the shape control of a bivariate rational interpolating surface. Literature [8] combines the methods in [4,6], and gets a new spline interpolation with more parameters, so it is more effective to modify the shape of curves and surfaces. Also literatures [9,10] discusses the rational interpolation with quadratic denominators and properties of it. The bivariate rational interpolation in literatures [4,6–8], in order to maintain smoothness of the interpolating surface, the parameters of y-direction must be limited, hence it is inconvenient to modify the shape of surfaces. For this reason, this paper presents a new weighted bivariate rational bicubic spline interpolation based on function values, which is C1-continuous for any positive parameters. Besides, this kind of interpolation only depends on the values of the function to be ⇑ Corresponding author at: Institute of Computer Images and Graphics, Zhejiang University, Hangzhou 310027, China. E-mail address:
[email protected] (G.-J. Wang). 0096-3003/$ - see front matter Ó 2010 Elsevier Inc. All rights reserved. doi:10.1016/j.amc.2010.11.016
4645
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
interpolated, without partial derivatives, which is convenient in the engineering. In addition, the paper discusses some properties of this type of interpolation, such as the bases, the matrix representation, the bounded property, the error, and gives some examples. The rest of this paper is organized as follows. In the next section we present a new weighted bivariate rational interpolation based on function values. In Section 3, we present the bases of the interpolation and the matrix representation of the interpolating function. In Section 4, we discuss the bounded property of this type of interpolation. In Section 5, we present two methods for the error estimation of the interpolation. Finally, in Section 6, we present some examples of this kind of interpolation and conclude in Section 7. 2. Interpolation Let X : ½a; b; c; d be the plane region, and fðxi ; yi ; fi;j Þ; i ¼ 1; 2; . . . ; n; n þ 1; j ¼ 1; 2; . . . ; m; m þ 1g be a given set of data points, where a ¼ x1 < x2 < < xn < xnþ1 ¼ b; c ¼ y1 < y2 < < ym < ymþ1 ¼ d are the knot spacings. Let yy i hi ¼ xiþ1 xi ; lj ¼ yjþ1 yj , and for any point ðx; yÞ 2 ½xi ; xiþ1 ; yj ; yjþ1 in the xy-plane, let h ¼ xx ; g ¼ lj j . First, for each hi y ¼ yj ; j ¼ 1; 2; . . . ; m þ 1, we construct the x-direct interpolating curve [4]:
pi;j ðxÞ ¼ p ðxÞj½xi ;xiþ1 ¼
ð1 hÞ3 ai;j fi;j þ hð1 hÞ2 V i;j þ h2 ð1 hÞW i;j þ h3 fiþ1;j ; ð1 hÞai;j þ h
i ¼ 1; 2; . . . ; n 1;
where V i;j ¼ ðai;j þ 1Þfi;j þ ai;j fiþ1;j ; W i;j ¼ ðai;j þ 2Þfiþ1;j hi Diþ1;j , with ai;j > 0 and Di;j ¼ rational cubic interpolation based on function values which satisfies:
pi;j ðxi Þ ¼ fi;j ;
pi;j ðxiþ1 Þ ¼ fiþ1;j ;
@pi;j ðxi Þ @x
¼ Di;j ;
@pi;j ðxiþ1 Þ @x
fiþ1;j fi;j . hi
ð1Þ
This interpolation is called the
¼ Diþ1;j :
Clearly, the interpolating function pi;j ðxÞ on ½xi ; xiþ1 is unique for the given data ðxr ; f ðxr ; yj ÞÞ; r ¼ i; i þ 1; i þ 2 and the parameter ai;j . Based on the x-direct interpolation pi;j ðxÞ, we now construct the bivariate C1-continuous rational interpolation as follows:
p1i;j ðx; yÞ ¼ p1 ðx; yÞj½xi ;xiþ1 ;yj ;yjþ1 ¼ ð1 gÞ3 pi;j ðxÞ þ gð1 gÞ2 V i;j þ g2 ð1 gÞW i;j þ g3 pi;jþ1 ðxÞ; i ¼ 1; 2; . . . ; n 1;
j ¼ 1; 2; . . . ; m 1;
ð2Þ
where V i;j ¼ 2pi;j ðxÞ þ pi;jþ1 ðxÞ; W i;j ¼ 3pi;jþ1 ðxÞ lj Di;jþ1 ðxÞ, with Di;j ðxÞ ¼ ðpi;jþ1 ðxÞ pi;j ðxÞÞ=lj . p1i;j ðx; yÞ is called the bivariate rational bicubic interpolation based on function values which satisfies:
p1i;j ðxr ; ys Þ ¼ f ðxr ; ys Þ ¼ fr;s ;
@p1i;j ðx; ys Þ @y
¼ Di;s ðxÞ;
r ¼ i; i þ 1;
s ¼ j; j þ 1:
Obviously, the interpolating function p1i;j ðx; yÞ on ½xi ; xiþ1 ; yj ; yjþ1 is unique ðxr ; ys ; f ðxr ; ys ÞÞ; r ¼ i; i þ 1; i þ 2; s ¼ j; j þ 1; j þ 2 and the parameters ai;s ; s ¼ j; j þ 1; j þ 2.
for
the
given
data
Lemma 1. For any positive parameter ai;j ; ði ¼ 1; 2 . . . ; n 1; j ¼ 1; 2; . . . ; m 1Þ, the interpolating function p1 ðx; yÞ on the region ½x1 ; xn ; y1 ; yn is C1-continuous. Proof. By Theorem 1 in literature [4], obviously, Lemma 1 holds. h In the above we firstly construct x-direct interpolation. Similarly, for each x ¼ xi ; i ¼ 1; 2; . . . ; n þ 1, we can construct ydirect interpolation as follows:
qi;j ðyÞ ¼ q ðyÞj½yj ;yjþ1 ¼
ð1 gÞ3 bi;j fi;j þ gð1 gÞ2 T i;j þ g2 ð1 gÞZ i;j þ g3 fi;jþ1 ; ð1 gÞbi;j þ g
j ¼ 1; 2; . . . ; m 1:
ð3Þ
where T i;j ¼ ðbi;j þ 1Þfi;j þ bi;j fi;jþ1 ; Z i;j ¼ ðbi;j þ 2Þfi;jþ1 lj di;jþ1 , with di;j ¼ ðfi;jþ1 fi;j Þ=lj , and bi;j > 0. Then we construct the bivariate rational interpolation as follows:
p2i;j ðx; yÞ ¼ p2 ðx; yÞj½xi ;xiþ1 ;yj ;yjþ1 ¼ ð1 hÞ3 qi;j ðyÞ þ hð1 hÞ2 T i;j þ h2 ð1 hÞZ i;j þ h3 qiþ1;j ðyÞ; i ¼ 1; 2; . . . ; n 1;
j ¼ 1; 2; . . . ; m 1;
ð4Þ
where T i;j ¼ 2qi;j ðyÞ þ qiþ1;j ðyÞ; Z i;j ¼ 3qiþ1;j ðyÞ hidiþ1;j ðyÞ, with di;j ðyÞ ¼ ðqiþ1;j ðyÞ qi;j ðyÞÞ=hi . p2i;j ðx; yÞ is also called the bivariate rational bicubic interpolation based on function values which satisfies:
p2i;j ðxr ; ys Þ ¼ f ðxr ; ys Þ ¼ fr;s ;
@p2i;j ðxr ; yÞ @x
¼ dr;j ðyÞ;
r ¼ i; i þ 1;
s ¼ j; j þ 1:
4646
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
Clearly, the interpolating function p2i;j ðx; yÞ on ½xi ; xiþ1 ; yj ; yjþ1 is unique ys ÞÞ; r ¼ i; i þ 1; i þ 2; s ¼ j; j þ 1; j þ 2 and the parameters br;j ; r ¼ i; i þ 1; i þ 2. For k 2 ½0; 1, let
pi;j ðx; yÞ ¼ pðx; yÞj½xi ;xiþ1 ;yj ;yjþ1 ¼ kp1i;j ðx; yÞ þ ð1 kÞp2i;j ðx; yÞ;
i ¼ 1; 2; . . . ; n 1;
for
the
given
j ¼ 1; 2; . . . ; m 1:
data
ðxr ; ys ; f ðxr ;
ð5Þ
Then p(x, y) is called the weighted bivariate rational bicubic interpolation which satisfies:
@pi;j ðxr ; ys Þ ¼ Dr;s ¼ di;j ðys Þ; @x
pi;j ðxr ; ys Þ ¼ fr;s ;
@pi;j ðxr ; ys Þ ¼ dr;s ¼ Di;j ðxr Þ; @y
ð6Þ
for r ¼ i; i þ 1; s ¼ j; j þ 1. Hence by Lemma 1 we can obtain the following theorem: Theorem 1. For given ai;j > 0; bi;j > 0; i ¼ 1; 2; . . . ; n; n 1; j ¼ 1; 2; . . . ; m; m 1, and ðxi ; yi ; fi;j Þ; i ¼ 1; 2; . . . ; n; n þ 1; j ¼ 1; 2; . . . ; m; m þ 1, and k 2 ½0; 1, the interpolating function p(x, y) which is defined by Eq. (5) and satisfies Eq. (6) is C1-continuous on the region ½x1 ; xn ; y1 ; yn . 3. The bases and the matrix representation Let 1;1
1;1
bi;j ¼ bi;j ðai;j ; hÞ ¼
1;2
ð1 hÞ2 ðai;j þ hÞ ; ð1 hÞai;j þ h
1;2
biþ1;j ¼ biþ1;j ðai;j ; hÞ ¼
1;3
1;3
biþ2;j ¼ biþ2;j ðai;j ; hÞ ¼
hð1 hÞai;j þ 2h2 h3 þ h2 ð1 hÞ hhiþ1i ð1 hÞai;j þ h h2 ð1 hÞ hhiþ1i ð1 hÞai;j þ h
;
;
2 1;1 c1;1 i;j ¼ ci;j ðgÞ ¼ ð1 gÞ ð1 þ gÞ;
1;2 2 2 c1;2 i;jþ1 ¼ c i;jþ1 ðgÞ ¼ gð1 þ g g Þ þ g ð1 gÞ
1;3 2 c1;3 i;jþ2 ¼ c i;jþ2 ðgÞ ¼ g ð1 gÞ
2;1
2;1
bi;j ¼ bi;j ðbi;j ; gÞ ¼
2;2
2;2
2;3
2;3
lj ; ljþ1
lj ; ljþ1
ð1 gÞ2 ðbi;j þ gÞ ; ð1 gÞbi;j þ g
bi;jþ1 ¼ bi;jþ1 ðbi;j ; gÞ ¼
bi;jþ2 ¼ bi;jþ2 ðbi;j ; gÞ ¼
lj gð1 gÞbi;j þ 2g2 g3 þ g2 ð1 gÞ ljþ1 ; ð1 gÞbi;j þ g
g2 ð1 gÞ l
lj jþ1
ð1 gÞbi;j þ g
;
2 2;1 c2;1 i;j ¼ ci;j ðhÞ ¼ ð1 hÞ ð1 þ hÞ;
2;2 2 2 c2;2 iþ1;j ¼ c iþ1;j ðhÞ ¼ hð1 þ h h Þ þ h ð1 hÞ
2;3 2 c2;3 iþ2;j ¼ c iþ2;j ðhÞ ¼ h ð1 hÞ
hi ; hiþ1
hi ; hiþ1
then the interpolating function pi;j ðx; yÞ defined by Eq. (5) can be written as follows:
pi;j ðx; yÞ ¼
jþ2 iþ2 X X r¼i
s¼j
ar;s fr;s ;
i ¼ 1; 2; . . . ; n 1;
j ¼ 1; 2; . . . ; m 1;
ð7Þ
4647
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
where 1;1
2;1
2;1 ai;j ¼ bi;j ðai;j ; hÞc1;1 i;j ðgÞk þ bi;j ðbi;j ; gÞc i;j ðhÞð1 kÞ; 1;2
2;1
1;3
2;1
2;2 aiþ1;j ¼ biþ1;j ðai;j ; hÞc1;1 i;j ðgÞk þ biþ1;j ðbiþ1;j ; gÞc iþ1;j ðhÞð1 kÞ; 2;3 aiþ2;j ¼ biþ2;j ðai;j ; hÞc1;1 i;j ðgÞk þ biþ2;j ðbiþ2;j ; gÞc iþ2;j ðhÞð1 kÞ; 1;1
2;2
1;2 2;1 ðgÞk þ bi;jþ1 ðbi;j ; gÞci;j ðhÞð1 kÞ; ai;jþ1 ¼ bi;jþ1 ðai;jþ1 ; hÞci;jþ1 1;2
2;2
1;3
2;2
2;2 aiþ1;jþ1 ¼ biþ1;jþ1 ðai;jþ1 ; hÞc1;2 i;jþ1 ðgÞk þ biþ1;jþ1 ðbiþ1;j ; gÞc iþ1;j ðhÞð1 kÞ; 2;3 aiþ2;jþ1 ¼ biþ2;jþ1 ðai;jþ1 ; hÞc1;2 i;jþ1 ðgÞk þ biþ2;jþ1 ðbiþ2;j ; gÞc iþ2;j ðhÞð1 kÞ; 1;1
2;3
1;3 2;1 ðgÞk þ bi;jþ2 ðbi;j ; gÞci;j ðhÞð1 kÞ; ai;jþ2 ¼ bi;jþ2 ðai;jþ2 ; hÞci;jþ2 1;2
2;3
1;3
2;3
2;2 aiþ1;jþ2 ¼ biþ1;jþ2 ðai;jþ2 ; hÞc1;3 i;jþ2 ðgÞk þ biþ1;jþ2 ðbiþ1;j ; gÞc iþ1;j ðhÞð1 kÞ; 2;3 aiþ2;jþ2 ¼ biþ2;jþ2 ðai;jþ2 ; hÞc1;3 i;jþ2 ðgÞk þ biþ2;jþ2 ðbiþ2;j ; gÞc iþ2;j ðhÞð1 kÞ:
far;s ; r ¼ i; i þ 1; i þ 2; s ¼ j; j þ 1; j þ 2g is called the bases of the interpolating function which satisfies: jþ2 iþ2 X X r¼i
ar;s ¼ 1:
ð8Þ
s¼j
Here suppose that ai;j ¼ ai;jþ1 ¼ ai;jþ2 ; bi;j ¼ biþ1;j ¼ biþ2;j , then it is easy to obtain the matrix representation for the interpolating function. Let
0
fi;j B F i;j ¼ @ fiþ1;j fiþ2;j B1i;j ¼ B2i;j ¼
fi;jþ1 fiþ1;jþ1 fiþ2;jþ1
1;1
biþ1;j
2;1
bi;jþ1
bi;j bi;j
fi;jþ2
1;2
2;2
1
C fiþ1;jþ2 A; fiþ2;jþ2
1;3 biþ2;j ; 2;3
bi;jþ2
T
C 1i;j ¼ ;
C 2i;j ¼
c1;1 i;j
c1;2 i;jþ1
2;1 ci;j
c2;2 iþ1;j
c1;3 i;jþ2
T
;
2;3 ciþ2;j ;
then the interpolation pi;j ðx; yÞ defined by Eq. (5) can be represented by matrices as follows:
pi;j ðx; yÞ ¼ kB1i;j F i;j C 1i;j þ ð1 kÞC 2i;j F i;j B2i;j ;
i ¼ 1; 2; . . . ; n 1;
j ¼ 1; 2; . . . ; m 1:
ð9Þ
4. The bounded property of the interpolation r¼iþ2;s¼jþ2 Let M ¼ maxr¼i;s¼j jfr;s j, we can obtain the following theorem:
Theorem 2. For any positive parameters ai;s ; br;j ; s ¼ j; j þ 1; j þ 2; r ¼ i; i þ 1; i þ 2, the upper boundary of the interpolating function defined by Eq. (5) on ½xi ; xiþ1 ; yj ; yjþ1 is determined by the following inequality:
2 2 2 2 !1=2 8lj hi 8hi lj jpi;j ðx; yÞj 6 M 1þ 1þ þ 1þ 1þ ; 27ljþ1 2hiþ1 27hiþ1 2ljþ1 ¼ 1; 2; . . . ; m 1:
i ¼ 1; 2; . . . ; n 1;
j ð10Þ
Especially when hi ¼ hiþ1 ; lj ¼ ljþ1 , Then
pffiffiffi 35 2M : jpi;j ðx; yÞj 6 18
ð11Þ
4648
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
Proof. By the bases representation (Eq. (7)) of the interpolation pi;j ðx; yÞ and k 2 ½0; 1; h 2 ½0; 1; g 2 ½0; 1, it follows that:
jþ2 jþ2 jþ2 X X iþ2 X iþ2 X iþ2 X X jpi;j ðx; yÞj ¼ ar;s fr;s 6 jar;s fr;s j 6 M jar;s j; r¼i s¼j r¼i s¼j r¼i s¼j 1;riþ1
By the definitions of c1;sjþ1 and br;s i;s
; r ¼ i; i þ 1; i þ 2; s ¼ j; j þ 1; j þ 2, it is easy to obtain that:
jþ2 X lj 8lj 1;sjþ1 61þ ; ci;s ¼ 1 þ 2g2 ð1 gÞ l 27l jþ1 jþ1 s¼j
iþ2 X 1;riþ1 br;s ¼1þ r¼i
2h2 ð1 hÞhi hi 61þ : ½ð1 hÞai;s þ hhiþ1 2hiþ1
Similarly, iþ2 X hi 8hi 2;riþ1 61þ ; cr;j ¼ 1 þ 2h2 ð1 hÞ h 27h iþ1 iþ1 r¼i
jþ2 X 2;sjþ1 br;s ¼1þ s¼j
2g2 ð1 gÞlj lj 61þ ; ½ð1 gÞbr;j þ gljþ1 2ljþ1
for r ¼ i; i þ 1; i þ 2; s ¼ j; j þ 1; j þ 2. Then jþ2 iþ2 X X r¼i
s¼j
8lj hi 8hi lj 1þ þ ð1 kÞ 1 þ 1þ jar;s j 6 k 1 þ 27ljþ1 2hiþ1 27hiþ1 2ljþ1 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2ffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 8lj hi 8hi lj 2 2 1þ þ 1þ 1þ 1þ 6 k þ ð1 kÞ 27ljþ1 2hiþ1 27hiþ1 2ljþ1 s ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 2 2ffi 8lj hi 8hi lj : 1þ þ 1þ 1þ 6 1þ 27ljþ1 2hiþ1 27hiþ1 2ljþ1
Hence Eq. (10) holds, and Eq. (11) can be derived directly from Eq. (10).
h
5. The error of the interpolation Here we give two methods for the error estimation of the interpolation. The first is based on the Taylor series, and the second is based on the Peano–Kernel. 5.1. The error estimation based on Taylor series Suppose ðx; yÞ 2 ½xi ; xiþ1 ; yj ; yjþ1 , and f ðx; yÞ 2 C 1 is the function to be interpolated, by Taylor series we have:
@ @ f ðn1 ; n2 Þ; f ðx; yÞ ¼ f ðxr ; ys Þ þ ðx xr Þ þ ðy ys Þ @x @y
r ¼ i; i þ 1; i þ 2;
s ¼ j; j þ 1; j þ 2;
where n1 2 ½xi ; xiþ2 ; n2 2 ½yj ; yjþ2 . Let
@f ¼ max @f ðx; yÞ; @x x2½xi ;xiþ2 @x Then
@f ¼ max @f ðx; yÞ: @y y2½yj ;yjþ2 @y
@f @f jf ðx; yÞ f ðxr ; ys Þj 6 ðhi þ hiþ1 Þ þ ðlj þ ljþ1 Þ ; @x @y
r ¼ i; i þ 1; i þ 2;
s ¼ j; j þ 1; j þ 2:
Suppose pi;j ðx; yÞ is defined by Eq. (5), and it can be rewritten by Eq. (7). By Eq. (8) it follows that:
f ðx; yÞ pi;j ðx; yÞ ¼
jþ2 iþ2 X X r¼i
ar;s ðf ðx; yÞ fr;s Þ;
s¼j
hence
iþ2 jþ2 @f @f X X jf ðx; yÞ pi;j ðx; yÞj 6 ðhi þ hiþ1 Þ þ ðlj þ ljþ1 Þ jar;s j: @x @y r¼i s¼j According to Theorem 2 and the proof of it, the following theorem can be obtained:
ð12Þ
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
4649
Theorem 3. Suppose pi;j ðx; yÞ is defined by Eq. (5) , and f ðx; yÞ 2 C 1 is the function to be interpolated, for any positive parameters ai;s ; br;j ; s ¼ j; j þ 1; j þ 2; r ¼ i; i þ 1; i þ 2, the error of the interpolation on ½xi ; xiþ1 ; yj ; yjþ1 is determined by the following inequality:
@f @f jf ðx; yÞ pi;j ðx; yÞj 6 ðhi þ hiþ1 Þ þ ðlj þ ljþ1 Þ @x @y " 2 2 2 2 #1=2 8lj hi 8hi lj 1þ 1þ þ 1þ 1þ ; 27ljþ1 2hiþ1 27hiþ1 2ljþ1 i ¼ 1; 2; . . . ; n 1;
j ¼ 1; 2; . . . ; m 1:
ð13Þ
Especially when hi ¼ hiþ1 ; lj ¼ ljþ1 ,
jf ðx; yÞ pi;j ðx; yÞj 6
pffiffiffi @f @f 35 2 hi þ lj : @x @y 9
ð14Þ
5.2. The error estimation based on the Peano–Kernel Suppose ðx; yÞ 2 ½xi ; xiþ1 ; yj ; yjþ1 , and f ðx; yÞ 2 C 2 is the function to be interpolated, pi;j ðx; yÞ is defined by Eq. (5). Here let
@p
@f ðx;yÞ
¼ maxx2½x ;x @pi;j ðx;yÞ;
@f ¼ maxx2½x ;x . @x
i
iþ1
@x
@x
i
@x
iþ1
According to Theorem 4.1.1 in literature [11], we have
jf ðx; yÞ pi;j ðx; yÞj 6 jf ðx; yÞ f ðxi ; yÞj þ jf ðxi ; yÞ pi;j ðxi ; yÞj þ jpi;j ðxi ; yÞ pi;j ðx; yÞj Z yjþ2 2 @f @p @f @p @ f ðxi ; sÞ R ½ðy s Þ d s 6 hi þ þ jf ðxi ; yÞ pi;j ðxi ; yÞj ¼ hi þ þ y þ 2 yj @x @x @x @x @y Z 2 y jþ2 @f @p @ f ðxi ; yÞ 6 hi þ þ jRy ½ðy sÞþ jds; @y2 yj @x @x
2
2
where Ry ½ðy sÞþ ¼ kR1y ½ðy sÞþ þ ð1 kÞR2y ½ðy sÞþ ; @ f@yðx2i ;yÞ ¼ maxy2½yj ;yjþ2 @ f@yðx2i ;yÞ, with
8 8 1 > > ðy sÞ c1;2 ðyjþ1 sÞ þ c1;3 ðyjþ2 sÞ ; yj < s < y; > i;jþ1 i;jþ2 > > < r 1 ðsÞ; yj < s < y; < 1 1;2 1;3 Ry ½ðy sÞþ ¼ c ðyjþ1 sÞ þ c ðyjþ2 sÞ ; ¼ r 12 ðsÞ; y < s < yjþ1 ; y < s < yjþ1 ; i;jþ1 i;jþ2 > > : 1 > > > r 3 ðsÞ; yjþ1 < s < yjþ2 ; : c1;3 ðy sÞ; y < s < y ; jþ1 jþ2 i;jþ2 jþ2 8 2;2 2;3 8 2 > > > ðy sÞ bi;jþ1 ðyjþ1 sÞ þ bi;jþ2 ðyjþ2 sÞ ; yj < s < y; > > < < r 1 ðsÞ; yj < s < y; 2;3 r 2 ðsÞ; y < s < yjþ1 ; R2y ½ðy sÞþ ¼ b2;2 ¼ y < s < yjþ1 ; i;jþ1 ðyjþ1 sÞ þ bi;jþ2 ðyjþ2 sÞ ; > > 2 > : > > r 23 ðsÞ; yjþ1 < s < yjþ2 : : b2;3 ðy sÞ; yjþ1 < s < yjþ2 ; i;jþ2 jþ2 2g½ð1gÞbi;j þg 2 i;j þ2gg
Let wðg; bi;j Þ ¼ ð1gÞb
Z
y
yi
þ
gðg2 2g1Þ½ð1gÞbi;j þ2gg2 þ 4g. ð1gÞbi;j þg
By simple integral calculation, we can obtain that
h
i h i 1 o kr ðsÞ þ ð1 kÞr 2 ðsÞds ¼ g þ kc1;2 þ ð1 kÞb2;2 ð2 gÞ þ kc1;3 þ ð1 kÞb2;3 ð2 gÞ þ 2ljþ1 =lj gl2 =2 i;jþ1 i;jþ2 j 1 1 i;jþ1 i;jþ2 n
2
¼ B1 ðg; k; bi;j ; ljþ1 =lj Þlj : ð15Þ Z y
yjþ1
1 kr ðsÞ þ ð1 kÞr 2 ðsÞds 6 2 2
Z y
yjþ1
jkr 12 ðsÞjds þ
Z
yjþ1 y
jð1 kÞr22 ðsÞjds
2 gð1 6g2 þ 6g3 þ 2g4 4g5 þ g6 Þ 1;3 2ð1 2g þ g3 Þljþ1 lj ¼k þ c þ ð1 i;jþ2 ð1 þ g g2 Þ ð1 þ g g2 Þlj 2 " # 2 2 2 2;3 2ð1 gÞ bi;j 2g ð1 gÞ ljþ1 lj kÞ wðg; bi;j Þ þ bi;jþ2 lj 2 ð1 gÞbi;j þ 2g g2 2
¼ B2 ðg; k; bi;j ; ljþ1 =lj Þlj ;
ð16Þ
4650
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
Z
yjþ2
yjþ1
2;3 2 2 jkr 13 ðsÞ þ ð1 kÞr 23 ðsÞjds ¼ kc1;3 i;jþ2 ð1 kÞbi;jþ2 ljþ1 =2 ¼ B3 ðg; k; bi;j ; ljþ1 =lj Þlj ;
ð17Þ
Then the following theorem can be derived: Theorem 4. Suppose pi;j ðx; yÞ is defined by Eq. (5) , and f ðx; yÞ 2 C 2 is the function to be interpolated, for any positive parameters ai;s ; br;j ; s ¼ j; j þ 1; j þ 2; r ¼ i; i þ 1; i þ 2, the error of the interpolation on the subregion ½xi ; xiþ1 ; yj ; yjþ1 is determined by the following inequality:
2 @f @p @ f ðxi ; yÞ x 2 jf ðx; yÞ pi;j ðx; yÞj 6 hi þ þ Bi;j lj ; 2 @y @x @x
i ¼ 1; 2; . . . ; n 1;
j ¼ 1; 2; . . . ; m 1;
ð18Þ
where Bxi;j ¼ max0
2
2
@ f ðx ;yÞ
@ f ðx ;yÞ Similarly, let @yiþ1
¼ maxy2½yj ;yjþ2 @yiþ1 , we have 2 2 Theorem 5. Suppose pi;j ðx; yÞ is defined by Eq. (5) , and f ðx; yÞ 2 C 2 is the function to be interpolated, for any positive parameters ai;s ; br;j ; s ¼ j; j þ 1; j þ 2; r ¼ i; i þ 1; i þ 2, the error of the interpolation on the subregion ½xi ; xiþ1 ; yj ; yjþ1 is determined by the following inequality:
2 @f @p @ f ðxiþ1 ; yÞ x 2 jf ðx; yÞ pi;j ðx; yÞj 6 hi þ þ Biþ1;j lj ; 2 @x @x @y
i ¼ 1; 2; . . . ; n 1;
j ¼ 1; 2; . . . ; m 1;
ð19Þ
where Bxiþ1;j ¼ max0
2
2
@p ðx;yÞ @f
@f ðx;yÞ
@ f ðx;y Þ
@ f ðx;y Þ ¼ maxy2½yj ;yjþ1 i;j@y ; @y Let @p
¼ maxy2½yj ;yjþ1 @y and @x2 s ¼ maxx2½xi ;xiþ2 @x2 s ; s ¼ j; j þ 1. Then by the symme@y
try, we can obtain the following theorem: Theorem 6. Suppose pi;j ðx; yÞ is defined by Eq. (5) , and f ðx; yÞ 2 C 2 is the function to be interpolated, for any positive parameters
ai;s ; br;j ; s ¼ j; j þ 1; j þ 2; r ¼ i; i þ 1; i þ 2, the error of the interpolation on the subregion ½xi ; xiþ1 ; yj ; yjþ1 is determined by the following inequality:
2 @f @p @ f ðx; ys Þ y 2 jf ðx; yÞ pi;j ðx; yÞj 6 lj þ þ B h ; @x2 i;s i @y @y
s ¼ j; j þ 1;
i ¼ 1; 2; . . . ; n 1;
¼ 1; 2; . . . ; m 1;
j ð20Þ
Byi;s
where ¼ max0
B3 ðg; k; bi;j ; 1Þ 6 0:14124;
ð21Þ
B2 ðg; k; bi;j ; 1Þ 6 0:27386;
ð22Þ
B1 ðg; k; bi;j ; 1Þ 6 0:33799;
ð23Þ
Proof
" B3 ðg; k; bi;j ; 1Þ ¼ 1=2 kðg2 g3 Þ þ ð1 kÞ
#
g2 g3 g2 g3 6 1=2 kðg2 g3 Þ þ ð1 kÞ ð1 gÞbi;j þ g g
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiqffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffi ¼ 1=2ðg g2 Þðgk þ ð1 kÞÞ 6 1=2ðg g2 Þ 1 þ g2 k2 þ ð1 kÞ2 6 1=2ðg g2 Þ 1 þ g2 6 0:14124; wðg; bi;j Þ 6
4g 7g2 þ 4g3 g4 ; 2g
therefore it is easy to know that:
4651
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
gð1 6g2 þ 6g3 þ 2g4 4g5 þ g6 Þ 2ð1 2g þ g3 Þðg2 þ g3 Þ þ 1 þ g g2 1 þ g g2 ! 2 2 3 4 2 4g 7g þ 4g g 2g ð1 gÞ þ 1=2ð1 kÞ þ 2g 2g 3 5 6 1 2g þ 2g 2g þ g 4 5g þ g3 ¼ 1=2 kg þ ð1 kÞ g 1 þ g g2 2g
B2 ðg; k; bi;j ; 1Þ 6 1=2k
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 2g þ 2g3 2g5 þ g6 Þ2 ð4 5g þ g3 Þ2 2 2 6 1=2 k þ ð1 kÞ g þ ð1 þ g g2 Þ2 ð2 gÞ2 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ð1 2g þ 2g3 2g5 þ g6 Þ2 ð4 5g þ g3 Þ2 6 0:27386; þ 6 ð1=2Þg ð1 þ g g2 Þ2 ð2 gÞ2
B1 ðg; k; bi;j ; 1Þ 6 1=2½g2 þ kgððg þ 2g2 2g3 Þð2 gÞ þ ðg2 þ g3 Þð4 gÞÞ þ ð1 kÞgð3g 2g2 Þð2 gÞ ¼ 1=2½g2 þ kgð2g g2 g3 þ g4 Þ þ ð1 kÞgð3g 2g2 Þð2 gÞ qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 6 1=2 g2 þ k2 þ ð1 kÞ2 g2 ð2 g g2 þ g3 Þ2 þ ð2g2 7g þ 6Þ2 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 6 1=2 g2 þ g2 ð2 g g2 þ g3 Þ2 þ ð2g2 7g þ 6Þ2 6 0:33799:
Remark 1. All the maximum values are calculated in Maple 12. Remark 2. According to the result of Theorem 7, the error estimation of the interpolation defined by Eq. (18) is stable. Similarly, the error estimations defined by Eqs. (19) and (20) are also stable.
6. Numerical examples qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Suppose the interpolated function is f ðx; yÞ ¼ 1 ð1 xÞ2 ð1 yÞ2 ; ðx; yÞ 2 ½0:5; 1:5; 0:5; 1:5, let hi ¼ lj ¼ 1=5; i ¼ 1; 2; 3; 4; 5; j ¼ 1; 2; 3; 4; 5, then xi ¼ 0:5 þ ði 1Þ=5; yj ¼ 0:5 þ ðj 1Þ=5; i ¼ 1; 2; 3; 4; 5; 6; j ¼ 1; 2; 3; 4; 5; 6. Also let ai;j ¼ 0:5þ ði 1Þ=2 þ ðj 1Þ=2; bi;j ¼ 0:5 þ ði 1Þ=2 þ ðj 1Þ=2; k ¼ 0:5. Then the surface f(x, y) is shown by the left picture of Fig. 1, the surface of interpolating function is shown by the central picture of Fig. 1, and the surface of the error is shown by the right picture of Fig. 1. Changing the parameters, let ai;j ¼ 1555 þ 2i þ 2j; bi;j ¼ 3555 þ 3i þ 3j, then the surface of the interpolation, written as p(x, y), and the surface of the error are shown by Fig. 2. From Fig. 1, it can be seen that the error of the interpolation is smaller than ±8 103. This means that under the condition of C1-continuous, a good approximation can be achieved. From Fig. 2, we can know that the error of the interpolation is smaller than ±15 103. This means that with large parameters, the error changes little, hence the interpolation is stable with parameters, therefore it can be used to modify the shape of the interpolating surface appropriately. In the following, we present a simple example to demonstrate how to construct the interpolating surface by solving the parameters ai;j ; bi;j ; k. For simplicity, we assume that ai;j ¼ bi;j a; k ¼ 0:5. Suppose the interpolated function is f ðx; yÞ ¼
−3
x 10 0.95
0.95
2
0.9
0.9
0
0.85
0.85
−2
0.8
0.8
0.75
0.75
0.7 1.5 1
y
0.5
0.4
0.6
0.8
x
1
1.2
1.4
z
4
z
1
z
1
−4 −6 −8 1.5
0.7 1.5 1
y
0.5
0.4
0.6
0.8
x
1
1.2
1.4 1
y
0.5
0.4
0.6
0.8
1
1.2
x
Fig. 1. The surfaces when ai;j ¼ bi;j ¼ 0:5 þ ði 1Þ=2 þ ðj 1Þ=2: (left) function’s surface, (center) interpolating surface, (right) error’s surface.
1.4
4652
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
−3
x 10 5
1 0.95
0
z
z
0.9 0.85 0.8
−5 −10
0.75 0.7 1.5
y
1 0.5
0.8
0.6
0.4
1
1.2
1.4
−15 1.5
y
x
1 0.5
0.4
0.6
0.8
1
1.2
1.4
x
Fig. 2. The surfaces when ai;j ¼ 1555 þ 2i þ 2j; bi;j ¼ 3555 þ 3i þ 3j: (left) interpolating surface, (right) error’s surface.
0.08
z
0.06 0.04 0.02 0 0.7
0.65
0.6
x
0.55
0.5 0.7
0.65
0.6
0.55
0.5
y
Fig. 3. Central point interpolating surface when s = 0.0245.
ð1 xÞ2 ð1 yÞ2 ; ðx; yÞ 2 ½0:5; 0:9; 0:5; 0:9, and let hi ¼ lj ¼ 1=5; i ¼ 1; 2; j ¼ 1; 2, then xi ¼ 0:5 þ ði 1Þ=5; yj ¼ 0:5þ ðj 1Þ=5; i ¼ 1; 2; 3; j ¼ 1; 2; 3. Suppose we need to determine the parameter a such that the interpolating surface interpolates the central point (0.6, 0.6, s), s is some real number. The matrices in Eq. (9) at the domain point (0.6, 0.6) can be expressed as follows:
0
F 1;1
1 f ð0:5; 0:5Þ f ð0:5; 0:7Þ f ð0:5; 0:9Þ B C ¼ @ f ð0:7; 0:5Þ f ð0:7; 0:7Þ f ð0:7; 0:9Þ A; f ð0:9; 0:5Þ f ð0:9; 0:7Þ f ð0:9; 0:9Þ
B11;1 ¼ ðB21;1 ÞT ¼ ðC 11;1 ÞT ¼ C 21;1 ¼
0:5aþ0:25 aþ1
3 8
3 4
1 8
0:5aþ1 aþ1
0:25 aþ1
;
:
Hence, we can obtain the following equation:
p1;1 ð0:6; 0:6Þ ¼ B11;1 F 1;1 C 11;1 þ C 21;1 F 1;1 B21;1 ¼ 2s: If we take s = f(0.6, 0.6), then the parameter can be solved as a = 1; If we take s = 0.0245, then the parameter can be solved as a = 0.1 (see Fig. 3). 7. Conclusion This paper presents a new weighted bivariate rational bicubic interpolation based on function values which is defined by Theorem 1, and discusses some properties of it under some conditions, see Theorems 2–7. In the numerical examples, we can see that the error of the interpolation is much small and the interpolation is stable with parameters. Besides, we can modify the shape of the interpolating surface appropriately, hence it is useful in computer aided geometric design. The other properties of the interpolating function (such as the relationship between the shape and parameters) will be presented in the future research. Moreover, the upper boundary in Theorem 7 would not be the best, so it will also be presented in the future research. Acknowledgements This work was supported by the National Nature Science Foundations of China under Grant Nos. 60873111 and 60933007. We thank to all referees for their comments which helped us to improve the paper.
W.-X. Huang, G.-J. Wang / Applied Mathematics and Computation 217 (2011) 4644–4653
4653
References [1] B. Barsky, D.P. Greenberg, Determining a set of B-spline control vertices to generate an interpolating surface, Comput. Graph. Image Proc. 14 (1980) 203–226. [2] B.K. Choi, H.Y. Shin, W.S. Yoo, Visually smooth composite surfaces for an unevenly spaced:3D data array, Comput. Aided Geom. Des. 10 (1993) 157– 171. [3] J. Peters, Local smooth surface interpolation: a classification, Comput. Aided Geom. Des. 7 (1990) 191–195. [4] Q. Duan, L. Wang, E.H. Twizell, A new bivariate rational interpolation based on function values, Inform. Sci. 166 (2004) 181–191. [5] Q. Duan, K. Djidjeli, W.G. Price, E.H. Twizell, A rational cubic spline based on funcion values, Comput. Graph. 22 (1998) 479–486. [6] Q. Duan, Y. Zhang, E.H. Twizell, A bivariate rational interpolation and the properties, Appl. Math. Comput. 179 (2006) 190–199. [7] Q. Duan, H. Zhang, Y. Zhang, Bounded Property and point control of a bavariate rational interpolating surface, Comput. Math. Appl. 52 (2006) 975–984. [8] Q. Duan, K. Djidjeli, W.G. Price, E.H. Twizell, Weighted rational cubic spline interpolation and its application, J. Comput. Appl. Math. 117 (2000) 121– 135. [9] Q. Duan, L. Wang, E.H. Twizell, A new weighted rational cubic interpolation and its approximation, Appl. Math. Comput. 168 (2005) 990–1003. [10] Q. Duan, L. Wang, E.H. Twizell, A new C2 rational interpolation based on function values and constrained control of the interpolant curves, Appl. Math. Comput. 161 (2005) 311–322. [11] M.P. George, Interpolation and Approximation by Polynomials, Springer-Verlag, 2003.