第八章边界条件

1、第八章边界条件第八章边界条件第八章 边界条件任何数值模拟都可以认为仅仅是在物理区 域或系统的一部分中进行的。 区域的断层产生了 人工边界，在这个断层中有我们处理的物理量。 此外，还有暴露在流体中的自然边界。 边界条件 的数值处理需要特别注意。 在实际的系统中处理 不当模拟就会出现偏差。 与此同时， 稳定性和求 解方案中的合成速度同样对数值模拟有消极的 影响。下边边界条件的类型是我们在欧拉方程和 N-S 方程中数值计算最常见的几种：固体壁面外表面的远场和流体内部流出或流出的表面对称面平整切割和周期性边界。平板间的边界这些边界条件的处理问题在以后几节中会进行详细的介绍。 对于文献中进一步涉及的边界

2、 条件，比如壁面上的热辐射或者是自由表面上的 （热辐射），读者可以在 3.4 节中了解。8.1 虚拟单元的概念在我们讨论边界条件时， 我们需要提到虚拟 单元（也可以被称作虚拟点）这个概念。在规则 的网格中这种方法非常的流行。 然而，在不规则 网格中， 虚拟单元仍然有很多的优点。 虚拟单元 是在物理区域外部附加层上的一些网格点。 这个 可以由图 8.1 中的二维规则网格中看到。正如我 们看到的，整个计算区域被两层虚拟单元包围着 （由虚线标出），虚拟单元（点）通常不会像区 域内的网格一样产生（除过多平板的网格） 。尽 管它仍然有几何形状，比如体积或者表面的矢 量，但是它仅仅是虚拟的。利用虚拟单元可

3、以简化计算沿边界的通量， 梯度， 散度等等。 这是由于在边界上可以将空间 离散的模型进行扩展。 正如图 8.1中我们看到的， 在物理区域内同样可以进行离散。 因此，我们可 以在所有的“物理”网格点中求解控制方程。这 种方法可以使离散工作非常简单。 此外，所有规 则的网格点可以存在在一个单独的区域内， 这在 矢量计算中非常很有用。 虚拟单元不但包含有守 恒变量，同时也有几何量。很明显的是，虚拟单 元层必须完全覆盖物理区域外。 几何量通常由边 界的控制体积来求得。 在多网格平板中 （3.1 节）， 所有的流体变量和几何变量可以从相邻的平板 求得。图 8.1 中灰色的虚拟单元产生了一个问题， 由于不

4、清楚怎么设置它的值 （如果没有相邻网格 板），通过环绕形的离散模型，我们不需要知道 它的值。但是，在计算梯度（粘性通量见 4.4 节）或者在多网格中计算转移量时就很有比必要 知道它的值。通常，如图 8.1 中用箭头表示的一 样，相邻的“规则的”虚拟单元的平均值是很有 用的。iIIiFigure B,1: Two layers of dummy cells (dashed line) around the coEiiputatiaiial daniam (thick lino) in 2D Filled circles represent the standard stencil of a 2n

5、d- t)rder cell-ver tex (dual) scheme, rpft；inglcs outline the steiu 11 of a 生 iiIymiIit cell-centiel scheme (see Section 1.3).& 2固体壁面非粘性流动在非粘性流动中，流体流过表面，由于没有 摩擦力，速度矢量与壁面相切，与就是说在壁面 的法线方向没有作用力。即：V ；=o在壁面上，（81）n表示壁面上的单位法相矢量。因此，相对 应的速度V等于0。即，对流通量的矢量简化成 只有压力项：0n xPw(FJ=n yPwn zpw 一0 一(&2) 其中Pw为壁面压力。规则的单元

6、中心方案在以单元中心的方案中，在单元的质心可以 求出压力。但是，在等式(8.2)中，边界表面 单元的Pw需要求解出。我们可以很容易的通过 内部区域推导出壁面的压力。考虑到图8.2,我们可以简单的设Pw=P2。通过任一个两两相连 点，我们可以取得较高的精确度。Pw =(3P2 - P3)(63) 或者一个三点相连的推导公式是1Pw (15P2-1OP3 3P4)8(64)为了解释网格的延伸，到壁面的距离可以替代常 数系数。【1】上式的推导公式8.3和8.4没有涉及到网格和表面的几何形状。一个可供选择的方法一Rizzi 发明的被称作为法向动量的关系。它是基于在非粘性流动中，壁面是流线型的。在式 流

7、线型壁面的法线流动方向的导数为 量方程可以替换成如下的公式：8.1中沿着0。并且动fw r v （v 人）n = n 人 p（8.5）等式（8.5）包括密度，速度和壁面法线方向的 压力。法向动量在【1】中给出了精确的结果。 但是在复杂的几何表面，法向动量的数值求解存 在着问题。推导中详细的描述和精确的比较可以 在文献【1】中找到。虚拟单元中守恒变量的值在壁面内部能够 被线性的推导出：W* =2N2 -w3Wo =3W2 - 2W3（8.6）公式（8.6）中的参数和图8.2中的相对应。 如果虚拟单元在空间上离散后可以被利用，那么 对流通量的计算就可以和公式（8.2） 致。Figure 8.2:

8、Solid wall ijouiidai y )ndition for the cvll-centrol s 1). Locations, where the convective fluxes Eq. (8.2) are evaluated, are marked by diamonds Compare also to the sketch iri Fig. 4.6.规则的单元顶点的方案通过重叠的控制体积（见422）,对单元顶 点的离散可以直接由方程（& 1）在边界条件上 应用。壁面的对流通量可以通过公式（ 8.2）来 求解。壁面的压力Pw可以通过计算平均节点的 值来表示，二维公式节点的值可

9、以通过（4.26） 计算，三维公式节点的值可以通过（4.27）计算 分类公式（二维公式中的（4.30）现在仅仅解 释了两个单元（在三维中是4个）。8.2）因此，通过在重叠的控制体积内，单元顶点的值可以通过许 多不同的方法求解。一种方法是通过公式（ 将壁面上每个表面的控制体积分离 图8.3,我们可以写出：010(nx )i J,2( Pw )i J/4,2(nx)i,2( Pw)i*/4,2(ny )i A2( Pw) i 4/4,2+(ny) i,2 ( Pw )i 出/4,2(nz)i 4,2( Pw)i 4/4,2(nz) i,2 ( Pw)i */4,2- 0 一_ 0 一(Fcw)i,

10、27）公式（8.7）中的压力（Pw）ig和（Pw）i*2可以通过线性的插值表示: 对应的三维公式可以表示未划分的网格。(Pw)i 1/4,2= 43(Pw)i，2 (Pw)i1,2(8.8)另外一种可能的方法是在各自的壁面节点 利用公式（&2）来进行计算。壁面的压力Pw可 以简单的设它等于 P2。单位法向矢量可以通过 节点2上所有的面的法向矢量的平均值来表示。 这种方法要求对速度矢量进行修正。通过时间步 进的求解，壁面上的速度矢量被投影到切平面。【3】【4】（Vj,2）corr =Vi,2 - Vi,2 （nav） （nav）（ 89）nav是平均单位法向矢量。通过这种方法，流 动将会变成与壁

11、面相切。为了对虚拟点进行赋值，在内部流场中，利 用与公式（&6）相关的公式推导守恒变量是非 常有效的。不规则的以单元为中心的方案公式（8.1 ）中壁面的边界条件也同样适用 于未构建以单元为中心的方案。如果边界单元是 四面体、六面体或者是棱形（壁面上有三角形的 表面），压力可以通过公式（8.3）来推导。相邻 的单元（图8.2中的3号）可以通过中表 面的数值来代替。对于三角形或者四面体单元，在文献【5】【6】中通过一层虚拟单元来表示。 壁面上边界单元的速度矢量可以表示虚拟单元 中速度分量。比如：在图8.2中的虚拟单元1速 度可以如下来表示：W =v2 _2V2 n(8.10)在这里V2 =业兀 v

12、2ny w2nz是推导的速度，n = nx,ny,nzT 表示的是壁面的单位法向矢量。我们假设虚拟单 元中的压力和密度与边界单元的值相等。 （即就 是 Pw = p2）。规则的双重中线的方案公式（XI）中的边界条件需要更多的注意 离散化的双重中线问题。图8.4和图8.5分别表 示了二维和三维的情况。公式（&2）中的对流 通量可以在壁面上的每一个控制体积上的表面 分开来计算。这种方法与第一种构建单元顶点的 方案相同。对于一个四面体单元（就像图8.4中的1-34-5）,压力可以通过公式8.8来代替。对 于六面体，棱形或者棱锥形，控制体积的表面是 四边形（就像图8.5的面1-4-5-6），推导的公式

13、 可以如下来表示：(8.11Pint(9口 3P4 3P6 P5)1）如果边界元素是四面体（或者二维中的三角 形），壁面的压力应该通过有限体积法来计算。【7】在壁面的1-2部分，例如，1-2*表面的压力 可以通过如下的公式计算：Pini =（5Pi P2）（ 8.12）6对于四面体，例如图8.5中的壁面1-2-3,压 力可以如下表示：Pint =g（6Pi P2 P3）（ 813）8Figure 8.4: Solid wall boundary condition for the 2-D unstructured, dualcontrol volume mixed-grid scheme. L

14、ocations, where the convective fluxes Eq. (8.2) are evaluated are marked by diamondsFigure 8.5: Solid wall boundary condition for the 3-D unstructured, mixed- grid scheme Locations, where the convective fluxes Eq. (82) are evaluated, are marked by dianionds.& 2.2粘性流动对于经过固体壁面的粘性流体，介于表 面和流体的中垂直于表面的为零。

15、因此，我 们我们可以称之为非滑动的边界条件。在一 个固定的壁面，速度分量在表面变为：u=v = w= O（ 8.14）对于非滑动的边界条件，有两个基本的 结果。第一点，我们不需要解壁面上的动量 方程，这个在单元顶点方案中已经应用。 第 二点，对流通量通过非滑动的壁面的公式已 经在公式（&2）中给出，并且在公式（2.24） 中的项已经简化成6 = kK，因此，对流通量 中的壁面压力同样可以通过非粘性流动来 求得。但是，虚拟的单元处理的方法不同。以单元为中心的方案公式（8.14）中非滑动边界条件的处理 可以通过利用虚拟单元来简化。在一个绝热 的壁面（没有热通量通过壁面），我们可以 取（如图8.2）

16、= 2，E厂 E215）同样的对于单元0和3。方法同样适用 于构建的以单元为中心的方案的和非构建 的以单元为中心的方案（参考文献【6】）。如果壁面的温度已知，速度分量仍然像 方程（8.15）中一样是相反的。利用壁面的 温度可以将周围的温度线性的表示。由于压 力的梯度垂直于壁面并且为零，边界元素的 压力同样可以用虚拟单元来表示（例如 p = p, = p2）。虚拟单元中的密度和总能量可以 通过推导的值来计算。单元顶点的方案由于动量方程不需要求解，壁面上的对 流通量没有作用。公式（2.23）中的粘性通 量仅仅对能量方程中垂直于壁面的温度梯 度有作用。对于一个绝热的壁面，工n为零。 因此，我们没有必

17、要求解任何壁面上的对流 或者粘性通量。为了防止避免节点上非零速 度分量的产生，动量方程中剩余的项都为零。在已知壁面温度的情况下，我们可以直 接设壁面（例如图8.3中的节点（，2）的总能量为（假设为理想气体）（8.16）在这里，Tw表示给定的壁面温度。剩余 的动量和能量方程都设为零。对于未建立一 单元为中心的方案同样适用。另外一种方法，对于一些应用来说比较 简略，并且根本没有解决壁面的控制方程。 同样的，密度和能量可以直接简化为4,2 二P,3Pi,3TwRi,217)方程（8.仃）中假设垂直于壁面没有压 力梯度（因此P2 = P3）。由于所有的守恒变量 都是一定的，剩余的所有的方程都等于零。

18、这个方法同样适用于未构建的网格。但是， 压力的推导在三角形或四边形网格中需要 一些附加条件。如果壁面是绝缘的，虚拟单元中的值如下所示:(8.18)3，巳,1 二 Ei,3Ui,l=Ui,3，v,1=7,3，Wi,i=Wi,3以上同样适用于节点0和4.如果壁面的 温度已知，虚拟单兀中的温度可以从内部推 导出：Tm=2Tw-Ti,3 和 Tie=3Tw-2Tu(8.19)按照图8.3所示，速度分量和公式 (8.18)中的同样是相反的。密度和能量可 以通过内部的温度值和压力Pi,3来算出。8.3远场外部流体通过机翼，翅膀，汽车和其它另 外的结构数值模拟通常是在一个有限的区域内 传导。正是如此，人为加

19、上的无限远的边界条 件就很有必要。无限远边界条件的数值表示通 常要满足两个基本的要求。第一，与无限的区 域相比较，无限远处的横截面应该对流体的求 解没有显著的影响。第二，任何流出的扰动不 能对原来的流场有影响。由于它们具有抛物线 的特点，低于或接近于音速的流体通常对远场 边界条件比较敏感。处理不当的话就会降低稳 定状态流动下的集中。此外，解的精确度也会 受到影响。在人为规定的边界上求解吸收流出 波浪的方法有很多。 【10】【 15】。我们在文 献【 16】中可以找到以前所学的非反映条件的 边界。在以下两个部分中，我们将讨论由Whitfield和Janus提出的特征变量的概念，我 们将在远场边界

20、条件下讨论升力体。特征变量的概念通过雅克比对流通量的特征值， 在特征变 量的条件下，信息通常在计算区域内外进行传 递。比如，在亚音速的来流条件下有四个流进的 特性（在三维情况下）和一个流出的特性。对于 亚音速去流时情况正好相反。根据Kress的一维 理论，边界外边的特性应该等于内部的内部的特 性。剩余的条件应该决定于区域内部的解。Whitfield和Janus12的方法是基于垂直于壁 面的一维欧拉方程的特性。 这种方法可以很好的 处理大量已构建的网格和未构建的网格问题。 它 不但适用于远场条件， 也同样适用于非粘性固体 壁面。在远场条件下两个基本的流动情形可以在 图 8.6 中看到。流体可以自

21、由的进出所限定的区域。因此，根据当地的马赫数，可以得到四种不 同的远场条件：1. 超音速来流2. 超音速去流3. 亚音速来流4. 亚音速去流(a)Boundary surfaceBoundary surfaceiFigure &6： Far field boundary: inflow (a) and outflow (b) situation. Position a is out side, b on the boundary, arid positioii d is ittsidc the physical domain. Thr unit tiorniul vector u = nTJ

22、n；T points out of the domain.超音速来流对于超音速来流，所有的特征值都有相同的 符号。由于流体进入一个物理区域，边界上的守 恒变量仅仅取决于自由液面。因此，Wb =Wa(8.20)Wa的值主要与已知的马赫数M：和两个流体的角度有关。(攻角和边切角)超音速去流在这种情况下，所有的特征值同样有相同的 符号。但是，流体流出物理区域并且所有的守恒 变量决定于区域内的解。这里可以设：Wb =Wd(8.21)亚音速来流这里，物理区域内有四个特征变量区域外有 一个。因此，根据自由液面的值来计算四个特征 变量。区域外的一个特征变量可以通过物理区域 的内部四个变量进行推导。边界条件如

23、下所示：1Pb =?Pa Pd 一 ： 0Co nx(Ua -Ud) n y(Va - Vd )(Wa - Wd)卜 a (PPa)/C2Ub 二虫 一 nx(Pa - Pb)/( oCo)(8.22)Vb =Va - ny(Pa - Pb)/( ;?oCo)Wb =Wa - nz(Pa - Pb)/( ；?oCo)这里，订和co代表初始状态，初始状态通常等 于内部的点(如图8.6中的d点)。A点的值取 决于自由液面的状态。亚音速去流在亚音速去流的情况下，四个流体变量(密 度和三个速度分量)可以通过物理区域内部进行推导。剩余的第五个变量(压力)可以通过外部 得出。远场边界的初始变量可以通过文献

24、【12】 得到：Pb 二 pa几二订( Pb - Pd)/c0?ub = ud nx( pd - pb) /( : 0C0)(8.2Vb =Vd ny( Pd - Pb)/(，oCo)Wb =Wd nz(Pd - Pb)/( ?oCo)3)Pa是静态压力。虚拟单元中的物理特性可以从状态b和d线性的推导。& 3.2升力体的修正上面在远场边界条件下假定是没有循环的， 这个对于升力体在亚音速和接近声速的情况下 是不正确的。正式由于这样的原因，远场边界将 不得不远离物体。否则，流体的求解将不太准确。 如果自由液面的流动包括一个单独的点的影响 (三维中的马蹄形的顶点)，那么到远场的距离 就会大大的缩短(

25、一个数量级)。这个顶点可以 设为升力体的中心。接下来，我们将会在二维和 三维中对顶点进行修正。二维中顶点的修正我们在这里讲的方法是参考的文献【18】自由液面的速度修正分量可以用以下的表达式表示（假设流体可以压缩）：u： =u：/lil22si n 日1 - M；s in 2(。-a)2、cd2 二d 1 -M；si n2(日1 .cos-)24）:是环量，（d,r）是极坐标中远场的点：是 攻角，M:是当地马赫数，特别的，环量可以 有下式求得：?（8.25) 通过利用库塔条件，在方程（8.25），a表示机翼的弦长，Cl是升力系数，可以通过表面压力 求得。公式（8.24）中的极坐标可以表示如下：d

26、 =(X -冷 )2 W 一 yref )2(8.2f、a y - yref廿=ta n (XXref y这里Xr e和Yr e是坐标的初始点（落在顶点上， 比如在弦的1/4处）。(8.27)修正的自由液面的压力表示为:?/（丫1）- 1*（i）/ y -仃 一2其中随也2心2，修正的自由液面的密度可以通过状态方程来表示：(8.28)修正量ugVqPr和匸带替公式（8.22）或者（8.23） 中的 Ua,Va, Pa,和:a O以上的顶点修正公式（8.24） （ 8.25）只 有在亚音速的条件下才能使用。 但是，对于自由 液面条件下的修正同样适用于接近音速的情况 下。这种情况在图8.7中可以说

27、明。在图中，升 力系数到远场边界的距离是已知的。远场的半径 可以设为5%，20%，50%和99%的弦长。正如 我们看到的，没有用顶点修正的仿真和远场的距 离相对独立。相反的是，在距离为 20时，计算 的结果比较的精确。这样导致了网格单元/点的减少。这一点在文献【19】中利用高阶项来进行顶点修正来说明，远场边界能够仅仅允许 5%弦 长的精度。0360.340320.30 with vortex correction 4- without vortex cxrrectk)n0204050801000.23distance to farffetdFigure 8*7: Effects of dist

28、ance to the far fie d boundary and of single vortex on the lift coefficient. NACA 0012 airfoil, = 0.& a = 1.25.三维中的顶点修正机翼中关于远场边界的影响可以近似的认 为是马蹄形的顶点。修订的自由液面的速度分量 可以从文献【20】【21】得到：I |.:.2u = u +A2兀*r z+lz-l * xP2 I畑-齐(z+i)2+y2B-(z_i)2 + y2Cx2 + y2B2A_W：二 W：-yyI2“z+l)2+y2(z-l)2+y2 一(8.29)其中：表示环量，（x,y,z）表

29、示笛卡尔坐标系中 的远场点，I表示一半的跨度。此外，在公式（8.29） 中，假设流体在x轴的正方向，机翼沿着z轴放 置。公式（8.29）中的A，B和C如下所示：Cz l Z - I(8.30)简写成如下的式子：J ,X22z|2.y22-x22 zi2 厂（8.3：=.1 -M ：1）其中Ms是当地的马赫数。环量可以通过式（&25）来计算。A表示实际的弦线。修正值压 力和密度可以分别从公式（8.27）和（8.28）来 求得。公式（8.22）或（&23）中的量Ua,Va,Wa,Pa，和:a 可以被它们的修正值代替。公式（8.29）中修正 的速度分量的表达式变成了无限的，其中顶点的线通过流出的边界

30、。这些点可以表示为：zH, y=0, z-l, y =0,并且X=Xfaf。为了避免出现数值奇点，在文献I 21】 中对值进行了限制(z l)2 y2和(z 一1)2 y2在公式(&29)中是对于四分之一的翼展， 例如1/2，这种测量在沿着顶点线Z = l和z = -l时降低 了距离1/2的修正速度v：和w的精度。文献【21】中的数值结果表示的是按照公式(8.29) 中提供的方法，升力和拉力系数精度都 会降低。通过实验发现距离远场边界7*1时两个系数会得到精确的值。8.4入口和出口的边界基于N-S方程下有很多的方法处理出口和 入口的边界条件。在这里，我们主要关注在涡轮 机械应用中发展起来的方法

31、。适合于非反映的进 口和入口边界条件在文献【27】-【30】中有所 描述。Giles的文章【31】和Hirsch和Verhoff32 提出了欧拉方程形式的非反映边界条件，这个边 界条件是介于主体和进口或者出口平面的一段 小小的距离在确定的条件下，进口和出口的边界是周期 性的，它与速度、压力以及温度的梯度有关。这 种类型的流体偶尔可以碰到，比如，在热交换的 仿真中。周期性的进口和出口的边界条件的问题 在文献【34】-【36】里关于不同频道的LES (发 射脱离系统)中可以找到。亚音速的进口一个普通的公式包括总压力，总温度和两个 流动角。一个特征变量从流体区域内可以推导出 来。一个可以通过黎曼不变

32、式解出，其定义如下：r+十(8.32)参数d表示区域内部的状态(如图 &6a)。 黎曼不变式用来决定是绝对速度还是边界上的 速度。通过实验可以发现，选择声音的速度可以 得到一个较稳定结果，特别是在低马赫数下，因 此，我们设=j (1)cos2e+2c2 丫_1 ；6-(Y-1)cos2e+2、cos $ (Y_1)(R_)22 :(8.33) ：是与边界相关的流动角，co是声音的滞止速度， 因此：COST =vd nvd并且Co =C2 2Td(8.34)(8.35)这里|vd 2表示的是在内部点d的合速度（图8.6） 公式（8.34）中的单位法向矢量n垂直区域向外。RTbbv=f 2cp(T

33、o _Tb)一些静态的值比如温度，压力，密度或者边 界上的绝对速度如下所示：(8.36)To和Po是已知的总的温度和压力，R和Cp分别 代表特定的气体常数和在常压下的热系数。进口 的速度分量可以将|vd|2沿着两个流体角度进行分 解。（在二维中是一个）亚音速的出口在涡轮机械中，静态压力通常用来描述出口。亚音速下的出口边界可以由公式（8.23）中流出的情况做相同的处理。仅仅是环境压力Pa换成了静态出口压力。在虚拟单元中的变量可以 通过边界线性推导以及内部的 d点8.5对称平面如果流体是线性或者平面对称，第一种情况 下肯定是没有通量穿过边界。也就是说垂直于对 称边界的速度为零。此外，下面的梯度消失

34、：垂直于边界的标量的梯度垂直于边界的切向速度的梯度沿着边界的法向速度的梯度(由于 0 ) 我们可以把上面的条件写成如下式所示：*P-n 人 U =0n I (v t) =0(8.37)+ f fe-t = (v n) = 0式中，U表示的是标量的变量，t表示的切 于对称边界的矢量。以单元为中心的方案对称边界条件的处理可以通过虚拟单元来 进行大大的简化。虚拟单元中的流体变量可以通 过反映的单元来表示。这个就意味着虚拟单元中 例如密度和压力的标量和相反的内部单元是相 同的。例如：对应的符号参考图 8.2。速度分量如公式 （8.10）中的边界上表示的一样。虚拟单元中的 法向速度的法向梯度与对称的单元

35、的法向梯度 相等，但是符号相反。单元顶点的方案（双重控制体） 两种不同的方法可以解决这个问题， 一种可 能性是通过反射边界的网格来建立缺失的半个 控制体积。通量和梯度通过在内部利用反映出的 流体变量来计算（同上） 。第二种方法是计算半 个控制体积的通量（但是不是沿着边界） 。剩余 的垂直于对称平面的分量都用零来表示。 同样对 于贴于边界的控制体积其表面的法向矢量可以 进行修正。（如图 8.4 中的 2*点）修正包括去点所 有表面矢量的分量， 其分量垂直于对称平面。 梯 度同样可以由公式（ 8.37）进行修正。8.6 坐标切割这种情况仅仅是在建成的网格中遇到。 坐标 切割表示的是人为的， 非物理

36、的边界， 它是一条 加在网格点上的线 （三维中是个平面） 它有不同 的计算坐标， 但是有相同的物理方位， 这就意味 着这个网格与自己重叠。 正如我们将在 中 看到的，坐标切割在被称作C_(Fig.11 5)或者 0一拓扑网格 中出现。流体变量和它们的梯度将通过切割连 续。最好的处理切割边界条件的方法是利用虚 拟单元。这种情况在图8.8进行了描述。正如我 们看到的，这里的虚拟层不是虚拟的，但是它们 与切面的反面的网格是一致的，或者在虚拟的点 可以直接从反向的单元得到。在以单元为中心的 方案中，通过边界单元表面的通量可以像内部流 场一样精确的计算。切面边界由两种不同的以单元为中心的方 案来解决。一

37、种是在切面产生一个完整的控制体 积。(第二部分可以通过图8.8的虚线来表示)。 利用虚拟点，通量可以当作区域内的点进行求 解。第二种方法是组合每一个分开的半个控制体 积。在图8.8b点2和点5剩余的控制体是加上 的。图中点2和点5中很好的组合起来非常重要。(a)13 2 Z111 0=5 jrFigure ft-fi： Coonluiatc cut (thick line): cdU-ccntrtd scheme (&), dual controldump srtwnip fb). Dummy cells (points) are numbered as 0 and 18.7周期性边界有一些关

38、于流场对一个或者多个坐标方向 是周期性变化的实际应用。在这种情况下，只要 仅仅模拟出重复地区中的一个就可以很好的解 决问题。剩余的物理区域通过修正使其符合周期 性的边界条件。我们可以区分两种基本的周期性 的边界。第一个是周期性的移动。这个意味着通 过坐标的移动一个周期性的边界可以转移到另 一个边界。第二种类型是周期性的边界，它是由 于坐标的旋转。因此我们可以称之为旋转的周期 性。在下面，我们将会求解以单元为中心的周期 性的边界条件以及和单元顶点方案的周期性边 界条件。我们同样需要考虑旋转的周期性。 进一 步关于周期性边界的细节可以在文献【37】【38】 中找到。以单元为中心的方案利用虚拟单元的

39、概念可以对周期性边界条件进行简化。我们可以从表示涡轮机械的图8.9为例。在垂直方向，涡轮机械的结构是周期性的。 阴影中的单元1和2分别放在较低和较高的周期 性边界。由于周期性的原因，第一个虚拟单元层 与相反的周期边界的边界单元相一致。第二层虚 拟单元与第二层的物理单元相一致等等等等。因 此，所有在虚拟单元中的标量值（密度，压力） 都对应了相应的物理单元。例如：Ui =6 和 U2 二U2（8.39）矢量值（速度，梯度）同样在移动性的周期 性变化中具有相同的关系。旋转周期边界需要一 个矢量变量的修正。这个将在后面进行讨论。单元顶点的方案（双重控制体积）这种情况在图&10中进行了描述。一种处 理周

40、期性边界的方法是由围在阴影控制体积的表面组成的合成通量来解决。剩余的在点1和点 2合起来组成完整的净通量，因此：Rl,sum = Rl R2和 R2,sum = R2 Rl(840)点1和点2的部分控制体积也应该加上。在 移动性的周期变化边界下，反相边界剩余的控制 体没有变化。例如：R1 = R1和R2 = R2于是就有 Ram =R2am，旋转性的周期变化边界在利用公式(8.40)之前需要对动量方程进行变化。Figure 8.9: Periodic boundaries (thick lines) in the case of 2D ini /structured, cell-cenlnxl

41、 scheme. Dumiiiy cells (dashod lin(9 are denoted by the (piinieci) numbers of the cornspotiding physical cellsFigure 8.10: Pei iodic boundaries (thick lines) in the case of 2-D ui卜 /structured, cell-vertex scheme with dual control volumes. The 4 dtuninv parts of the control volumes (dashed line) are

42、 denoted by llw (primed) nunibei-s of the concrsponding control volumes at the opposite boundary. The same holds also for the1 duininy points 31 and 4Figure 8.11: Rotationally periodic boutidaries (A and 11). The rotaticual axis- is assumed to coincide with the i-axis.旋转性周期旋转性周期的情况是基于旋转的坐标系统。 因此，所有的

43、矢量项，例如速度或者标量的梯度 都会随着一同转换。标量项例如压力或者密度对 于坐标的旋转是不变的。如果我们假设旋转坐标 和x轴平行，旋转矩阵变为：_ 一1 0 0 R= 0 cos-sin (80 sin cos伞 _41)介于周期性边界 A和B的角顺时针是正 的。因此，比如，速度矢量从边界 A移动到B 可以写成(图8.9中的单元1,2和图&10中的点 1 -4)VB = Rva（8.42） 很容易看出X轴分量VA（i.e,UA）通过旋转没有变 化。因此Ub=Ua。所有的流体量的梯度以相同的方 式转移。正如上面陈述的，在以单元为顶点的方案 中，剩余的动量方程在使用公式（&40）必须要 进行修正

44、。式（8.41）中的旋转矩阵的应用可以 变成如下的公式：U,V,WU, V,wU, V,wRB,sum = RbRRa(8.43)图（8.43）的上标表示的是三个动量方程。8.8介于平板网格的界面在讨论3.1节中构建的网格进行空间离散 时，很明显通常它可能在复杂的几何区域内产生 一个单独的网格。（如图3.4）。我们提出两种可 能的方案来解决这个问题。第一种方法是多平板 方法，第二种是Chimera技术。在下文，我们将主 要阐述多平板方法。关于多平板进一步的讨论， 读者可以参考文献【39】-【45】。文献【46】对 多平板方法很好的进行了介绍。更多的关于 Chimera 技术的可以在文献【 47

45、】-【 51】，这里我们 不做讨论。在多平板方法中， 物理区域被分割成确定数 量的虚拟部分。 因此，计算区域同样分成了相同 的平板数。总体上，求解一个特定的平板的物理 量需要流场在相邻的一个或多个平板。 因此，我 们需要提供一个可以有效的在平板之间进行信 息交换的数据结构。 如果不同的处理器来解决平 板的控制方程的话，这个结构也要求其可以交 流。第一部分的数据结构包括平板边界的编号。 图 8.12 表示了一种特殊的编号方法。这种编号 的思想总结起来如下所示：boundary 1: i =IBEG boundary 2: i =IEND boundary 3: j =JBEG boundary

46、4: j =JEND boundary 5: k =KBEG boundary 6: k =KEND 所有的平板用相同的编号方案是非常重要的。 计 算区域中的网格点指数 i， j， k 定义如下：IBEG i IEND JBEGW j JEND KBEG k KEND单元的指数 IJK 需要以相同的方法定义为 以单元为中心的方案。 由于多平板方法经常通过 虚拟的单元来解决， 因此物理单元在开始或者结 束时将会有一定的偏移（见图 8.1）每一个平板的边界被分成一些不重叠的小 块。在相同的边界条件下， 这种分类可以规范不 同的边界条件。图 8.13 对此进行了描述。唯一 区别每一个小块的方法是相应的平板