CN103136431B - 一种自由曲面的保面积参数特性优化方法 - Google Patents

Abstract

本发明公开了一种自由曲面的保面积参数特性优化方法，具体步骤为：步骤一：输入空间自由曲面；步骤二：根据自由曲面类型的不同，分别计算空间自由曲面优化方程的系数；步骤三：通过计算两条平面参数曲线的交点获取自由曲面最优变换；步骤四：对空间自由曲面进行变换。在实际应用中，本发明的方法简单有效，在不改变用户给定的自由曲面形状和次数的前提下，通过对自由曲面进行线性重新参数化，只改变自由曲面的控制顶点、权因子以及节点向量，实现了对自由曲面的保面积特性的优化，并且保证了最优保面积曲面的唯一性。

Description

一种自由曲面的保面积参数特性优化方法

技术领域

本发明涉及一种自由曲面的保面积参数特性优化方法。

背景技术

国际标准化组织(ISO)于1991年正式确定把NURBS方法作为定义产品形状的唯一数学方法，越来越多的CAD/CAM系统采用NURBS作为其模型，NURBS在工业产品的外形设计中获得广泛的应用。基于NURBS曲面的几何算法大多是参数相关的，算法的结果和自由曲面当前的参数特性息息相关。

现有关于保面积参数特性的研究都针对离散三角面片模型，到目前为止还没有其他方法对NURBS曲面保面积参数特性进行优化。计算机辅助设计中（CAD），机械以及人体的模型都是基于自由曲面表示的，基于NURBS曲面的分析和可视化算法都是参数相关的，曲面的保面积特性直接影响自由曲面的分析和可视化的结果，特别是对于人脑、关节以及具备丰富凹凸特征的机械曲面。设计人员通常先对自由曲面的保面积特性进行优化，以满足后续分析和可视化算法的要求。因此基于NURBS曲面保面积参数特性优化方法具有很强的实用价值。

发明内容

本发明的目的就是为了解决上述问题，提供一种自由曲面的保面积参数特性优化方法，它具有算法简单有效、不破坏曲面形状、保持曲面次数和节点向量的优点。

为了实现上述目的，本发明采用如下技术方案：

一种自由曲面的保面积参数特性优化方法，具体步骤为：

步骤一：输入空间自由曲面；

步骤二：根据自由曲面类型的不同，分别计算空间自由曲面优化方程的系数；

步骤三：通过计算两条平面参数曲线的交点获取自由曲面最优变换；

步骤四：对空间自由曲面进行变换。

所述步骤一中空间自由曲面的表达式为：

$S(u,v)=\frac{{\mathrm{\Σ}}_{i=0}^{n_u-1}{\mathrm{\Σ}}_{j=0}^{n_v-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\ω}}_{i,j}{P}_{i,j}}{{\mathrm{\Σ}}_{i=0}^{n_u-1}{\mathrm{\Σ}}_{j=0}^{n_v-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\ω}}_{i,j}},u,v\∈\left[\mathrm{0,1}\right]$

其中，u和v是曲面S(u,v)的参数，n_u和n_v是曲面S(u,v)所含的控制顶点在u方向和v方向上的个数，P_i,j是曲面S(u,v)中序号为(i,j)的控制顶点的坐标值，ω_i,j是曲面S(u,v)中序号为(i，j)的控制顶点的权重，是定义在节点向量序列U之上的序号为i的p阶B-spline基函数；是定义在节点向量序列V之上的序号为j的q阶B-spline基函数，i＝0，1，…，n_u-1；j＝0，1，…，n_v-1，i和j均为大于等于0的整数，

其中p，q是曲面S(u,v)在u方向和v方向上的次数，u₀,u₁,…,u_n+p为节点向量序列U中的节点，其中u₀＝u₁＝…＝u_p＝0，u_n＝u_n+1＝…＝u_n+p＝1，0<u_p+1，u_n+2，…,u_n-1＜1；v₀,v₁,…,v_n+q为节点向量序列V中的节点，其中v₀=v₁＝…＝v_q＝0，v_n＝v_n+1＝…＝v_n+q＝1，0<v_q+1,v_q+2，…,v_n-1<1,所述的B-spline基函数采用Matlab系统样条工具库中B-spline基函数的定义方式。

所述步骤二的具体步骤为：

如果输入曲面是Bézier曲面，确定Bézier曲面的表示形式；

确定Bézier曲面在向量序列u/v方向上的偏微分，并确定Bézier曲面的微分第一基本形式；

根据Bézier曲面的微分第一基本形式确定优化方程的系数C_i,j。

如果输入曲面是B-spline曲面，通过节点插入法把B-spline曲面细分成n个Bézier曲面，用表示第r个Bézier曲面的局部参数域,确定第r个子曲面的优化方程系数；利用所有子曲面的优化方程系数确定B-spline曲面优化方程系数c_i,j。

如果输入曲面是NURBS曲面，取样2m个点在u区间内单调，2n个点在v区间内单调，

$\left\{\begin{array}{c}0=u_0<u_1<u_2<...<u_{2m}=1\\ 0=v_0<v_1<v_2<...<v_{2n}=1\end{array}\right.$

其中，u₀,u₁,…,u_2m为u方向采样点，v₁,…,v_2n为v方向采样点；使得

$I\&ap;T=\frac{h_1{g}_1}{4}\underset{i_1=0}{\overset{2m}{\mathrm{\&Sigma;}}}\underset{j_1=0}{\overset{2n}{\mathrm{\&Sigma;}}}{A}_{i_1}{A}_{j_1}{\left|\right|X_u(u_{i_1},v_{j_1})\&times;X_v(u_{i_1},v_{j_1})\left|\right|}^2{B}_i^2\left(u_{i_1}\right)B_j^2\left(v_{j_1}\right)---\left(1\right)$

满足：

$\left\{\begin{array}{cc}{h}_1=\frac{1}{2m}& \\ g_1=\frac{1}{2n},& A_{i_1}=\left\{\begin{array}{cc}1& i_1=\mathrm{0,2}m\\ 2& i_1=\mathrm{1,2},...2m-1\end{array}\right.\\ u_{i_1}=i_1{h}_1& \\ v_{j_1=j_1{g}_1}& \end{array}\right.$

$A_{j_1}=\left\{\begin{array}{cc}1& j_1=\mathrm{0,2}n\\ 2& j_1=\mathrm{1,2},...2n-1\end{array}\right.$

使用双重积分的算法计算C_i，j：

（2-1）设置式（1）中u，v两个方向采样密度m，n的值m₀=n₀=10;

（2-2）确定式（1）中数值积分T的值T₀=T,其中T表示数值双重积分；

（2-3）令m_r=2m_r-1，n_r=2n_r-1，其中m_r,n_r分别为u/v方向的采样密度；

（2-4）确定T的新近似值T_r;

（2-5）如果|T_r-T_r-1|＜∈，结束；其中∈是整个计算过程的终止条件；否则迭代次数r=r+1，回到步骤（2-3）。

所述步骤三的具体步骤为：

（3-1）确定当前曲面参数化与保面积参数化的偏差;

（3-2）利用公式等价法确定当前曲面参数化与保面积参数化偏差的最小值;

（3-3）通过计算两条平面曲线的交点确定最优参数变换的系数。

所述步骤四的具体步骤为：

（4-1）对节点向量参数u，v分别进行变换：

设s，t是对u，v变换后的参数；分别把参数s，t映射到单位区间上,使之满足：s和t为0时，u和v也为0；s和t为1时，u和v也为1；s和t为0.5时，u和v分别为1-α和1-β，由变换方程，确定参数s和t的表达式，α和β表示变换的两个自由参数；

（4-2）把s和t的表达式代入自由空间曲面的方程形式，确定重新参数化曲面X(s，t)，确定重新参数化曲面新的节点向量和

本发明的有益效果：在实际应用中，本发明的方法简单可行，不破坏曲面形状，在不改变用户给定的自由曲面形状和次数的前提下，通过对自由曲面进行线性重新参数化，只改变自由曲面的控制顶点、权因子以及节点向量，实现了对自由曲面的保面积特性的优化，并且保证了最优保面积曲面的唯一性。

附图说明

图1是整个方法的步骤流程图；

图2为输入的一个人体关节NURBS曲面S(u，v);

图3为对NURBS曲面S(u，v)进行保面积参数优化的最优结果。

具体实施方式

下面结合附图与实施案例对本发明作进一步说明。

一种自由曲面的保面积参数特性优化方法，具体步骤为：

步骤一：输入空间自由曲面；

步骤二：根据自由曲面类型的不同，分别计算空间自由曲面优化方程的系数；

步骤三：通过计算两条平面参数曲线的交点获取自由曲面最优变换；

步骤四：对空间自由曲面进行变换。

所述步骤一中空间自由曲面的表达式为：

$S(u,v)=\frac{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\&omega;}}_{i,j}{P}_{i,j}}{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\&omega;}}_{i,j}},u,v\&Element;\left[\mathrm{0,1}\right]$

其中，其中u和v是曲面S(u,v)的参数，n_u和n_v是曲面S(u，v)所含的控制顶点在u方向和v方向上的个数，P_i，j是曲面S(u，v)中序号为(i，j)的控制顶点的坐标值，ω_i,j是曲面S(u，v)中序号为(i，j)的控制顶点的权重，是定义在节点向量序列U之上的序号为i的p阶B-spline基函数；是定义在节点向量序列V之上的序号为j的q阶B-spline基函数，i＝0，1，…，n_u-1；j＝0，1，…，n_v-1，i和j均为大于等于0的整数，

其中p，q是曲面S(u,v)在u方向和v方向上的次数，u₀，u₁，…，u_n+p为节点向量序列U中的节点，其中u₀＝u₁＝…＝u_p＝0，u_n＝u_n+1＝…＝u_n+p＝1，0<u_p+1，u_n+2,…，u_n-1<1；v₀,v₁,…,v_n+q为节点向量序列V中的节点，其中v₀=v₁＝…＝v_q＝0，v_n＝v_n+1＝…＝v_n+q＝1，0<v_q+1，v_q+2，…,v_n-1<1，所述的B-spline基函数采用Matlab系统样条工具库中B-spline基函数的定义方式。

所述步骤二的具体步骤为：

（2-1）如果输入曲面是Bézier曲面，给定Bézier曲面的表示形式X(u，v)：

$X(u,v)=\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{n}{\mathrm{\&Sigma;}}}{B}_i^m\left(u\right)B_j^n\left(v\right)P_{i,j},u,v\&Element;\left[\mathrm{0,1}\right]$

其中，P_i,j是控制顶点，和是关于参数u和v的Bernstein基函数，m和n分别是曲面的u方向和v方向的次数，实际工程应用中，m和n一般取值为3。

确定X(u,v)在u/v方向上的偏微分X_u(u,v)和X_v(u。v)，

$X_u(u,v)=m\underset{i=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{n}{\mathrm{\&Sigma;}}}{B}_i^{m-1}\left(u\right)B_j^n\left(v\right)Q_{i,j},u,v\&Element;\left[\mathrm{0,1}\right]$

$X_v(u,v)=n\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{n-1}{\mathrm{\&Sigma;}}}{B}_i^m\left(u\right)B_j^{n-1}\left(v\right)R_{i,j},u,v\&Element;\left[\mathrm{0,1}\right]$

确定Q_i,j和R_i,j：

$\left\{\begin{array}{c}{Q}_{i,j}=P_{i+1,j}-P_{i,j}\\ R_{i,j}=P_{i,j+1}-P_{i,j}\end{array}\right.$

其中P_i+1,j，P_i,j+1，P_i,j是Bézier曲面控制顶点。

确定自由曲面X(u,v)的微分第一基本形式E，G，F：

$E=X_u\&CenterDot;X_u$

$=\frac{m!m!n!n!}{(2m-2)!\left(2n\right)!}\underset{i=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{k=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{l=0}{\overset{n}{\mathrm{\&Sigma;}}}\frac{(2m-2-i-k)!(i+k)!(2n-j-l)!(j+l)!}{(m-1-i)!i!(m-1-k)!k!(n-j)!j!(n-l)!l!}{B}_{i+k}^{2m-2}\left(u\right)B_{j+l}^{2n}\left(v\right)Q_{i,j}\&CenterDot;Q_{k,l}$

其中i,j,k,l,是循环变量，和分别是u和v方向上的Bernstein基函数，Q_k,l的计算方法是

Q_k,l=P_k+1,l-P_k,l

$G=X_v\&CenterDot;X_v=\frac{m!m!n!n!}{(2n-2)!\left(2m\right)!}\underset{i=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\underset{k=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{l=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\frac{(2m-i-k)!(i+k)!(2n-2-j-l)!(j+l)!}{(n-1-j)!j!(n-1-l)!l!(m-i)!i!(m-k)!k!}{B}_{i+k}^{2m}\left(u\right)B_{j+l}^{2n-2}\left(v\right)R_{i,j}$

$\&CenterDot;R_{k,l}$

R_k,l的计算方法是

R_k,l=P_k,l+1-P_k,l

$F=X_u\&CenterDot;X_v$

$=\frac{m!m!n!n!}{(2m-1)!(2n-1)!}\underset{i=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{k=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{l=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\frac{(2m-1-i-k)!(i+k)!(2n-1-j-l)!(j+l)!}{(m-1-i)!i!(m-k)!k!(n-j)!j!(n-1-l)!l!}{B}_{i+k}^{2m-1}\left(u\right)B_{j+l}^{2n-1}\left(v\right)Q_{i,j}$

$\&CenterDot;R_{k,l}$

最后得到优化方程的系数C_i，j:

$c_{i,j}=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_u\&times;X_v\left|\right|}^2{B}_i^2\left(u\right)B_j^2\left(v\right)\mathrm{dudv}$

$=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}(\mathrm{FG}-F^2)B_i^2\left(u\right)B_j^2\left(v\right)\mathrm{dudv}$

$=\left(\begin{array}{c}2\\ i\end{array}\right)\left(\begin{array}{c}2\\ j\end{array}\right)\frac{m!m!n!n!m!m!n!n!}{(4m+1)(4n+1)}\underset{i_1=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{j_1=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{k_1=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{l_1=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{i_2=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j_2=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\underset{k_2=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{l_2=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\frac{(4m-i_1-i_2-k_1-k_2-i)!(i_1+i_2+k_1+k_2+i)!}{(m-1-i_1)!i_1!(m-1-k_1)!k_1!(n-j_1)!j_1!(n-l_1)!l_1!}$

$\frac{(4n-j_1-j_2-l_1-l_2-j)!(j_1+j_2+l_1+l_2+j)!}{(n-1-j_2)!j_2!(n-1-l_1)!l_1!(m-i_2)!l_2!(m-k_2)!k_2!}(Q_{i_1,j_1}\&CenterDot;Q_{k_1,l_1})(R_{i_2,j_2}\&CenterDot;R_{k_2,l_2})$

$-\left(\begin{array}{c}2\\ i\end{array}\right)\left(\begin{array}{c}2\\ j\end{array}\right)\frac{m!m!n!n!m!m!n!n!}{(4m+1)(4n+1)}\underset{i_1=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{j_1=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{k_1=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{l_1=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{i_2=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{j_2=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\underset{k_2=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{l_2=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\frac{(4m-i_1-i_2-k_1-k_2-i)!(i_1+i_2+k_1+k_2+i)!}{(m-1-i_1)!i_1!(m-k_1)!k_1!(n-j_1)!j_1!(n-1-l_1)!l_1!}$

$\frac{(4n-j_1-j_2-l_1-l_2-j)!(j_1+j_2+l_1+l_2+j)!}{(m-1-i_2)!i_2!(m-k_2)!k_2!(n-j_2)!j_2!(n-1-l_2)!l_2!}(Q_{i_1,j_1}\&CenterDot;R_{k_1,l_1})(Q_{i_2,j_2}\&CenterDot;R_{k_2,l_2})$

$=\left(\begin{array}{c}2\\ i\end{array}\right)\left(\begin{array}{c}2\\ j\end{array}\right)\frac{m!m!n!n!m!m!n!n!}{(4m+1)(4n+1)}\underset{i_1=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{i_2=0}{\overset{m-1}{\mathrm{\&Sigma;}}}\underset{j_1=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\underset{j_2=0}{\overset{n-1}{\mathrm{\&Sigma;}}}\underset{k_1=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{k_2=0}{\overset{m}{\mathrm{\&Sigma;}}}\underset{l_1=0}{\overset{n}{\mathrm{\&Sigma;}}}\underset{l_2=0}{\overset{n}{\mathrm{\&Sigma;}}}\frac{(4m-i_1-i_2-k_1-k_2-i)!(i_1+i_2+k_1+k_2+i)!}{(m-1-i_1)!i_1!(m-1-i_2)!i_2!(n-l_2)!l_2!(n-l_1)!l_1!}$

$\frac{(4n-j_1-j_2-l_1-l_2-j)!(j_1+j_2+l_1+l_2+j)!}{(n-1-j_2)!j_2!(n-1-j_1)!j_1!(m-k_1)!k_1!(m-k_2)!k_2!}[(Q_{i_1,l_1}\&CenterDot;Q_{i_2,l_2})(R_{k_1,j_1}\&CenterDot;R_{k_2,j_2})-(Q_{i_1,l_1}$

$\&CenterDot;R_{k_1,j_1}\left)(Q_{i_2,l_2}\&CenterDot;R_{k_2,j_2})\right]$

其中，i₁，i₂，k₁，k₂，j₁，j₂，l₁，l₂是循环变量，和分别是关于u和v的Bernstein基函数， 计算方法是

$\left\{\begin{array}{c}{Q}_{i,j}=P_{i+1,j}-P_{i,j}\\ R_{i,j}=P_{i,j+1}-P_{i,j}\end{array}\right.$

（2-2）如果输入曲面是B-spline曲面，通过节点插入法把B-spline曲面细分成n个Bézier曲面，用表示第r个Bézier曲面的局部参数域,确定第r个子曲面的优化方程系数C_r,i，j：

其中，表示s对u的偏导数，表示t对v的偏导数，分别表示子曲面参数u的区间取值， 分别表示子曲面参数v的区间取值，和τ分别表示u和v的线性变换参数，参数u和v为：

计算B-spline曲面优化方程系数c_i,j：

$c_{i,j}=\underset{r=1}{\overset{n}{\mathrm{\&Sigma;}}}{c}_{r,i,j}$

其中，n表示划分成Bézier曲面的个数。

（2-3）如果输入曲面是NURBS曲面，取样2m个点在u区间内单调，2n个点在v区间内单调，

$\left\{\begin{array}{c}0=u_0<u_1<u_2<...<u_{2m}=1\\ 0=v_0<v_1<v_2<...<v_{2n}=1\end{array}\right.$

其中，u₀，u₁，…,u_2m为节点向量序列U中的节点，v₁,...,v_2n为节点向量序列V中的节点；使得

$I\&ap;T=\frac{h_1{g}_1}{4}\underset{i_1=0}{\overset{2m}{\mathrm{\&Sigma;}}}\underset{j_1=0}{\overset{2n}{\mathrm{\&Sigma;}}}{A}_{i_1}{A}_{j_1}{\left|\right|X_u(u_{i_1},v_{j_1})\&times;X_v(u_{i_1},v_{j_1})\left|\right|}^2{B}_i^2\left(u_{i_1}\right)B_j^2\left(v_{j_1}\right)---\left(1\right)$

其中： $\left\{\begin{array}{c}{h}_1=\frac{1}{2m}\\ g_1=\frac{1}{2n}\\ u_{i_1}=i_1{h}_1\\ v_{j_1}=j_1{g}_1\end{array}\right.,$ $A_{i_1}=\left\{\begin{array}{cc}1& i_1=\mathrm{0,2}m\\ 2& i_1=\mathrm{1,2},...2m-1\end{array}\right.,$

$A_{j_1}=\left\{\begin{array}{cc}1& j_1=\mathrm{0,2}n\\ 2& j_1=\mathrm{1,2},...2n-1\end{array}\right.,$

（2-4）使用双重积分的算法计算C_i,j。

具体步骤为：

（2-4-1）设置式（1）中u，v两个方向采样密度m，n的值m₀=n₀=10;

（2-4-2）计算式（1）中数值积分T的值T₀=T，其中T表示数值双重积分；

（2-4-3）令m_r=2m_r-1，n_r=2n_r-1，其中m_r，n_r分别为u/v方向的采样密度；

（2-4-4）计算T的新近似值T_r;

（2-4-5）如果|T_r-T_r-1|＜∈，结束，其中∈是整个计算过程的终止条件；否则迭代次数r=r+1，回到步骤（2-4-3）。

所述步骤三的具体步骤为：

（3-1）确定当前曲面参数化X(S,t)与保面积参数化的偏差：

$\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left(\right||X_s\&times;X_t||-A)}^2\mathrm{dsdt}---\left(2\right)$

其中，X_s表示X(s,t)在s方向上的偏导数，X_t表示X(s,t)在t方向上的偏导数，A表示给定曲面的面积，||·||表示三维空间的向量长度；

（3-2）为计算公式（2）的最小值，将公式等价为：

$J(\mathrm{\&alpha;},\mathrm{\&beta;})=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_s\&times;X_t\left|\right|}^2\mathrm{dsdt}$

（3-3）根据链式法则确定J（α，β），

$J(\mathrm{\&alpha;},\mathrm{\&beta;})=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_s\&times;X_t\left|\right|}^2\mathrm{dsdt}$

$=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_u\&times;X_v\left|\right|}^2\&CenterDot;{\left(\frac{{\&PartialD;}_u}{{\&PartialD;}_s}\right)}^2\&CenterDot;{\left(\frac{{\&PartialD;}_v}{{\&PartialD;}_t}\right)}^2\mathrm{dsdt}$

$=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_u\&times;X_v\left|\right|}^2\&CenterDot;\frac{{\&PartialD;}_u}{{\&PartialD;}_s}\&CenterDot;\frac{{\&PartialD;}_v}{{\&PartialD;}_t}\mathrm{dudv}$

$=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_u\&times;X_v\left|\right|}^2$

$\&CenterDot;\frac{{(1-\mathrm{\&alpha;}+2\mathrm{u\&alpha;}-u)}^2}{\mathrm{\&alpha;}(1-\mathrm{\&alpha;})}\&CenterDot;\frac{{(1-\mathrm{\&beta;}+2\mathrm{v\&beta;}-v)}^2}{\mathrm{\&beta;}(1-\mathrm{\&beta;})}\mathrm{dudv}$

$=\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{2}{\mathrm{\&Sigma;}}}\frac{c_{i,j}{B}_i^2\left(\mathrm{\&alpha;}\right)B_j^2\left(\mathrm{\&beta;}\right)}{\left(\begin{array}{c}2\\ i\end{array}\right)\left(\begin{array}{c}2\\ j\end{array}\right)\mathrm{\&alpha;}(1-\mathrm{\&alpha;})\mathrm{\&beta;}(1-\mathrm{\&beta;})}$

其中，X_u和X_v表示X(u,v)在u/v方向上的偏微分X_u(u,v)和X_v(u,v)，c_i,j是步骤二中计算得到的优化方程的系数，其中，J(α,β)的解满足：

$0=\frac{\&PartialD;J(\mathrm{\&alpha;},\mathrm{\&beta;})}{\&PartialD;\mathrm{\&alpha;}},$ $0=\frac{\&PartialD;J(\mathrm{\&alpha;},\mathrm{\&beta;})}{\&PartialD;\mathrm{\&beta;}}---\left(3\right)$

（3-4）确定J(α,β)的等价条件：

$0=\frac{\&PartialD;J(\mathrm{\&alpha;},\mathrm{\&beta;})}{\&PartialD;\mathrm{\&alpha;}}$

$=\frac{B_2^2\left(\mathrm{\&alpha;}\right)(c_{\mathrm{2,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{2,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{2,2}}{B}_2^2\left(\mathrm{\&beta;}\right))-B_0^2\left(\mathrm{\&alpha;}\right)(c_{\mathrm{0,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{0,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{0,2}}{B}_2^2\left(\mathrm{\&beta;}\right))}{{\mathrm{\&alpha;}}^2{(1-\mathrm{\&alpha;})}^2\mathrm{\&beta;}(1-\mathrm{\&beta;})}---\left(4\right);$

（3-5）确定公式（4）等价条件：

$\frac{B_2^2\left(\mathrm{\&alpha;}\right)}{B_0^2\left(\mathrm{\&alpha;}\right)}=\frac{c_{\mathrm{0,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{0,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{0,2}}{B}_2^2\left(\mathrm{\&beta;}\right)}{c_{\mathrm{2,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{2,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{2,2}}{B}_2^2\left(\mathrm{\&beta;}\right)};$

（3-5）计算出J(α,β)对β偏微分为0的等价条件：

$\frac{B_2^2\left(\mathrm{\&beta;}\right)}{B_0^2\left(\mathrm{\&beta;}\right)}=\frac{c_{\mathrm{0,0}}{B}_0^2\left(\mathrm{\&alpha;}\right)+\frac{c_{\mathrm{1,0}}{B}_1^2\left(\mathrm{\&alpha;}\right)}{2}+c_{\mathrm{2,0}}{B}_2^2\left(\mathrm{\&alpha;}\right)}{c_{\mathrm{0,2}}{B}_0^2\left(\mathrm{\&alpha;}\right)+\frac{c_{\mathrm{1,2}}{B}_1^2\left(\mathrm{\&alpha;}\right)}{2}+c_{\mathrm{2,2}}{B}_2^2\left(\mathrm{\&alpha;}\right)}$

对两个偏导公式（3）表示的两条平面代数曲线进行求交，确定交点坐标即为最后所求的最优线性变换的系数α和β。

所述步骤四的具体步骤为：

（4-1）对参数u，v分别进行变换：

$u=f_1\left(s\right)=\frac{(\mathrm{\&alpha;}-1)s}{2\mathrm{\&alpha;s}-s-\mathrm{\&alpha;}},\mathrm{\&alpha;}\&Element;\left(\mathrm{0,1}\right)$

$v=f_2\left(t\right)=\frac{(\mathrm{\&beta;}-1)t}{2\mathrm{\&beta;t}-t-\mathrm{\&beta;}},\mathrm{\&beta;}\&Element;\left(\mathrm{0,1}\right)$

其中，s,t是对u，v变换后的参数；分别把参数s,t映射到单位区间上,使之满足：

f₁(0)=0,，f₁(1)=1,f₁(0.5)=1-α

f₂(0)=0,f₂(1)=1,f₂(0.5)=1-β

由变换方程，确定参数s和t:

$s=f_3\left(u\right)=\frac{\mathrm{\&alpha;u}}{1-\mathrm{\&alpha;}+2\mathrm{\&alpha;u}-u}$

$t=f_4\left(v\right)=\frac{\mathrm{\&beta;v}}{1-\mathrm{\&beta;}+2\mathrm{\&beta;v}-v};$

其中，α和β表示变换的两个自由参数。

（4-2）把s和t的表达式代入NURBS曲面的方程形式，确定重新参数化曲面X(s,t):

$X(s,t)=\frac{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(s\right)N_j^q\left(t\right){\widetilde{\mathrm{\&omega;}}}_{i,j}{P}_{i,j}}{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(s\right)N_j^q\left(t\right){\widetilde{\mathrm{\&omega;}}}_{i,j}},s,t\&Element;\left[\mathrm{0,1}\right]$

其中，是曲面X(s,t)中序号为(i,j)的控制顶点新的权重，

${\widetilde{\mathrm{\&omega;}}}_{i,j}=\frac{{\mathrm{\&omega;}}_{i,j}}{{\mathrm{\&Pi;}}_{k=1}^m(1-\mathrm{\&alpha;}+2{\mathrm{\&alpha;u}}_{i+k}-u_{i+k}){\mathrm{\&Pi;}}_{l=1}^n(1-\mathrm{\&beta;}+2\mathrm{\&beta;}{u}_{j+l}-u_{j+l})},$

确定新的节点向量和

本发明的方法以Windows7中的Visual Studio2008为开发平台，其具体实施方式如下：

图1是输入的一个原始NURBS曲面S，对于曲面S(u,v)其形式为：

$S(u,v)=\frac{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\&omega;}}_{i,j}{P}_{i,j}}{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\&omega;}}_{i,j}},u,v\&Element;\left[\mathrm{0,1}\right]$

其中u和v是曲面S(u,v)的参数，n_u和n_v，是曲面S(u,v)所含的控制顶点在u方向和v方向上的个数，P_i,j是曲面S(u,v)中序号为(i，j)的控制顶点的坐标值，由用户输入，i＝0，1，…，n_u-1；j＝0，1，…，n_v-1，采用三维直角坐标表示，ω_i,j是曲面S(u，v)中序号为(i，j)的控制顶点的权重，由用户输入，i＝0，1，…，n_u-1；j＝0，1，…，n_v-1，是定义在节点向量序列U之上的序号为i的p阶B-spline基函数，i＝0，1，…，n_u-1；是定义在节点向量序列V之上的序号为j的q阶B-spline基函数，j＝0，1，…，n_v-1，

其中p，q是曲面S(u,v)在u方向和v方向上的次数，由用户输入，u₀，u₁,…,u_n+p等为节点向量序列U中的节点，其中u₀＝u₁＝…＝u_p＝0，u_n＝u_n+1＝…＝u_n+p＝1，u_p+1，u_n+2，…,u_n-1由用户输入；v₀,v₁,...,v_n+q等为节点向量序列V中的节点，其中v₀＝v₁＝…＝v_q＝0，v_n＝v_n+1＝…＝v_n+q＝1，v_q+1，v_q+2，…,v_n-1由用户输入，上述的B-spline基函数采用Matlab系统样条工具库中B-spline基函数的定义方式。

图2为对NURBS曲面S(u,v)进行保面积参数优化的最优结果；

（1）计算NURBS曲面的优化方程的系数，其步骤如下：

取样2m个点在u区间内单调，2n个点在v区间内单调，计算如下：

$\left\{\begin{array}{c}0=u_0<u_1<u_2<...<u_{2m}=1\\ 0=v_0<v_1<v_2<...<v_{2n}=1\end{array}\right.$

使得

$I\&ap;T=\frac{h_1{g}_1}{4}\underset{i_1=0}{\overset{2m}{\mathrm{\&Sigma;}}}\underset{j_1=0}{\overset{2n}{\mathrm{\&Sigma;}}}{A}_{i_1}{A}_{j_1}{\left|\right|X_u(u_{i_1},v_{j_1})\&times;X_v(u_{i_1},v_{j_1})\left|\right|}^2{B}_i^2\left(u_{i_1}\right)B_j^2\left(v_{j_1}\right)---(1-1)$

满足：

$\left\{\begin{array}{c}{h}_1=\frac{1}{2m}\\ g_1=\frac{1}{2n}\\ u_{i_1}=i_1{h}_1\\ v_{j_1}=j_1{g}_1\end{array}\right.,$ $A_{i_1}=\left\{\begin{array}{cc}1& i_1=\mathrm{0,2}m\\ 2& i_1=\mathrm{1,2},...2m-1\end{array}\right.,$

$A_{j_1}=\left\{\begin{array}{cc}1& j_1=\mathrm{0,2}n\\ 2& j_1=\mathrm{1,2},...2n-1\end{array}\right.,$

使用双重积分的算法计算C_i，j,其计算步骤如下：

第一步：设置式（1-1）中u,v两个方向采样密度m,n的值m₀=n₀=10;

第二步：计算式（1-1）中数值积分T的值T₀=T,其中T表示数值双重积分；

第三步：令m_r=2m_r-1,n_r=2n_r-1,其中m_r，n_r分别为u/v方向的采样密度；

第四步：计算T的新近似值T_r;

第五步：如果|T_r-T_r-1|＜∈，结束，其中∈是整个计算过程的终止条件；否则迭代次数r=r+1，回到第三步。

（2）根据公式变换，计算当前曲面参数化X(s,t)与保面积参数化的偏差，偏差计算方法如下：

$\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left(\right||X_s\&times;X_t||-A)}^2\mathrm{dsdt}---(2-1)$

其中，X_s表示X(s,t)在s方向上的偏导数，X_t表示X(s,t)在t方向上的偏导数，A表示给定曲面的面积，||·||表示三维空间的向量长度。

计算偏差值的最小值，等价于计算J(α,β)的最小值，计算方法如下：

$J(\mathrm{\&alpha;},\mathrm{\&beta;})=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_s\&times;X_t\left|\right|}^2\mathrm{dsdt}$

根据链式法则，进一步计算J(α,β)，计算方法如下：

$J(\mathrm{\&alpha;},\mathrm{\&beta;})=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_s\&times;X_t\left|\right|}^2\mathrm{dsdt}$

$=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_u\&times;X_v\left|\right|}^2\&CenterDot;{\left(\frac{{\&PartialD;}_u}{{\&PartialD;}_s}\right)}^2\&CenterDot;{\left(\frac{{\&PartialD;}_v}{{\&PartialD;}_t}\right)}^2\mathrm{dsdt}$ $=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_u\&times;X_v\left|\right|}^2\&CenterDot;\frac{{\&PartialD;}_u}{{\&PartialD;}_s}\&CenterDot;\frac{{\&PartialD;}_v}{{\&PartialD;}_t}\mathrm{dudv}$

$=\underset{0}{\overset{1}{\&Integral;}}\underset{0}{\overset{1}{\&Integral;}}{\left|\right|X_u\&times;X_v\left|\right|}^2$ $\&CenterDot;\frac{{(1-\mathrm{\&alpha;}+2\mathrm{u\&alpha;}-u)}^2}{\mathrm{\&alpha;}(1-\mathrm{\&alpha;})}\&CenterDot;\frac{{(1-\mathrm{\&beta;}+2\mathrm{v\&beta;}-v)}^2}{\mathrm{\&beta;}(1-\mathrm{\&beta;})}\mathrm{dudv}$

$=\underset{i=0}{\overset{2}{\mathrm{\&Sigma;}}}\underset{j=0}{\overset{2}{\mathrm{\&Sigma;}}}\frac{c_{i,j}{B}_i^2\left(\mathrm{\&alpha;}\right)B_j^2\left(\mathrm{\&beta;}\right)}{\left(\begin{array}{c}2\\ i\end{array}\right)\left(\begin{array}{c}2\\ j\end{array}\right)\mathrm{\&alpha;}(1-\mathrm{\&alpha;})\mathrm{\&beta;}(1-\mathrm{\&beta;})}$

其中，J(α,β)的解满足：

$0=\frac{\&PartialD;J(\mathrm{\&alpha;},\mathrm{\&beta;})}{\&PartialD;\mathrm{\&alpha;}},$

$0=\frac{\&PartialD;J(\mathrm{\&alpha;},\mathrm{\&beta;})}{\&PartialD;\mathrm{\&beta;}},---(2-2)$

计算出其等价条件，计算方法如下：

$0=\frac{\&PartialD;J(\mathrm{\&alpha;},\mathrm{\&beta;})}{\&PartialD;\mathrm{\&alpha;}}$

$=\frac{B_2^2\left(\mathrm{\&alpha;}\right)(c_{\mathrm{2,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{2,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{2,2}}{B}_2^2\left(\mathrm{\&beta;}\right))-B_0^2\left(\mathrm{\&alpha;}\right)(c_{\mathrm{0,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{0,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{0,2}}{B}_2^2\left(\mathrm{\&beta;}\right))}{{\mathrm{\&alpha;}}^2{(1-\mathrm{\&alpha;})}^2\mathrm{\&beta;}(1-\mathrm{\&beta;})}---(2-3)$

计算出（2-3）式的等价条件，计算方法如下：

$\frac{B_2^2\left(\mathrm{\&alpha;}\right)}{B_0^2\left(\mathrm{\&alpha;}\right)}=\frac{c_{\mathrm{0,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{0,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{0,2}}{B}_2^2\left(\mathrm{\&beta;}\right)}{c_{\mathrm{2,0}}{B}_0^2\left(\mathrm{\&beta;}\right)+\frac{c_{\mathrm{2,1}}{B}_1^2\left(\mathrm{\&beta;}\right)}{2}+c_{\mathrm{2,2}}{B}_2^2\left(\mathrm{\&beta;}\right)}$

计算出J(α,β)对β偏微分为0的等价条件，计算方法如下：

$\frac{B_2^2\left(\mathrm{\&beta;}\right)}{B_0^2\left(\mathrm{\&beta;}\right)}=\frac{c_{\mathrm{0,0}}{B}_0^2\left(\mathrm{\&alpha;}\right)+\frac{c_{\mathrm{1,0}}{B}_1^2\left(\mathrm{\&alpha;}\right)}{2}+c_{\mathrm{2,0}}{B}_2^2\left(\mathrm{\&alpha;}\right)}{c_{\mathrm{0,2}}{B}_0^2\left(\mathrm{\&alpha;}\right)+\frac{c_{\mathrm{1,2}}{B}_1^2\left(\mathrm{\&alpha;}\right)}{2}+c_{\mathrm{2,2}}{B}_2^2\left(\mathrm{\&alpha;}\right)}$

对两个偏导公式（2-2）表示的两条平面代数曲线进行求交，求得的交点坐标即为最后所求的最优线性变换的系数α和β。

（3）对NURBS空间曲面进行变换，其步骤如下：

对参数u，v分别进行变换，计算方法如下：

$u=f_1\left(s\right)=\frac{(\mathrm{\&alpha;}-1)s}{2\mathrm{\&alpha;s}-s-\mathrm{\&alpha;}},\mathrm{\&alpha;}\&Element;\left(\mathrm{0,1}\right)$

$v=f_2\left(t\right)=\frac{(\mathrm{\&beta;}-1)t}{2\mathrm{\&beta;t}-t-\mathrm{\&beta;}},\mathrm{\&beta;}\&Element;\left(\mathrm{0,1}\right)$

其中，s，t是对u，v变换后的参数。

分别把参数s，t映射到单位区间上,使之满足：

$\left\{\begin{array}{c}{f}_1\left(0\right)=0,f_1\left(1\right)=1,f_1\left(0.5\right)=1-\mathrm{\&alpha;}\\ f_2\left(0\right)=0,f_2\left(1\right)=1,f_2\left(0.5\right)=1-\mathrm{\&beta;}\end{array}\right.$

由变换方程，反求参数s和t,得到表达式：

$s=f_3\left(u\right)=\frac{\mathrm{\&alpha;u}}{1-\mathrm{\&alpha;}+2\mathrm{\&alpha;u}-u}$

$t=f_4\left(v\right)=\frac{\mathrm{\&beta;v}}{1-\mathrm{\&beta;}+2\mathrm{\&beta;v}-v}$

把s和t的表达式代入NURBS曲面的方程形式，计算重新参数化曲面X(s,t),计算方法如下：

$X(s,t)=\frac{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(s\right)N_j^q\left(t\right){\widetilde{\mathrm{\&omega;}}}_{i,j}{P}_{i,j}}{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(s\right)N_j^q\left(t\right){\widetilde{\mathrm{\&omega;}}}_{i,j}},s,t\&Element;\left[\mathrm{0,1}\right]$

其中是曲面X(s,t)中序号为(i,j)的控制顶点新的权重，计算方法如下：

${\widetilde{\mathrm{\&omega;}}}_{i,j}=\frac{{\mathrm{\&omega;}}_{i,j}}{{\mathrm{\&Pi;}}_{k=1}^m(1-\mathrm{\&alpha;}+2{\mathrm{\&alpha;u}}_{i+k}-u_{i+k}){\mathrm{\&Pi;}}_{l=1}^n(1-\mathrm{\&beta;}+2\mathrm{\&beta;}{u}_{j+l}-u_{j+l})}$

产生新的节点向量和计算方法如下：

上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。

Claims (7)

1.一种自由曲面的保面积参数特性优化方法，其特征是简单有效、在保持曲面次数不变，不破坏曲面形状的前提下优化曲面的保面积参数特性，具体步骤为：

步骤一：输入空间自由曲面；

步骤二：根据自由曲面类型的不同，分别计算空间自由曲面优化方程的系数；

步骤三：通过计算两条平面参数曲线的交点获取自由曲面最优变换；

步骤四：对空间自由曲面进行变换；

所述步骤三的具体步骤为：

(3-1)确定当前曲面参数化与保面积参数化的偏差；

(3-2)利用公式等价法确定当前曲面参数化与保面积参数化偏差的最小值；

(3-3)通过计算两条平面曲线的交点确定最优参数变换的系数。

2.如权利要求1所述一种自由曲面的保面积参数特性优化方法，所述步骤一中空间自由曲面的表达式为：

$S(u,v)=\frac{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_v-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\&omega;}}_{i,j}{P}_{i,j}}{{\mathrm{\&Sigma;}}_{i=0}^{n_u-1}{\mathrm{\&Sigma;}}_{j=0}^{n_u-1}{N}_i^p\left(u\right)N_j^q\left(v\right){\mathrm{\&omega;}}_{i,j}},u,v\&Element;\left[\mathrm{0,1}\right]$

其中，u和v是曲面S(u，v)的参数，n_u和n_v是曲面S(u，v)所含的控制顶点在u方向和v方向上的个数，P_i，j是曲面S(u，v)中序号为(i，j)的控制顶点的坐标值，ω_i，j是曲面S(u，v)中序号为(i，j)的控制顶点的权重，是定义在节点向量序列U之上的序号为i的p阶B-spline基函数；是定义在节点向量序列V之上的序号为j的q阶B-spline基函数，i＝0，1，…，n_u-1；j＝0，1，…，n_v-1，i和j均为大于等于0的整数，

其中p，q是曲面S(u，v)在u方向和v方向上的次数，u₀，u₁，…，u_n+p为节点向量序列U中的节点，其中u₀＝u₁＝…＝u_p＝0，u_n＝u_n+1＝…＝u_n+p＝1，0＜u_p+1，u_n+2，…，u_n-1＜1；v₀，v₁，…，v_n+q为节点向量序列V中的节点，其中v₀＝v₁＝…＝v_q＝0，v_n＝v_n+1＝…＝v_n+q＝1，0＜v_q+1，v_q+2，…，v_n-1＜1，所述的B-spline基函数采用Matlab系统样条工具库中B-spline基函数的定义方式。

3.如权利要求l所述一种自由曲面的保面积参数特性优化方法，所述步骤二的具体步骤为：

如果输入的曲面是Bézier曲面，首先确定Bézier曲面的表示形式及确定Bézier曲面在向量序列u/v方向上的偏微分，然后确定Bézier曲面的微分第一基本形式；

根据Bézier曲面的微分第一基本形式确定优化方程的系数c_i，j。

4.如权利要求l所述一种自由曲面的保面积参数特性优化方法，所述步骤二的具体步骤为：如果输入曲面是B-spline曲面，通过节点插入法把B-spline曲面细分成n个Bézier曲面，用表示第r个Bézier曲面的局部参数域，确定第r个子曲面的优化方程系数；利用所有子曲面的优化方程系数确定B-spline曲面优化方程系数c_i，j。

5.如权利要求l所述一种自由曲面的保面积参数特性优化方法，所述步骤二的具体步骤为：

(2-1)如果输入曲面是NURBS曲面，取样2m个点在u区间内单调，2n个点在v区间内单调，

$\left\{\begin{array}{c}0=u_0<u_1<u_2<...<u_{2m}=1\\ 0=v_0<v_1<v_2<...<v_{2n}=1\end{array}\right.$

其中，u₀，u₁，…，u_2m为u方向采样点，v₁，…，V_2n为v方向采样点；使得

$I\&ap;T=\frac{h_1{g}_1}{4}\underset{i_2=0}{\overset{2m}{\mathrm{\&Sigma;}}}\underset{j_2=0}{\overset{2n}{\mathrm{\&Sigma;}}}{A}_{i_2}{A}_{j_2}{\left|\right|X_u(u_{i_1},v_{j_1})\&times;X_v(u_{i_1},v_{j_1})\left|\right|}^2{B}_i^2\left(u_{i_1}\right)B_j^2\left(v_{j_2}\right)---\left(1\right)$

满足： $\left\{\begin{array}{c}{h}_1=\frac{1}{2m}\\ g_1=\frac{1}{2n}\\ u_{i_1}=i_1{h}_1\\ v_{j_1}=j_1{g}_1\end{array}\right.,A_{i_1}=\left\{\begin{array}{cc}1& i_1=\mathrm{0,2}m\\ 2& i_1=\mathrm{1,2},...2m-1\end{array}\right.,$

$A_{j_1}=\left\{\begin{array}{cc}1& j_1=\mathrm{0,2}n\\ 2& j_1=\mathrm{01,2},...2n-1\end{array}\right.;$

(2-2)使用双重积分的算法计算c_i，j。

6.如权利要求l所述一种自由曲面的保面积参数特性优化方法，所述步骤四的具体步骤为：

(4-1)对参数u，v分别进行变换：

设s，t是对u，v变换后的参数；分别把参数s，t映射到单位区间上，使之满足：s和t为0时，u和v也为0；s和t为1时，u和v也为1；s和t为0.5时，u和v分别为1-α和1-β，由变换方程，确定参数s和t的表达式，α和β表示变换的两个自由参数；

(4-2)把s和t的表达式代入自由空间曲面的方程形式，确定重新参数化曲面X(s，t)，

确定重新参数化曲面的节点向量和。

7.如权利要求5所述一种自由曲面的保面积参数特性优化方法，所述步骤(2-2)的具体步骤为：

(2-2-1)设置式(1)中u，v两个方向采样密度m，n的值m₀＝n₀＝10；

(2-2-2)确定式(1)中数值积分T的值T₀＝T，其中T表示数值双重积分；

(2-2-3)令m_r＝2m_r-1，n_r＝2n_r-l，其中m_r，n_r分别为u/v方向的采样密度；

(2-2-4)确定T的新近似值T_r；

(2-2-5)如果|T_r-T_r-l|＜∈，结束；其中∈是整个计算过程的终止条件；否则迭代次数r＝r+1，回到步骤(2-2-3)。