CN111079326B - 二维各向异性网格单元度量张量场光滑化方法 - Google Patents

Abstract

本发明公开了二维各向异性网格单元度量张量场光滑化方法。目前尚无较好的各向异性网格度量张量的光滑化方法。本发明在飞机翼形几何模型区域的三角形背景网格上设定飞机翼形的初始各向异性网格度量张量场，将每个网格点上的各向异性网格度量张量进行取对数运算，然后建立限制各向异性网格度量张量变化的梯度约束，根据梯度约束对飞机翼形的初始各向异性网格度量张量场进行光滑化处理，从而建立各向异性尺寸度量张量场的光滑过渡数学模型，得到光滑各向异性网格度量张量场。本发明将变化剧烈的飞机翼形度量张量场处理成光滑的度量张量场，从而提高各向异性网格质量及数值模拟精度，使得最终生成的飞机翼形网格单元质量高且网格单元数最少。

Description

二维各向异性网格单元度量张量场光滑化方法

技术领域

本发明涉及数值模拟领域的前处理网格生成过程中的单元尺寸设置方法，具体涉及面向各向异性网格生成的二维各向异性网格单元度量张量场光滑化方法。

背景技术

网格生成是有限元法、有限体积法和有限差分法等数值模拟技术中的前处理过程，该过程将连续的几何区域划分成有限个基本几何形体的组合，这些基本几何形体被称为网格单元，常用的网格单元类型有三角形单元、四边形单元、四面体单元和六面体单元等。网格单元的质量和数量直接影响数值计算的精度和效率，好的网格需要用尽量少的自由度获得尽量高的数值计算精度。另一方面，很多数值模拟问题往往在某些局部区域存在很强的各向异性特征，如流体里的激波和边界层，及燃烧时的火焰面等。在这些局部区域，物理量沿某一个方向的变化远大于沿其它方向的变化。对这类具有强各向异性特征的数值模拟问题，最佳的网格配置方式应该是在这些存在各向异性特征的局部区域生成长宽比很大的各向异性网格单元，而在其它区域生成各向同性网格单元。这样的网格配置可以在计算精度和计算效率上获得最佳平衡。

区别于各向同性网格单元，各向异性网格单元外观主要由尺寸、形状和方向三个因素决定。在二维情形下，这些决定因素可由2×2的对称正定矩阵来表示，该矩阵也称为各向异性网格度量张量。要获得数值仿真几何模型的各向异性网格，通常需要先获得覆几何模型的度量张量场。这一张量场可以根据模型的几何特征计算得到，如沿高曲率方向的网格尺寸要小于其它方向的网格尺寸。多年来，求解自适应技术被成功应用于各向异性网格生成中，这类方法中的各向异性网格度量张量场根据求解误差进行重构得到，进而调整初始网格单元的尺寸、形状和方向并得到符合物理解特征的各向异性网格。

然而，求解自适应各向异性网格生成方法仍面临着一个较大的技术挑战，即通过数值解重构得到的各向异性网格度量张量场通常不太光滑，在某些区域从一个位置到另一位置内度量张量变化过大。不光滑的度量张量场会导致生成低质量的各向异性网格单元，从而降低数值解的精度。为此，在将根据数值解重构得到的度量张量场应用于各向异性网格生成前需要先对其进行光滑处理。相对于各向同性网格尺寸场光滑化处理过程中只需要考虑网格尺寸的变化，各向异性网格度量张量场的光滑化更加复杂，既需要考虑张量的主方向变化还需要考虑主方向上的网格尺寸变化。当前，尚无较好的各向异性网格度量张量的光滑化方法提出，而在工业方面有着相关诉求，例如NACA0012飞机翼形几何模型生成的初始各向异性网格度量张量场在用于网格生成时质量不足以使用，需要做优化处理。

发明内容

本发明的目的是针对现有构建网格单元张量场方法存在的缺点，提供一种面向自适应曲面网格生成的二维各向异性网格单元度量张量场光滑化方法，将变化剧烈的飞机翼形度量张量场处理成光滑的度量张量场，从而提高各向异性网格质量及数值模拟精度；不仅可以减少运算内存和时间消耗，同时还可构建飞机翼形单元尺寸值过渡合理、生成网格单元数目最少的尺寸场，使得最终生成的飞机翼形网格单元质量高且网格单元数最少。

本发明采用的技术方案是：

本发明二维各向异性网格单元度量张量场光滑化方法，具体如下：

步骤1、生成覆盖飞机翼形几何模型区域的三角形背景网格，三角形背景网格包含网格点和网格单元，设定飞机翼形的初始各向异性网格度量张量场在网格点B_q上的各向异性网格度量张量为M_q，q＝1,2,...,n，n为网格点总数，各向异性网格度量张量M_q为二阶张量。

步骤2、将飞机翼形的初始各向异性网格度量张量场中每个网格点上的各向异性网格度量张量进行取对数运算，从而将飞机翼形的初始各向异性网格度量张量场转换到对数空间。取对数运算具体如下：

步骤3、建立限制各向异性网格度量张量变化的梯度约束，具体如下：设三角形背景网格的网格单元T覆盖的参数平面区域为Ω_T，

i＝0,1,2为网格单元T的三个网格点上各向异性网格度量张量在对数空间中的值，通过线性插值计算出参数平面区域Ω_T内任意点(x,y)的各向异性网格度量张量：

式中，exp代表指数运算，t_i,i＝0,1,2为点(x,y)在参数平面区域Ω_T内的面积坐标，即：

式中，a_i＝x_jy_k-x_ky_j，b_i＝y_j-y_k，c_i＝x_k-x_j，i＝0,1,2,j＝(i+1)mod3，k＝(i+2)mod3，mod代表求余运算，(x_i,y_i),i＝0,1,2分别为参数平面区域Ω_T的三个顶点坐标，A为参数平面区域Ω_T的面积。

点(x,y)的各向异性网格度量张量进行取对数运算后在参数平面区域Ω_T内变化的梯度如下：

式中，

为梯度符号，

为偏导符号，T为转置符号。

使三角形背景网格的网格单元T内的尺寸值满足如下梯度约束：

式中，

为

的模，参数β为大于0的一个设定值。

步骤4、将三角形背景网格记为M，设三角形背景网格的网格单元数为m，三角形背景网格M表达如下：

M＝{E＝{e_p|p＝1,2,…,m},B＝{B_q|q＝1,2,…,n}},

其中，E和B分别为网格单元集合和网格点集合。

计算每一网格单元e_p区域内点(x,y)的各向异性网格度量张量进行取对数运算后的变化梯度

若对于所有网格单元均满足

则三角形背景网格M上获得光滑过渡的各向异性张量场，否则按照尺寸张量修改模型修改不满足

的网格单元网格点上的尺寸张量，直到所有网格单元均满足

得到光滑各向异性网格度量张量场。

进一步，建立尺寸张量修改模型如下：

式中，

为网格点B_q上修改后的尺寸张量在对数空间中的值，

为网格点B_q上修改前的尺寸张量在对数空间中的值。

本发明具有的有益效果是：

本发明采用非结构化曲面背景网格存储飞机翼形网格单元尺寸值，从而构建单元尺寸合理过渡的张量场；非结构化网格具有拓扑结构灵活、局部修改简单等特点，提高了张量场构建的存储性能和时间性能；另外由于只需要生成曲面背景网格，所需单元数少；本发明将单元尺寸梯度约束构建为一个非线性的凸优化问题，该优化问题存在全局最优解，该解可以在整个计算域内满足单元尺寸合理过渡要求，同时使得背景网格点上的尺寸值改变最小，从而使得最终生成的飞机翼形网格单元质量高且网格单元数最少，为飞机翼形下游数值模拟过程提供最佳输入。

附图说明

图1为NACA 0012飞机翼形的几何模型。

图2为NACA 0012飞机翼形的初始各向异性网格度量张量场及局部放大图。

图3为本发明处理后的NACA 0012飞机翼形光滑各向异性网格度量张量场及局部放大图。

图4为NACA 0012飞机翼形基于初始各向异性尺寸张量场生成的网格及局部放大图。

图5为NACA 0012飞机翼形基于光滑各向异性网格度量张量场生成的网格及局部放大图。

图6为基于光滑各向异性网格度量张量场及其网格的仿真结果图。

具体实施方式

下面结合附图对本发明作进一步说明。

本发明二维各向异性网格单元度量张量场光滑化方法，具体如下：

1、生成覆盖飞机翼形几何模型(如图1所示)区域的三角形背景网格，如图4所示，三角形背景网格包含网格点和网格单元，设定飞机翼形的初始各向异性网格度量张量场在网格点B_q上的各向异性网格度量张量为M_q，q＝1,2,...,n，n为网格点总数，各向异性网格度量张量M_q为二阶张量，其主方向及主方向上的半径定义了各向异性网格的主方向及沿主方向上的尺寸值。如图2所示为飞机翼形定义在正方形区域内的初始各向异性网格度量张量场。

2、将飞机翼形的初始各向异性网格度量张量场中每个网格点上的各向异性网格度量张量进行取对数运算，从而将飞机翼形的初始各向异性网格度量张量场转换到对数空间。取对数运算具体如下：

3、建立限制各向异性网格度量张量变化的梯度约束，具体如下：设三角形背景网格的网格单元T覆盖的参数平面区域为Ω_T，

i＝0,1,2为网格单元T的三个网格点上各向异性网格度量张量在对数空间中的值(进行取对数运算后的值)，通过线性插值计算出参数平面区域Ω_T内任意点(x,y)的各向异性网格度量张量：

式中，exp代表指数运算，t_i,i＝0,1,2为点(x,y)在参数平面区域Ω_T内的面积坐标，即：

式中，a_i＝x_jy_k-x_ky_j，b_i＝y_j-y_k，c_i＝x_k-x_j，i＝0,1,2,j＝(i+1)mod3，k＝(i+2)mod3，mod代表求余运算，(x_i,y_i),i＝0,1,2分别为参数平面区域Ω_T的三个顶点坐标，A为参数平面区域Ω_T的面积。

点(x,y)的各向异性网格度量张量进行取对数运算后在参数平面区域Ω_T内变化的梯度(用于表示各向异性网格度量张量变化的剧烈程度)如下：

式中，

为梯度符号，

为偏导符号，T为转置符号。

为使点(x,y)的各向异性网格度量张量在参数平面区域Ω_T内缓慢过渡，保证三角形背景网格的网格单元T内的尺寸值满足如下梯度约束：

式中，

为

的模，参数β为大于0的一个设定值。

4、根据步骤3中建立的梯度约束对飞机翼形的初始各向异性网格度量张量场进行光滑化处理，从而建立各向异性尺寸度量张量场的光滑过渡数学模型，得到光滑各向异性网格度量张量场，具体如下：将三角形背景网格记为M，设三角形背景网格的网格单元数为m，三角形背景网格M表达如下：

M＝{E＝{e_p|p＝1,2,…,m},B＝{B_q|q＝1,2,…,n}},

其中，E和B分别为网格单元集合和网格点集合。

计算每一网格单元e_p区域内点(x,y)的各向异性网格度量张量进行取对数运算后的变化梯度

若对于所有网格单元均满足

则三角形背景网格M上获得光滑过渡的各向异性张量场，否则按照尺寸张量修改模型修改不满足

的网格单元网格点上的尺寸张量，直到所有网格单元均满足

得到光滑各向异性网格度量张量场。为了使网格点上的尺寸值改变最小，建立尺寸张量修改模型如下：

式中，

为网格点B_q上修改后的尺寸张量在对数空间中的值，

为网格点B_q上修改前的尺寸张量在对数空间中的值，min代表“在后面的式子取最小值的前提下”，s.t.代表“使得”。

本实施例中采用NACA 0012飞机翼形，图1为NACA 0012飞机翼形的几何模型，图2中的NACA 0012飞机翼形初始各向异性网格度量张量场按照本发明步骤进行光滑化后的各向异性尺寸张量场如图3所示。图4和5分别为基于初始各向异性网格度量张量场和光滑各向异性网格度量张量场生成的网格，可以看到基于光滑各向异性网格度量张量场生成的网格尺寸过渡很光滑，该网格更利于仿真分析，图6为仿真分析结果，图6中e为以10为底的指数函数。

Claims (2)

1.二维各向异性网格单元度量张量场光滑化方法，其特征在于：该方法具体如下：

步骤1、生成覆盖飞机翼形几何模型区域的三角形背景网格，三角形背景网格包含网格点和网格单元，设定飞机翼形的初始各向异性网格度量张量场在网格点B_q上的各向异性网格度量张量为M_q，q＝1,2,...,n，n为网格点总数，各向异性网格度量张量M_q为二阶张量；

步骤2、将飞机翼形的初始各向异性网格度量张量场中每个网格点上的各向异性网格度量张量进行取对数运算，从而将飞机翼形的初始各向异性网格度量张量场转换到对数空间；取对数运算具体如下：

步骤3、建立限制各向异性网格度量张量变化的梯度约束，具体如下：设三角形背景网格的网格单元T覆盖的参数平面区域为Ω_T，

为网格单元T的三个网格点上各向异性网格度量张量在对数空间中的值，通过线性插值计算出参数平面区域Ω_T内任意点(x,y)的各向异性网格度量张量：

式中，exp代表指数运算，t_i,i＝0,1,2为点(x,y)在参数平面区域Ω_T内的面积坐标，即：

式中，a_i＝x_jy_k-x_ky_j，b_i＝y_j-y_k，c_i＝x_k-x_j，i＝0,1,2,j＝(i+1)mod3，k＝(i+2)mod3，mod代表求余运算，(x_i,y_i),i＝0,1,2分别为参数平面区域Ω_T的三个顶点坐标，A为参数平面区域Ω_T的面积；

点(x,y)的各向异性网格度量张量进行取对数运算后在参数平面区域Ω_T内变化的梯度如下：

式中，

为梯度符号，

为偏导符号，T为转置符号；

使三角形背景网格的网格单元T内的尺寸值满足如下梯度约束：

式中，

为

的模，参数β为大于0的一个设定值；

步骤4、将三角形背景网格记为M，设三角形背景网格的网格单元数为m，三角形背景网格M表达如下：

M＝{E＝{e_p|p＝1,2,…,m},B＝{B_q|q＝1,2,…,n}},

其中，E和B分别为网格单元集合和网格点集合；

计算每一网格单元e_p区域内点(x,y)的各向异性网格度量张量进行取对数运算后的变化梯度

若对于所有网格单元均满足

则三角形背景网格M上获得光滑过渡的各向异性张量场，否则按照尺寸张量修改模型修改不满足

的网格单元网格点上的尺寸张量，直到所有网格单元均满足

得到光滑各向异性网格度量张量场。

2.根据权利要求1所述的二维各向异性网格单元度量张量场光滑化方法，其特征在于：建立尺寸张量修改模型如下：

式中，

为网格点B_q上修改后的尺寸张量在对数空间中的值，

为网格点B_q上修改前的尺寸张量在对数空间中的值。

Images