張 莉，趙 林，檀結(jié)慶,2

(1.合肥工業(yè)大學(xué) 數(shù)學(xué)學(xué)院， 安徽 合肥 230009； 2.合肥工業(yè)大學(xué) 計(jì)算機(jī)與信息學(xué)院， 安徽 合肥 230009)

帶互異權(quán)值的漸進(jìn)迭代逼近算法及其應(yīng)用

張 莉1，趙 林1，檀結(jié)慶1,2

(1.合肥工業(yè)大學(xué) 數(shù)學(xué)學(xué)院， 安徽 合肥 230009； 2.合肥工業(yè)大學(xué) 計(jì)算機(jī)與信息學(xué)院， 安徽 合肥 230009)

在計(jì)算機(jī)輔助幾何設(shè)計(jì)(CAGD)領(lǐng)域,漸進(jìn)迭代逼近(PIA)算法因其具有很好的自適應(yīng)性和收斂穩(wěn)定性，被廣泛應(yīng)用于插值與逼近問(wèn)題.其中帶權(quán)漸進(jìn)迭代逼近(WPIA)算法通過(guò)調(diào)整向量加權(quán)明顯加快了收斂速度.提出了一種帶互異權(quán)值的漸進(jìn)迭代逼近算法,不僅操作靈活,還可根據(jù)需要對(duì)各控制頂點(diǎn)進(jìn)行調(diào)整，實(shí)現(xiàn)不同的迭代效果;同時(shí)通過(guò)引入一個(gè)參數(shù),給出了可調(diào)權(quán)值迭代算法,當(dāng)參數(shù)取合適值時(shí)，該算法的收斂速度比帶權(quán)PIA算法更快,且權(quán)值取法不依賴于配置矩陣的特征值.最后用數(shù)值實(shí)例,通過(guò)對(duì)Bézier曲線、張量積Bézier曲面,以及三角Bézier曲面進(jìn)行迭代,展示了該算法的有效性.

漸進(jìn)迭代逼近；帶權(quán)漸進(jìn)迭代逼近；插值與逼近；Bézier曲線曲面；三角Bézier曲面

Progressive iterative approximation with different weights and its application. Journal of Zhejiang University(Science Edition), 2017,44(1):022-027

在計(jì)算機(jī)輔助幾何設(shè)計(jì)與逆向工程中,構(gòu)造一組不同精度的擬合曲線(曲面)進(jìn)行插值或逼近于一個(gè)給定的有序點(diǎn)集是一類很重要的課題.為了更好地解決這一問(wèn)題,早在1975年,齊東旭等[1]就提出了均勻三次B-spline曲線的盈虧修正算法;1979年DE BOOR[2]在報(bào)告中也闡述了類似的思想.

2004年,藺宏偉等[3]用逐步迭代非均勻B-spline曲線曲面的方法擬合給定點(diǎn)集,證明了非均勻三次B-spline曲線(曲面)具有盈虧修正性質(zhì);2005年,藺宏偉等[4]進(jìn)一步證明了所有的標(biāo)準(zhǔn)全正基(NTP bases)混合曲線(曲面)都具有這一性質(zhì),并將此迭代方法命名為漸進(jìn)迭代逼近(progressive iterative approximation, PIA).DELGADO等[5]發(fā)現(xiàn)規(guī)范B基在所有NTP基函數(shù)中收斂速度最快.2011年,陳杰等[6]給出了三角域上的PIA方法,證明了三角Bézier曲面和有理三角Bézier曲面在均勻參數(shù)化下的PIA性質(zhì).

為研究PIA方法與傳統(tǒng)的代數(shù)插值方法之間的關(guān)系,鄧少輝等[7]推導(dǎo)了代數(shù)插值和PIA方法之間的等價(jià)關(guān)系.藺宏偉[8-9]給出了PIA方法的局部性質(zhì)和收斂性證明,并將局部性質(zhì)應(yīng)用于自適應(yīng)數(shù)據(jù)擬合技術(shù)中,取得了很好的迭代效果.

陸立正[10]通過(guò)對(duì)調(diào)整向量加權(quán)的方式對(duì)PIA方法進(jìn)行加速(簡(jiǎn)稱為WPIA),并給出了最優(yōu)的權(quán)值取法.陳杰等[11]給出了WPIA方法的2種推廣:帶權(quán)的局部PIA方法和均勻參數(shù)化下的三角Bézier曲面的帶權(quán)PIA方法.2015年,劉曉艷等[12]結(jié)合Jacobi迭代法,提出了非均勻三次B樣條曲線插值的Jacobi-PIA算法,對(duì)PIA迭代進(jìn)行加速.

2011年,藺宏偉等[13]給出了逼近型PIA迭代格式,實(shí)現(xiàn)了PIA對(duì)大量數(shù)據(jù)的擬合.2014年,鄧重陽(yáng)等[14]提出了一種新的逼近型PIA方法,迭代的極限曲線或曲面是給定數(shù)據(jù)點(diǎn)的最小二乘擬合.基于諸多學(xué)者對(duì)PIA的大量研究與應(yīng)用,藺宏偉[15]從理論和應(yīng)用兩方面對(duì)PIA方法進(jìn)行了綜述,展示了PIA方法在多個(gè)領(lǐng)域成功應(yīng)用的實(shí)例.

本文提出一種帶互異權(quán)值的PIA方法,在迭代過(guò)程中,每一個(gè)控制頂點(diǎn)所對(duì)應(yīng)的調(diào)整向量都乘以不同的權(quán)值,增加了實(shí)際操作的靈活性.當(dāng)基函數(shù)是三次非均勻B樣條基、權(quán)值取某些特殊值時(shí),文獻(xiàn)[12]中的方法即為本文方法的特例;當(dāng)所有權(quán)值都相等時(shí),本文方法就是WPIA方法,數(shù)值實(shí)例表明，當(dāng)取一組合適的權(quán)值時(shí),本文方法比WPIA方法收斂速度更快;此外,還給出一種權(quán)值的取法，在減少計(jì)算量的同時(shí)(不需要求矩陣特征值),收斂速度更快.

1 參數(shù)曲線帶互異權(quán)值的PIA方法

1.1 迭代格式的導(dǎo)出

是一個(gè)全正矩陣.

類似地,第k+1步的迭代為

為加速收斂的同時(shí)增加操作的靈活性,對(duì)上述迭代步驟做如下調(diào)整：

Δk=(I-BW)Δk-1=(I-BW)kΔ0,

(1)

其中,I是n+1階單位矩陣,

B為基函數(shù){Bi(t)}在參數(shù){ti}(i=0,1,…,n)下的配置矩陣,W=diag(μ0,μ1,…,μn)稱為權(quán)矩陣.

1.2 迭代格式(1)的收斂性

定理1 迭代格式(1)是收斂的,即迭代的極限曲線插值于給定的數(shù)據(jù)點(diǎn).

證明 由基函數(shù)的非負(fù)性及歸一化性質(zhì)可得

0<λi(BW)<2,

所以迭代矩陣D=I-BW的特征值滿足

-1<λi(D)=1-λi(BW)<1,

D的譜半徑ρ(D)<1,則迭代格式(1)收斂.那么

即極限曲線插值于初始數(shù)據(jù)點(diǎn).

定理2 當(dāng)?shù)袷?1)的權(quán)值取

時(shí)，其中a∈(0,2),迭代格式(1)亦收斂.

證明 矩陣BW的列范數(shù)

矩陣BW的特征值滿足

-1<λi(D)=1-λi(BW)<1,

所以ρ(D)<1，迭代格式(1)收斂.

2 張量積曲面PIA方法的改進(jìn)

2.1 迭代格式的導(dǎo)出

設(shè)第k+1次迭代的曲面為

其中，

由此得到如下迭代格式：

Δk=DΔk-1=DkΔ0,

(2)

其中，

2.2 迭代格式(2)的收斂性

類似于曲線的情況,對(duì)于張量積曲面,同樣給出了關(guān)于收斂性分析的2個(gè)定理.并且定理的證明與曲線的情況類似，注意矩陣B為2個(gè)配置矩陣的Kronecker積即可.

定理3 迭代格式(2)是收斂的,即極限曲面插值于給定的數(shù)據(jù)點(diǎn).

定理4 當(dāng)?shù)袷?2)中的權(quán)值取

時(shí)(p=0,1,…,m,q=0,1,…,n),迭代格式(2)收斂.

3 三角Bézier曲面PIA方法的改進(jìn)

3.1 迭代格式的導(dǎo)出

三角域T∶={(u,v,w)∶u,v,w≥0,u+v+w=1}上的一個(gè)n階三角Bézier曲面為

隨著迭代的進(jìn)行,可得到第l+1次的曲面為

為便于表達(dá),引入字典順序的定義.

定義[16]給定2個(gè)d維向量α,β,若α-β=(α1-β1,α2-β2,…,αd-βd)的第1個(gè)非零元為正數(shù)，則向量α排在向量β之前,記作α>β.那么,定義在三角域T上的n階三角Bernstein基函數(shù),可表示為(n+1)(n+2)/2維的行向量

參數(shù)序列為

則得到如下迭代格式

Δl+1=(I-BW)Δl=(I-BW)l+1Δ0，

(3)

3.2 迭代格式(3)的收斂性

定理5 迭代格式(3)是收斂的,即極限曲面插值于給定的數(shù)據(jù)點(diǎn).

證明 雙變量Bernstein基函數(shù)Bn在參數(shù)tn下的配置矩陣為

由于雙變量Bernstein基函數(shù)的非負(fù)性和歸一性，

則矩陣B的無(wú)窮范數(shù)為

所以矩陣B的特征值λi,j,k(B)滿足

0<λi,j,k(B)≤1, i+j+k=n.

又由于μi,j,k∈(0,2),那么矩陣W的特征值和無(wú)窮范數(shù)滿足

0<λi,j,k(W)=μi,j,k<2，

0<λi,j,k(BW)<2,

-1<λi,j,k(D)=1-λi,j,k(BW)<1,

矩陣D的譜半徑ρ(D)<1,則迭代格式(3)收斂.即

那么，極限曲線插值于初始數(shù)據(jù)點(diǎn).

定理6 當(dāng)?shù)袷?3)中的權(quán)值取

證明 矩陣BW的列范數(shù)為

此時(shí),矩陣BW和D的特征值分別滿足

-1<λi,j,k(D)=1-λi,j,k(BW)<1,

所以ρ(D)<1,迭代格式(3)收斂.

4 數(shù)值實(shí)例

用Bézier曲線、張量積Bézier曲面,以及三角Bézier曲面進(jìn)行迭代,逐步插值于給定數(shù)據(jù)點(diǎn),并將本文方法與經(jīng)典PIA[4]和取最佳權(quán)值的WPIA方法[10-11]進(jìn)行比較.

數(shù)值實(shí)例均采用均勻參數(shù)化方法;曲線，張量積曲面和三角曲面第k步迭代的誤差分別取

例1 在圓弧p(t)=(sin t,cos t)上取樣本p(ti)(ti=πi/6,i=0,1,…,10)作為數(shù)據(jù)點(diǎn)，并作Bézier曲線.圖1中黑色的是使用本文方法取a=1.4迭代后的曲線.表1展示了圓弧在相同迭代次數(shù)下不同PIA方法的誤差比較.

圖1 使用本文方法迭代例1中的數(shù)據(jù)點(diǎn)Fig.1 Iterate the data in example 1 using our method

所用方法迭代次數(shù)03510PIA方法0．29300．02140．00840．0023WPIA方法0．29300．07850．05060．0176本文方法0．29300．01370．00420．0009

由表1知,在相同迭代次數(shù)下,本文方法的誤差明顯較小.需要注意的是,WPIA算法的加速收斂作用會(huì)隨著迭代次數(shù)的增加而漸趨明顯,而在迭代之初,其加速作用不明顯.

例2 圖2展示的是以類peaks曲面上的點(diǎn)(1,1,1),(1,2,4),(1,3,2),(1,4,2),(1,5,2),(2,1,2),(2,2,2),(2,3,3),(2,4,6),(2,5,4),(3,1,3),(3,2,2),(3,3,4),(3,4,4),(3,5,3),(4,1,4),(4,2,6),(4,3,1),(4,4,1),(4,5,2)為數(shù)據(jù)點(diǎn),取a=1.9,采用本文方法迭代后的張量積Bézier曲面.表2列出了各方法的誤差.

表2 例2中相同迭代次數(shù)下不同方法的誤差比較

例3 圖3展示的是以字典順序排列的點(diǎn)列{(0,3,0.3),(-0.5,2,3),(0.5,2,4),(-1,1.5,1),(0,1.4,3),(1,1.5,0.5)}為數(shù)據(jù)點(diǎn),取a=1.2，使用本文方法迭代后的三角Bézier曲面.表3列出了相應(yīng)的誤差,可以看出迭代效果有顯著改善.

圖2 使用本文方法迭代類例2中的數(shù)據(jù)點(diǎn)Fig. 2 Iterate the data in example 2 using our method

圖3 使用本文方法迭代例3中的6個(gè)數(shù)據(jù)點(diǎn)Fig. 3 Iterate the data in example 3 using our method

所用方法迭代次數(shù)0159PIA方法1．80430．90220．05640．0035WPIA方法1．80430．60140．00749．2×10-5本文方法1．80430．36095．8×10-49．2×10-7

5 結(jié) 論

給出了一種帶互異權(quán)值的PIA方法.在調(diào)整過(guò)程中，對(duì)不同的數(shù)據(jù)點(diǎn)分配不同的權(quán)值,操作較靈活；同時(shí)給出了一組含參數(shù)的權(quán)值,在迭代次數(shù)不變的情況下通過(guò)調(diào)控參數(shù)可以得到更小的誤差；并且本文方法不必依賴于矩陣特征值;此外，WPIA的加速收斂作用在迭代之初不明顯，需多次迭代才顯現(xiàn)，而本文方法在迭代前期誤差就很小.

數(shù)值實(shí)驗(yàn)證明,本文方法不僅適用于Bézier曲線和張量積Bézier曲面,亦適用于三角Bézier曲面.

[1] 齊東旭,田自賢,張玉心，等.曲線擬合的數(shù)值磨光方法[J].數(shù)學(xué)學(xué)報(bào),1975,18(3):173-184. QI D X, TIAN Z X, ZHANG Y X, et al. The method of numeric polish in curve fitting [J]. Acta Mathematica Sinica,1975,18(3):173-184.

[2] DE BOOR C. How does Agee’s smoothing method work[J/OL]. [2016-7-10]. http://ftp.cs.wisc.edu/debooron/ agee.pdf.

[3] 藺宏偉,王國(guó)瑾,董辰世.用迭代非均勻B-spline曲線(曲面)擬合給定點(diǎn)集[J].中國(guó)科學(xué):E輯,2004,33:315-331. LIN H W, WANG G J, DONG C S. Constructing iterative non-uniform B-spline curve and surface to fit data points [J]. Science in China: Ser E,2004,33:315-331.

[4] LIN H W, BAO H J, WANG G J. Totally positive bases and progressive iteration approximation [J]. Computers and Mathematics with Applications,2005,50(3/4):575-586.

[5] DELGADO J, PENA J M. Progressive iterative approximation and bases with the fastest convergence rates [J]. Computer Aided Geometric Design,2007,24(1):10-18.

[6] CHEN J, WANG G J. Progressive iterative approximation for triangular Bézier surfaces [J]. Computer-Aided Design,2011,43(8):889-895.

[7] 鄧少輝,汪國(guó)昭.漸進(jìn)迭代逼近方法的數(shù)值分析[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào),2012,24(7):879-884. DENG S H, WANG G Z. Numerical analysis of the progressive iterative approximation method [J]. Journal of Computer-Aided Design and Graphics,2012,24(7):879-884.

[8] LIN H W. Local progressive-iterative approxim-ation format for blending curves and patches [J]. Computer Aided Geometric Design,2010，27(4):322-339.

[9] LIN H W. Adaptive fitting by the progressive-iterative approximation [J]. Computer Aided Geometric Design,2012,29(7):463-473.

[10] LU L Z. Weighted progressive iteration approximation and convergence analysis[J]. Computer Aided Geometric Design,2010,27:129-137.

[11] 陳杰,王國(guó)瑾,金聰健.兩類推廣的漸近迭代逼近[J].自動(dòng)化學(xué)報(bào),2012,38(1):135-139. CHEN J, WANG G J, JIN C J. Two kinds of generalized progressive iterative approximations [J]. Acta Automatica Sinica,2012,38(1):135-139.

[12] 劉曉艷,鄧重陽(yáng).非均勻三次B樣條曲線插值的Jacobi-PIA算法[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào),2015,27(3):484-491. LIU X Y, DENG C Y. Jacobi-PIA algorithm for non-uniform cubic B-spline curve interpola-tion[J]. Journal of Computer-Aided Design& Graphics,2015,27(3):484-491.

[13] LIN H W, ZHANG Z Y. An extended iterative format for the progressive-iteration approximation [J]. Computers & Graphics,2011,35:967-975.

[14] DENG C Y, LIN H W. Progressive and iterative approximation for least squares B-spline curve and surface fitting[J]. Computer-Aided Design,2014,47(1):32-44.

[15] 藺宏偉.幾何迭代法及其應(yīng)用綜述[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào),2015,27(4):582-589. LIN H W. Survey on geometric iterative methods with applications[J]. Journal of Computer-Aided Design& Graphics,2015,27(4):582-589.

[16] DUNLL C, XU Y. Orthogonal Polynomials of Several Variables[M]. Cambridge: Cambridge University Press,2001.

ZHANG Li1, ZHAO Lin1, TAN Jieqing1,2

(1.SchoolofMathematics,HefeiUniversityofTechnology,Hefei230009,China; 2.SchoolofComputerandInformation,HefeiUniversityofTechnology,Hefei230009,China)

In CAGD, progressive iterative approximation (PIA) method is widely used to solve interpolation and approximation problems due to its perfect adaptability and convergence stability. Weighted progressive iterative approximation (WPIA) can accelerate the convergence rate by assigning an appropriate weight for each adjusting vectors. One new PIA method with mutually different weights is presented. It not only provides more flexibility in operation, but also achieves satisfactory iterative result for different control vertices. A set of weights with an adjustable parameter has also been put forward, which can be obtained without resorting to the eigenvalue of collocation matrices and can speed up the convergence rate compared with the WPIA method. Numerical examples of Bézier curves, tensor-product Bézier surfaces and triangular Bézier surfaces demonstrate the effectiveness of the method.

progressive iterative approximation(PIA); weighted progressive iterative approximation (WPIA); interpolation and approximation; Bézier curves and surfaces; triangular Bézier surfaces

2016-07-15.

國(guó)家自然科學(xué)基金重點(diǎn)資助項(xiàng)目(U1135003);國(guó)家自然科學(xué)基金資助項(xiàng)目(61472466,61100126);中國(guó)博士后科學(xué)基金面上資助項(xiàng)目(2015M571926);浙江大學(xué)CAGD&CG國(guó)家重點(diǎn)實(shí)驗(yàn)室開(kāi)放課題(A1607).

張 莉(1976-),ORCID:http://orcid.org/0000-0002-9208-1949，女，博士，教授，主要從事CAGD研究，E-mail：lizhang@hfut.edu.cn.

10.3785/j.issn.1008-9497.2017.01.003

TP 391.41

A

1008-9497(2017)01-022-06