EI、Scopus 收录
中文核心期刊

黏性边界层网格自动生成

甘洋科, 刘剑飞

甘洋科, 刘剑飞. 黏性边界层网格自动生成[J]. 力学学报, 2017, 49(5): 1029-1041. DOI: 10.6052/0459-1879-17-154
引用本文: 甘洋科, 刘剑飞. 黏性边界层网格自动生成[J]. 力学学报, 2017, 49(5): 1029-1041. DOI: 10.6052/0459-1879-17-154
Gan Yangke, Liu Jianfei. AUTOMATIC VISCOUS BOUNDARY LAYER MESH GENERATION[J]. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(5): 1029-1041. DOI: 10.6052/0459-1879-17-154
Citation: Gan Yangke, Liu Jianfei. AUTOMATIC VISCOUS BOUNDARY LAYER MESH GENERATION[J]. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(5): 1029-1041. DOI: 10.6052/0459-1879-17-154
甘洋科, 刘剑飞. 黏性边界层网格自动生成[J]. 力学学报, 2017, 49(5): 1029-1041. CSTR: 32045.14.0459-1879-17-154
引用本文: 甘洋科, 刘剑飞. 黏性边界层网格自动生成[J]. 力学学报, 2017, 49(5): 1029-1041. CSTR: 32045.14.0459-1879-17-154
Gan Yangke, Liu Jianfei. AUTOMATIC VISCOUS BOUNDARY LAYER MESH GENERATION[J]. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(5): 1029-1041. CSTR: 32045.14.0459-1879-17-154
Citation: Gan Yangke, Liu Jianfei. AUTOMATIC VISCOUS BOUNDARY LAYER MESH GENERATION[J]. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(5): 1029-1041. CSTR: 32045.14.0459-1879-17-154

黏性边界层网格自动生成

基金项目: 

国家自然科学基金 11172005

国家重点基础研究发展计划 2010CB832701

详细信息
    通讯作者:

    2) 刘剑飞, 副教授, 主要研究方向:网格生成, 计算几何, 计算机图形学.E-mail:liujianfei@pku.edu.cn

  • 中图分类号: O242.21

AUTOMATIC VISCOUS BOUNDARY LAYER MESH GENERATION

  • 摘要: 高雷诺数黏性流动在壁面附近存在边界层,在计算模拟中自动生成可靠且有效的计算网格仍然是计算流体力学存在的瓶颈.三棱柱/四面体混合网格技术在一定程度上缓解了这个困难.然而,对于复杂外形的情况,在边界层内自动高效生成高质量的三棱柱单元仍然十分困难.常用的层推进法在凹凸区域及角点处生成的边界层网格单元质量较差,边界层网格最外层尺寸不均匀.针对这些问题,发展了一种黏性边界层网格自动生成方法,通过顶点周围边的二面角识别物面网格特征确定多生长方向,预估并调整生长高度处理相交情况.同时提出一种三维前沿尺寸调节方式,提高了边界层网格单元的正交性,保证了边界层网格与远场网格尺寸的光滑过渡.通过ONERA M6翼型以及带发动机短舱的DLR-F6翼身组合体等外形的网格生成实例及绕流数值模拟,将计算值与标准实验值进行对比,结果表明:该方法能够自动高效地生成满足数值计算需求的混合网格.
    Abstract: High Reynolds number fluid flows have boundary layers at the wall. To automatically generate robust and valid boundary layer mesh for the simulations is still the bottleneck problem of computational fluid dynamics. Prisms/Tetrahedra hybrid mesh leads to significant savings in mesh size and solution costs. However, it's still difficult to generate prismatic elements of high quality within boundary layers of complex models. Previous advancing layers techniques sometimes lead to invalid meshes and poor quality elements at concave/convex ridges and sharp corners. To improve these situations, we present a strategy for automatically generating viscous boundary layer mesh. In this method, multiple growth directions at ridges and corners are well defined by the dihedral angles around the vertices and the growth heights are adjusted appropriately. Therefore, boundary layer mesh grows well at sharp corners, convex and concave ridges of the domain. We also decrease the number of global intersection checks by predefining the total growth heights before generating elements through one global check, which improves the efficiency of mesh generation. At the same time, we develop a 3D strategy of mesh size control to get a size uniform triangular mesh of the outer boundary of the boundary layer mesh, which is beneficial to generate far-field isotropic mesh of high quality. Finally, mesh examples and the viscous flow simulations including 2D and 3D are presented. In 3D, the hybrid mesh over the standard ONERA M6 and DLR-F6 configurations are generated with the present method. The numerical results agree very well with experiment data which indicates that the hybrid viscous meshes generated by the proposed method are effective and efficient.
  • 近年来,计算流体力学成为学术界研究的热点,在航空航天、汽车、热能动力等工程领域都有着重要的地位,广泛用于设计、分析和优化.在计算流体力学模拟中,需要借助计算机数值求解物理量的控制方程.当前主要采用的计算方法,包括有限差分法、有限体积法、有限元法等,都需要对计算区域进行网格划分,进而在网格单元节点上离散求解偏微分方程. 30多年来,尽管计算流体力学取得了很大的成就,网格生成依然是需要解决的关键问题之一.大批商业网格生成软件,包括Gambit,Gridgen以及ANSYS ICEM-CFD等,虽然具有功能全面、使用方便、技术服务好等优点,但在生成含边界层网格时的自动性及可靠性较差,对于复杂外形需要大量手工操作,有时甚至无法生成.高雷诺数黏性流动中不能自动生成可靠且有效的网格仍然是计算流体力学中的瓶颈问题.

    边界层是流体运动中贴近壁面的薄层,其厚度与物体特征尺寸的比值为一个小量,通常小于千分之一.如果采用各向同性网格,则生成的边界层网格单元数量巨大,大大提高了计算的代价,对于大规模的复杂问题,其计算代价是目前工程上不能接受的.而如果计算域都采用大长宽比的各向异性单元则会得到质量很差的数值解[1].因此,边界层网格生成的难点在于如何用尽可能少的网格单元达到捕捉层内流动特性的目的.

    三维黏性混合网格策略,即在边界层内生成高质量的半结构化三棱柱单元而在其他远场区域生成各向同性的四面体网格单元,是目前最常用的方法,这种方法兼具结构网格和非结构网格的优点. Pirzadeh[2]提出的层推进法(advancing layer method)是其中最具影响的一种.层推进方法从需要生长边界层的表面三角网格出发,沿通过表面顶点的某一个方向布置各向异性网格的节点,然后根据表面网格的拓扑关系连接这些布置的节点来生成一层一层的三棱柱单元,有时由于求解器的原因,三棱柱单元可以转化为四面体,如图 1所示.

    图  1  三棱柱单元生长示意图
    Figure  1.  Prism elements growing example

    围绕层推进法,很多学者做了不同程度的改进,Kallingderis等[3]改进了推进向量的求法,提出了可见区域的概念,利用该方法求得的向量可以保证生成的单元都有正的体积. Löhner[4]提出了连接模型边界上生成的半结构化网格与其他区域生成的各向同性网格的方法,在半结构化边界层网格生成过程中,检测形状差, 反转、相交、尺寸不适合的单元,将相应的单元删除. Aubry等[5]提出了一种求节点法向的方法,并将其用于表面各向异性网格的生成[6-7],最近又利用Voronoi图来确定非流形角点处的生长法向,进而更好地生成边界层网格,并集成到软件中[8]. Connell等[9]实现了全四面体黏性网格的生成,在确保相容性的前提下将三棱柱单元全部转化为四面体单元,通过删除单元的方式处理相交问题. Garimella等[10]引进多生长方向,提高了非流形角点处的单元正交性及单元质量,提出采用收缩、平滑、裁剪等操作处理层推进过程中的单元自交及前沿相交情况. Athanasiadis等[11-12]则提出了构造虚拟点来实现多法向生长及法向融合操作,提高了凹凸区域黏性网格的正交性及结构化特性. Sharov等[13]通过生成各向异性的表面网格来适应表面凹凸区域的半结构化黏性网格生成.王刚等[14-15]在这方面做了一些工作,将层推进方法应用到黏性网格生成.陈建军等[16]也在已有的研究基础上,集成各种处理方法,发展了一套混合网格生成程序.

    以上所述的方法,在出现相交情况的时候采用删除单元或停止推进的方式进行处理,这样会导致各向同性网格生成器的初始表面网格尺寸不均匀,存在长宽比很大的狭长单元,不利于各向同性网格的生成; 同时边界层网格的层状结构也会破坏.另外,该方法在凹凸脊及角点等几何特征处需要特殊处理,现有的处理方法都会产生一定的问题,尚未合理地解决通用性的问题.

    除了层推进方法,一些学者也提出了其他方法. Dyedov等[17]将偏移面(face offsetting)方法[18]应用到黏性混合网格生成中,考虑面追踪,并直接推进表面三角形单元,事先确定推进距离的上限,避免相交检测,并提出采用变分优化的方法对三棱柱网格进行质量改善,该方法通过求解拉格朗日演化方程的方式来生成三棱柱边界层网格,并将该混合网格生成方法应用到生物力学计算中,取得了一定效果.

    与层推进方法不同,Ito等[19]以及Alauzet等[20]发展了一种黏性混合网格的生成方法,先生成全域的非结构各向同性网格,然后采用挤压内部网格的方式生长边界层半结构化网格.前者不改变各向同性网格的拓扑结构,后者结合拓扑可变网格变形方法.该方法可以避免网格生长过程中的相交检测,但受限于网格变形技术,由于网格变形能力有限,在凹角点区域可以出现单元反转的现象,同时由于挤压,内部非结构网格的质量下降,影响计算的精度.

    张来平等[21-23]在生成壁面三棱柱网格方面发展了一种聚合的方法,将物面附近各向异性的四面体单元聚合生成三棱柱单元,从而得到半结构化的边界层网格,并将生成的混合网格应用到非定常数值计算中[24-25].

    本文在前人研究成果的基础上,基于前沿层推进方法,发展了一套黏性边界层网格自动生成程序,其中的细节处理方法,包括前沿多法向的确定以及生长高度的调节等,很好地解决了凹凸区域及角点附近的边界层网格生成稳定性问题,能够自动高效地生成可靠且有效的边界层网格.同时发展了一种三维尺寸调节方式,可以为各向同性网格生成器提供均匀尺寸的初始表面网格,进而得到高质量的远场各向同性网格,有利于提高流体计算解的质量.

    本文算法步骤如下:

    输入:计算区域物面三角形网格,边界层首层厚度$ h_{\rm f}$,厚度增长比例$h_{\rm r}$,层数$n_{\rm l}$.

    输出:边界层网格,及最外层三角形网格前沿.

    Step1读入物面三角形网格数据,建立网格数据结构及拓扑关系.

    Step2调整各节点总生长高度,此时生长高度为用户给定总生长高度

    $$ h_{\rm t}= \dfrac{h_{\rm f} (1 - h_{\rm r}^{n_l})}{1 - h_{\rm r}} $$ (1)

    Step 2.1提取前沿特征信息,包括特征边及角点等信息,确定节点生长线数量及方向.

    Step 2.2根据特征信息生长对应的单元,采用收缩高度的方式处理单元自交(反转)问题,得到各点的合法生长高度.

    Step 2.3对新前沿进行全局相交检测,收缩高度避免全局相交,得到各点新的生长高度.采用拉普拉斯方法进行高度光滑处理得到最终的总生长高度.

    Step3根据调整后的总生长厚度,得到每层生长高度,依次重复Step2.1——Step2.2,直到达到边界层数$n_l$,对最外层前沿进行全局相交检测,采用收缩高度的方式处理.

    边界层网格生成过程中的前沿为表面三角形网格.物面初始三角形网格剖分为初始前沿.前沿上节点$V$,其周围相邻节点记为$P_i $ ($i =1, 2, \cdots$),对应的邻边记为$e_i $ ($i =1, 2, \cdots$).

    $V$的邻边$e_i$两侧的三角面片之间的有向二面角记为$\theta _i$,$\theta _i$较大(大于240$^\circ$)的边定义为特征边.邻边中存在三条及以上特征边的节点称为角点.如图 2标出了边$e_3$及对应的二面角$\theta _3$.

    图  2  生长法向确定
    Figure  2.  Calculate growth directions

    计算每个节点周围的有向二面角,确定其邻边中的特征边,根据特征边将节点周围的三角面片分组,如图 2,$V$点周围存在三条特征边$VP_1 $, $VP_3 $, $VP_4 $,则将周围的三角面片分为3组.每一组确定一条生长方向${\pmb n}_{i} $ ($i=1, 2, 3$),由该组三角面片的外法向加权得到,$N_i $为第$i$组三角面片数,${\pmb m}_{j}$为第$j$个三角片外法向,权值$w_j $可简单取为1或者采用文献[26]中提出的角度权值

    $$ {\pmb n}_{i} = \sum\limits_{j = 1}^{N_i } (w_j {\pmb m}_{j} ) \Big / \sum\limits_{j = 1}^{N_i } {w_j } $$ (2)

    以上的法向计算策略能够很好地满足推进单元正交性及可视化条件,而且三条法向就足够了.因此本文设定节点的生长方向最多有三条.当节点周围存在三条以上特征边时,则从中挑出有向二面角最大的三条特征边作为最后的分组标准.

    节点生长法向确定后,周围三角片在该节点处有一个对应的生成法向;特征边在该点处有一到两个生长方向,可由边两侧的三角片的对应情况确定.当节点存在三个生长法向时,该节点会单独生长出一个四面体单元.

    初始生长高度为用户给定,即程序的输入中包含的第一层网格高度$h_{\rm f}$,生长高度增长比率$h_{\rm r}$及层数$n_{\rm l}$.由于节点的生长法向不一定都垂直于物面,因此有必要对初始生长高度进行适当的调整.按照调整后的高度信息逐层推进生长边界层,在推进过程中,由于边界形状的复杂性,可能出现单元自交及前沿全局相交的问题.因此也需要对当前的生长高度进行调节.

    在用户给定初始生长高度等信息作为输入后,由于上一节确定的节点生长方向不一定垂直于节点周围所有的三角面片,因此为了尽量保证推进前沿与当前前沿的距离符合用户给定的信息,需要对初始生长高度进行适当的调整.

    对于二维情况,由于生长法向是节点外角的$N$等分线($N$为生长法向数目),因此只需要考虑法向与物面垂线夹角即可. 图 3为二维示意图,前沿推进高度为$h$,节点$V_1 $有一个法向,则生长高度为$h / \sin \theta $,点$V_2 $处$\theta =90^\circ$,所以生长高度为$h$.

    图  3  二维初始高度调节
    Figure  3.  2D initial height setting

    对于三维情况,由于采用多生长方向,同时法向的确定为分组三角面片法向的加权平均,因此为了保证推进前沿与当前前沿尽量平行,则需要采用加权法调整每个节点的初始生长高度.如图 4,节点$V$周围的一个三角片$T_i $,其单元法向为${\pmb n}_{\boldsymbol T}$,它对应的分组节点生长法向为${\pmb n}_{\boldsymbol i} $(节点$V$可能有多个法向,图 4只显示了对应于$T_i $的法向${\pmb n}_{\boldsymbol i} $),则初始高度调整后为

    $$ h_i = \dfrac{h}{\cos \theta } = \dfrac{h}{{\pmb n}_{{\boldsymbol T}_i} {\pmb n}_{\boldsymbol i}} $$ (3)
    图  4  三维初始高度调节
    Figure  4.  3D initial height setting

    遍历节点$V$周围所有三角片加权平均得到节点$V$调整后的初始节点生长高度$h_V $,记$T(V)$为节点$V$周围的三角片集合,$N_{T(V)} $为集合元素个数,则有

    $$ h_V = \dfrac{1}{N_{T(V)} }\sum\limits_{T_i \in T(V)} {h_i } $$ (4)

    自交即单元出现反转,如图 5(a)所示.当出现自交时,收缩自交单元节点的生长高度,二维单元则收缩至交点处,即法向融合,这样可以避免前沿小尺寸出现,同时也可以减少之后处理自交的频次.

    图  5  自交处理二维示意图
    Figure  5.  Dealing with 2D self-intersection

    对于三维情况,当三棱柱单元出现反转时,如图 6,记三角形$ABC$推进到三角形$A'B'C'$,面积由$S$变为$S'$,则采取收缩高度的方式处理,收缩因子$\gamma $为

    $$ \gamma = \dfrac{\sqrt S }{3(\sqrt S + \sqrt {S'} )} $$ (5)
    图  6  自交处理三维示意图
    Figure  6.  Dealing with 3D self-intersection

    自交处理完以后,需要采用拉普拉斯方法对前沿各点的生长高度进行光滑处理,这样可以避免三棱柱单元邻边突变的情形.

    对于复杂的计算域,前沿层推进的过程中,不同物面推进的前沿之间可能存在相交的情况,对于前沿规模非常大的问题,全局相交检测是一个耗时的过程.为了减少相交检测的范围,可以建立八叉树(二维四叉树)结构使得全局检测局部化.尽管八叉树结构提高了全局相交检测的效率,但以往的层推进方法每一层都需要检测,使效率下降.为了提高边界层网格生成的效率,本文的方法先对全局预测一个总生长高度,即先采用总生长高度生长一层边界层,在单元生长的过程中进行自相交处理,得到每个表面节点的生长高度,然后进行拉普拉斯光滑处理,最后进行一次全局相交检测,采用收缩的方式得到每个节点的总生长高度,进而得到合适的每层生长高度.采用预测得到的生长高度进行边界层网格推进则可以减少全局相交检测的次数,提高效率.

    图 7所示的是一个二维全局收缩的示意图,收缩后的总生长高度能够较好地避免内部边界层网格前沿的全局相交,因此,在很大程度上能够避免多次的全局相交检测及收缩操作.

    图  7  全局相交处理二维示意图
    Figure  7.  Dealing with global intersection (2D example)

    根据第2节可知,前沿存在三种特征信息:角点,特征边以及三角片.因此,本文方法的单元生长方式共有如图 8四种情况,粗点及粗线在当前前沿上,箭头表示生长方向:角点生长出一个四面体单元,三角片生长出一个三棱柱单元,特征边则可以生长出一个三棱柱单元或者一个退化的三棱柱单元.

    图  8  单元生长示意图
    Figure  8.  Element growing types

    由于大多数流体求解器局限于纯四面体网格[7],为了适应大多数流体求解程序的要求,本文实现了将得到的边界层网格转化全四面体网格的过程:连接三棱柱或退化三棱柱的侧面对角线,将三棱柱细分为三个各向异性的四面体单元,如图 9.

    图  9  三棱柱侧面对角线连接方式
    Figure  9.  Connecting types of diagonals for lateral faces of a prism

    为了保持单元之间的相容性,即单元相邻当且仅当单元具有共同的节点、边或面.然而,不是所有的连接方式都能对三棱柱单元进行有效的剖分,如图 9(a)的S型连接方式则不能对三棱柱单元进行剖分得到三个四面体单元.因此在连接三棱柱侧边四边形对角线的时候需要进行特殊的处理.

    排除以上的S型连接方式,可以得到合法连接的充要条件是存在两条有公共端点的对角线.如图 9(b)连接方式将三棱柱剖分为$ABCF$, $ABFE$, $AEFD$三个四面体.

    因此本文处理连接对角线的方式为:对角线统一由编号小的点连向编号大的生长点,如图 10编号$i_1 > i_0 > i_2 $,则对角线由$P_{i_0 } $和$P_{i_2 } $连向编号大的点$P_{i_1 } $的生长点$P'_{i_1 } $,$P_{i_2 } $连向$P_{i_0 } $的生长点$P'_{i_0 } $.对于全部的三棱柱和退化的三棱柱的侧面都进行如此连接处理,则可以得到全局相容的四面体网格.

    图  10  节点拓扑连接关系示意图
    Figure  10.  Diagonal choice for prism tetrahedronization

    边界层网格生成结束之后,边界层最外层前沿将作为各向同性网格生成器的表面网格,因此最外层前沿三角网格的质量直接关系到生成的各向同性网格的质量.所以,很有必要在生成边界层网格的同时进行前沿尺寸调节,为之后的各向同性网格生成提供高质量的三角形网格输入.

    在边界层生长推进的过程中,其前沿的尺寸可能不断发生变化.二维情况下,在生成边界层网格时,表面线段在推进后的尺寸可能会变大(变小的情况按自交情况处理),如图 11中的$BC$推进到$AD$,$AD$的长度变大,记

    $$ AD = \alpha \cdot BC \ \quad \alpha > 1 $$ (6)
    图  11  二维前沿尺寸变大
    Figure  11.  2D front size control

    文献[27]在生成二维混合网格时对这种情况进行了处理,提出当$\alpha $大于给定值(本文取1.5) 时,则取$AD$中点$E$,将四边形单元$ABCD$细分为三个三角形.

    本文将其推广应用到三维尺寸控制中.由于三维边界前沿在推进过程中只出现四种情形,角点产生四面体,特征边和三角面片产生三棱柱,不管哪种情形,其表面只存在两种形状:三角形或四边形.首先标记尺寸大的边,三角形的处理方式为标记边的中点与对应的顶点相连,然后消除标记,依次进行,直至所有标记边处理结束; 而对于四边形,则采用图 11的细分方式,同时也加入三角形的操作方式; 参照图 12,红色粗线为标记的尺寸较大的边,蓝色线为新加入线.对应于三维中的调节操作可参见附录A, 另外有两种情况需要进行一次边置换操作,如图 13.

    图  12  三角形及四边形处理方式
    Figure  12.  Subdivision of triangles and quads
    图  13  边置换操作(蓝色线)
    Figure  13.  Edge swap operation (blue line)

    程序在输入物面三角形网格之后,能够根据给定的边界层网格第一层高度、层数及总高度,自动生成边界层网格.以下实例的远场网格均采用前沿推进法[28]生成.同时,为了验证本文生成的混合网格的计算有效性,采用本文方法生成的混合网格对二维NACA0012翼型、三维ONERA M6翼型、带发动机短舱F6翼身组合体进行黏性绕流计算.

    图 14是自交处理和尺寸调节结合的二维实例图.从图中可以看出,调整生长高度融合生长法向之后,单元自交得到很好的处理(凹角区域).

    图  14  自交处理和尺寸调节二维实例
    Figure  14.  2D example dealing with self-intersection and size controlling

    图 15图 16为一个凹槽黏性网格图,展示了凹角区域收缩生长高度的效果.

    图  15  凹槽黏性网格
    Figure  15.  Viscous mesh example of a cavity
    图  16  凹槽黏性网格内部放大图
    Figure  16.  Enlarged inside view of the cavity mesh

    图 17是一个二维NACA0012翼型混合网格生成实例,生成20层边界层网格,含8 394个四边形,102个三角形; 远场包括11 302个三角形. 图 18为翼型尾部放大图.从二维NACA0012实例可以看出,采用尺寸调节后,最外层前沿线段长度比较均匀,因此边界层网格与远场网格衔接较好.

    图  17  NACA0012二维混合网格
    Figure  17.  2D hybrid mesh for NACA0012
    图  18  NACA0012二维混合网格尾部放大图
    Figure  18.  Enlarged view of the 2D hybrid mesh of NACA0012

    图 19是一个三维NACA0012混合网格生成实例,生成10层边界层网格,包括36 819个四面体;远场包含22 787个四面体. 图 20为内部视图,可以看出多生长方向提高了边界层网格的正交性.

    图  19  NACA0012三维网格
    Figure  19.  3D hybrid mesh of NACA0012
    图  20  NACA0012三维网格内部视图
    Figure  20.  Inside view of the 3D hybrid mesh of NACA0012

    二维NACA0012翼型绕流的入口边界马赫数为0.8,攻角为0$^\circ$,采用Spalart-Allmaras湍流模型,得到的压力分布云图如图 21所示,其上表面壁面附近的速度矢量图如图 22所示,从结果可以看出本文采用的混合网格能够捕捉壁面附近及流场内部的信息,网格是有效的.

    图  21  二维NACA0012压力分布云图
    Figure  21.  Pressure contours of NACA0012 wing
    图  22  机翼上表面附近速度矢量图
    Figure  22.  Velocity vectors near upper wing surface

    采用本文方法生成M6的黏性混合网格如图 23所示,图 24为内部视图.

    图  23  ONERA M6混合网格
    Figure  23.  Hybrid mesh of ONERA M6 wing
    图  24  ONERA M6混合网格内部视图
    Figure  24.  Inside view of hybrid mesh of M6 wing

    为了与现有实验数据进行对比,三维M6翼型绕流采用标准实验工况,$Ma =0.839\, 5$,攻角为3.06$^\circ$,采用$k$-$\omega$湍流模型,机翼不同站位上的压力系数分布计算值与实验数据[29]对比如图 25~图 29.选取的站位$y/b$分别为0.2,0.44,0.65,0.8,0.9.

    图  25  $y/b=0.2$压力系数分布
    Figure  25.  Pressure coeffients on the wing surface at $y/b=0.2$
    图  26  $y/b=0.44$压力系数分布
    Figure  26.  Pressure coeffients on the wing surface at $y/b=0.44$
    图  27  $y/b=0.65$压力系数分布
    Figure  27.  Pressure coeffients on the wing surface at $y/b=0.65$
    图  28  $y/b=0.8$压力系数分布
    Figure  28.  Pressure coeffients on the wing surface at $y/b=0.8$
    图  29  $y/b=0.9$压力系数分布
    Figure  29.  Pressure coeffients on the wing surface at $y/b=0.9$

    从图中对比可以看出,压力系数分布总体上吻合得较好,各站位上的激波位置捕捉较为理想.总之,从对比结果可以看出本文生成的黏性混合网格能够有效应用到黏性绕流计算中.

    采用本文方法对F6翼身组合体进行黏性混合网格生成,图 30为内部网格剖分图.

    图  30  DLR-F6混和网格内部视图
    Figure  30.  Inside view of hybrid mesh of DLR-F6

    计算中,为了与实验数据进行对比,采用$Ma = 0.75$,$Re = 3\times 10^6$,$\alpha =1.738$的计算条件,文献[30]表明计算结果与湍流模型的选取有关,本文的湍流模型采用$k$-$\omega$模型.将得到的机翼不同站位上的压力系数与实验值[31]进行对比. 图 31~图 38展示了计算结果与所有8个实验站位上的对比图.从图中可以看出,计算结果与实验得到的物面压力系数分布较为符合,同时激波位置的捕捉也较为准确,因此本文生成的混合网格运用到F6翼身组合体的绕流计算是有效的.

    图  31  $y/b=0.15$压力系数分布
    Figure  31.  Pressure coeffients on the wing surface at $y/b=0.15$
    图  32  $y/b=0.239$压力系数分布
    Figure  32.  Pressure coeffients on the wing surface at $y/b=0.239$
    图  33  $y/b=0.331$压力系数分布
    Figure  33.  Pressure coeffients on the wing surface at $y/b=0.331$
    图  34  $y/b=0.377$压力系数分布
    Figure  34.  Pressure coeffients on the wing surface at $y/b=0.377$
    图  35  $y/b=0.411$压力系数分布
    Figure  35.  Pressure coeffients on the wing surface at $y/b=0.411$
    图  36  $y/b=0.514$压力系数分布
    Figure  36.  Pressure coeffients on the wing surface at $y/b=0.514$
    图  37  $y/b=0.638$压力系数分布
    Figure  37.  Pressure coeffients on the wing surface at $y/b=0.638$
    图  38  $y/b=0.847$压力系数分布
    Figure  38.  Pressure coeffients on the wing surface at $y/b=0.847$

    本文提出的黏性边界层网格自动生成策略可以很好地解决物面凹凸区域及角点处边界层网格的生长问题,避免产生非法单元,同时减少全局相交检测的次数,提高了网格生成的效率; 另外还发展了一种三维尺寸调节方法,保证了边界层网格最外层三角形网格的尺寸统一,有利于生成高质量的各向同性远场网格.最后将本文生成的二维及三维的黏性混合网格应用到高雷诺数黏性绕流计算中,得到的计算结果表明本方法生成的黏性混合网格能够满足实际计算的要求,得到合理的计算结果.结合网格变形方法,如弹簧类方法,代数插值类方法[32]等, 本文提出的混合网格生成策略也将很快应用到非定常流场计算[33-34]中.

    本文在处理自交或全局相交时只采用收缩因子的形式,在容易产生自交的凹角区域,其法向推进距离小于其他区域,前沿尺寸也相应减小,如图 39所示为凹角区域收缩效果,在该区域的未推进空间过小,生成各向同性网格会有一定困难.目前的处理方式是增加该区域表面网格的密度,减小表面网格尺寸,即生成的初始表面网格在凹凸脊及角点区域的单元会加密.同样,在一些狭缝区域,两侧壁面同时推进的情况下,由于空间狭小,会造成未推进区域过小,给各向同性网格生成造成困难.本文实例中尚未涉及此种特别狭窄的区域,全局的前沿合并策略是解决这类情况的一种好的方式,将在今后的工作中进一步考虑和完善.

    图  39  凹角网格收缩效果
    Figure  39.  Mesh shrinking at concave ridges

    另外,本文的程序将表面网格直接作为输入,计算域边界的表面三角剖分对于边界层网格的顺利进行具有重要意义[6, 35],因此,在生成边界层网格之前,如何获得适用于边界层推进的高质量的表面网格值得进一步研究.

  • 图  1   三棱柱单元生长示意图

    Figure  1.   Prism elements growing example

    图  2   生长法向确定

    Figure  2.   Calculate growth directions

    图  3   二维初始高度调节

    Figure  3.   2D initial height setting

    图  4   三维初始高度调节

    Figure  4.   3D initial height setting

    图  5   自交处理二维示意图

    Figure  5.   Dealing with 2D self-intersection

    图  6   自交处理三维示意图

    Figure  6.   Dealing with 3D self-intersection

    图  7   全局相交处理二维示意图

    Figure  7.   Dealing with global intersection (2D example)

    图  8   单元生长示意图

    Figure  8.   Element growing types

    图  9   三棱柱侧面对角线连接方式

    Figure  9.   Connecting types of diagonals for lateral faces of a prism

    图  10   节点拓扑连接关系示意图

    Figure  10.   Diagonal choice for prism tetrahedronization

    图  11   二维前沿尺寸变大

    Figure  11.   2D front size control

    图  12   三角形及四边形处理方式

    Figure  12.   Subdivision of triangles and quads

    图  13   边置换操作(蓝色线)

    Figure  13.   Edge swap operation (blue line)

    图  14   自交处理和尺寸调节二维实例

    Figure  14.   2D example dealing with self-intersection and size controlling

    图  15   凹槽黏性网格

    Figure  15.   Viscous mesh example of a cavity

    图  16   凹槽黏性网格内部放大图

    Figure  16.   Enlarged inside view of the cavity mesh

    图  17   NACA0012二维混合网格

    Figure  17.   2D hybrid mesh for NACA0012

    图  18   NACA0012二维混合网格尾部放大图

    Figure  18.   Enlarged view of the 2D hybrid mesh of NACA0012

    图  19   NACA0012三维网格

    Figure  19.   3D hybrid mesh of NACA0012

    图  20   NACA0012三维网格内部视图

    Figure  20.   Inside view of the 3D hybrid mesh of NACA0012

    图  21   二维NACA0012压力分布云图

    Figure  21.   Pressure contours of NACA0012 wing

    图  22   机翼上表面附近速度矢量图

    Figure  22.   Velocity vectors near upper wing surface

    图  23   ONERA M6混合网格

    Figure  23.   Hybrid mesh of ONERA M6 wing

    图  24   ONERA M6混合网格内部视图

    Figure  24.   Inside view of hybrid mesh of M6 wing

    图  25   $y/b=0.2$压力系数分布

    Figure  25.   Pressure coeffients on the wing surface at $y/b=0.2$

    图  26   $y/b=0.44$压力系数分布

    Figure  26.   Pressure coeffients on the wing surface at $y/b=0.44$

    图  27   $y/b=0.65$压力系数分布

    Figure  27.   Pressure coeffients on the wing surface at $y/b=0.65$

    图  28   $y/b=0.8$压力系数分布

    Figure  28.   Pressure coeffients on the wing surface at $y/b=0.8$

    图  29   $y/b=0.9$压力系数分布

    Figure  29.   Pressure coeffients on the wing surface at $y/b=0.9$

    图  30   DLR-F6混和网格内部视图

    Figure  30.   Inside view of hybrid mesh of DLR-F6

    图  31   $y/b=0.15$压力系数分布

    Figure  31.   Pressure coeffients on the wing surface at $y/b=0.15$

    图  32   $y/b=0.239$压力系数分布

    Figure  32.   Pressure coeffients on the wing surface at $y/b=0.239$

    图  33   $y/b=0.331$压力系数分布

    Figure  33.   Pressure coeffients on the wing surface at $y/b=0.331$

    图  34   $y/b=0.377$压力系数分布

    Figure  34.   Pressure coeffients on the wing surface at $y/b=0.377$

    图  35   $y/b=0.411$压力系数分布

    Figure  35.   Pressure coeffients on the wing surface at $y/b=0.411$

    图  36   $y/b=0.514$压力系数分布

    Figure  36.   Pressure coeffients on the wing surface at $y/b=0.514$

    图  37   $y/b=0.638$压力系数分布

    Figure  37.   Pressure coeffients on the wing surface at $y/b=0.638$

    图  38   $y/b=0.847$压力系数分布

    Figure  38.   Pressure coeffients on the wing surface at $y/b=0.847$

    图  39   凹角网格收缩效果

    Figure  39.   Mesh shrinking at concave ridges

  • [1]

    Sahni O, Jansen KE, Shephard MS, et al. Adaptive boundary layer meshing for viscous flow simulations. Engineering with Computers, 2008, 24(3):267-285 doi: 10.1007/s00366-008-0095-0

    [2]

    Pirzadeh S. Unstructured viscous grid generation by the advancinglayers method. AIAA Journal, 1994, 32(8):1735-1737 doi: 10.2514/3.12167

    [3]

    Kallinderis Y, Ward S. Prismatic grid generation for threedimensional complex geometries. AIAA Journal, 1993, 31(10):1850-1856 doi: 10.2514/3.11858

    [4]

    Lohner R. Matching semi-structured and unstructured grids for Navier-Stokes calculations. AIAA Paper 1993-3348, Orlando, FL, U.S.A. 1993 doi: 10.2514/6.1993-3348

    [5]

    Aubry R, Löhner R. On the 'most normal' normal. Communications in Numerical Methods in Engineering, 2008, 24(12):1641-1652 https://www.reddit.com/r/anime/comments/31klc6/what_are_the_most_normal_animes_you_know/

    [6]

    Aubry R, Dey S, Mestreau E, et al. An entropy satisfying boundary layer surface mesh generation. SIAM Journal on Scientific Computing, 2015, 37(4):A1957-A1974 doi: 10.1137/140963625

    [7]

    Aubry R, Dey S, Karamete K, et al. Smooth anisotropic sources with application to three-dimensional surface mesh generation. Engineering with Computers, 2016, 32(2):313-330 doi: 10.1007/s00366-015-0420-3

    [8]

    Dey S, Aubry R, Karamete B, et al. Capstone:A geometry-centric platform to enable physics-based simulation and design of systems. Computing in Science & Engineering, 2016, 18(1):32-39 http://ieeexplore.ieee.org/document/7274258/

    [9]

    Connell SD, Braaten ME. Semistructured mesh generation for threedimensional Navier-Stokes calculations. AIAA Journal, 1995, 33(6):1017-1024 doi: 10.2514/3.12522

    [10]

    Garimella RV, Shephard MS. Boundary layer mesh generation for viscous flow simulations. International Journal for Numerical Methods in Engineering, 2000, 49(1):193-218 doi: 10.1007/s00366-008-0095-0

    [11]

    Athanasiadis AN, Deconinck H. Object-oriented three-dimensional hybrid grid generation. International Journal for Numerical Methods in Engineering, 2003, 58(2):301-318 doi: 10.1002/(ISSN)1097-0207

    [12]

    Athanasiadis AN, Deconinck H. A folding/unfolding algorithm for the construction of semi-structured layers in hybrid grid generation. Computer Methods in Applied Mechanics and Engineering, 2005, 194(48):5051-5067 http://www.sciencedirect.com/science/article/pii/S0045782505000782

    [13]

    Sharov D, Luo H, Baum JD, et al. Unstructured Navier-Stokes grid generation at corners and ridges. International Journal for Numerical Methods in Fluids, 2003, 43(6-7):717-728 doi: 10.1002/fld.v43:6/7

    [14] 王刚, 叶正寅, 陈迎春.三维非结构黏性网格生成方法.计算物理, 2001, 18(5):402-406 http://www.cnki.com.cn/Article/CJFDTOTAL-JGXK200903015.htm

    Wang Gang, Ye Zhengyin, Chen Yingchun. Generation of three-dimensional unstructured viscous grids. Chinese Journal of Computational Physics, 2001, 18(5):402-406 (in Chinese) http://www.cnki.com.cn/Article/CJFDTOTAL-JGXK200903015.htm

    [15] 王刚, 叶正寅.三维非结构混合网格生成与N-S方程求解.航空学报, 2003, 24(5):385-390 http://www.cnki.com.cn/Article/CJFDTOTAL-HKXB200305001.htm

    Wang Gang, Ye Zhengyin. Generation of three dimensional mixed and unstructured grids and its application in solving Navier-Stokes equations. Acta Aeronautica et Astronautica Sinica, 2003, 24(5):385-390(in Chinese) http://www.cnki.com.cn/Article/CJFDTOTAL-HKXB200305001.htm

    [16] 陈建军, 曹建, 徐彦等.适应复杂外形的三维黏性混合网格生成算法.计算力学学报, 2014, 31(3):363-370 doi: 10.7511/jslx201403014

    Chen Jianjun, Cao Jian, Xu Yan, et al. Hybrid mesh generation algorithm for viscous computations of complex aerodynamics configurations. Chinese Journal of Computational Mechanics, 2014, 31(3):363-370(in Chinese) doi: 10.7511/jslx201403014

    [17]

    Dyedov V, Einstein DR, Jiao X, et al. Variational generation of prismatic boundary-layer meshes for biomedical computing. International Journal for Numerical Methods in Engineering, 2009, 79(8):907-945 doi: 10.1002/nme.v79:8

    [18]

    Jiao X. Face offsetting:A unified approach for explicit moving interfaces. Journal of Computational Physics, 2007, 220(2):612-625 doi: 10.1016/j.jcp.2006.05.021

    [19]

    Ito Y, Nakahashi K. Improvements in the reliability and quality of unstructured hybrid mesh generation. International Journal for Numerical Methods in Fluids, 2004, 45(1):79-108 doi: 10.1002/(ISSN)1097-0363

    [20]

    Alauzet F, Marcum D. A closed advancing-layer method with connectivity optimization-based mesh movement for viscous mesh generation. Engineering with Computers, 2015, 31(3):545-560 doi: 10.1007/s00366-014-0385-7

    [21]

    Zhang L, Zhao Z, Chang X, et al. A 3D hybrid grid generation technique and a multigrid/parallel algorithm based on anisotropic agglomeration approach. Chinese Journal of Aeronautics, 2013, 26(1):47-62 doi: 10.1016/j.cja.2012.12.002

    [22]

    Zhong Z, Zhang L P, Xin H E. Hybrid grid generation technique for complex geometries based on agglomeration of anisotropic tetrahedrons. Acta Aerodynamica Sinica, 2013, 31(1):34-39 http://www.en.cnki.com.cn/Article_en/CJFDTotal-KQDX201301007.htm

    [23]

    Zhang Laiping, Zhang Hanxin. A Cartesian/triangular hybrid grid solver for complex inviscid flow fields. Acta Mechanica Sinica, 1998, 30(1):104-108(in Chinese) http://en.cnki.com.cn/article_en/cjfdtotal-lxxb801.013.htm

    [24] 张来平, 王志坚, 张涵信.动态混合网格生成及隐式非定常计算方法.力学学报, 2004, 36(6):664-672 http://lxxb.cstam.org.cn/CN/abstract/abstract141195.shtml

    Zhang Laiping, Wang Zhijian, Zhang Hanxin. Hybrid moving grid generation and implicit algorithm for unsteady flows. Acta Mechanica Sinica, 2004, 36(6):664-672(in Chinese) http://lxxb.cstam.org.cn/CN/abstract/abstract141195.shtml

    [25] 张来平, 常兴华, 段旭鹏, 等.基于动态混合网格的不可压非定常流计算方法.力学学报, 2007, 39(5):577-586 http://lxxb.cstam.org.cn/CN/abstract/abstract141588.shtml

    Zhang Laiping, Chang Xinghua, Duan Xupeng, et al. Implicit algorithm based on dynamic hybrid mesh for incompressible unsteady flows. Acta Mechanica Sinica, 2007, 39(5):577-586(in Chinese) http://lxxb.cstam.org.cn/CN/abstract/abstract141588.shtml

    [26]

    Thürrner G, Wüthrich C A. Computing vertex normals from polygonal facets. Journal of Graphics Tools, 1998, 3(1):43-46 doi: 10.1080/10867651.1998.10487487

    [27]

    Dussin D, Fossati M, Guardone A, et al. Hybrid grid generation for two-dimensional high-Reynolds flows. Computers & Fluids, 2009, 38(10):1863-1875 http://www.sciencedirect.com/science/article/pii/S004579300900067X

    [28]

    Geuzaine C, Remacle JF. Gmsh:A 3-D finite element mesh generator with built-in pre-and post-processing facilities. International Journal for Numerical Methods in Engineering, 2009, 79(11):1309-1331 doi: 10.1002/nme.v79:11

    [29]

    Schmitt V, Charpin F. Pressure distributions on the ONERA-M6-Wing at transonic Mach numbers, experimental data base for computer program assessment. Report of the Fluid Dynamics Panel Working Group 04, AGARD AR-138, 1979. http://sandi.co.in/v2/home/projects/aerospace/onera_m6_wing/

    [30] 石磊, 杨云军, 周伟江等.两种湍流模型在高速旋转翼身组合弹箭中的对比研究.力学学报, 2017, 49(1):84-92 http://lxxb.cstam.org.cn/CN/abstract/abstract146220.shtml

    Shi Lei, Yang Yunjun, Zhou Weijiang. A comparative study of two turbulence models for Magnus effect in spinning projectile. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(1):84-92(in Chinese) http://lxxb.cstam.org.cn/CN/abstract/abstract146220.shtml

    [31]

    Laflin KR, Brodersen O. Summary of data from the Second AIAA CFD Drag Prediction Workshop. AIAA2004-0555, 2004 doi: 10.2514/6.2013-46

    [32] 刘中玉, 张明锋, 聂雪媛等.一种基于径向基函数的两步法网格变形策略.力学学报, 2015, 47(3):534-538 doi: 10.6052/0459-1879-14-280

    Liu Zhongyu, Zhang Mingfeng, Wei Xueyuan, et al. A two-step mesh deformation strategy based on radial basis function. Chinese Journal of Theoretical and Applied Mechanics, 2015, 47(3):534-538(in Chinese) doi: 10.6052/0459-1879-14-280

    [33] 常兴华, 马戎, 张来平等.基于计算流体力学的"虚拟飞行"技术及初步应用.力学学报, 2015, 47(4):596-604 doi: 10.6052/0459-1879-14-320

    Chang Xinghua, Ma Rong, Zhang Laiping, et al. Study on CFD-based numerical virtual flight technology and preliminary application. Chinese Journal of Theoretical and Applied Mechanics, 2015, 47(4):596-604 (in Chinese) doi: 10.6052/0459-1879-14-320

    [34] 童福林, 李新亮, 唐志共.激波与转捩边界层干扰非定常特性数值分析.力学学报, 2017, 49(1):93-104 http://lxxb.cstam.org.cn/CN/abstract/abstract146222.shtml

    Tong Fulin, Li Xinliang, Tang Zhigong. Numerical analysis of unsteady motion in shock wave/transitional boundary layer interaction. Chinese Journal of Theoretical and Applied Mechanics, 2017, 49(1):93-104(in Chinese) http://lxxb.cstam.org.cn/CN/abstract/abstract146222.shtml

    [35]

    Aubry R. Ensuring a smooth transition from semi-structured surface boundary layer mesh to fully unstructured anisotropic surface mesh//53rd AIAA Aerospace Sciences Meeting and Exhibit, AIAA, Reston, VA, 2015, AIAA-2015-1507

  • 期刊类型引用(9)

    1. 刘田田,董庆利,杨洋,冷珏琳,郑澎. 一种基于双前沿推进的边界层网格生成算法. 计算力学学报. 2022(01): 42-48 . 百度学术
    2. 王年华,鲁鹏,常兴华,张来平. 基于机器学习的非结构网格阵面推进生成技术初探. 力学学报. 2021(03): 740-751 . 本站查看
    3. 任一鹏,刘枫,林其,胡华伟. 网格生长动力学模型的建立和求解及其混合黏性网格生成中的应用. 力学与实践. 2021(03): 406-413 . 百度学术
    4. 李新平,黄明智,王刚,徐坤,刘婷婷. 基于神经网络和平直节理接触模型的细观参数标定方法. 力学与实践. 2021(03): 393-405 . 百度学术
    5. 陈坚强. 国家数值风洞(NNW)工程关键技术研究进展. 中国科学:技术科学. 2021(11): 1326-1347 . 百度学术
    6. 王年华,鲁鹏,常兴华,张来平,邓小刚. 基于人工神经网络的非结构网格尺度控制方法. 力学学报. 2021(10): 2682-2691 . 本站查看
    7. 王硕,庞宇飞,肖素梅,刘杨,齐龙,何雨阳,谢冬香. 基于网格框架的非结构附面层网格生成技术. 计算力学学报. 2021(06): 819-824 . 百度学术
    8. 靳宏宇,吴壮志,王奇,贾贺,荣伟. 面向降落伞稳态CFD计算的网格生成方法研究. 航天返回与遥感. 2019(04): 30-37 . 百度学术
    9. 杜敏韬,许善燎,吕霄,肖盼. 边界层网格对汽车风阻系数仿真计算的影响研究. 数字制造科学. 2019(03): 196-200+240 . 百度学术

    其他类型引用(9)

图(39)
计量
  • 文章访问数:  1770
  • HTML全文浏览量:  389
  • PDF下载量:  572
  • 被引次数: 18
出版历程
  • 收稿日期:  2017-05-03
  • 网络出版日期:  2017-07-06
  • 刊出日期:  2017-09-17

目录

/

返回文章
返回