一种半径可控的参数化三维前缘钝化设计方法与流程
本发明涉及高超声速飞行器技术领域，尤其涉及一种半径可控的参数化三维前缘钝化设计方法。
背景技术：
在高超声速飞行时，空气受到剧烈压缩而产生的气动加热使得飞行器面临的巨大的热防护压力，因此高超声速飞行器在工程应用中必须要对存在的尖锐前缘进行钝化处理。在旋成体、升力体、翼身融合体等常见高超声速飞行器气动构型之中，乘波体具有升阻比大、设计过程简单等优点，是一种应用前景广阔的高超声速飞行器气动构型方案。但乘波体的设计原理决定了其外形必然带有尖锐前缘，给热防护和加工制造带来一定困难，在设计时有必要对其前缘进行钝化处理。
目前，虽然已有相应三维钝化方法可使钝化前缘与飞行器表面光滑连接，但尚不能指定各前缘点处的钝化半径，故在设计飞行器气动外形时无法根据热防护需求所提出的前缘约束对钝化前缘进行针对性设计。
技术实现要素：
针对现有技术存在的上述缺陷，本发明目的在于提供一种半径可控的参数化三维前缘钝化设计方法，其能够解决现有的前缘钝化方法不能兼顾钝化前缘与飞行器上下表面光滑连接与指定各前缘点钝化半径的问题。本方法使钝化前缘能在满足与飞行器上下表面光滑连接且钝化前缘形状沿横向连续变化的同时，还可根据热防护需求来设计各前缘点处不同的钝化半径，有效地实现了热防护与气动外形的多学科耦合设计。
为实现本发明之目的，采用以下技术方案予以实现：
一种半径可控的参数化三维前缘钝化设计方法，首先给定初始乘波体外形，然后对初始乘波体外形中的尖锐前缘进行钝化。给定的初始乘波体外形，既可是吻切锥乘波体，但不局限于该乘波体或锥导乘波体等其他类型的乘波体外形，还可用于其他带有尖锐前缘的飞行器外形。也就是本发明接下来重点介绍的钝化方法可以应用于带有尖锐前缘的乘波体或者是带有尖锐前缘的飞行器，实现乘波体或飞行器的三维前缘钝化。
具体地，对尖锐前缘进行钝化，包括以下步骤：
s1将乘波体上表面抬升一定高度；
给定初始乘波体外形，即已知乘波体未钝化时的原始前缘点为p0和未钝化前的乘波体纵向截面曲线；未钝化前的乘波体纵向截面曲线由前缘点、上表面曲线和下表面曲线连接而成。
设对乘波体的上表面沿法向抬升高度h，抬升后的乘波体上表面端点为p1，抬升后的乘波体下表面端点为乘波体未钝化前的原始前缘点p0。抬升后的乘波体纵向截面曲线由上表面曲线、下表面曲线以及钝化前缘曲线连接而成，其中上表面曲线和下表面曲线即是未钝化前的乘波体纵向截面曲线上的上表面曲线以及和下表面曲线，因此是已知的，钝化前缘曲线是待求解的；抬升后的乘波体上表面端点p1是上表面曲线上的点，抬升后的乘波体下表面端点p0是下表面曲线上的点。
s2获取抬升后的乘波体前缘点p2的插值条件；
抬升后的乘波体前缘点p2的插值条件为抬升后的乘波体下表面端点p0处的切矢v0、曲率k0和抬升后的乘波体上表面端点p1处切矢v1和曲率k1。
s3给出抬升后的钝化前缘曲线上不同横向站位处的钝化半径。
s4设置抬升后的乘波体前缘点p2以及乘波体的上表面沿法向抬升高度h。
s5求出所有钝化前缘曲线，然后放样得到钝化前缘曲面。
本发明s2中，对上表面曲线进行均匀离散，获得多个离散点，获得的多个离散点中包括抬升后的乘波体上表面端点p1，对多个离散点进行差分，能够获得抬升后的乘波体上表面端点p1的1阶和2阶导矢，进而获得抬升后的乘波体上表面端点p1处的切矢v1和曲率k1；
同样的，对下表面曲线进行均匀离散，获得多个离散点，获得的多个离散点中包括抬升后的乘波体上表面端点p1，对多个离散点进行差分，能够获得抬升后的乘波体下表面端点p0的1阶和2阶导矢，进而获得抬升后的乘波体下表面端点p0处的切矢v0和曲率k0。
本发明s3中，未钝化前的乘波体纵向截面曲线中的前缘曲线已知，故前缘曲线上各横向站位处的前缘后掠角已知，前缘曲线上各横向站位均一一对应钝化前缘曲线上的横向站位，钝化前缘曲线上各横向站位处的当地前缘后掠角与前缘曲线上与之对应的横向站位处的当地前缘后掠角看作相同，即钝化前缘曲线上各横向站位处的当地前缘后掠角λ已知。
设定抬升后的乘波体前缘点p2处的前缘钝化半径，由预先设定的抬升后的乘波体前缘点p2处的前缘钝化半径和钝化前缘曲线上各横向站位处的当地前缘后掠角λ来确定抬升后的钝化前缘曲线上各横向站位处的前缘钝化半径，如式(3)所示：
rswept＝rsymmetry(cosλ)2.2(3)
其中，rswept为钝化前缘曲线上各横向站位处的前缘钝化半径；rsymmerty为预先设定的抬升后的乘波体前缘点p2处的前缘钝化半径。
本发明s4中，抬升后的乘波体前缘点p2和乘波体的上表面沿法向抬升高度h的设置应满足以下两个要求：
1)抬升后的乘波体前缘点p2和乘波体的上表面沿法向抬升高度h的设置能使乘波体钝化前缘曲线能够成功生成；
2)抬升后的乘波体前缘点p2的曲率大于乘波体钝化前缘曲线上各横向站位处的曲率。
进一步地，本发明将乘波体的上表面沿法向抬升高度h设置为当地前缘钝化半径的1～2倍。
抬升后的乘波体前缘点p2设置在上表面曲线以及下表面曲线的延长线所构成夹角的角平分线上，抬升后的乘波体前缘点p2距p0点的水平距离l可由式(4)得出：
l＝rswept/cosλ(4)
其中，rswept为当地前缘钝化半径。
本发明s5中，根据点p0和p1以及s2获取的抬升后的乘波体前缘点p2的插值条件求出插值曲线即钝化前缘曲线；将求解钝化前缘曲线的问题看作插值于两个端点的ghi问题，利用3次bézier曲线构造钝化前缘曲线，求出所有的钝化前缘曲线，选取其中同时满足s4中要求的钝化前缘曲线进行放样，从而得到所需的钝化前缘曲面。具体地，在本发明中是以p0、p2为端点利用一段3次bézier曲线构造点p0和p2间的钝化前缘曲线，以p2、p1为端点利用一段3次bézier曲线构造点p2和p1间的钝化前缘曲线，完整的钝化前缘曲线即由点p0和p2间的钝化前缘曲线和点p2和p1间的钝化前缘曲线连接而成，求出所有的钝化前缘曲线，选取其中同时满足s4中要求的钝化前缘曲线进行放样，从而得到所需的钝化前缘曲面。
3次bézier曲线共需要4个控制点来进行曲线的求解，因此首先根据已知的两个控制点(即点p0和p2或者点p2和p1)求得插值曲线中的另外两个控制点，然后根据已知的两个控制点以及求解得到的另外两个控制点，即4个控制点来进行插值曲线的求解。
下面介绍如何以p0、p2为端点利用一段3次bézier曲线构造点p0和p2间的钝化前缘曲线，方法如下：
已知的两个控制点p0和p2(即两个端点)与待求的两个控制点p3和p4，这4个控制点间的关系可由下式表示
其中：δ0和δ1为设定的中间参数，为待求值；v0和v2分别为端点p0和p2处的切矢。
两个端点p0和p2处的曲率为
由s2获取的抬升后的乘波体前缘点p2的插值条件可以得到关于δ0和δ1的方程组
其中d是两个端点p0和p2之间的位移矢量。
通过求解方程组(7)可以得到未知参数δ0和δ1的值，从而计算得到其余两个控制点p3和p4的坐标，进而得到一条曲率连续的插值曲线，即利用3次bézier曲线构造的点p0和p2间的初步钝化前缘曲线。
方程组(7)的解可能有0至3个解。由于乘波体钝化前缘曲线为三维空间曲线，钝化前缘曲线上各个横向站位处的端点插值条件不一，可能导致方程无解或者多解。如要将该方法应用于三维外形上，则对于方程组解的讨论将是主要障碍之一。此处对插值条件进行简化，即不要求钝化前缘曲线与上表面曲线、下表面曲线间的曲率连续，将钝化前缘曲线与上表面曲线、下表面曲线连接处的曲率值设置为0。
方程组(7)可简化为
原方程组(7)变为一个一元一次方程与一个二元二次方程的组合，不难看出，方程组(8)必然有且只有一个解。通过对抬升后的乘波体前缘点p2的插值条件的调整，以降低钝化前缘曲线与上表面曲线、下表面曲线连接的光滑性为代价，简化了bézier曲线的求解过程，有利于保持不同横向站位处钝化前缘曲线之间的连续性。
以p2、p1为端点利用一段3次bézier曲线构造点p2和p1间的钝化前缘曲线的方法与以p0、p2为端点利用一段3次bézier曲线构造点p0和p2间的钝化前缘曲线的方法完全相同，在此不再赘述。
本发明的有益技术效果是：
本发明基于增加材料的前缘钝化思路和基于bézier曲线和有理bézier曲线的前缘钝化方法的半径可控的参数化三维前缘钝化设计方法，在保证钝化前缘曲线与上、下表面曲线光滑连接的同时，还可以根据热防护需求指定各前缘点处的钝化半径，有效地实现了热防护与气动外形的多学科耦合设计。
综上所述，该方法主要有以下几个特点：
1)乘波体前缘钝化过程实现了参数化，有利于程序实现和优化设计，减少了飞行器设计的工作量；
2)钝化前缘曲面与乘波体表面光滑连接，可以避免飞行器表面不光滑对性能所带来的不利影响；
3)前缘点以及前缘点处的钝化半径可以直接指定，使得气动外形设计与热防护系统设计可以实现有效衔接。
附图说明
图1是实施例给定的乘波体外形形体示意图；
图2是移除材料和增加材料两种前缘钝化思路对比示意图，上部分为乘波体尖锐前缘的侧视图和后视图，中间部分为移除材料钝化前缘示意图，下面部分是增加材料钝化前缘示意图；
图3是乘波体上表面“抬升”前后的乘波体纵向截面示意图；
图4是乘波体上表面抬升后的新前缘点p2设置示意图，其中l为钝化前后前缘点的水平距离；
图5是三次bézier曲线示意图，其中：p0、p1、p2、p3为设计插值曲线的四个控制点，v0、k0分别是p0点处的切矢和曲率，相应地，v2、k2分别是p2点处的切矢和曲率，d为p0和p2之间的的位移矢量；
图6是利用2段3次bézier曲线构造钝化前缘曲线示意图，p0、p2和p2、p1分别是2段3次bézier曲线的端点，v、k则分别是各点处的切矢和曲率；
图7是脊形吻切锥乘波体等轴测视图及头部钝化细节图；
图8是有粘和无粘升力系数随攻角变化曲线；
图9是有粘和无粘阻力系数随攻角变化曲线；
图10是有粘和无粘俯仰力矩系数随攻角变化曲线；
图11是有粘和无粘升阻比数随攻角变化曲线。
具体实施方式
下面将结合本发明实施例图中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，做进一步详细说明，但本发明的实施方式不仅限于此。
一种半径可控的参数化三维前缘钝化设计方法，其首先给定乘波体外形。本实施例以脊形吻切锥乘波体为例进行前缘钝化设计，参照图1，为本实施例给定的乘波体外形形体示意图。特别强调，本发明所提出的钝化方法，既可用于吻切锥乘波体，但不局限于该乘波体或锥导乘波体等其他类型的乘波体，还可用于其他带有尖锐前缘的飞行器外形，本实施例仅以吻切锥乘波体为例进行说明。
现有的各种乘波体前缘钝化方法在材料的增减方面存在差异：steven提出了通过移除材料和增加材料两种途径进行前缘钝化的思路，参见stevensdr.practicalconsiderationsinwaveriderapplications[m].aiaaaircraftdesignsystemsmeeting.1992.。如图2所示，其中移除材料的方法被一些学者在研究中所采用。考虑到乘波体的前缘较薄，移除材料的方法可能对升力面的几何形状会产生较大影响，进而使得升阻比和容积减小。然而实际上采用增加材料的方法对乘波体容积、升力面和流场结构的影响较小；本发明采用增加材料的钝化思路，沿流向设计各个横向站位上钝化前缘曲线。
具体的，对给定的乘波体外形中的尖锐前缘进行钝化，包括以下步骤：
s1将乘波体上表面抬升一定高度
给定初始乘波体外形，即已知乘波体未钝化时的原始前缘点为p0和未钝化前的乘波体纵向截面曲线；未钝化前的乘波体纵向截面曲线由前缘点、上表面曲线和下表面曲线连接而成。
设对乘波体的上表面沿法向抬升高度h，抬升后的乘波体上表面端点为p1，抬升后的乘波体下表面端点为乘波体未钝化前的原始前缘点p0。抬升后的乘波体纵向截面曲线由上表面曲线、下表面曲线以及钝化前缘曲线连接而成，其中上表面曲线和下表面曲线即是未钝化前的乘波体纵向截面曲线上的上表面曲线以及和下表面曲线，因此是已知的，钝化前缘曲线是待求解的；抬升后的乘波体上表面端点p1是上表面曲线上的点，抬升后的乘波体下表面端点p0是下表面曲线上的点。参照图3是乘波体上表面“抬升”前后的乘波体纵向截面示意图。
s2获取抬升后的乘波体前缘点p2的插值条件
抬升后的乘波体前缘点p2的插值条件为抬升后的乘波体下表面端点p0处的切矢v0、曲率k0和抬升后的乘波体上表面端点p1处切矢v1和曲率k1。
对上表面曲线进行均匀离散，获得多个离散点，获得的多个离散点中包括抬升后的乘波体上表面端点p1，对多个离散点进行差分，能够获得抬升后的乘波体上表面端点p1的1阶和2阶导矢，进而获得抬升后的乘波体上表面端点p1处的切矢v1和曲率k1；
同样的，对下表面曲线进行均匀离散，获得多个离散点，获得的多个离散点中包括抬升后的乘波体上表面端点p1，对多个离散点进行差分，能够获得抬升后的乘波体下表面端点p0的1阶和2阶导矢，进而获得抬升后的乘波体下表面端点p0处的切矢v0、曲率k0；
其中：具体的差分方法参见杨炯,宁涛,席平.前缘点曲率可控的曲率连续前缘几何设计[j].计算机辅助设计与图形学学报,2016,28(7):1195-200.以抬升后的乘波体上表面端点p1为例，具体的差分方法为：由于上表面曲线已知，所以可以将上表面曲线离散为多个点p1、p2、p3、p4、p5…pn，pi＝(xi，yi)，根据导矢的定义，点pi对应的一阶导矢可由点位置的差分计算得到:
一阶差分：δf(xi+1)＝f(xi+2)-f(xi)
一阶导矢：ui+1＝(2h，yi+2-yi)，其中h为相邻离散点间的x轴向距离
二阶差分：δ2f(xi+1)＝f(xi+2)+f(xi)-2f(xi+1)
二阶导矢：u′(xi+1)＝(2h，yi+2+yi-2yi+1)
即可得到插值条件。
s3给出抬升后的钝化前缘曲线上不同横向站位处的钝化半径；
抬升后的乘波体前缘点p2(乘波体纵向对称面处的驻点)是热防护压力最大的部位，该处的前缘钝化半径应由热防护材料和飞行工况等条件确定，是预先设定好的。
未钝化前的乘波体纵向截面曲线中的前缘曲线已知，故前缘曲线上各横向站位处的前缘后掠角已知，前缘曲线上各横向站位均一一对应钝化前缘曲线上的横向站位，钝化前缘曲线上各横向站位处的当地前缘后掠角与前缘曲线上与之对应的横向站位处的当地前缘后掠角相同，即钝化前缘曲线上各横向站位处的当地前缘后掠角λ已知。
由预先设定的抬升后的乘波体前缘点p2处的前缘钝化半径和钝化前缘曲线上各横向站位处的当地前缘后掠角λ来确定钝化前缘曲线上各横向站位处(即不同横向坐标下前缘曲线上的点)的前缘钝化半径，如式(3)所示。
rswept＝rsymmetry(cosλ)2.2(3)
其中，rswept为钝化前缘曲线上各横向站位处的前缘钝化半径(由公式(3)得到)，rsymmetry为抬升后的乘波体前缘点p2处的前缘钝化半径。
s4设置抬升后的乘波体前缘点p2以及乘波体的上表面沿法向抬升高度h；
抬升后的乘波体前缘点p2和乘波体的上表面沿法向抬升高度h的设置应满足以下两个要求：
1)参数设置能使钝化前缘曲线能够成功生成；2)抬升后的乘波体前缘点p2的曲率大于乘波体钝化前缘曲线上各横向站位处的曲率。
进一步地，本发明优选将乘波体的上表面沿法向抬升高度h设置为当地前缘钝化半径的1～2倍即可。
将抬升后的乘波体前缘点p2设置在上表面曲线以及下表面曲线的延长线所构成夹角的角平分线上，如图4所示，其距p0点的水平距离l可由式(4)得出，
l＝rswept/cosλ(4)
s5求出所有钝化前缘曲线，然后放样得到钝化前缘曲面；
根据点p0和p1以及s2获取的抬升后的乘波体前缘点p2的插值条件求出插值曲线即钝化前缘曲线，将求解钝化前缘曲线的问题看作插值于两个端点的ghi问题，以p0、p2为端点利用一段3次bézier曲线构造点p0和p2间的钝化前缘曲线，以p2、p1为端点利用一段3次bézier曲线构造点p2和p1间的钝化前缘曲线，完整的钝化前缘曲线即由点p0和p2间的钝化前缘曲线和点p2和p1间的钝化前缘曲线连接而成，求出所有的钝化前缘曲线，选取其中同时满足s4中要求的钝化前缘曲线进行放样，从而得到所需的钝化前缘曲面。
根据端点条件构造插值曲线的问题，可以看作插值于两个端点的ghi(geometrichermiteinterpolation，几何厄米特插值)问题。该问题由deboor首次提出，参见boorcd,
3次bézier曲线共需要4个控制点来进行曲线的求解，因此首先根据已知的两个控制点(即点p0和p2或者点p2和p1)求得插值曲线中的另外两个控制点，然后根据已知的两个控制点以及求解得到的另外两个控制点，即4个控制点来进行插值曲线的求解,如图5所示。
下面介绍如何以p0、p2为端点利用一段3次bézier曲线构造点p0和p2间的钝化前缘曲线，方法如下：
已知的两个控制点p0和p2(即两个端点)与待求的两个控制点p3和p4，这4个控制点间的关系可由下式表示
其中：δ0和δ1为设定的中间参数，为待求值；v0和v2分别为端点p0和p2处的切矢。
两个端点p0和p2处的曲率为
由s2获取的抬升后的乘波体前缘点p2的插值条件可以得到关于δ0和δ1的方程组
其中d是两个端点p0和p2之间的位移矢量。
通过求解方程组(7)可以得到未知参数δ0和δ1的值，从而计算得到其余两个控制点p3和p4的坐标，进而得到一条曲率连续的插值曲线，即利用3次bézier曲线构造的点p0和p2间的初步钝化前缘曲线。
方程组(7)的解可能有0至3个解。由于乘波体钝化前缘曲线为三维空间曲线，钝化前缘曲线上各个横向站位处的端点插值条件不一，可能导致方程无解或者多解。如要将该方法应用于三维外形上，则对于方程组解的讨论将是主要障碍之一。此处对插值条件进行简化，即不要求钝化前缘曲线与上表面曲线、下表面曲线间的曲率连续，将钝化前缘曲线与上表面曲线、下表面曲线连接处的曲率值设置为0。
方程组(7)可简化为
原方程组(7)变为一个一元一次方程与一个二元二次方程的组合，不难看出，方程组(8)必然有且只有一个解。通过对抬升后的乘波体前缘点p2的插值条件的调整，以降低钝化前缘曲线与上表面曲线、下表面曲线连接的光滑性为代价，简化了bézier曲线的求解过程，有利于保持不同横向站位处钝化前缘曲线之间的连续性。
以p2、p1为端点利用一段3次bézier曲线构造点p2和p1间的钝化前缘曲线的方法与以p0、p2为端点利用一段3次bézier曲线构造点p0和p2间的钝化前缘曲线的方法完全相同，在此不再赘述。
本实施例中，利用2段3次bézier曲线构造钝化前缘曲线，具体如图6所示。图6是利用2段3次bézier曲线构造钝化前缘曲线示意图，p0、p2和p2、p1分别是2段3次bézier曲线的端点，v、k则分别是各点处的切矢和曲率。钝化前缘曲线由两段曲线连接而成，即点p0和p2间的钝化前缘曲线和点p2和p1间的钝化前缘曲线。同样的将求解钝化前缘曲线的问题看作插值于两个端点的ghi问题，以p0、p2为端点利用一段3次bézier曲线构造点p0和p2间的钝化前缘曲线，以p2、p1为端点利用一段3次bézier曲线构造点p2和p1间的钝化前缘曲线，求出所有的钝化前缘曲线。实际生成前缘钝化曲面时，较难保证所有生成的钝化前缘曲线均能满足前述两点要求，因此选取其中同时满足s4中所有要求的钝化前缘曲线进行放样，从而得到所需的钝化前缘曲面。本实施例最终得到钝化后的吻切锥乘波体外形如图7所示。
对钝化后的乘波体体进行气动性能分析，由于乘波体为面对称外形，在只研究俯仰特性时可取一半流场进行数值模型，计算域包括了乘波体的底部。计算网格为非结构网格，为提高对边界层的预测精度，在壁面附近划分了棱柱网格以提高边界层的预测精度。流动控制方程为可压流的三维雷诺平均navier-stokes方程(rans)，在计算中采用二阶迎风格式进行空间离散，时间项为隐格式。采用rng-两方程湍流模型和增强壁面函数模拟湍流，壁面边界条件设置为绝热无滑移壁面。在马赫数ma＝8.0和高度h＝31km时，气动性能随攻角变化规律如图8～11所示。
不难看出，升力系数和俯仰力矩系数与攻角基本上呈线性关系，且值随着攻角的增大而增大。粘性对升力系数和俯仰力矩系数的影响不大，在飞行器设计初期可通过无粘计算快速获取升力特性和俯仰特性。由于高超声速流动具有强粘性的特点，乘波体阻力系数的有粘和无粘计算结果存在较大的差异。乘波体的最大升阻比出现在8°攻角，与一般研究认为的最大升阻比出现在小攻角的结论不一致。可能的原因是本文采用了脊形的设计，上表面与一般乘波体所用的自由流面设计有较大区别，进而导致了升阻比和俯仰特性存在一定的差异。
以上包含了本发明优选实施例的说明，这是为了详细说明本发明的技术特征，并不是想要将发明内容限制在实施例所描述的具体形式中，依据本发明内容主旨进行的其他修改和变型也受本专利保护。本发明内容的主旨是由权利要求书所界定，而非由实施例的具体描述所界定。
- 还没有人留言评论。精彩留言会获得点赞！