基于非结构网格规则化重构技术的高精度热流计算方法

摘要

本发明涉及飞行器流场气动热求解技术领域，公开了一种基于非结构网格规则化重构技术的高精度热流计算方法。该方法从利用初始流场进行特征面提取思想出发，结合激波探测技术和三维激波曲面拟合技术，多次迭代逐步逼近找到符合真实物理解的空间激波面，对物面驻点热流网格进行规则化处理以及对空间激波面附近网格进行二次贴体重构，从而消除因网格随机性排列和扰动引起的数值误差。本发明提出的方法在保证整体非结构网格快速生成的前提下，通过半自动化网格修复技术，反复迭代重构形成符合物理特征的规则化网格，排除因网格扰动导致的数值计算误差，为高超声速飞行器气动热模拟提供了一种新的快速而精细的模拟方法。

说明书

技术领域

本发明涉及飞行器流场气动热求解技术领域，具体涉及一种高精度求解高超声速飞行器驻点附近热流的求解方法。

背景技术

精确预测气动热环境是高超声速飞行器热防护设计的核心技术。当前普遍采用的方法是基于计算流体力学（CFD）技术进行N-S方程离散求解得到空间流场信息，通过分析流场信息进一步得到飞行器的气动热相关的特性参数。当前的CFD计算方法主要包括基于结构化的网格和非结构化网格两大类，由于热流计算要求对飞行器表面温度梯度进行精准预测，因而基于结构化网格的方法由于网格的规则化排布特性较好通常可以得到较高精度的模拟结果，但结构化网格存在的一个显著问题是网格生成耗时过长，且在复杂外形局部区域拓扑设计困难，这给复杂外形的热流计算带来了很大的挑战；而非结构网格因具有自动化、灵活性和高效率的优点，可以较好地解决飞行器复杂外形的网格生成效率问题，因而基于非结构网格的CFD技术在气动热的数值模拟计算领域得到了越来越多的关注。但非结构网格的应用当前还存在如下问题：（1）由于非结构网格排列的不规则性，会导致热流预测精度不够;(2)当空间网格没有严格沿激波排列时，激波后流场会出现较大的扰动误差，这种扰动误差会传播到驻点附近，在极端情况下，误差可能会支配流场，造成激波向外突出，形成所谓的“红玉”现象，导致出现非物理的流场解等。而基于非结构网格自适应的方法虽然能通过自动化的手段加密激波附近的网格，使激波模拟更精细，但依然无法实现网格严格沿激波排列，无法消除网格带来的扰动误差。同样激波装配方法虽然可以实现网格严格沿激波排列，但是计算复杂，稳定性差。因此为了实现基于非结构网格的热流模拟，亟需解决非结构网格的规则化排布尤其是激波间断区域的网格正交性问题。

发明内容

本发明所要解决的技术问题是：针对上述存在的问题，提供了一种基于非结构网格规则化重构修正技术的高精度热流计算新方法。该方法从利用初始流场进行特征面提取思想出发，结合激波探测技术和三维曲面拟合技术，多次迭代逐步逼近找到符合真实物理解的空间激波面，对物面驻点热流网格进行规则化处理以及对空间激波面附近网格进行二次贴体重构，从而消除因网格随机性排列和扰动引起的数值误差。

本发明采用的技术方案如下：

一种基于非结构网格规则化重构技术的高精度热流计算方法，包括：

步骤1：针对飞行器外形，通过网格软件生成初始的非结构混合网格，基于雷诺平均N-S方程的流场求解器进行计算，得到收敛的流场解，获取流场中的压力分布信息。

步骤2：根据压力分布信息，在非结构混合网格中判断出重构驻点位置；采用交叉剖分技术规则化重构驻点附近的网格，具体为：首先根据上一轮计算得到的流场中的压力分布信息，获取压力极值区域，即为物面驻点区域；然后在物面驻点区域附近生成各向同性的规则型四边形网格，在四边形网格的基础上，任选一个单元，沿对角线剖分成两个三角形，然后沿已剖分单元的剖分顶点剖分相邻单元，循环进行，直到将所有单元剖分完毕，得到规则排列的三角形网格。

步骤3：使用“激波探测”方法进行激波面提取。

具体步骤包括：根据步骤1初始模拟流场的空间压力分布，采用基于压力间断的方法进行激波探测，遍历所有网格单元进行激波面网格标记；

针对网格单元

将所有被标记为激波面网格单元的体心坐标输出，形成一个

步骤4：根据步骤3基于压力间断判断得到的激波面网格单元的附近单元可能存在少量被错误标记，通过采用孤立森林模型算法将这些异常点进行过滤删除。

其中，孤立森林模型算法包括：

训练阶段，对数据矩阵中的

评估阶段，对数据矩阵中的

其中，

删除数据矩阵中的

步骤5：采用局部加权回归散点平滑法将经过异常点筛选后剩下的

步骤6：将三维激波曲面导入步骤1的网格软件，按照在三维激波曲面上生成各向同性三角形网格以及在三维激波曲面上沿其法向方向外推生成三棱柱网格这两个要求重新生成网格，得到新网格。

步骤7：利用步骤1的流场计算结果，插值到新网格上作为初始流场，将初始流场导入到流场解算器进行续算。观察初始流场计算结果是否满足热流分布光滑要求，若满足要求，则停止计算；否则重复步骤2到步骤6，直到满足要求为止。

进一步的，所述最大压力突变判据的计算公式如下：

其中，

进一步的，所述训练阶段，在每一步分割过程中，

（1）孤立树达到了限制的深度，取为

（2）节点上只有一个样本点；

（3）节点上样本点的特征

最终，构建得到一个孤立树的数据结构。

与现有技术相比，采用上述技术方案的有益效果为：该方法基于激波探测、机器学习异常点筛选、三维激波曲面拟合等技术对激波面附近的网格和表面驻点附近网格进行规则化二次修复重构，从而能有效解决利用非结构网格的CFD计算进行精确热流模拟的问题。

附图说明

图1是本发明一实施例提供的整体计算方法流程图。

图2是驻点附近网格剖分示意图，其中，图2（a）、图2（b）、图2（c）皆表示在规则型四边形网格上生成三角形网格，图2（d）为规则排列的三角形网格。

图3是标记为激波面的网格单元位置示意图。

图4是孤立森林算法中的数据分割示意图。

图5是孤立森林算法中的孤立树结构示意图。

具体实施方式

下面结合附图对本发明做进一步描述。

为了解决非结构网格的物面和空间随机分布的缺点引起的数值计算误差问题，如图1所示，本实施例提出了一种基于激波探测和曲面拟合相结合的非结构网格规则化重构技术来实现高精度的热流计算。考虑到热流计算误差主要由驻点附近网格排列不规则及激波附近空间网格没有沿激波排列造成的，因此，为了消除网格排列引起的误差，根据初始流场计算结果，利用激波探测技术，获取物面驻点和空间激波的位置信息，对驻点附近表面网格进行规则化重构，对激波面进行拟合，进而进行空间网格规则化修正优化，从而较好地模拟飞行器驻点附近滞止区的热流，获取更可靠精细的热流计算结果。

具体的步骤如下：

步骤1：针对飞行器外形，通过网格软件生成初始的非结构混合网格，基于雷诺平均N-S方程的流场求解器进行计算，得到收敛的流场解，获取流场中的压力分布信息。其中，网格软件可选用Pointwise或Gridstar。

步骤2：根据飞行器物面的压力分布信息，判断重构驻点位置。在本实施例中，采用交叉剖分技术规则化重构驻点附近的物面网格，具体操作如图2所示。首先根据上一轮计算得到的流场中的压力分布信息，获取压力极值区域，即为物面驻点区域；然后在物面驻点区域附近生成各向同性的规则型四边形网格，在四边形网格的基础上，任选一个单元，沿对角线剖分成两个三角形，然后沿已剖分单元的剖分顶点剖分相邻单元，循环进行，直到将所有单元剖分完毕，得到规则排列的三角形网格。

步骤3：使用“激波探测”方法进行激波面提取，输出激波点。具体步骤如下：

根据步骤1初始模拟流场的空间压力分布，采用基于压力间断的方法进行激波探测，遍历所有网格单元进行激波面网格标记。

对于网格单元

其中，

将所有被标记为激波面网格单元的体心坐标输出，形成一个

步骤4：去除异常点。

根据步骤3基于压力间断判断得到的激波面附近单元可能存在少量被错误标记，需要将这些异常点进行过滤删除。具体采用“孤立森林（Isolation Forest）”模型算法：

i.训练阶段，对数据矩阵中的

ii.评估阶段，对数据矩阵中的

其中，

iii.删除数据矩阵中的

步骤5：采用局部加权回归散点平滑法将经过异常点筛选后剩下的

步骤6：将三维激波曲面导入步骤1的网格软件，按照在三维激波曲面上生成各向同性三角形网格以及在三维激波曲面上沿其法向方向外推生成三棱柱网格这两个要求重新生成网格，得到新网格。

步骤7：利用步骤1的流场计算结果，插值到新网格上作为初始流场，将初始流场导入到流场解算器进行续算。观察初始流场计算结果是否满足热流分布光滑要求，若满足要求，则停止计算；否则重复步骤2到步骤6，直到满足要求为止。

下面提供一个具体的实施实例。

针对某钝锥头部开展初始热流计算，该算例来流条件为来流马赫数

根据步骤1得到的流场压力值，进行物面驻点附近规则化网格生成，如图2所示。以驻点为中心到压力平滑过渡区域生成规则型四边形网格，然后按照步骤2的规则剖分成三角形网格。

根据步骤3，利用压力梯度，进行激波面网格标记，并将所有被标记为激波面网格单元的体心坐标值输出。按照步骤4，对这些数据进行异常点检查，利用孤立森林模型进行训练识别，评估器数量(n_estimator)为100，最大样本(max_samples)为256，找到远离激波面的异常点并删除。

按照步骤5，利用局部加权回归散点平滑法对经过异常点筛选后剩下的数据进行三维激波曲面拟合。得到拟合后的三维激波曲面，并将STL的模型导入到原网格软件，在激波面上生成表面网格并按照在三维激波曲面上生成各向同性三角形网格以及在三维激波曲面上沿其法向方向外推生成三棱柱网格这两个原则，生成新网格。

在新网格的基础上进行重新计算，收敛后判断结果是否达到要求，否则回到步骤1重新迭代。此算例最终迭代三次后满足要求停止。

本发明并不局限于前述的具体实施方式。本发明扩展到任何在本说明书中披露的新特征或任何新的组合，以及披露的任一新的方法或过程的步骤或任何新的组合。如果本领域技术人员，在不脱离本发明的精神所做的非实质性改变或改进，都应该属于本发明权利要求保护的范围。