基于物理量梯度修正的紧支基函数多场耦合数据传递方法与流程

本发明涉及多场耦合中的数据传递方法领域，具体涉及一种基于物理量梯度修正的紧支基函数多场耦合数据传递方法。

背景技术：

当前飞行器气动力/热/结构多物理场耦合问题的求解主要是松耦合方法，其基本思路为按照设定的顺序在各自的求解器中分别求解流体控制方程和结构控制方程，通过流固耦合交界面把两者的计算结果互相传递，同时能够充分的利用现有的cfd、csd求解方法和程序模块，只需要增加少量的耦合界面数据交互模块即可完成整个分析过程。径向基函数(radialbasisfunction,rbf)插值方法形式简单，不依赖于求解器离散格式及网格拓扑结构，易于编程且精度较高，因而在多场耦合数据插值领域得到广泛应用。比较常用的径向基函数包括薄板样条(tps)等样条函数、多重二次曲面双协调(mq)曲面拟合及紧支c2基函数。tps、mq等全局域基函数其插值精度较高，但插值时需要所有的网格点，当实际工程中网格量较大时，所需的内存和时间都是非常巨大的。而紧支c2基函数在插值过程中只选择较少的网格点，能够有效的减少内存和时间的开销，因而已成为当前多场耦合计算领域重要的一种插值方法，但是紧支c2基函数存在的瓶颈是在有效节省了插值计算量的情况下如何进一步提升其精度。

技术实现要素：

针对现有技术中的上述不足，本发明提供的一种基于物理量梯度修正的紧支基函数多场耦合数据传递方法解决了紧支c2基函数在有效节省了插值计算量的情况下如何进一步提升其精度的问题。

为了达到上述发明目的，本发明采用的技术方案为：一种基于物理量梯度修正的紧支基函数多场耦合数据传递方法，包括以下步骤：

s1、对于所有已知物理量网格点进行物理量梯度的计算；

s2、根据梯度值对紧支基函数进行各向异性修正，并根据修正结果对紧支半径计算公式进行更新；

s3、根据更新后的紧支半径计算公式得到修正后的径向基的取值，并按照紧支基函数标准插值流程进行数据传递。

进一步地：所述步骤s1中的网格点包括内部点、网格块角点和边界点。

进一步地：所述步骤s1的具体步骤为：

s11、遍历所有内部点对应的网格块角点和边界点；

s12、根据内部点的相邻网格块角点的物理量和该内部点的物理量计算内部点的梯度；

s13、将边界点的相邻内部点梯度赋值为该边界点的梯度；

s14、将网格块角点最近的边界点梯度赋值为该网格块角点的梯度。

进一步地：所述步骤s12中内部点的梯度的计算公式为：

上式中，ui为网格点i上物理量的值，uj为网格点j上物理量的值，δ(·)ij＝(·)j-(·)i，i，j均为网格点的序号，(x,y,z)为网格点的坐标，n为网格点个数，θj为权系数。

进一步地：所述步骤s2的具体步骤为：

s21、计算所有网格点在x,y,z方向上的梯度总和；

s22、根据网格点在x,y,z方向上的梯度总和计算网格点在x,y,z方向上的平均梯度；

s23、对网格点在x,y,z方向上的平均梯度实现几何外形上的尺度归一；

s24、根据尺度归一后的平均梯度对所有网格点的坐标进行缩放；

s25、根据径向基函数公式和缩放后的坐标得到更新后的紧支半径计算公式。

进一步地：所述步骤s22中网格点在x,y,z方向上的平均梯度的计算公式为：

上式中，当i＝1时，ave_gradx1＝ave_gradx，all_gradx1＝all_gradx；当i＝2时，ave_gradx2＝ave_grady，all_gradx2＝all_grady；当i＝3时，ave_gradx3＝ave_gradz，all_gradx3＝all_gradz；ave_gradx、ave_grady、ave_gradz分别为网格点在x、y、z方向上的平均梯度，all_gradx、all_grady、all_gradz分别为网格点在x、y、z方向上的梯度总和，i，j均为网格点的序号，n为网格点个数。

进一步地：所述步骤s23中尺度归一的计算公式为：

上式中，ave_gradxi'为尺度归一后的平均梯度，ave_gradx'、ave_grady'、ave_gradz'为尺度归一后x、y、z方向上的平均梯度。

进一步地：所述步骤s24中缩放的计算公式为：

x'＝x·ave_gradx'

y'＝y·ave_grady'

z'＝z·ave_gradz'

上式中，x'、y'、z'为缩放后的网格点坐标。

进一步地：所述步骤s25中更新后的紧支半径计算公式为：

上式中，r为修正后的径向基。

本发明的有益效果为：本发明通过在原始紧支基函数计算方法中，根据物理量梯度大小，对径向基计算中的xyz三方向进行系数缩比调控，以求在相同的紧支半径范围内，选取更多物理量变化剧烈的点进行插值，从而使用来插值的点能够更具有代表性和聚集性，具有实际物理量分布的各向异性特征，从而使插值的结果更能表征物理实际分布特征，从而提高精度。因此本发明可以为飞行器气动力/热/结构多场耦合计算提供一种可行的具备更高精度的数据传递方法。

附图说明

图1为本发明流程图；

图2为本发明中结构化网格点相邻点的示意图；

图3为本发明中y型剖分网格点相邻点的示意图；

图4为本发明实施例中的误差改进效果示意图。

具体实施方式

下面对本发明的具体实施方式进行描述，以便于本技术领域的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。

如图1所示，一种基于物理量梯度修正的紧支基函数多场耦合数据传递方法，其特征在于，包括以下步骤：

s1、对于所有已知物理量网格点进行物理量梯度的计算；网格点包括内部点、网格块角点和边界点。分别如附图2所示的q1、q2、q3。对于内部点一般有四个单元的顶点与其相连，如附图2中与q1点相连的四个节点p1、p2、p3、p4，对于存在y型剖分的网格某一内部节点只有三个相连的节点，如附图3中与q4点相连的三个节点p8、p9、p10；对于网格块角点q2找到其所在单元e7，然后找出q2在该单元内相邻点p5、p6；对于边界上的点q3，找到与之相邻的内部点p7。

表面网格节点梯度求解的的具体步骤为：

s11、遍历所有内部点对应的网格块角点和边界点；

s12、根据内部点的相邻网格块角点的物理量和该内部点的物理量计算内部点的梯度；内部点的梯度的计算公式为：

上式中，ui为网格点i上物理量的值，uj为网格点j上物理量的值，δ(·)ij＝(·)j-(·)i，i、j均为网格点的序号，(x,y,z)为网格点的坐标，n为网格点个数，θj为权系数，如面积、体积或距离等，一般也可以取1，本申请中的θj值均为1。

s13、将边界点的相邻内部点梯度赋值为该边界点的梯度；

s14、将网格块角点最近的边界点梯度赋值为该网格块角点的梯度。

v.对于存在其他形式的内部点相连形式(如：某一内部点连接5个相邻点)，梯度求解方法是找到最近的相邻节点，将该相邻节点梯度赋值给该点。

s2、根据梯度值对紧支基函数进行各向异性修正，并根据修正结果对紧支半径计算公式进行更新；

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

s21、计算所有网格点在x,y,z方向上的梯度总和；

s22、根据网格点在x,y,z方向上的梯度总和计算网格点在x,y,z方向上的平均梯度；

网格点在x,y,z方向上的平均梯度的计算公式为：

上式中，当i＝1时，ave_gradx1＝ave_gradx，all_gradx1＝all_gradx；当i＝2时，ave_gradx2＝ave_grady，all_gradx2＝all_grady；当i＝3时，ave_gradx3＝ave_gradz，all_gradx3＝all_gradz；ave_gradx、ave_grady、ave_gradz分别为网格点在x、y、z方向上的平均梯度，all_gradx、all_grady、all_gradz分别为网格点在x、y、z方向上的梯度总和，i，j均为网格点的序号，n为网格点个数。

s23、对网格点在x,y,z方向上的平均梯度实现几何外形上的尺度归一；

尺度归一的计算公式为：

上式中，ave_gradxi'为尺度归一后的平均梯度，ave_gradx'、ave_grady'、ave_gradz'为尺度归一后x、y、z方向上的平均梯度。

s24、根据尺度归一后的平均梯度对所有网格点的坐标进行缩放；

缩放的计算公式为：

x'＝x·ave_gradx'

y'＝y·ave_grady'

z'＝z·ave_gradz'

上式中，x'、y'、z'为缩放后的网格点坐标。

s25、根据径向基函数公式和缩放后的坐标得到更新后的紧支半径计算公式。

更新后的紧支半径计算公式为：

上式中，r为修正后的径向基。

s3、根据更新后的紧支半径计算公式得到修正后的径向基的取值，并按照紧支基函数标准插值流程进行数据传递。

在本发明的一个实施例中，针对某一翼面开展飞行器插值，该翼面流场表面网格节点数4073个，单元数为4016个；结构表面网格节点数5317个，单元数为10520个，已知流场结构网格点上每一点的热流值。首先根据步骤s1计算流场表面4016个单元对应的4073个点的热流梯度。

根据步骤s1得到的每一点上的热流梯度，得到xyz各个方向的梯度总和all_gradx：1374.2、all_grady:1920.2、all_gradz：448.99，求各个方向的平均梯度ave_gradx:1.06、ave_grady：2.64、ave_gradz:16.38。

已知该外形的min_x:8.91e^-7、min_y:0min_z:-6.88e^-2、max_x:2.716max_y:1.27max_z:6.88e^-2，根据步骤s2中尺度归一公式计算得到每一点在不同方向的缩放因子。根据缩放因子，计算缩放后的径向基的值。

用缩放后径向基的值，代替标准紧支径向基函数中的径向基，完成数据的插值。附图4给出了改进后的实际效果，插值的结果更能表征物理实际分布特征，从而提高精度。