笛卡尔网格生成技术概述（节选）

x99

......（省略了）
§1．2 计算流体动力学中网格生成技术的发展
计算流体动力学作为计算机科学、流体力学、偏微分方程数学理论、计算几何、数值分析等学科的交叉融合，它的发展除依赖于这些学科的发展外，更直接表现于对网格生成技术、数值计算方法发展的依赖。
在计算流体动力学中，按照一定规律分布于流场中的离散点的集合叫网格（Grid），分布这些网格节点的过程叫网格生成（Grid Generation）。网格生成对CFD至关重要，直接关系到CFD计算问题的成败。1974年Thompson等提出采用求解椭圆型方程方法生成贴体网格，在网格生成技术的发展中起到了开创作用。随后Steger等又提出采用求解双曲型方程方法生成贴体网格。但直到二十世纪八十年代中期，相比于计算格式和方法的飞跃发展，网格生成技术未能与之保持同步发展。因而从二十世纪八十年代开始，各国计算流体和工业界都十分重视网格生成技术的研究。二十世纪九十年代以来迅速发展的非结构网格和自适应笛卡尔网格等方法，使复杂外形的网格生成技术呈现出了更加繁荣发展的局面。现在网格生成技术已经发展成为CFD的一个重要分支，它也是计算流体动力学近二十年来一个取得较大进展的领域。也正是网格生成技术的迅速发展，才实现了流场解的高质量，使工业界能够将CFD的研究成果——求解Euler/NS方程方法应用于型号设计中。
随着CFD在实际工程设计中的深入应用，所面临的几何外形和流场变得越来越复杂，网格生成作为整个计算分析过程中的首要部分，也变得越来越困难，它所需的人力时间已达到一个计算任务全部人力时间的60%左右。在网格生成这一“瓶颈”没有消除之前，快速地对新外形进行流体力学分析，和对新模型的实验结果进行比较分析还无法实现。尽管现在已有一些比较先进的网格生成软件，如ICEM、Gridgen、Gambit等等，但是对一个复杂的新外形要生成一套比较合适的网格，其需要的时间还是比较长，而对于设计新外形的工程人员来说，一两天是他们可以接受的对新外形进行一次分析的最大周期。CFD已经成功地缩短了新外形设计中所需要的风洞实验时间，但在CFD对任意外形成为一种适时的分析工具以前，新外形设计中所需要的风洞实验时间依然非常巨大。要将CFD从专业的研究团体中脱离出来，并且能让工程设计人员应用到实际的设计中去，就必须首先解决网格生成的自动化、即时性问题，R.Consner等人在他们的一篇文章中，详细地讨论了这些方面的问题，并提出：CFD研究人员的关键问题是“你能把整个设计周期缩短多少天？”。而缩短设计周期的主要途径就是缩短网格生成时间和流场计算时间。因此，生成复杂外形网格的自动化和及时性已成为应用空气动力学、计算流体力学最具挑战性的任务之一。
当今，有众多研究人员对复杂外形的网格生成技术从分区结构网格、非结构网格和笛卡尔网格三个不同的方向展开研究。
1．2．1分区结构网格方法
分区结构网格方法将原始的物理区域按不同的空间拓扑结构分成若干区域块，每个单块网格的拓扑结构简单，易于生成贴体网格，然后合并这些单域贴体网格来形成复杂外形的空间网格。常用的传统单域贴体网格生成方法可分为代数网格生成方法、求解椭圆微分方程生成方法和求解双曲微分方程生成方法，以及从求解椭圆微分方程生成方法发展而成的求解抛物微分方程生成方法。近十多年来发展了不少新的分区结构网格。比较成熟的分区结构网格方法有：
1）分区对接网格方法。  
2）分区重叠网格方法。
它们都属于分区网格生成方法，且与流场的分区计算方法以及并行计算方法密切相关，并构成这些算法的基础。对于分区对接网格方法，虽然网格块之间的关系和数据交换比较简单，但相连网格块之间必须共用网格线（二维）或网格面（三维），为了提高计算的精度，要求块与块之间的连接应尽量光滑，这给网格生成带来诸多限制和不便。随着分区重叠网格方法的引入，网格线（面）的公用限制被打破，网格生成的难度得以降低，但随之而来的是网格块之间数据交换变得复杂，对于三维复杂外形绕流，划分的子块可能上百块，合理处理块与块之间的边界就变得非常繁琐，往往需要付出大量的手工劳动。
采用结构网格方法的优势在于它易于生成物面附近的边界层网格、有许多成熟的计算方法和比较好的湍流计算模型，因此它仍然是目前复杂外形飞行器气动力数值模拟的主要方法，计算技术最成熟。但是比较长的物面离散时间、单块网格边界条件的确定以及网格块之间各种相关信息的传递，又增加了快速计算分析的难度，而且对于不同的复杂外形，它得构造不同的网格拓扑结构，因而无法实现网格生成的“自动”，生成网格费时费力。其发展方向是朝着减少工作量，实现网格的自动生成和自适应加密，具有良好的人机对话及可视化，具有与CAD良好的接口，并强调更有效的数据结构等方面进一步发展。
1．2．2非结构三角形网格方法
复杂外形网格生成的第二方向是最近应用比较广泛的非结构三角形网格方法，它利用三角形（二维）或四面体（三维）在定义复杂外形时的灵活性，以Delaunay法或推进波阵面法【27】为基础，全部采用三角形（四面体）来填充二维（三维）空间，它消除了结构网格中节点的结构性限制，节点和单元的分布可控性好，因而能较好地处理边界，适用于模拟真实复杂外型。非结构网格生成方法在其生成过程中采用一定的准则进行优化判断，因而能生成高质量的网格，很容易控制网格的大小和节点的密度，它采用随机的数据结构有利于进行网格自适应。一旦在边界上指定网格的分布，在边界之间可以自动生成网格，无需分块或用户的干预，而且不需要在子域之间传递信息。因而，近年来非结构网格方法受到了高度的重视，有了很大发展。
非结构网格方法的一个不利之处就是不能很好地处理粘性问题，在附面层内只采用三角形或四面体网格，其网格数量将极其巨大。现在比较好的方法就是采用混合网格技术【45】，即先贴体生成能用于粘性计算的四边型或三棱柱网格，然后以此为物面边界，生成三角形非结构网格，但是生成复杂外型的四边形或三棱柱网格难度很大。
非结构网格方法的另一个不利之处就是对于相同的物理空间，网格填充效率不高，在满足同样流场计算条件的情况下，它产生的网格数量要比结构网格的数量大得多（一个长方体要划分为5个四面体）。随机的数据结构也增加了流场参数交换的时间，因此此方法要求较大的计算机内存，计算时间长。在物面附近，非结构网格方法，特别是对于复杂外形如凹槽、细缝等处比较难以处理。
非结构网格与结构网格一样都属于贴体网格，模型表面网格的好坏直接关系到空间网格的质量，因而它们的模型表面网格必须同时与网格拓扑结构和当地的几何外形特性相适应，为了更好地适应其中一方面，有时不得不在另一方面作出让步，因而往往顾此失彼。因此，在生成非结构网格和结构网格时，处理模型表面又成为一个关键而费时的工作。
1．2．3笛卡尔网格方法
笛卡尔网格是CFD计算中最早使用，也是最易生成的一种网格，它不同于传统的贴体网格，笛卡尔网格中的单元基本按照笛卡尔坐标方向（X，Y，Z）排列，流场可以采用有限体积法进行模拟计算，在与模型表面相交的单元处需要给出特殊的处理，为此，必须准确计算和判断网格单元与模型表面的相交情况。笛卡尔网格可以通过简单的再划分来达到准确拟合几何物面的目的，当前，采用笛卡尔网格的计算方法虽不如采用前两种网格的计算方法那样众多，但该方法具有网格建立简单、快速、数据结构简单、网格自适应容易等特点【84】，而且可以实现网格生成的自动化而在近几年受到人们更多的关注。笛卡尔网格技术的发展应用为CFD的非专业化、工程实际应用化提供了一种新的思路，使为工程设计人员提供操作简单、计算快速的CFD分析软件成为可能。
近年来人们开始采用自适应（AMR：Adaptive Mesh Refinement）的笛卡尔网格来计算复杂几何形状的流场，即在原始的均匀笛卡尔网格基础上根据几何外形特点或流场特点在局部区域内不断进行网格细化，得到精度符合要求、分布又是最理想的非均匀笛卡尔网格，达到准确模拟外形和捕捉激波等目的。相比于结构网格和非结构网格，采用笛卡尔网格和实现网格自适应具有以下优点：
1．由于笛卡尔网格的生成不是从模型表面出发，而是采用先空间后物面的方式，模型表面网格仅仅用于物理外型的描述，因此对模型表面网格的要求不如结构网格和非结构网格那样严格，对于多部件模型，可以采用模型部件分开描述的方式，容易重新移动、旋转部件，而且不用考虑部件之间的相互关系，可以一次性生成计算所需的计算网格，使网格生成过程简单、省时。
2．相比于贴体结构网格，不需要从物理空间向计算空间的转换，不需要在分块网格之间交换复杂的流场信息，使流场计算简单，节约计算时间。流场计算中实现自适应也比较容易，较简单。
3．笛卡尔网格不存在分区结构网格中不同外形有不同的网格拓扑结构的要求，网格生成过程容易统一，对模型表面处理的依赖程度较低，因而容易写出通用的网格生成程序，网格生成过程中不需要人为干预，因而可以实现网格生成的自动化。
4．笛卡尔网格对流场空间的填充效率高，能够缩短流场计算时间。
相对于贴体结构网格和非结构网格，笛卡尔网格虽然在数据交换、数据结构、空间网格生成等方面有优势，但它在生成贴近物面的一层网格（物面层网格）时却需要做大量的工作，换句话说，在编制笛卡尔网格生成程序时，以处理物面层网格的复杂性、多样性，替代了分区结构网格中网格拓扑结构的复杂性、多样性。但一个非常重要的结果是它可以实现网格生成的自动化。因此，近年来笛卡尔网格生成方法受到人们的重视，有了快速的发展。
§1．3 笛卡尔网格生成方法的发展
1．3．1 对笛卡尔网格生成方法发展的简单回顾
在1976年，Reyhner率先采用非贴体的笛卡尔网格模拟绕进口的跨音速流动，在他发表的论文中，他预见性地深入讨论了发展笛卡尔网格生成方法所要面临的许多问题，包括准确地确定物面边界条件、相邻网格单元间光滑过渡和网格自适应等问题【1】。1978年，Purvis和Burkhalter将笛卡尔网格与有限体积法结合起来求解全速势方程，成功地求解了轴对称外形问题【2】，这种思想被Wedan和South所发展，并应用于多段翼型和内流问题【6】。在1985年，Clarke、Hassan、Salas将Wedan和South的工作再向前推进了一步：采用Runge-Kutta时间推进法求解Euler方程【13】。在1986年，Grossman and Whitaker采用笛卡尔网格和全速势方程求解三维圆锥体的超音速问题【14】，其中值得注意的是，他们给出了解决薄体问题的一种方法。在1987年，Gaffney、Hassan和Salas首次采用笛卡尔网格和Euler方程求解三维机翼问题【15】。在1988年，Choi和Grossman成功地应用二维矢通量分裂法求解Euler方程，解决了喷口处超音速流场模拟问题【20】。
在这前十年的后一阶段，采用笛卡尔网格进行跨音速数值模拟不如采用贴体结构网格和非结构网格那样广泛。一个例外就是波音公司开发的TRANAIR软件，该软件首次对任意外形采用笛卡尔网格求解三维全速势方程。
除了TRANAIR外，早期针对笛卡尔网格的各种应用都不包含网格自适应，因此，一个很大的问题就是为了捕捉几何外形细节和流场特征，在实际应用中必然存在一些网格单元效率低的特点，这增加了计算过程中对CPU和内存的要求。当前，随着自适应笛卡尔网格的应用，人们在精确的边界条件、稳定性、矢量化、高效、内存需求和高效的数据结构等方面展开了更广泛的研究。
在LeVeque、Berger的一系列文章中【21，24，29，30】，详细地阐述了采用自适应笛卡尔网格时的各种物面边界条件，并给出了包含强激波的二维无粘非定常流情况下的处理结果。在Bell、Colella和Pembe的文章中也有关于非定常数值模拟的论述【32，43】。在1989年，Epstein、Luntz和Nachson在其文章中描述了采用重叠笛卡尔网格用多重网格法求解三维流场【26】。在1991年，De Zeeuw和Powell提出了一种在求解二维定常跨音速Euler方程情况下进行个别单元再划分的方法【33】。Morinishi提出了一种针对笛卡尔网格的求解二维Euler方程的有限差分方法，并用Runge-Kutta积分法计算了多段翼的跨音速流场【36】。在1992年，Quirk提出了对用Bezier曲线描述的外形，采用自适应笛卡尔网格求解二维Euler方程的方法【35，53】。Coirier和Powell采用笛卡尔网格方法计算了Ringleb流，并将结果与采用贴体结构网格的计算结果进行了对比【47】。Melton、Enomoto和Berger将CAD/CAM软件与采用笛卡尔网格的计算软件结合应用，使不需生成模型表面网格而直接生成空间网格成为可能【48】。
现在已有采用笛卡尔网格方法，而且形成比较成熟的商业软件，如CART3D、MGAERO，在其中包含了基于部件几何定义、多重网格Euler方程求解和用户指定参数网格自适应等内容【34，50，51，52】。Landsber、Young和Boris已将笛卡尔网格方法应用到三维外形的并行计算【72】。近期，Z .J.WANG、R.F.CHEN等采用直角/三棱柱网格和粘性笛卡尔网格求解Navier-Stokes方程【78，79，80，81，82】。 由于笛卡尔网格生成可以做到自动化和实时性，象美国NASA等航天航空机构已对笛卡尔网格生成方法有了相当的重视。在国外，有如下团体在从事笛卡尔网格生成方法的研究和开发：
学院：
Courant Institute group led by Prof. Marsha Berger
Florida State University group led by Dr. Ching Jen Chen
UC Berkeley group led by Prof. Phil Colella
U of Michigan group led by Prof. Ken Powell
James Quirk at GALCIT
Prof. Frank Evans at the University of Colorado, Boulder
商业团体：
MGAERO from Analytical Methods, Inc.
TRANAIR Development Group at Boeing - NAS 92/93
Dr. William Coirier at CFDRC
Dr. Steve Karman, Jr., at Lockheed - NAS 94/95
政府部门：
Michael Aftosmis at NASA Ames
Dr. Jay Boris at NRL
Dr. John Melton at NASA Ames
1．3．2 直角/三棱柱混合网格、“粘性”笛卡尔网格技术
生成笛卡尔网格的难点在于物面层网格的处理，在早期，人们在处理物面层网格时，直接通过笛卡尔网格的细分，用相对较细的锯齿型台阶来模拟物面，此时就不存在笛卡尔网格单元与物面的相交问题，但为了更进一步精确地模拟物面和保证网格单元的平滑过度，物面附近的网格单元就会很小，因而整个网格数量就很大。这种早期方法的另一个不足就是在计算时会有“台阶”效应，有时导致计算不稳定。
后来，人们开始采用求交切割的方法【24，44，48，58，68，69】，切割法首先在整个流场中生成笛卡尔网格，第二步去除物面内部的网格单元，第三步对与物面相交的网格单元进行切割处理，然后针对物面几何特征，如相邻两个表面三角形法向矢量的夹角的大小等，进行网格自适应划分，重复第二、三、四步，直到达到要求。采用切割的方法产生的物面层网格具有任意的形状，而且容易产生非凸的、很小的网格单元，这给网格标记、法向计算、求中心点以及对控制体采用有限体积法积分求解时，带来很大的不便，很小的网格单元也会带来求解的不稳定，而此方法最大的缺点就是不适合解决粘性问题，处理几何重叠、薄体等问题时比较困难。
为了解决粘性问题，Karman【60，62】 和Wang【74】分别提出了直角/三棱柱混合网格解决方法，Karman提出的是自适应的笛卡尔网格/固定的三棱柱网格方法，此方法的不足点在于固定的三棱柱网格部分抵消了自适应笛卡尔网格的有效性，为此Wang提出了自适应的笛卡尔网格/自适应的三棱柱网格方法，实现了用于二维粘性计算网格的自动生成【78】。生成直角/三棱柱混合网格时，同样采用切割的方法，因此同样具有切割法的缺点，即对物面描述的封闭性要求，产生不规则、非凸、小体积网格单元等等。
现在一种很有发展前景的方法是采用“投影”的方法，生成可用于粘性计算的“粘性”笛卡尔网格【75，76，77】。该方法中表面描述仅仅用于笛卡尔网格单元的划分和投影，物面层网格在笛卡尔网格前锋面向表面进行投影时自动产生，所生成的物面网格具有与前锋面相同的网格分布结构，只要将物面层网格在法向进行加密划分，就可得到用于粘性计算的“粘性”笛卡尔网格。而且对几何重叠、薄体等难题有很好的处理能力。现在Wang等人已经将此方法推广到采用2N数据结构、各向异性自适应笛卡尔网格等方面【79】。
1．3．3 笛卡尔网格生成方法的发展方向
笛卡尔网格生成方法经过二十多年的发展，在工程设计过程中得到了成功的应用。目前笛卡尔网格生成方法正向以下几个方面发展：
;发展各向异性的笛卡尔网格，以更利于网格自适应、更真实地模拟外形和提高网格填充效率，减少网格数量。
;发展“粘性”笛卡尔网格或笛卡尔网格与三棱柱网格的混合网格，用于模拟粘性流场。
;发展应用效率更高的数据结构，将笛卡尔网格应用到并行计算中。
;将生成笛卡尔网格的软件直接与CAD软件对接，进一步减少模型表面处理时间。
;将笛卡尔网格方法应用到有相对运动物体的复杂流场计算中，如外挂投放。
;将现有成熟的各种方法融合一体，实现该技术的软件化、实用化、商品化。