六面体网格生成办法

2 α = xη + yη2 ,
ห้องสมุดไป่ตู้
β = xξ xη + yξ yη ,
γ = xξ2 + yξ2
和代数法相比，偏微分方程法的计算较为复杂，需要求解偏微分方程，右端强迫函数的选取 也不方便，可以通过强迫函数的选取得到希望的映射网格。
・21・
上海交通大学博士学位论文
3.2.3 超单元映射法
等 参 超 单 元 网 格 生 成 法 是 由 Zienkiewicz 和 Phillips 最早提出的[31,42]，其基本思想是将划分区域分 成更简单的子块（超单元） ，子块被映射为自然坐标系 中的单位直角形（正方形） ，根据每个方向给定的分级 权系数，将自然坐标系中的直角形（正方形）离散并 反向变换回超单元。在三维情况使用超单元生成法时， 要划分的实体用 20 节点六面体单元分成数个子区域 （图3-3） 。在超单元内一点的坐标 x，y 和 z 和自然坐 标之间的关系由下式确定
・20・
第 3 章 六面体单元生成方法及相关技术
x (ξ , η , ζ ) F (ξ ,η , ζ ) = y (ξ , η , ζ ) z (ξ ,η , ζ )
(0 ≤ ξ ≤ 1,
0 ≤ η ≤ 1, 0 ≤ ζ ≤ 1)
(3-1)
多变量的坐标变换式(3-1)使用无限插值(Transfinite Interpolation)来完成。在三维情况下，使 用混合函数和与之相联系的参数（特定点的位置及偏导数）来显式地决定式(3-1)，然后通过对每 个单变量的循环完成无限插值。一般来说，这些混合函数和参数都选取为区域边界处的函数和参 数。文献[40]给出了几种无限插值函数及其参数的选取。 代数法无限插值网格生成法的特点是计算简单，可以采用中间变量的方法方便地控制网格的 密度，对边界简单的区域，可以生成质量较高的网格，但缺点是不适应复杂边界的划分，边界不 规则时生成的网格的质量很差，并可能产生奇异性。可以通过将划分区域分解为子区域的方法， 在子区域上应用，可以在一定程度上克服这些缺点，但不易实现自动划分。
扫描线
目标
源
图3-4 扫描法生成六面体单元 Fig. 3-4 Hex elements generated by sweeping
第3章 六面体单元生成方法及相关技术
3.1 引言
有限元法是求解偏微分方程描述的连续体问题的一种近似工程方法。为克服实际连续体问题 难于处理的问题，它将分析区域离散化，将偏微分方程转化为线性方程组，采用数值计算方法求 出连续体问题的近似解。用有限元法进行工程分析的主要过程包括三个阶段： （1）有限元模型的 建立和数据准备； （2）用软件分析计算； （3）分析结果的判断和评定。迅速而合理地划分有限元 网格是完成有限元分析的前提和保证。 有限元网格生成技术发展到现在，已经出现了大量的、不同的实现方法 [22,23,24,25,26,27,28,29,30,31, 32,33,34,35,36,37,38,39] 。网格生成方法一般可以分为两类：结构化网格生成方法和非结构网格生成方法。 所谓结构化网格，严格的讲，是指网格内部节点具有相同的邻接单元的网格，即网格内部节点具 有相同的拓扑结构。结构化网格生成通常使用复合的迭代平滑技术，用边界或物理区域来排列单 元，划分区域的边界不能太复杂，以便能将区域分解为拓扑结构的块。所谓非结构化网格，放松 了对节点的连接要求，允许一个节点属于任意数量的单元。虽然四边形和六面体单元也可以是非 结构网格，但通常非结构网格是指三角形和四面体网格。目前二维区域网格自动划分发展较为成 熟，许多商业软件都提供了成熟的工具。而三维区域的网格自动划分，除个别商业软件提供了四 面体网格的自动划分外，还不成熟，特别是形状较规范的六面体网格的自动生成。 一般来说，由程序自动划分的网格总有一部分的单元的形状不是很好，需要对网格的质量进 行优化。节点序号的的标注方法直接影响刚度矩阵存储容量的大小从而影响计算的效率，因而需 要对网格的节点重新编号，以得到具有较小带宽的网格系统。 本章主要讨论了三维六面体网格生成的方法及其相关技术。
N i (ξ , η , ζ ) =
图3-3 20 节点六面体单元 Fig. 3-3 20-node hexahedral element
1 (1 + ξξi )(1 +ηηi )(1 + ζζi ), i = 1, 2, L, 8 8
(3-6)
在各棱边的中间节点定义为
1 1 − ξ 2 (1 + ηη i )(1 + ζζi ), i = 9, 11, 17, 19 4 1 N i (ξ , η , ζ ) = (1 + ξξi ) 1 − η 2 (1 + ζζi ), i = 10 , 12, 18, 20 4 1 N i (ξ , η , ζ ) = (1 + ξξi )(1 + ηηi ) 1 − ζ 2 , i = 13, 14, 15, 16 4 N i (ξ ,η , ζ ) =
3.2 结构化网格生成
结构化网格生成方法主要分为两类：代数法和 PDE（Partial Differential Equation，偏微分方 程） 方法。 代数法包括超限插值法、 等参映射法和保角映射法； PDE 方法包括椭圆系统方程、 Poisson 系统方程、抛物线系统方程和双曲线系统方程方法。应用最多的是等参映射法。结构化网格的优 点是易于实现，在每个子区域内网格可以得到很好的控制，生成规则的结构化网格，而且能生成 曲面网格。缺点是对不规则的形体有时生成的网格质量很差，需要事先根据所要产生的网格类型 将目标区域分割为一系列可映射的子区域，这一步通常需要手工完成，因而自动化程度低，不适 合网格全自动生成。提供结构化网格生成的软件有[55]3DMAGGS（NASA，Langley & Lockhead ， 四边形和六面体单元） 、 CAGI （ERC，Miss. State，四边形和六面体单元） 、CSCMDO（CSC/NASA LaRC，四边形单元） 、EMC2（INRIA，四边形单元） 、GENIE++（ERC，Miss. State，四边形和六
(3-2)
右端的强迫函数为
N ξ − ξ m −c m ξ −ξ m e − ∑ bi ξ − ξm m=1 i =1 M N η − η m − cm η−ηm Q (ξ ,η ) = − ∑ am e − ∑ bi η − ηm m =1 i =1
P (ξ , η ) = − ∑ a m
M
ξ − ξ m − d i [(ξ −ξ i )2 +(η−ηi )2 ]2 e ξ − ξm 1 η − η m −d i [(ξ −ξi )2 +(η−ηi )2 ]2 e η −ηm
3.2.1 代数法网格生成
代数法网格生成技术[40]是通过插值函数将理想的直角坐标系表示的计算区域变换为实际物 理区域来实现的（图3-1） 。它的具体做法是将计算区域的直角坐标用均匀的间隔划分为计算区域 的网格，通过变换，计算区域上均匀分布的坐标被映射为物理区域上的坐标（图 3-2） 。网格点数 目的多寡并不影响变换的特性。从计算区域到物理区域的变换是通过一个向量函数完成的
计算区域 物理区域
ζ η ξ
z y x
图3-1计算区域和物理区域 Fig. 3-1 Computational doma in and physical domain 计算网格 物理网格
x = x(ξ, η, ζ) y = y(ξ, η, ζ) z = z(ξ, η, ζ)
ζ η
z
ξ
y x
图3-2计算区域网格和物理区域网格 Fig. 3-2 Grids in the computational domain and in the physical domain
3.2.2 偏微分方程法
类似于代数法，偏微分方程法也是通过映射将曲的边界映射为网格坐标线，只不过这里使用 的映射函数满足某一类微分方程。以二维椭圆偏微分方程为例[41]，满足最大条件的椭圆型微分方 程组可以为
ξ xx + ξ yy = P(ξ , η ) η xx + η yy = Q(ξ , η )
・19・
上海交通大学博士学位论文 面体单元） 、ICEM CFD（ICEM CFD Engineering，四边形和六面体单元） 、QulkGrid （Perspective Edge，四边形单元） 、TrueGrid（XYZ Scientific Applications ，Inc.，四边形和六面体单元） 、VGM （NASA，Langley & Lockhead ，四边形和六面体单元）等。
x = ∑ N i xi , y = ∑ N i y i , z = ∑ N i z i
i =1 i=1 i =1 20 20 20
20 5 16 13 12 1
8
19 18
ζ η ξ
7
17
6 15 14 3
4 9
11 10 2
(3-5)
式中，Ni(ξ, η, ζ)是与按曲线坐标系ξ，η和ζ定义的与每 个节点相联系的形函数， ξ， η和ζ的取值范围为−1 到 1。 20 节点六面体单元形函数在各角点定义为
(
)
(
)
(3-7)
(
)
如节点的坐标已知，则与特定坐标ξ，η和ζ相对应的点的坐标可由式(3-5)计算出。 沿每个局部自然坐标，给定划分的数目和划分的间距，这些子区域（超单元）被划分为更小 的单元。以 ξ为例，划分时自然坐标按下式增加
ξ i = ξ i −1 + 2 (Wξ )i WξT
(3-8)
这里WξT = ∑ (Wξ )i ，ξ0 = −1。 超单元映射法计算简单，通过对自然坐标的不同划分，可以得到密度变化的网格，但超单元 描述的边界和划分区域的实际边界存在误差，不能适应复杂的边界。