基于NURBS曲面插值的船体曲面重构

基于NURBS曲面插值的船体曲面重构钱宏1,刘敏2,贺庆1,刘朕明1,荣焕宗1（1.中国船舶及海洋工程设计研究院，上海200011；2.中国科学技术大学，合肥230000）摘要本文提出了一种船体曲面重构方法。首先，用特征线把船体曲面划分为曲面片，便于表达复杂船体曲面和保留船体特征。其次，对船体线框模型的网格均匀化处理，得到型值点的位置，并给它们加上特征，生成NURBS网格。最后，用带重节点的B样条曲面插值带切矢和重点的NURBS网格，并进行曲面片拼接，生成无缝隙的NURBS曲面。因此，可以用较少的曲面片（几至十几块）重构带有球首、球尾、平边线、平底线、部分折角线和方尾的复杂船体曲面。这给船舶结构设计、CAM和CFD计算提供了良好的船体NURBS曲面模型。关键词：船舶设计；NURBS曲面插值；船体曲面重构中图分类号：U66.22文献标识码：A0前言1在船体型线设计中，船体线框和曲面设计方法[1]是以B样条曲线为基础，先用B样条曲线插值算法，根据给定的型值点生成横剖线、水线、纵剖线、首尾轮廓线和折角线等，创建船体线框模型。然后用NURBS曲面片覆盖此线框模型重构船体曲面，最后生成船体曲面模型。由于该方法符合设计船舶型线的习惯，船舶设计师容易掌握，在国内还有很大市场。重构船体曲面的方法主要有孔斯（Coons）曲面、蒙皮曲面和插值曲面方法。孔斯曲面方法[2]采用超限插值法构造满足给定边界和跨界切矢条件的曲面，使构造出来的曲面片光滑拼合。方法简洁、易于实现，且适应各种船型。但是，要得到比较光顺的船体曲面，必须对线框模型进行仔细处理和光顺，而且生成的曲面片太多（几十至几百块），曲面片之间会有缝隙，这给进一步交互修改船体、生成光顺船体曲面带来困难。瑞典TribonM3[3]和芬兰NAPA[4]等软件属于这类方法，它们占领着中国市场。蒙皮曲面方法[2]是先在关键位置处布置若干条描述曲面形状的截面曲线，然后生成一张通过这些曲线的曲面。这种方法直观、易于实现，生成的曲面片少（几块），但是对于复杂船型很难实现。闫秋莲和杨启[5]采用一种改进的蒙皮法重构船体曲面，利用广义站改进了船体曲面的首尾控制网格，并利用节点删除优化了相容性处理后的节点数目，具有一定的工程应用价值，但是对于复杂船型仍然很难实现。插值曲面方法[6-7]采用插值法构造通过给定型值点或点云数据的曲面，常用NURBS曲面片插值或三角曲面片插值等。这种方法需对型值点进行分块，对各分块的型值点分别进行插值，得到分片插值曲面，然后采用曲面拼接技术构造复杂曲面。该方法生成的曲面片太多（几十至几百块），而且曲面片之间会有缝隙。对于船型设计来说，常常需对船体曲面进行修改以满足设计要求，显然，太多的曲面片使得交互修改很困难，而且工作量大。虽然在逆1收稿日期：2015-02-04；修改稿收稿日期：2015-2向工程的曲面重构中有广泛的应用，但是还没有见到用于重构船体曲面的文献。本文提出了一种NURBS曲面插值方法，直接插值船体线框模型的网格型值点，生成船体NURBS曲面，成功地处理了带切矢的点和重点等复杂情况，实现了用较少的曲面片（几至十几块）重构带有球首、球尾、平边线、平底线、部分折角线和方尾的复杂船体曲面。由于NURBS曲面由较少的曲面片组成，所以可以比较方便地对它进行修改，得到满足设计要求的、光顺的船体曲面。1B样条（B-spline）曲线曲面定义1.1B样条曲线定义给定n+1个空间矢量(0,,)iinP…，节点矢量{}0mu,,uU…，1iiuu，0,,-1im…，称p次（degree）或p+1阶（order）参数曲线10()()11nii,ppniuNuuuu,pnCP(1)为p次B样条曲线。式中，m=n+p+1，)(uNi,p为p次B样条基函数，(0,,)iinP…为控制顶点。p次B样条曲线的定义域为][1npu,u。为便于统一处理，通常将B样条曲线的定义域取成规范参数域10u。B样条基函数用deBoor-Cox-Mansfield递推公式计算。本文只讨论非周期B样条曲线，简称为B样条曲线[2]。其节点矢量形式如下：111{0,,0,,,,1,,1}pnppuuU(2)1.2B样条曲面定义给定(n+1)(m+1)个空间矢量(0,,;0,,)i,jinjmP……，节点矢量1{}0npu,,uU…与1{}0mqv,,vV…，称pq次参数曲面1100(,)()(),nmi,ji,pj,qpnqmijuvNuNvuuuvvvSP(3)为pq次B样条曲面。式中，)(uNi,p和)(vNj,q为p和q次B样条基函数，(0,,;0,,)i,jinjmP……为控制顶点。本文只讨论非周期B样条曲面，简称为B样条曲面[2]。其节点矢量形式如下：3111111{0,,0,,,,1,,1}{0,,0,,,,1,,1}pnppqmqquuuuUV(4)2B样条曲线插值为了表示一组数据点（型值点）的位置和切矢，将三维空间的的点)(x,y,zQ附加上切矢控制条件，用ida表示左导数条件，idf表示右导数条件，根据它们可以确定该点的左右切矢。称它们为型值点的特征，可取如下值，切矢长度取数据点的弦长。ORDINARY：无给定的左右切矢，称为普通点。DIRLINE：以两点所在直线方向为切方向。若切方向为0，切矢条件退化为重顶点条件。DIRCURVE：以此点出发的连续三个数据点插值出的抛物线在此点处的切方向为切方向。若不够三个数据点，退化成DIRLINE；若切方向为0，切矢条件退化为重顶点条件。作为端点切矢的默认值。DIRX：以DIRLINE在X轴的投影为切方向。DIRY：以DIRLINE在Y轴的投影为切方向。DIRZ：以DIRLINE在Z轴的投影为切方向。DIRMIDCURVE：以此点和左右两个点插值出的抛物线在此点处的切方向为切方向。用于光滑拼接。若三个数据点中有且仅有两个相等，退化成DIRLINE；若切方向为0，切矢条件退化为重顶点条件。INVERSE：表示该方向的切矢与其相反方向切矢相等。用于180º的光滑拼接。ZERO：切矢取为0。当型值点中某点有左或右切矢时，用带有2重节点的B样条曲线插值；当型值点中某点有左右切矢时，用带有3重节点的B样条曲线插值，解法详见文献[8]。当型值点中有重型值点时，把重型值点处理成重控制顶点，解法与文献[8]类似，由于篇幅限制，不再详述。一般，两重点会产生明显的凸起，但是，如果在该处加上方向相反的左右切矢条件，就可使之比较平滑；三重点及更高重会产生折线。3B样条曲面插值为了表示一组网格形式的数据点（型值点）的位置和切矢，将三维空间的一个点z)y,(,xQ附加上切矢控制条件，用ida_u和idf_u分别表示沿着u方向的左、右导数条件，ida_v和idf_v分别表示沿着v方向的左、右导数条件，根据它们可以确定该点的u/v方向左右切矢，称它们为网格型值点的特征，它们的取值见第2节，切矢长度取子网格内对应等参数线的弦长。4设切矢条件总是一行行或一列列地出现，这种有切矢控制的线称之为特征线（边界也可以认为是特征线）。如果某一行/列只有部分数据有切矢控制，根据DIRMIDCURVE规则将该行缺失的特征补全。当网格型值点中某列/行有左或右切矢时，在u/v方向节点向量中用2重节点；当有左右切矢时，用3重节点；当型值点中有重型值点时，把重型值点处理成重控制顶点。下面将详细讨论边界有切矢，内部u、v方向特征线都有切矢的曲面插值算法。类似地，可以得到其它情况的算法，由于篇幅限制，不再详述。给定数据点,(0,,;0,,)ijinjmQ……，四条边界切矢,00nviiT（左边界）、,0nvimiT（右边界）、0,0mujjT（上边界）、,0munjjT（下边界）以及四个角点的扭矢0,0uvD、,0uvnD、0,uvmD、,uvnmD。再给出两条内部特征线的左右切矢：第3列网格线,20niiQ的左右切矢,20nviiT及扭矢0,2uvD、0,2uvD、,2uvnD、,2uvnD；第3行网格线2,0mjjQ的左右切矢2,0mujjT及扭矢2,0uvD、2,0uvD、2,uvmD、2,uvmD；以及内部交点处扭矢2,2uvD、2,2uvD、2,2uvD、2,2uvD。所有切矢可以根据型值点的特征ida_u、idf_u、ida_v和idf_v估算它们，详见第2节，切矢的长度取其所在等参数线的总弦长，即等参数线上所有型值点连接线的长度。为了方便起见，将扭矢都取为0，因而曲面上相应的地方可能会比较平坦，也可用其它方法确定扭矢，如文献[8]的简单扭矢等。网格数据点、切矢和扭矢如图1所示。图1u、v方向都带左右切矢的网格数据点、切矢和扭矢用带重节点的双三次B样条曲面插值这些条件5443300nmi,ji,j,ij(u,v)B(u)N(v)PP(5)这里，控制顶点为)5()5(mn的网格，u向节点数为9n，v向节点数为9m。对于曲面，既可以采用对每条网格线的参数化做平均来得到u、v的参数，也可采用与曲线参数无关的均匀参数化方法。经过试验，当数据分布不均匀时（网格线不是几乎平行的），给予曲线参数化然后取平均的方法效果并不好，因为此时的参数不再有几何意义。这里我们采用均匀参数化确定参数值和节点向量。u向参数：01111[01]inninu,uuu,u,,nnnU…………u向节点向量（含3重节点）：45674122210,0,0,0,,,,1,1,1,10,0,0,0,,1,1,1,1nnu,u,uuuu,u,u,u,uU……v向参数：01111[01]jmmjmv,vvv,v,,mmmV…………v向节点向量（含3重节点）：456741222m-10,0,0,0,,,1,1,1,10,0,0,0,v,,1,1,1,1mv,v,v,v,v,v,v,v,vV……插值条件为}m,{j},n,{i,)v,u(},{j},n,{i,)v,u(}m,{j},,{i,)v,u(},{j},,{i,)v,u(}m,,{j},n,,{i,)v,u(}m,,{j},n,,{i,)v,u(}m,,{j},n,,{i,)v,u(vuj,ijivuvuj,ijivuvuj,ijivuvuj,ijivuvj,ijivuj,ijiuj,iji22202220202020002000DPDPDPDPTPTPQP(6)将其写成矩阵形式SM(7)式中，6)v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,u()v,