霍彦妏，蔡占川



基于多结点样条磨光函数的几何迭代法

霍彦妏，蔡占川

(澳门科技大学资讯科技学院，澳门 999078)

几何迭代法在计算机辅助几何设计(CAGD)中有广泛地应用，为了提高传统的B-样条曲线插值在几何迭代中的收敛速度和迭代精度，提出了基于多结点样条磨光函数的几何迭代法，引入多结点样条磨光函数，在曲线拟合时把多结点样条磨光方法和几何迭代方法结合，经过磨光和迭代，在L-BFGS迭代算法的最优解下构造具有高逼近性的曲线拟合方法。实验结果表明，在相同精度下，该方法不仅减少了迭代次数，且提高了迭代速度，可以用于飞机、汽车等外形设计上，亦可用于文物、房屋等外形重构和重建，以及卫星图形图像的处理中。

几何迭代法；多结点样条磨光；L-BFGS算法；B-样条

1 研究背景

几何迭代法[1-2](geometric iteration method，GIM)，早期被称为渐进迭代逼近(progressive- iterative approximation，PIA)，其是迭代逼近算法在几何设计中的应用。该方法是将原始的离散数据点集用迭代的方法调整控制点，形成逐渐逼近控制点的新样条函数，最终形成一条光滑且具有更高逼近性的曲线或曲面的一种图形数据拟合方法[3]。

几何迭代法在计算机辅助几何设计(computer- aided geometric design，CAGD)中，解决的主要问题是工业产品几何形状的数字描述，例如汽车车轮、涡轮发动机、飞机外形的设计。CAGD技术起源于航空工业，由于飞机外形流线型的特殊结构是由多个自由的曲线及曲面组成的，所以CAGD技术的研究与曲线、曲面的发展息息相关。

1975年，齐东旭等[4]提出均匀三次B-样条的盈亏修正曲线拟合数值磨光方法，之后又在文献[5-6]中提出多结点样条的磨光方法，即用高阶的多结点样条磨光来达到更高的逼近效果，使“盈亏修正”方法得到更好地拓展和延伸。DE BOOR[7]于1979年也研究盈亏修正并且证明了该方法的收敛速度。2004年，LIN等[8]在非均匀三次B-样条的研究上证明了其同样具有盈亏修正的性质，随后在2006年提出了渐进迭代逼近(progressive iterative approximation，PIA)[9-13]的概念，并证明了几何迭代的性质同样适用于归一化的全正基混合曲线和曲面中。

2007年DELGADO和PEÑA[14]对比了在迭代过程中当取不同基函数的样条函数时，各自的收敛速度不同，且证明了B-样条的基函数迭代具有最高的收敛速度。同年，MAEKAWA等[15]用相邻数据点的参数值代替原始点的参数值并提出了几何插值(geometric interpolation，GI[16])的概念。熟知，由于PIA和GI[17]方法的相似性，在CAGD研究领域，被统称为几何迭代法。

为了提高函数的迭代精度以及收敛性，本文提出基于多结点样条磨光函数方法的几何迭代法，采用L-BFGS算法求型值点和调整后控制点之间距离的最小值即最优解，选取最恰当的点，以得到更加光滑且更加逼近型值点的几何图形。实验结果表明，该方法在相同的精度下，比通用的三次B-样条曲线插值迭代次数减少，且收敛速度成倍加快，不仅有效地保证了几何外形的精确性，还保持了其光滑的特性。

2 相关工作(GIM基本思想)

按此规律，生成的样条曲线的表达式如下：经过+1次迭代后，向量的差值为

图1 迭代示意图

3 多结点样条磨光函数

自由曲线或曲面的数据拟合是几何图形设计中不可或缺的部分。为了提高三次B-样条曲线或曲面的收敛性以及精确性，本文探讨多结点样条磨光方法并研究其迭代的几何意义。

当前，在CAD/CAM系统中，B-样条曲线曲面已经发展为被广泛使用、很成熟的一类样条方法。在齐东旭等[4]提出均匀的三次B-样条曲线拟合数值磨光方法；之后，为了推广B-样条函数，在保证样条曲线的精确性和光滑性同时兼顾的情况下，又引入级磨光算子，推出了“多结点” B-样条函数的磨光法[5-6]。事实上，B-样条函数是将次数为0的磨光算子在宽度=1的情况下，不断地磨光计算得出的函数[5]。

在研究曲线或曲面的拟合逼近问题时，齐东旭[18]研究了单位算子，并推出含有微分算子的样条函数磨光公式，使微分算子的逼近运算精度更高[19]。文献[20]研究了层次多结点样条算法曲线曲面的拟合和应用。

多结点样条磨光函数的定义为

图2 边界延拓示意图

样条函数的次数可以有0次，1次，2次， 3次······，由次数分别为0，1，2，3的表达式形成的基函数的图像，如图3所示。

多结点样条磨光法的优点：

(1) 一般性。该磨光方法只储存型值点，即不论几何图形是猫、狗亦或是松鼠这类型值点个数不同的几何图形，但其应用曲线的表达式都是统一的，在处理结果的误差上相比常用到的最小二乘法更小。

(2) 整体性。不论型值点如何增加，迭代的次数如何改变，几何图形的形状均能保持和基本初始的数据形状一致。

(3) 局部性。改变个别的数据点后，只会随之影响这些点对应的个别部分，而整体的几何图形形状不会改变。

图3 K=0,1,2,3时多结点样条基函数图像

4 基于多结点样条磨光函数的几何迭代法

熟知，牛顿法在进行Hessian的逆运算时，不仅费时而且占用大量内存空间，给计算带来很大的不便。相比牛顿法，BFGS算法由于Hessian的正定性即参数化的数据点唯一有且只对应一个数据点，并不需要求解Hessian，也就不需要占用很大的内存空间来存储每次生成的矩阵结果，既节省存储空间又提高了运行速度。而L-BFGS算法更是被称为“节省内存的BFGS算法”。

1989年LIU和NOCEDAL[21]提出了L-BFGS (Limited-memory BFGS)算法，该算法常应用于规模较大的数据处理的数值计算中。

其中，

上式是=3的情况，由于实验研究的是三次多结点样条磨光的基函数，所以其他3种情况的表达式不做说明。

在使用均匀三次B-样条函数进行数据拟合时，需求解一组线性方程的根来确定该几何形状控制顶点的位置。反向求解方式为确定样条函数的系数带来了很大的不便，既增大运算量，又延长了运算时间。本文采用多结点样条磨光方法来克服B-样条函数的不足。由图4可知，使用三次多结点样条磨光函数和三次B-样条函数做对比，多结点样条磨光函数的逼近效果要好于B-样条函数，提高了几何轮廓的精确性，而且更加逼近原始型值点。

图4 多结点样条磨光函数与B-样条对比

基于多结点样条磨光函数几何迭代方法算法如下：

输出：找到最小化函数值对应的坐标、步长、梯度、精度。

步骤1. 初始化阶段，选择起始结点，原始型值点参数化。

步骤2. 形成初始磨光样条函数曲线。

步骤3. 迭代次数=1；

步骤3.1. 阈值判定；

步骤3.2. 迭代正向或逆向求解。

步骤4. 阈值判定；

步骤4.1.若符合条件，画出曲线；

步骤4.2. 若不符合，=+1，循环步骤3，4。

本文算法流程图如图5所示。

5 实 验

本文的实验操作是基于16 G内存，i7-6700 CPUintel (R)处理器，在MATLAB编译器上完成本文的实验内容。

图5 本文实现过程流程图

5.1 实验结果

研究几何迭代法离不开研究拟合曲线迭代的速度和精度，即算法收敛性。文献[28]中在三次B-样条函数的基础上对各种几何图形做出了相应的迭代拟合。图6~20是本文算法与文献[28]的实验结果的对比图。当原始型值点个数不同且接近100个数据点时，最大的迭代次数在5次以内就可以达到所需的精度值。为了验证本文算法的收敛性，通过以下5种实例来验证算法的迭代效果以及精度分别与迭代次数和时间的关系。

在精度与迭代次数和时间的对比关系图中，显示相同的精度设置下的对比结果，蓝色线条为三次B-样条函数方法实现的迭代结果，红色线条为多结点样条磨光函数方法实现的迭代实验结果。由图6~20可知，在相同精度下，多结点样条磨光函数方法可比三次B-样条用更少的迭代次数和运行时间更加迅速地达到阈值条件。

5.2 实验分析

由表1和图21不同几何图形的精度对比以及表2的5种迭代详细情况对比可知，5种不同的几何图形用本文的算法迭代后在5次之内均可达到阈值迭代精度。

图6 牛迭代实验结果

图7 牛误差与迭代次数关系

图8 牛误差与时间关系

图9 松鼠迭代实验结果

图10 松鼠误差与迭代次数关系

图11 松鼠误差与时间关系

图12 大象迭代实验结果

图13 象误差与迭代次数关系

图14 象误差与时间关系

图15 米老鼠迭代实验结果

图16 米老鼠误差与迭代次数关系

图17 米老鼠误差与时间关系

图18 螃蟹迭代实验结果

图19 螃蟹误差与迭代次数关系

图20 螃蟹误差与时间关系

表1 不同几何图形误差对比

图21 5种几何图像误差对比

表2 采用多结点样条磨光方法与三次B-样条函数方法的几何图形迭代详细情况对比

6 结 论

本文基于多结点样条磨光函数，提出了相应的优化几何迭代方法。在获取几何图形的原始型值点后，利用多结点样条磨光函数的高逼近性和局部性，采用L-BFGS算法寻找最佳迭代逼近数据点的位置，通过数次磨光和迭代产生数据的最优解，在迭代达到指定精度值后，得到所需具有“较高精确性”的光滑曲线。实验结果表明，本文提出的磨光方法，比三次B-样条函数曲线在几何迭代的曲线拟合上逼近性更高，收敛性更快，在精度和速度上都有所提升。

今后，将继续探讨基于多结点样条磨光函数的几何迭代方法在三维的立体几何图形以及曲面中的应用。

[1] OKANIWA S, NASRI A, LIN H, et al. Uniform B-spline curve interpolation with prescribed tangent and curvature vectors [J]. IEEE transactions on visualization and computer graphics, 2012, 18(9): 1474-1487.

[2] LIN H, MAEKAWA T, DENG C. Survey on geometric iterative methods and their applications [J]. Computer-Aided Design, 2018, 95: 40-51.

[3] LIN H, ZHANG Z. An efficient method for fitting large data sets using T-splines [J]. SIAM Journal on Scientific Computing, 2013, 35(6): A3052-A3068.

[4] 齐东旭, 田自贤, 张玉心, 等. 曲线拟合的数值磨光方法[J]. 数学学报, 1975, 18(3): 173-184.

[5] 齐东旭, 梁振珊. 多结点样条磨光(Ⅰ)[J]. 高等学校计算数学学报, 1979, (2): 196-209.

[6] 齐东旭, 梁振珊. 多结点样条磨光(Ⅱ)[J]. 高等学校计算数学学报, 1981, (1): 68-81.

[7] DE BOOR C. How does Agee’s smoothing method work? [R]. Washington D C: Army Research Office 79-3, 1979: 299-302.

[8] 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: Series F, 2004, 47(3): 315-331.

[9] 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.

[10] LIN H, WANG G, DONG C. Constructing iterative non-uniform B-spline curve and surface to fit data points [J]. Science in China Series: Information Sciences, 2004, 47(3): 315-331.

[11] 邓少辉, 汪国昭. 渐进迭代逼近方法的数值分析[J]. 计算辅助设计与图形学学报, 2012, 24(7): 879-884.

[12] LIU M, LI B, GUO Q, et al. Progressive iterative approximation for regularized least square bivariate B-spline surface fitting [J]. Journal of Computational and Applied Mathematics, 2018, 327: 175-187.

[13] 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: 32-44.

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

[15] MAEKAWA T, MATSUMOTO Y, NAMIKI K. Interpolation by geometric algorithm [J]. Computer-Aided Design, 2007, 39(4): 313-323.

[16] LIN H. The convergence of the geometric interpolation algorithm [J]. Computer-Aided Design, 2010, 42(6): 505-508.

[17] XIONG Y, LI G, MAO A. Convergence analysis for B-spline geometric interpolation [J]. Computers & Graphics, 2012, 36(7): 884-891.

[18] 齐东旭. 关于计算机辅助几何造型数学方法的若干注记[J]. 北方工业大学学报, 1991, 3(1): 1-8.

[19] 李岳生, 齐东旭. 样条函数方法[M]. 北京: 科学出版社, 1979: 148-176.

[20] CAI Z C, LAN T, ZHENG C. Hierarchical MK splines: Algorithm and applications to data fitting [J]. IEEE Transactions on Multimedia, 2017, 19(5): 921-934.

[21] LIU D C, NOCEDAL J. On the limited memory BFGS method for large scale optimization [J]. Mathematical Programming, 1989, 45(1): 503-528.

[22] QIN Z W. The relationships between CG, BFGS, and two limited-memory algorithms [J]. Furman University Electronic Journal of Undergraduate Mathematics, 2016, 12(1): 5-20.

[23] NASH S G, NOCEDAL J. A numerical study of the limited memory BFGS method and the truncated- newton method for large scale optimization [J]. SIAM Journal on Optimization, 1991, 1(3): 358-372.

[24] SHANNO D F. On the convergence of a new conjugate gradient algorithm [J]. SIAM Journal on Numerical Analysis, 1978, 15(6): 1247-1257.

[25] LIU Y, WANG W P, LÉVY B, et al. On centroidal voronoi tessellation-energy smoothness and fast computation [J]. ACM Transactions on Graphics, 2009, 28(4): 1-17.

[26] WANG Y C, DONG L, XUE X, et al. On the design of constant modulus sequences with low correlation sidelobes levels [J]. IEEE Communications Letters, 2012, 16(4): 462-465.

[27] WANG Y C, WANG X, LIU H, et al. On the design of constant modulus probing signals for MIMO radar [J]. IEEE Transactions on Signal Processing, 2012, 60(8): 4432-4438.

[28] 刘超, 辛士庆, 舒振宇, 等. 几何迭代法的加速[J]. 计算机辅助设计与图形学学报, 2016, 28(11): 1838-1843.

Geometric Iteration Method Based on Many-Knot Spline Polishing Functions

HUO Yan-wen, CAI Zhan-chuan

(Faculty of Information Technology, Macau University of Science and Technology, Macau 999078, China)

Geometric iteration method has been widely used in computer aided geometric design (CAGD). In order to improve the convergence speed and iterative accuracy of the traditional B-spline curve interpolation in geometric iterations, this study proposes the geometric iteration method based on many-knot spline polishing functions, which introduces many-knot spline polishing functions, and combines many-knot spline polishing functions method and geometric iteration method in curve fitting. After polishing operator and iterating, the curve fitting method with high approximation under the optimal solution of L-BFGS iterative algorithm is constructed. Experimental results show that the proposed method not only reduces the times of iterations, but also improves the iterative speed under the same accuracy. The proposed geometric iteration method can be used in the shape design of airplanes, automobiles, etc. It can also be used to reconstruct and rebuild the shape of cultural relic houses and satellite image processing.

geometric iteration method; many-knot spline polishing functions; limited-memory BFGS algorithm; B-splines

TP 391

10.11996/JG.j.2095-302X.2019010015

A

2095-302X(2019)01-0015-09

2018-07-02；

2018-07-21

国家基础研究计划“973”项目(2011CB302400)；澳门科技发展基金项目(048/2016/A2，0012/2018/A1，0069/2018/A2)；国家自然科学基金面上项目(61272364)；浙江大学CAD&CG国家重点实验室开放课题(A1910)；北京理工大学珠海学院科研发展基金项目(XK-2018-04)

霍彦妏(1992-)，女，山西太原人，硕士。主要研究方向为计算机图形学等。E-mail：huohyw@163.com

蔡占川(1973-)，男，河北馆陶人，教授，博士，博士生导师。主要研究方向为计算机图形图像处理、数值分析。E-mail：zccai@must.edu.mo