第二十二章 解算器
本章描述了FLUENT解算器的结构以及
使用方法
消防栓的使用方法指针万用表的使用方法84消毒液使用方法消防灭火器使用方法铁材计算器使用方法
。
目录
数值格式概况
离散
分离解
耦合解
解算器的使用概况
离散差分格式的选择
选择压力速度耦合方法
设定松弛因子
改变Courant数
Turning On FAS Multigrid
设定解的限制
解的初始化
计算
监视解的收敛性
计算期间命令的执行
收敛性与稳定性
数值格式概况
FLUENT为我们提供了两种数值方法:分离解法("FLUENT/UNS")和耦合解法("RAMPANT")。
Fluent的两种解法都可以解决守恒型积分方程,其中包括动量、能量、质量以及其他标量如湍流和化学组分的守恒。在两种情况下都应用了控制体技术,它包括:
使用计算网格对流体区域进行划分
对控制方程在控制区域内进行积分以建立代数方程,这些代数方程中包括各种相关的离散变量如:速度、压力、温度以及其他的守恒标量
离散方程的线化以及获取线性方程结果以更新相关变量的值
两种数值方法采用相似的离散过程——有限体积,但线化的方法以及离散方程的解法是不同的。首先我们在离散解法与耦合解法中讨论一般的解法,然后讨论一下线性显式与隐式中的线化方法
分离解方法
分离解算器原来是FLUENT 4和FLUENT/UNS所用的算法。使用该方法,控制方程是分离解出的(即:一个一个的解)。因为控制方程是非线性的(还是耦合的),所以在得到收敛解之前,必须进行迭代。下面是对每步迭代的介绍:
1.
在当前解的基础上,更新流体属性(如果计算刚刚开始,流体的属性用初始解来更新)
2.
为了更新流场,u,v和w的动量方程用当前压力和
表
关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf
面质量流量按顺序解出。
3.
因为第一步得到的速度可能在局部不满足连续性方程,所以从连续性方程和线化动量方程推导出压力校正的泊松方程。然后解出压力校正方程获取压力和速度场以及表面质量流量的必要校正从而满足连续性方程。
4.
在适当的地方,用前面更新的其它变量的数值解出湍流、能量、组分与及辐射等标量。5.
当包含相间耦合时,可以用离散相轨迹计算来更新连续相的源项。
6.
检查设定的方程的收敛性。
直到满足收敛判据才会结束上述步骤。
Figure 1: 分离解算器方法概述
耦合解方法
耦合解算器原来用于RAMPANT。该方法同时解连续性、动量、能量以及组分输运的控制方程(即:耦合在一起)。然后分离解方法中的分离解算器程序解附加的标量控制方程(即:和耦合方程是分离的)。因为控制方程式非线性的和耦合的,所以在获取收敛解之前需要进行适当的解循环的迭代。组成每一步迭代的步骤见上图,现概括如下:
1.
在当前解的基础上更新流体属性(如果刚刚开始计算则用初始解来更新)。
2.
同时解连续性、动量、能量和组分输运方程。
3.
在适当的地方,用前面更新的其它变量的数值解出如湍流和及辐射等标量。
4.
当包含相间耦合时,可以用离散相轨迹计算来更新连续相的源项。
5.
检查设定的方程的收敛性。
直到满足收敛判据才会结束上述步骤。
Figure 1: 耦合解方法概述
线化:隐式和显式的比较
在分离和耦合解方法中,离散,非线性控制方程被线化为每一个计算单元中相关变量的方程组。然后用线化方程组的解来更新流场。
控制方程的线化形式可能包括关于相关变量的隐式或显式形式。隐式和显式的意义如下:
· 隐式:对于给定变量,单元内的未知值用邻近单元的已知和未知值计算得出。因此,每一个未知值会在不止一个方程中出现,这些方程必须同时解来给出未知量。
· 显示:对于给定变量,每一个单元内的未知量用只包含已知量的关系式计算得到。因此未知量只在一个方程中出现,而且每一个单元内的未知量的方程只需解一次就可以给出未知量的值。
在分离解算器中,每一个离散控制方程都是该方程的相关变量的隐式线化。从而区域内每一个单元只有一个方程,这些方程组成一个方程组。因为每一个单元只有一个方程,所以常常会被称为标量系统方程。点隐式(高斯-塞德尔)线化方程解算器和代数多重网格方法 (AMG)一起被用于解单元内相关变量的标量系统方程。例如,x向动量方程被线化得到速度u未知的方程系统。方程系统的共时解(用标量AMG解算器)更新了u方向上的速度场。
总而言之,分离解方法同时考虑所有单元来解出单个变量的场(如:p)。然后再同时考虑所有单元来解出下一个变量的场,直至全部解出。分离解算器没有什么明确的选项。
在耦合解方法中你可以选择控制方程的隐式或者显式线化形式。这一选项只用于耦合控制方程组。与耦合方程组分开解的附加标量,如湍流、辐射等,的控制方程是采用和分离解方法中介绍的相同程序来线化和解出的。不管你选择的是显式还是隐式格式,解的过程都要遵循上图中的耦合解方法。
如果你选择耦合解算器的隐式选项,耦合控制方程组的每一个方程都是关于方程组中所有相关变量的隐式线化。这样我们便得到了区域内每一个单元的具有N个方程的线化方程系统,其中N是方程组中耦合方程的数量。因为每一个单元中有N个方程,所以这通常被称为方程的块系统。因为每个单元有N个方程,所以它通常被称为方程的块系统。点隐式(块结构高斯-塞德尔)线化方程解算器和代数多重网格方法 (AMG)一起被用于解单元内N个相关变量的块系统方程。例如,连续性方程和x,y,z方向动量方程以及能量方程的耦合会产生一个方程系统,在这个方程系统中,p,u,v,w和T都是未知的。用块AMG解算器同时解这些方程就会马上更新压力、三个坐标轴方向上的速度以及温度场。
总而言之,耦合隐式解算器同时在所有单元内解出所有变量(p, u, v, w, T)。
如果你选择耦合解算器的显式选项,耦合的一组控制方程都用显式的方式线化。和隐式选项一样,通过这种方法也会得到区域内每一个单元的具有N个方程的方程系统。同样地,方程系统中的所有相关变量都同时更新。然而,方程系统中都是未知的因变量。例如,x向动量方程写成的形式是为了保证更新后的x速度为流场变量已知值的函数。正因为如此,我们不需要线化方程解算器。取而代之的是,解的更新是使用多步(Runge-Kutta)解算器来完成的。在这里你可以选择全近似存储(FAS)多重网格格式来加速多步解算器。FAS多重网格的耦合显示解算器原来是用于RAMPANT中的。
总而言之,耦合显式方法同时解一个单元内的所有变量(p,u,v,w,T)。
注意FAS多重网格是显式解算器方法的一个可选部分,而AMG方法是分离和耦合隐式方法需要的部分。
离散
FLUENT使用基于控制体的方法将控制方程转换为可以用数值方法解出的代数方程。该方法,在每一个控制体内积分控制方程,从而产生基于控制体的每一个变量都守恒的离散方程。
考虑标量
输运的定常状态守恒方程可以很容易的说明控制方程的离散。下面就是写成对于控制体积V的积分形式的方程:
其中
r
=
密度
v
=
速度矢量(= u ,i(hat) + v ,j(hat) in 2D)
A
=
曲面面积矢量
G_f
=
f的扩散系数
?f
=
f的梯度 (= (秄/?x) ,i(hat) + (秄/?y) ,j(hat) in 2D)
S_f
=
每一单位体积f的源项
上面的方程被应用于区域内每一个控制体积或者单元。下面图中的二维三角单元就是控制体积的一个例子。在给定单元内离散上面的方程有:
where
N_faces
=
封闭单元的面的个数
f_f
=
通过表面f的对流量
v_f
=
通过表面的质量流量
A_f
=
表面
的面积,|A| (= | A_x i(hat) + A_y j(hat) | in 2D)
(?f)_n
=
magnitude of ?f normal to face f
V
=
单元体积
由FLUENT所解的方程和上面所给出的一般形式相同,而且很容易扩展到多维情况和友人以多边形、多面体组成的非结构网格。
Figure 1: 用于显示标量输运方程离散的控制体积
FLUENT在单元的中心(上图的c0和c1)存贮标量f的离散值。然而,方程3的对流项中需要表面值f_f,因此必须从单元中心插值。这个任务由迎风格式完成。
迎风的意思就是,表面值f_f是从单元上游或者说迎风处的量推导出来的,这个上游是指相对于方程3法向速度v_n的方向而言的,FLUENT允许你选择几种迎风格式:一阶迎风,二阶迎风,幂率和QUICK格式。这些格式在一阶迎风格式一节中介绍[95]。
方程3中的扩散项是中心差分而且一般具有二阶精度。
一阶迎风格式
当需要一阶精度时,我们假定描述单元内变量平均值的单元中心变量就是整个单元内各个变量的值,而且单元表面的量等于单元内的量。因此,当选择一阶迎风格式时,表面值f_f被设定等于迎风单元的单元中心值。
幂率格式
幂率离散格式使用一维对流扩散方程的精确解来插值变量f在表面处的值。
其中G和r u是通过间隔?x的常值。积分方程1可得如下f随x的变化关系:
其中:
f_0
=
f|_x=0
f_L
=
f|_x=L
Pe是Peclet数。
下图所示为不同Pelect数下f(x)在x=0和x=L之间的变化关系。该图表明对于较大的Pe,f在x=L/2处的值近似等于迎风值。这就意味着当流动由对流项主导时,只需要让变量表面处的值等于迎风或者上游值就可以完成插值。这是FLUENT的标准一阶格式。
Figure 1: 变量f在x=0和x=L之间的变化(方程1)
如果选择幂率格式,FLUENT用方程3等价的幂率格式[118]作为插值格式。
如一阶迎风格式所述,上图表明,对于较大的Pe,f在x=L/2处的值近似等于迎风值。当Pe=0(无流动或者纯扩散)图1表明f可以用x=0到x=1之间简单的线性平均来实现插值。当Peclet数的值适中时,f在x=L/2处的插值必须使用方程3等价的幂率插值格式来得到。
二阶迎风格式
当需要二阶精度时,使用多维线性重建方法[5]来计算单元表面处的值。在这种方法中,通过单元中心解在单元中心处的泰勒展开来实现单元表面的二阶精度值。因此,当使用二阶迎风格式时,用下面的方程来计算表面值f_f:
其中f和?f分别是单元中心值和迎风单元的梯度值,Ds是从迎风单元中心到表面中心的位移矢量。在这种情况下需要确定每个单元内的梯度?f。我们使用散度定理来计算这个梯度,其离散格式如下:
在这里,表面处的值f _f由邻近表面的两个单元的f的平均值来计算。最后,限制梯度?f以保证不会引进新的最大值和最小值。
QUICK格式
对于四边形和六面体网格,我们可以确定它们唯一的上游和下游表面以及单元。
FLUENT还提供了计算对流变量
在表面处高阶值的QUICK格式。QUICK类型的格式[95]是通过变量的二阶迎风与中心插值加上适当的权因子得到的,具体可以写成:
Figure 1: 一维控制体
在上式中
q=1就是中心二阶插值,
q=0就是二阶迎风值。传统的QUICK格式对应的
q= 1/8。FLUENT中使用一个变量,解相关的q值,以避免引进新的解的极值。
当结构网格和流动方向一致时,QUICK格式明显具有较高精度。需要注意的是FLUENT也允许对非结构网格或者混合网格使用QUICK格式,在这种情况下,常用的二阶迎风离散格式将被用于非六面体单元表面或者非四边形单元表面。当使用并行解算器时,二阶迎风格式还被用于划分的边界处。
离散方程的线化形式
离散标量输运方程(离散一节中的方程3)包括了单元中心的标量f的未知值,还包括周围相邻单元出的未知值。一般说来这些方程关于这些变量是非线性的。离散一节中方程3的线化形式为:
其中下标nb是指相邻单元,a_p和a_nb分别是f和f_nb.的线化系数。
每一个单元的邻近单元的数量取决于网格的拓扑结构,但是一般说来都等于围成该单元的表面的数量(边界单元除外)。
在网格中每一个单元都可以写出相似的方程。这样就产生了具有稀疏系数矩阵的代数方程。对于标量方程,FLUENT是用点隐式(Gauss-Seidel)线化方程解算器与代数多重网格方法(AMG)连接起来从而解出这个线性系统。
亚松驰
由于FLUENT所解方程组的非线性,我们有必要控制
的变化。一般用亚松驰方法来实现控制,该方法在每一部迭代中减少了f的变化量。亚松驰最简单的形式为:单元内变量f等于原来的值f_old加上亚松驰因子a与f变化的积:
分离解算器
在本节中,我们主要介绍了几个练习,这几个练习和动量方程和连续性方程的离散以及通过分离解算器计算出相应的解有关。通过考虑定常状态的连续性和动量方程一般形式可以和容易介绍这几个练习:
其中I是单位矩阵,t是应力张量,F是力矢。
动量方程的离散
标量输运方程的离散中所讨论的离散格式也用于离散动量方程。例如,设定f=u可以很容易得到x向动量方程:
如果已知压力场和表面质量通量,我们就可以用离散一节所介绍的方法来解方程1并获取一个速度场。但是,压力场和表面质量流率预先是未知的,因此必须作为解的一部分。关于压力的存储和压力梯度项的离散有重要的问题需要讨论,下一步我们就会讨论这些问题。
FLUENT使用共同定位格式,在这种情况下,压力和速度共同存储在单元中心。然而,如离散一节的图1所示,方程1需要单元c0和c1之间的表面处的压力值。因此,我们需要使用插值格式来从单元处压力值计算表面压力值。
压力插值格式
在默认的情况下,FLUENT使用动量方程系数[132]来处理表面处的压力插值。只要单元中心的压力变化比较光滑,这个程序都可以算得很好。当控制体之间的动量源项有跳跃或者较大的梯度时,在单元表面处的压力轮廓会有较大的梯度,此时不可以用这种格式来插值。。如果还是用这种格式,就会出现单元速度过大和过小的矛盾。
标准压力插值格式有问题的流动包括:具有大体积力的流动,如强烈的漩流、高Rayleigh数自然对流等。在这种情况下,有必要在大梯度区域处理网格从而能够较好的粪便压力的变化。
FLUENT的另一个误差源就是假定壁面法向压力梯度为零。这种假定对于边界层是有效的,但是对于出现体积力或者曲率的情况就不合适了。此外,压力梯度计算最明显的错误就是速度矢量指向或者指出壁面。
当标准压力插值格式不可用时,FLUENT提供了几种可选的方法:
· 线性格式用相邻单元压力的平均值来计算表面压力。
· 二阶格式用二阶精度对流项(见二阶迎风格式一节)中使用的方法来重建表面压力。这一格式可能会比标准和线性格式有所改善,但是如果将其用于初始计算和/或较差网格计算就可能会出现一些麻烦。
· 体积力加权格式假定压力和体积力之间差值梯度为常数来计算表面压力。如果动量方程中的体积力预先知道(如:浮力和轴对称涡流的计算),那么这一格式会处理得很好。
· PRESTO! (PREssure STaggering Option)格式使用离散连续性平衡来计算交错压力的表面处的交错体积。这一程序和非结构网格中[118]使用的交错网格格式的思想相似。
注意:PRESTO!格式只用于四边形或者六面体网格。
当使用这些可选格式时,推荐大家参阅一下选择压力插值格式一节。
连续性方程的离散
分离解算器中的方程1可以在离散一节中的图1所示的控制体上积分,从而得到如下离散方程:
其中J_f是通过表面f, r v_n A_f的质量流速。因此:
如数值格式概述一节所述,动量方程和连续性方程式按顺序解出的。在这个顺序格式中,连续性方程是作为压力的方程使用的。但是对于不可压流动,压力并没有明确出现在方程3中,这是因为密度并不是直接和压力相关。SIMPLE (Semi-Implicit Method for Pressure-Linked Equations)的一组算法[118]将压力引进连续性方程。具体可以参阅SIMPLE一节。
为了将来能够继续下去,有必要将速度的表面值v_n和速度在单元的中心的存储值关联起来。单元中心速度到表面的线性插值会导致压力的非物理校验板(checker-boarding)。FLUENT用Rhie和Chow[132]概述的相似程序来避免校验板。速度v_n的表面值不是线性平均而是动量权重平均,这个权重因子是基于动量方程离散一节中的方程1的a_P系数的。采用这种程序,表面流动速度J_f可以写成:
其中,p_c0和p_c1是表面两边单元的压力,J(hat)_f包含了这些单元中速度的影响(见离散一节中的图1)。d_f向可以写成:
a(bar)_P项是表面两边单元中动量方程a_P系数的平均。
压力-速度耦合
使用连续性方程的离散一节中的方程5来实现压力速度耦合,从而从离散连续性方程(连续性方程离散一节中的方程3)推导出一个压力方程。FLUENT提供了三种可选的压力速度耦合算法:SIMPLE,SIMPLEC和PISO。关于这些算法的选择请参阅选择压力速度耦合方法一节。
SIMPLE
SIMPLE算法使用压力和速度之间的相互校正关系来强制质量守恒并获取压力场。
如果用猜测压力场p^*来解动量方程,从连续性方程离散一节中的方程5所得到的表面流量J^*_f为:
它并不满足连续性方程。因此将校正项J^'_f加入到表面流速J^*_f中来校正质量流速J_f:
此时满足了连续性方程。SIMPLE假定J^'_f写成如下形式:
其中p^'是单元压力校正。
SIMPLE算法将流量校正方程(方程3和5)代入到离散连续性方程(连续性方程的离散一节中的方程3)从而得到单元内压力校正p^'的离散方程。
其中,源项b是流入单元的净流速。
压力校正方程(方程7)可以用代数多重网格一节中所介绍的代数多重网格方法来解。一旦得到解,使用下面的方程校正单元压力和表面流动速度:
在这里,a_p是压力亚松驰因子(请参阅亚松驰方面的介绍)。校正后的表面流速J_f在每一部迭代中同一地满足离散连续性方程。
SIMPLEC
基本SIMPLE算法的很多变量都可以在相关文献资料中查到。除了SIMPLE算法之外,FLUENT还提供了SIMPLEC (SIMPLE-Consistent)算法[164]。SIMPLE算法是默认的,但是对于很多问题如果使用SIMPLEC可能会更好一些,具体可以参阅SIMPLE与SIMPLEC的比较一节。
SIMPLEC程序和SIMPLE程序相似。两种算法所使用的表达式唯一的区别就是表面流动速度校正项。和SIMPLE中一样,校正方程可写为:
但是系数d_f重新定义为:
可以看出,在压力速度耦合是得到解的主要因素时,使用修改后的校正方程可以加速收敛。
PISO
压力隐式分裂算子(PISO)的压力速度耦合格式是SIMPLE算法族的一部分,它是基于压力速度校正之间的高度近似关系的一种算法。SIMPLE和SIMPLEC算法的一个限制就是在压力校正方程解出之后新的速度值和相应的流量不满足动量平衡。因此必须重复计算直至平衡得到满足。为了提高该计算的效率,PISO算法执行了两个附加的校正:相邻校正和偏斜校正。
PISO算法的主要思想就是将压力校正方程[69]中解的阶段中的SIMPLE和SIMPLEC算法所需的重复计算移除。经过一个或更多的附加PISO循环,校正的速度会更接近满足连续性和动量方程。这一迭代过程被称为动量校正或者邻近校正。PISO算法在每个迭代中要花费稍多的CPU时间但是极大的减少了达到收敛所需要的迭代次数,尤其是对于过渡问题,这一优点更为明显。
对于具有一些倾斜度的网格,单元表面质量流量校正和邻近单元压力校正差值之间的关系是相当简略的。因为沿着单元表面的压力校正梯度的分量开始是未知的,所以需要进行一个和上面所述的PISO邻近校正中相似的迭代步骤[51]。初始化压力校正方程的解之后,重新计算压力校正梯度然后用重新计算出来的值更新质量流量校正。这个被称为偏斜矫正的过程极大的减少了计算高度扭曲网格所遇到的收敛性困难。PISO偏斜校正可以使我们在基本相同的迭代步中,从高度偏斜的网格上得到和更为正交的网格上不相上下的解。
多相流中强体积力的特定处理
当多项流中存在较大的体积力(如:重力或者表面张力),动量方程中的体积力项和压力梯度项几乎是平衡的,相比较来说,对流项和粘性项的贡献就较小了。除非考虑压力梯度和体积力的局部平衡,否则分离算法的收敛性会很差。FLUENT提供了一种可选的隐式体积力处理,这种处理考虑了上面所说的影响从而使得解更具有鲁棒性。
基本的程序是将包含体积力的校正项增加到表面流动校正方程中(SIMPLE中的方程13)。这样,SIMPLE中的方程9就多了一个额外的体积力校正项,从而使得流动在迭代过程中提早得到真实的压力场。
这一选项只在多项流计算中使用,但是在默认情况下是关闭的。设定多相流计算的说明中包括了打开隐式体积力处理的说明,具体可以在以下几节中找到相关说明:在VOF计算中包括体积力,在气穴计算中包括体积力,在代数滑移混合计算中包括体积力。
除此之外,通过使用体积力的亚松驰因子,FLUENT允许你控制体积力中的变化。
耦合解算器
FLUENT中的耦合解算器同时解连续性、动量、(合适的话)能量和组分输运,并将它们作为一组控制方程或者方程的矢量来处理。随后会按顺序解附加标量的控制方程(也就是说这些附加标量方程相互之间是分离的而且和耦合方程组之间是分离的)。
矢量形式的控制方程
将曲面面积微分dA作为控制体积,积分描述平均流动属性的单组分流体的控制方程系统,相应的笛卡尔坐标形式为:
其中,矢量W,F和G分别定义为:
,
,
矢量H包含了诸如体积力、能量等源项。
在这里,r、v、E和p分别是流体的密度、速度、单位质量的总能量以及压力。T是粘性应力张量,q是热流量。
总能量E和总焓H的关系为:
其中:
在方程1中所表示的Navier-Stokes方程在低马赫数下会非常具有刚性,这主要是因为流体速度v和声速c相差太大。对于不可压流来说这种情况尤其真实,不管流动速度有多大,不可压流体中的声速都是无穷大。在这些情况下,方程的数值刚性会导致较差的收敛速度。在FLUENT中,我们采用耦合解算器中的一种被称为(时间导数)预处理[175]的方法克服了这种困难。
预处理
时间导数预处理方法,是用预处理矩阵先乘以矢量形式控制方程中方程1的时间导数项。这一步就重新标度了所解方程系统的声速(特征值),从而减轻了低马赫数和不可压流动中会遇到的数值刚性的影响。
推导预处理矩阵,首先是用微分学中的一个法则(具体看方程就知道了)将守恒量W在控制方程的方程1的因变量形式变形为原始变量Q的形式,结果如下:
其中,Q是矢量{p, u, v, w, T }^T,雅克比(Jacobian)矩阵?W/?Q为:
其中:
,
对于理想气体d = 1,对于不可压流体d = 0。
选择原始变量Q作为因变量有几个原因。首先,解不可压流动时它是自然的选择。其次,当我们使用二阶精度时,为了得到粘性流动中更高精度的速度和温度梯度以及无粘流动的压力梯度,我们需要重建Q而不是W。最后,选择压力作为因变量允许系统中的声波的传播被挑选出来[165]。
我们用预处理矩阵G来替换雅克比矩阵?W/?Q(方程3)来实现方程的预处理,这样,预处理系统的守恒形式为:
其中:
参数
为:
方程9中的参考速度U_r在当地选取,从而保证系统的特征值关于对流和耗散时间尺度能够调节的很好[175]。
预处理之后的方程系统(方程5)的特征值为:
其中:
在这里,d Q是空间差分Q_R - Q_L。流量F_R = F(Q_R)和F_L = F(Q_L)是用表面“右边”和“左边”的(重建解)矢量Q_R和Q_L计算出来的。矩阵|A(hat)|定义为:
其中L是特征值的对角阵,M是G^-1A对角化的模式矩阵,其中A是无粘流量Jacobian矩阵?F/?Q。
对于非预处理系统(和理想气体),方程13化简为在Roe平均用于估计G,|A(hat)|时的Roe流量差分分裂[133]。目前所使用的是Q_R和Q_L的平均。
在当前形式中,方程13可以看成二阶中心差分加上附加的矩阵耗散。附加的矩阵耗散项不仅对对流变量迎风的产生和超音速流动的压力和流量速度有作用,而且还提供了定常状态的压力速度耦合,并对低速和不可压流的稳定性和高效收敛有作用。
定常流的时间步进法
FLUENT中耦合控制方程组(预处理中的方程5)在时间上的离散既可以用于定常计算也可以用于非定常计算。在定常情况下,进行时间步进计算直至达到定常解。耦合方程可以用隐式也可以用显式时间步进来进行时间离散。这两种算法在下面的显式和隐式格式两节中介绍。
显式格式
在显式格式中,预处理的方程5是用多步,时间步进算法[73]来离散时间导数的。从第n步迭代到第n+1步迭代是用m步Runge-Kutta格式来提高解的,m步Runge-Kutta格式为:
时间步D t是从CFL(Courant-Friedrichs-Lewy)条件计算得到的
其中l_max是预处理一节中方程11定义的局部特征值的最大值。
使用全近似多重网格方法(FAS)可以加快显示格式的收敛速度,关于FAS请参阅后面将要介绍的内容。
通过与相邻残差的隐式平均增加格式的支持可以在后面的计算中增加最大时间步长。残差通过Laplacian 光滑算子进行过滤,其中Laplacian 光滑算子为:
该方程可以用下面的Jacobi迭代解出:
(译者注:分母为什么这样写?)
对于e=0.5,两次Jacobi迭代就足可以允许时间步增加一倍。
隐式格式
在显式格式中,控制方程(预处理一节中的方程5)的时间Euler隐式离散和流量的牛顿类型线化界和得到下面delta形式的线化系统[174]:
中心D和对角线下面的系数矩阵S_j,k分别为:
残差矢量R^n和时间步D分别在显式格式中的方程1和方程3中定义。
方程1是用点Gauss-Seidel格式和适应耦合的一组方程的代数多重网格方法解出的。
非定常流的时间离散
对于时间精度计算,可以使用显式和隐式时间步进格式。(隐式方法也被称为双重时间步进法)
显式时间步进
在显式时间步进方法中,使用显式格式一节所介绍的方法,在区域的每一个单元内使用相同的时间步,取消预处理选项。
双重时间步进
为了得到预处理方程的高效的,时间精度解,我们采用双重时间步进法和多步格式。在这里我们在矢量形式的控制方程的方程1中引入伪时间导数:
其中t是物理时间,t是时间步进程序中的伪时间。注意:方程1左手边的第二项t??消失了,而重新获得了矢量形式控制方程中的方程1。
方程1中的时间相关项是用一阶或者二阶精度后向时间差分进行隐式离散的。写成半隐式形式如下:
其中f_e是精确解。在解收敛之前可能会有一个和近似解f有关的误差d:
我们寻找一个校正f的y,这样,精确解由下式给出:
将方程7代入到方程3中有:
现在,使用方程5和方程9可以得到:
这是根据初始细化层面算子A和误差d的校正方程。假如局部(高频)误差在精细层面已经充分消除,就可以在下一个较粗糙的层面光滑校正y并更有效的解决问题。
限制和延拓
在粗化层面解校正方程需要将精细层面(限制)的误差转移过来,还需要计算校正方程,然后再将校正结果从粗化网格层面(延拓)传回来。粗化层面校正y^H可以写为:
其中A^H是粗化层面算子,R是限制算字,负责将精细层面的误差传到粗化层面。根据更新之后的精细层面解可以给出方程1的解为:
其中P是延拓算子,用于将粗化层面的校正传给细化层面。
非结构多重网格
在非结构网格中使用多重网格的主要困难在于粗化网格的创建和使用。在结构网格中,可以通过简单的从精细网格中移去每隔一个的网格线来得到粗化网格而且限制和延拓算子很容易用公式表示出来(如:注入和双线性插值)。
FLUENT中使用两种多重网格方法分开的方式克服了在非结构网格上应用多重网格的困难。迄今为止所讨论的基本原则和多重网格循环一节中讨论的循环策略是相同的,限制、延拓和粗化网格算子的建构技巧是不同的,具体在代数多重网格(AMG)和全近似存储(FAS)多重网格两节中分别讨论了AMG和FAS方法。
译者注:
多重网格法基本原理
微分方程的误差分量可以分为两大类,一类是频率变化较缓慢的低频分量;另一类是频率高,摆动快的高频分量。一般的迭代方法可以迅速地将摆动误差衰减,但对那些低频分量,迭代法的效果不是很显著。高频分量和低频分量是相对的,与网格尺度有关,在细网格上被视为低频的分量,在粗网格上可能为高频分量。
多重网格方法作为一种快速计算方法,迭代求解由偏微分方程组离散以后组成的代数方程组,其基本原理在于一定的网格最容易消除波长与网格步长相对应的误差分量。该方法采用不同尺度的网格,不同疏密的网格消除不同波长的误差分量,首先在细网格上采用迭代法,当收敛速度变缓慢时暗示误差已经光滑,则转移到较粗的网格上消除与该层网格上相对应的较易消除的那些误差分量,这样逐层进行下去直到消除各种误差分量,再逐层返回到细网格上。
目前两层网格方法从理论上已证明是收敛的,并且其收敛速度与网格尺度无关[哈克布思,1988]。 多重网格法是迭代法与粗网格修正的组合,经过证明,迭代法可迅速地将那些高频分量去掉,粗网格修正则可以帮助消除那些光滑了的低频分量,而对那些高频分量基本不起作用。
在多重网格计算中,需要一些媒介把细网格上的信息传递到粗网格上去,同时还需要一些媒介把粗网格上的信息传递到细网格上去。限制算子Iih(i-1)h是把细网格i-1层上的残余限制到粗网格i层上的算子,最简单的算子是平凡单射,另外还有特殊加权限制;插值算子Iih(i-1)h是把粗网格i层上的结果插值到细网格i-1层上的算子,一般采用线性插值或完全加权限制算子。
需要说明的是在多重网格迭代方法中,粗网格修正之前,细网格必须进行光滑迭代,以消除高频误差,使粗网格修正最有效地发挥其作用;在粗网格修正之后,不可避免的引入高频误差,所以也必须进行光滑迭代,不过高频误差能很快的通过光滑迭代消除。
Multigrid Cycles
多重网格循环可以定义为在每一个网格层面通过网格层次时在网格层面内应用的递归程序。FLUENT中有四种多重网格循环:V,W,F以及灵活("flex")循环。V和W循环可以用在AMG和FAS中,F和灵活循环只限用于AMG方法。(W和灵活AMG循环由于要花费大量的计算而不可用于解耦合方程组。)
The V and W Cycles
下面两个图分别是V和W多重网格循环(定义如下)。在每一个图中,多重网格循环被描述为正方形,然后被扩展来表示循环内执行的个别步骤。当你读下面的步骤时,你可能想要遵循图中所描述的步骤。
Figure 1: V循环多重网格
Figure 2: W循环多重网格
对于V和W循环,各个层面的转换由三个参数控制:b_1,b_2和b_3:
1.
b_1“滤波”(有时称为预松弛扫描),在当前网格层面进行滤波减少误差(当地误差)的高频部分。
在图1和2中,这一步被描述为一个循环,并标记了多重网格循环的起点。误差的高频分量必须减小直至所剩下的误差在下一个较粗糙的网格中没有明显的混淆现象。
如果这是最粗糙的网格层面,在这个层面的多重网格循环就结束了(在图1和2中有3个粗糙网格层,所以描述层面3多重网格循环的四边形等价于每一个图中最后一个图表的循环)
注意:在AMG方法中b_1的默认值是0(即:没有预松弛)。
2.
接着就是用适当的限制算子将问题限制映射到下一个粗糙网格层面。
在图1和2中,从精细网格限制映射到较粗糙网格层面用向下的斜线来表示。
3.
通过执行b_2多重网格循环可以减少粗化网格的误差(在图1和2中用四边形表示)。通常说来,对于固定的多重网格策略,b_2是1或者2,分别对应V循环和W循环多重网格。
4.
然后,使用适当的延拓算子,在粗化网格上积累的校正用内插替换返回到精细网格并加到精细网格解中。
在图1和2中,延拓算子用斜向上的线来表示。
现在在精细网格层面上出现了高频误差,这个误差是由于延拓程序使用传递矫正而造成的。
5.
在最后一步,执行b_3光滑(后松弛)来去掉b_2网格循环在粗糙网格中引进的高频误差。
在图1和2中,松弛程序用一个三角形表示。
注意:在FAS方法中,在这一步中使用显式格式中的拉普拉斯光滑算子,校正又会被进行光滑处理,b_3的默认值是零(也就是说,不进行这种光滑处理);在最精细的网格层面的循环结束处从来不进行b_3光滑。
代数多重网格(AMG)
我们将会看到,这一算法之所以被称为代数多重网格格式,主要是因为粗糙层面方程的生成不需要在粗糙层使用任何几何图形或者重新离散;这一功能使得在使用非结构网格时AMG尤其有用。优点在于:不必创建或者存储粗糙网格,而且不需要在粗糙层面估计任何流量或者源项。这一方法和FAS (有时被称为几何学的)多重网格方法形成鲜明对比,在FAS方法中需要网格的层级,和每一个层面的离散方程。。从理论上说,FAS优于AMG方法的地方在于,对于非线性问题前者可以做得更好,这是因为系统的非线性可以通过重新离散传到粗糙层面;当使用AMG时,一旦系统被线化,直到细化层面算子被更新,解算器才会“感觉到”非线性。
AMG限制和延拓算子
这里所使用的限制和延拓算子是基于附加校正(AC)策略的,结构网格的附加校正是由Hutchinson和Raithby[68]描述的。层间的传递是用分段线性插值和延拓完成的。任何粗化层面单元的缺点由细化层面包含的缺点的总和给出的,通过注入粗化层面的值来获取细化层面的校正。在这种方式中,延拓算子由限制算子的转置给出。
限制算子是用细化层面单元粗化或者分组到粗化层面单元来定义。在这个过程中,每一个细化层面的单元用一个或者更多的最强相邻网格来分组,对当前未分组的相邻网格优先选择。这个算法尝试将单元集中到成组的固定尺寸中,一般是二或者四组,但是可以指定任何数。在分组的前后关系中,最强是指当前单元i的相邻单元j的系数A_ij是最大的。对于成组的耦合方程,A_ij是一个块矩阵,它的大小的度量被简单的看成第一个单元的大小。除此之外,对于给定单元的一组耦合方程被一起处理,而且在不同的粗糙单元中不再划分。这样,就对系统中的每一个方程产生了相同的粗化。
AMG粗糙层面算子
使用Galerkin方法来建构粗化层面算子A^H。在这里,当转移到粗化层面时,与细化层面解有关的缺陷必须消除。因此我们可以写出:
在双重时间步中的方程5和限制和延拓中的方程3的基础上替换d^new 和f^new可以得到:
现在重新整理和使用双重时间步一节中的方程5有:
比较限制和延拓一节中的方程1和方程5可以得到下面的粗化层面算子的表达式:
因此,粗化层面算子的建构减少到所有精细层面单元的对角线与组中的相应对角线下面的块之和来形成该组粗化单元的对角线块。
F循环
多重网格F循环本质上是V和W循环的结合。
多重网格循环是一个递归程序,该程序通过在当前层面完成单一网格循环来扩展到下一个粗糙网格层面。参阅V和W循环一节中的图1中的V和W循环以及图2中的V和W循环我们可以看出,F循环就是用0-1层循环(每一个图中的第二个图表)程序来替换当前层(描述单一循环)的方块。可以看出,V循环包括:
W循环包括:
F循环则是用W循环之后进行V循环来构成:
正如预料的那样,F循环比V循环需要更多的计算,但是比W循环花费要少一些。但是它的收敛性比V循环要好,大致和W循环的收敛性差不多。对于耦合解算器设置来说,F循环是默认的AMG循环类型。
灵活循环
对于灵活循环,使用粗化网格的计算由下图中所示的逻辑控制多重网格程序来实现。这一逻辑控制可以保证当前层面的网格的残差减小速度足够慢时就会调用较粗糙网格的计算。除此之外,当当前粗糙网格层上的校正迭代解充分收敛而且因此应该转到下一层精细网格时,多重网格控制会做出适当的处理。这两个判断由下图中所示的两个参数a和b来控制。详细请参阅下面内容。注意多重网格程序的逻辑是这样的:在某一方程上的单一全局迭代过程中,可能会重复处理网格层面。例如:对于4层多重网格,其中包括0,1,2,3四层,解给定输运方程的灵活多重网格程序可能会按这样的顺序处理网格层面:0-1-2-3-2-3-2-1-0-1-2-1-0:
Figure 1灵活多重网格的逻辑控制
灵活循环和V,W循环之间的主要区别是:灵活循环会通过残差减小的公差和终止判据的满足情况来确定什么时候,按什么样的频率来处理每一层网格,而V和W循环则明确定义了各个层面之间的转换模式。
当当前层面的误差减小速度不够快时,多重网格程序就会调用下一个网格层面的计算,具体判断公式如下:
其中,R_i是在当前层面经过第i步迭代之后从当前网格层面计算出的残差和的绝对值。上面的公式表明,如果经过i步迭代后迭代解的残差比某一分数b(在0和1之间)大,就需要处理下一个粗化网格层面了。因此b被称为残差减小公差,用于确定什么时候放弃当前网格层面的迭代解,并转移到下一个更粗糙的网格上解校正方程。B的值控制了处理的粗化网格层面的频率。默认值是0.7。如果b的值较大,就会处理较小的频率,反之亦然。
假定残差减小速度足够块,校正方程会在当前网格层收敛,这样就可以转到下一个更精细的网格层面计算。
当校正解的误差减小到该网格层初始误差的某一分数a(在0和1之间)时,当前网格层上的校正方程就可以被认为是充分收敛了。
其中,R_i是在当前层面经过第i步迭代之后从当前网格层面计算出的残差。R_0时在当前全局计算中该网格层面最初获得的残差。参数a被称为终止判据,默认值是0.1。注意,上面的方程也用于多重网格程序中最低(最精细)网格层面的终止计算。因此,继续在每一个网格层面(包括最精细的网格层面)进行松弛,直到满足该方程的判据(或者直到完成最大的迭代步,而所指定的判据还没有实现)。
全近似格式(FAS)多重网格
FLUENT中FAS形成多重网格层级的方法是将精细网格的成组单元结合起来形成粗化网格单元。如下图所示,粗化网格单元是用节点周围的单元聚束形成的。根据网格的拓扑结构的不同,这一方法可以形成不规则网格和不定的表面数。但是,网格层面的简单创建和植入也就对应了简单的延拓和松弛算子。
Figure 1: 节点聚束形成粗化网格单元
有趣的是,虽然粗化网格单元看起来很不规则,但是在单元表面上还是看不到锯齿状的结构。离散只使用单元表面的区域映射,因此将两个不规则外形单元分开的每组锯齿状单元表面,在二维问题中等价于连接锯齿单元终点的一条直线,在三维问题中等价于从不规则但连续的几何形状映射而来的区域。这种优化方法减少了所需的内存和计算时间。
FAS限制和延拓算子
FAS既需要精细网格解f的限制又需要其残差d(误差)的限制。这个用于将解传到下一个较粗糙的网格层面的算子是用全近似格式构成的[22]。也就是说,取精细网格单元上解的体积平均值来得到粗糙单元的解。粗糙网格上的残差等于植入的精细网格单元内的残差之和。
通过简单的设定精细网格校正为与其相关的粗糙网格的值,我们构建了用于将校正传到细化层面的延拓算子P。
粗糙网格校正y^H是从粗糙层面f^H计算得到的解和限制到粗糙层面Rf的初始解之间的差值计算得到的,它产生于粗糙层面并被应用于精细层面解。因此精细层面解的校正变成:
FAS粗化层面算子
FAS粗化网格算子是很简单的,它源于粗糙层面网格上控制方程的再离散。因为在离散一节中出现的离散方程和耦合解算器对组成单元的很多表面不加限制,所以在由不规则外形单元组成的粗糙网格上进行的再离散不存在问题。
当在不规则粗糙网格单元上使用有限体积格式时会丢失一些精度,但是因为多重网格解的精度主要有最精细的网格决定,所以解的精度不会受到粗糙网格离散的影响。
为了保持精细网格解的精度,我们修改了粗糙层面方程来包括源项[72],这就保证了如果精细网格d^h上的残差为零,那么在粗糙网格f^H上计算得到的校正也会是零。因此,粗糙网格方程用公式表示为:
其中,d^H是从当前粗糙网格解f^H中计算得到的粗糙网格残差,d^H(R f)是从限制的精细层面解R f计算得到的粗糙网格残差。起初,这两项是相同的(因为在初始时我们有f^H = R f),删除这两项化简上面的方程可得:
所以当细化网格残差d^h为零时就不会有粗化层面校正了。
解算器使用概述
当你确定了模型和解算器(参阅使用解算器形式一节),你就可以运行解算器了。下面是你可以参阅的一般步骤:
1.
选择离散格式,如果使用分离解算器,你还可以选择压力插值格式(见选择离散格式一节)。
2.
(只用于分离解算器)选择压力速度耦合方法(参阅选择压力速度耦和方法一节)。
3.
设定亚松驰因子(见设定亚松驰因子一节)。
4.
(只用于耦合显式解算器)打开FAS多重网格(见打开FAS网格一节)。
5.
对解算器设定作任何附加的修改,具体的修改可以参阅介绍你所使用的模型的相关章节。
6.
初始化解(参阅初始化解一节)
7.
激活适当的解监视器(参阅监视解的收敛一节)。
8.
开始计算(参阅执行计算一节)。
9.
如果出现收敛性问题,请尝试收敛性和稳定性一节中所讨论的方法。
前三项的默认设定可以满足大多数问题而不必修改。下面将要介绍这三个以及其它解参数可以做什么样的修改,以及什么时候你希望做出修改。上面所列出的步骤应用于所有定常状态的计算。具体可以参阅非定常解程序的时间相关模拟。
选择离散格式
FLUENT允许你选择每一个控制方程的对流项的离散格式(对于粘性项自动使用二阶精度)当使用分离解算器时,在默认的情况下,所有的方程都是使用一阶迎风离散来解对流问题。当使用耦合解算器时,在默认的情况下,流动方程是二阶格式解出的,其它方程是用一阶格式解出的。关于FLUENT中离散格式的完全介绍请见离散一节。
除此之外,当你使用分离解算器时,你可以指定压力插值格式,有关压力插值格式的相关内容可以参阅前面的介绍。
一阶与二阶的比较
当流动和网格成一条线时(如:矩形网格或者六面体网格模拟矩形导管的层流流动),可以使用一阶迎风离散格式。但是,当流动和网格不在一条线上时(即:流动斜穿网格线)一阶对流离散增加了对流离散的误差(数值耗散)。。对于三角形和四面体网格,流动从来就不会和网格成一条线,此时一般要使用二阶离散来获取更高精度的结果。对于四边形或者六面体网格,如果使用二阶离散格式,尤其是对于复杂流动来说,你可以获取更好的结果。
总而言之,一阶离散一般会比二阶离散收敛得好,但是精度要差,尤其是对于三角形或者四面体网格精度更差。关于控制收敛性的信息请参阅收敛性和稳定性一节。
对于大多数情况,你可以在计算的开始使用二阶格式。对于有些情况,你应该以一阶离散开始计算,在进行了初步迭代之后再转到二阶格式。例如,如果你解高马赫数流动问题,初始解科所预期的解相差较大,你就应该先用一阶格式迭代几步然后打开二阶格式继续计算直至收敛。
对于与网格成一条线的简单流动(如:划分为矩形网格或者六面体网格的矩形导管的层流流动),数值耗散自然会很低,所以一般使用一阶格式替代二阶格式而不损失精度。
最后,如果你使用二阶格式遇到收敛性问题,你就应该尝试使用一阶格式。
幂率和QUICK格式
对于旋转和涡流来说,如果你使用四边形或者六面体网格,QUICK离散格式会比二阶格式产生更为精确的结果。但是,一般说来二阶格式已经足够了,QUICK格式也未必会提高多少精度。
我们还可以选择幂率格式,但是一般说来它不会比一阶格式精确多少。
选择压力插值格式
如压力插值格式所述,当使用分离解算器时我们可以采用很多压力插值格式。对于大多数情况,标准格式已经足够了,但是对于特定的某些模型使用其它格式可能会更好:
· 对于具有较大体积力的问题,推荐使用体积力加权格式。
· 对于具有高涡流数,高Rayleigh数自然对流,高速旋转流动,包含多孔介质的流动和高度扭曲区域的流动,使用PRESTO!格式。
注意:PRESTO!只能用于四边形或者六面体网格。
· 对于可压流动推荐使用二阶格式。
· 当其它格式不适用时,使用二阶格式来提高精度(如:对于流过具有非六面体或者非四边形网格的曲面边界的流动)。
用户输入
你可以在解控制面板(如下图)中指定离散格式,如果是使用分离解算器还可以选择压力插值格式。菜单:Solve/Controls/Solution...。
Figure 1: 分离解算器的解控制面板
在离散框下面列出的每一个标量方程(对于分离解算器,有动量、能量、湍流动能、能量等,对于耦合解算器,有湍流动能、湍流耗散速度)中,你可以在下拉列表中选择一阶迎风、二阶迎风、幂率或者QUICK格式。对于耦合解算器的流动方程(其中包括动量和能量),你既可以选择一阶迎风也可以选择二阶迎风格式。上面的面板是分离解算器对应的面板。
如果你使用分离解算器,在离散框中的压力后面的下拉列表中选择压力插值格式。你可以选择标准、线性、二阶、体积力权重或者(只用于四边形或六面体网格)PRESTO!。
如果你要改变离散的设定,但是你要回到FLUENT解算器默认的设定中修改,你可以点击默认按钮。此时所有的设定都改为默认设定,同时默认按钮就会变成重置按钮。要改回设定,你可以点击重置按钮。
选择压力速度耦合方法
在分离解算器中,FLUENT提供了压力速度耦和的三种方法:SIMPLE,SIMPLEC以及PISO。定常状态计算一般使用SIMPLE或者SIMPLEC方法,对于过渡计算推荐使用PISO方法。PISO方法还可以用于高度倾斜网格的定常状态计算和过渡计算。需要注意的是压力速度耦合只用于分离解算器,对于耦合解算器你不可以使用它。
SIMPLE与SIMPLEC比较
在FLUENT中,可以使用标准SIMPLE算法和SIMPLEC(SIMPLE-Consistent)算法,默认是SIMPLE算法,但是对于许多问题如果使用SIMPLEC可能会得到更好的结果,尤其是可以应用增加的亚松驰迭代时,具体介绍如下。
对于相对简单的问题(如:没有附加模型激活的层流流动),其收敛性已经被压力速度耦合所限制,你通常可以用SIMPLEC算法很快得到收敛解。在SIMPLEC中,压力校正亚松驰因子通常设为1.0,它有助于收敛。但是,在有些问题中,将压力校正松弛因子增加到1.0可能会导致不稳定。对于这种情况,你需要使用更为保守的亚松驰或者使用SIMPLE算法。对于包含湍流和/或附加物理模型的复杂流动,只要用压力速度耦合做限制,SIMPLEC会提高收敛性。它通常是一种限制收敛性的附加模拟参数,在这种情况下,SIMPLE和SIMPLEC会给出相似的收敛速度。
PISO
对于所有的过渡流动计算,强烈推荐使用PISO算法邻近校正。它允许你使用大的时间步,而且对于动量和压力都可以使用亚松驰因子1.0。对于定常状态问题,具有邻近校正的PISO并不会比具有较好的亚松驰因子的SIMPLE或SIMPLEC好。
对于具有较大扭曲网格上的定常状态和过渡计算推荐使用PISO倾斜校正。
当你使用PISO邻近校正时,对所有方程都推荐使用亚松驰因子为1.0或者接近1.0。如果你只对高度扭曲的网格使用PISO倾斜校正,请设定动量和压力的亚松驰因子之和为1.0(比如:压力亚松驰因子0.3,动量亚松驰因子0.7)。如果你同时使用PISO的两种校正方法,推荐参阅PISO邻近校正中所用的方法。
用户输入
你可以在解控制面板(如下图)中指定压力速度耦合方法。菜单:Solve/Controls/Solution...。
在离散框下面的压力速度耦合下拉菜单中选择SIMPLE,SIMPLEC或者PISO。
如果你选择PISO,面板就会扩展到PISO参数设定面板。在默认的情况下,倾斜校正和邻近校正都是打开的。如果你只想使用一种校正,你可以关闭适当的选项。默认的迭代数为1,你不必修改这个值。
用户输入
设定亚松驰因子
如亚松驰一节所述,分离解算器使用亚松驰来控制