近些年来,随着计算机性能的提高,计算流体力学(computational fluid dynamics,简称CFD)在很多工作中得到广泛运用。在实际的工程应用中,很多问题涉及到非定常流动,并且其计算量远大于定常流动,比如,气动弹性动力学问题、涡轮叶片绕流问题、分离涡脱落问题、动导数的求解等等。由此而促进了各种加速收敛技术和时间推进格式的发展。
为了提高计算效率,前人发展了很多加速收敛措施,如当地时间步长、隐式残值光顺方法和多重网格方法[1]。隐式时间推进方法可以在保证稳定性的同时大幅提高时间步长,进而提高计算效率,常用的有近似因子分解(AF)[2]、超松弛迭代方法(SSOR)[3]、直接求解方法和Newton-Krylov方法[4]。在近似因子分解方法中应用最广的是LUSGS算法[5]。后来Chen和Wang[6]对LUSGS进行了改进,在有限体积方法中发展了BLUSGS算法。Newton-Krylov方法中最受欢迎的是Newton-GMRES[7],其中无矩阵存储的GMRES方法[8]在CFD中应用更为广泛。
非定常流场的计算方法可分为2大类:双时间步方法(dual-time stepping)和物理时间迭代方法(physical time subiteration)。双时间步法引入了伪时间项,对于每一个物理时刻看做定常流动来求解。双时间步法的好处是原来用于定常计算的预处理、当地时间步、对角化、多重网格都可以应用到计算中;真实时间步长可以取得很大;双时间法可以降低通量线性化误差和近似因子分解误差,放宽了稳定性限制[9]。这些方法都能有效的提高收敛速度。
流场求解效率不仅仅与收敛速度有关,还与初值与收敛解间的差距有关。因此,在基于双时间方法的非定常流场求解中,为了提高计算效率,通常将上一个时刻的流场作为下一个时刻流场求解的初值。至于从初值角度进一步提高计算效率,相关研究鲜有报道。本文从该研究角度出发,将前2~3个物理时刻的流场信息外插求解出待求时刻的流场作为伪时间步迭代的初场,这样就进一步缩小了初值与收敛解的差距,以期减小迭代步数,提高计算效率。
基于以上的思想,本文通过采用内迭代初值外插,发展了一种提高非定常流场求解效率的有效策略。采用绕圆柱非定常流动的求解算例来验证本方法的计算效果,并研究了不同初值外插格式、空间离散格式、时间步长和收敛标准对初值外插方法效果的影响。
1 数值方法 1.1 双时间步方法为了提高非定常流场的时间计算精度,同时又要求具有较高的计算效率,Jameson 提出了一种双时间步方法,即在冻结的物理时间点上加入类似 Newton 迭代的虚拟时间迭代过程,通过增加内迭代过程提高 LU-SGS 等隐式线性化方法所损失的时间精度。双时间步法思想简单,且在定常流计算程序基础上进行改造的工作量小,得到了广泛应用,其实现形式如下。
非定常流动控制方程离散后的形式为
方程(2)直接求解比较困难,双时间步法采用虚时间迭代技术对方程(2)进行求解,引入一个虚拟时间τ,将控制方程(2)改写为
初值外插方法是通过记录前几个时刻的流场的基本信息进行外插得到下一个物理时刻的内迭代初值。
一般情况下,采用外插时必须要保证插值对象的连续性。若其出现比较强的间断,则外插出现的偏差往往会比较大,这样就会与我们的初衷背道而驰。因此,合理选择外插对象至关重要。相对于原始变量,守恒变量的连续性更好,因此我们将流场的守恒量进行外插。但对于流场本身存在强间断情况,比如有激波出现,此时用初值外插方法容易造成更大偏差,不宜使用。在时间方向上的泰勒级数展开如下:
本文还提出了一阶与三阶交替的外插格式。即在第n个物理时刻,取2m为周期有:
若mod(n,2m)<m,则下一时刻的初值为Qn+1=2Qn-Qn-1
若mod(n,2m)>m-1,则下一时刻的初值为Qn+1=4Qn-6Qn-1+4Qn-2-Qn-3
根据多次计算对比结果,取m=3时效果更好。因此后面的交替外插都采用的是每6个时间点为1个周期,三步一阶外插三步三阶外插的格式。
2 数值算例验证 2.1 算例介绍以非定常静止圆柱绕流为例,对于不同外插格式和空间离散格式下的初值外插方法进行效率对比并分析时间步长和收敛标准对初值外插方法效率的影响。来流马赫数为0.1,雷诺数为100,层流,采用非结构三角形网格,圆柱表面网格节点数为100,网格总数为29 788,计算域网格和圆柱表面局部网格示意图如图 1所示。
2.2 不同因素影响对比在初值外插方法介绍中给出了一阶、二阶、三阶和交替外插4种外插格式。本节针对采用不同的外插格式、空间离散格式对初值外插方法效率的影响进行分析并确定使该方法效率最好的外插格式和空间离散格式,并分析不同的时间步长和收敛标准对初值外插方法效率的影响。
在非定常圆柱绕流算例中,空间离散采用中心格式,取时间步长为1个周期120步,设1个周期的步数为Tn,则有Tn=120。设每个周期平均耗时为Tp。收敛要求为最大残值小于5.0×10-9。具体状态在2.1节中已经交代。对一阶、二阶、三阶和交替外插精度的外插效果进行对比。以下是计算结果:
外插格式 | 每周期迭代步数 | 平均每迭代一步所用时间/s | Tp/s | 效率提高/% |
无 | 8 869 | 0.100 | 884.21 | |
一阶 | 5 553 | 0.107 | 595.48 | 48.49 |
二阶 | 4 471 | 0.110 | 493.23 | 79.27 |
三阶 | 11 183 | 0.112 | 1 256.48 | -29.63 |
交替 | 4 099 | 0.110 | 452.27 | 95.50 |
由图 2可以看出,无论是初始迭代的最大误差还是平均误差,均随着外插精度的提高而降低。而迭代步数并不满足该规律,初始迭代误差最小的三阶精度外插的迭代步数反而是最多的,换而言之,采用三阶精度外插格式会出现收敛困难的现象。对于出现收敛困难现象的原因,在下一节中会具体分析。
从表 1中不难看出一阶精度、二阶精度和交替外插格式均可明显减少每周期迭代步数,但由于采用初值外插方法每一步迭代增加了计算量和存储因而每一步迭代所用的时间会有所增加。比较后发现平均每迭代一步所用的时间增加量很小,大约10%左右。因此整体的效率仍然是提高的。在4种格式中采用三步一阶三步三阶的交替外插格式效率最高。因此,本文后面的计算没有特殊说明均采用三
步一阶三步三阶的交替外插格式。从图 3中可以看出采用初值外插方法前后计算结果完全吻合,证明了初值外插方法的准确性。
采用不同的空间离散格式,分析对初值外插方法的效率是否会有所影响以及影响程度如何。本文对比了中心格式和迎风ROE格式2种空间离散格式,时间步长取1个周期120步,收敛标准为5.0×10-9。以下是对比结果:
从表 2中可以看出采用不同的空间离散格式对于初值外插方法效率会有一定的影响。空间离散采用中心格式初值外插方法的效果更好。对比中心格式和迎风格式不难发现,中心格式利用了较多的网格系统作为模板,因此对于模板的迭代初值更为敏感。当迭代初值更接近收敛解时,采用中心格式比采用迎风格式提高的效率就会更多。
时间步长和收敛标准选取的不同对初值外插方法的效率也会有所影响。本文通过改变时间步长和收敛标准观察初值外插方法的效率变化情况。
不改变其他的参数,仅对时间步长进行改变,分析对该方法效率的影响。以中心格式为例,对于时间步长分别取1个周期60步、120步、200步和400步的情况,收敛标准为最大残值小于1.0×10-8。图 4和表 3显示了初值外插方法的效率变化情况。
Tn | 无外插Tp/s | 外插法Tp/s | 效率提高/% |
60 | 630.92 | 340.16 | 85.48 |
120 | 734.42 | 381.11 | 92.71 |
200 | 830.92 | 424.64 | 95.68 |
400 | 1014.67 | 537.28 | 88.85 |
从以上图表中不难看出,时间步长的改变虽对外插方法的计算效率提升量有一定影响,但影响不大。因此,可以认为在很大的时间步数变化范围内初值外插方法均有非常明显的效果。
下面研究收敛标准对外插方法效率提升的影响。空间离散采用中心格式,固定时间步长为1个周期120步,收敛标准分别取最大误差小于5.0×10-9、1.0×10-8、5.0×10-8和1.0×10-7时观察初值外插方法的效果。
收敛标准 | 无外插Tp/s | 外插法Tp/s | 效率提高/% |
5.0×10-9 | 884.21 | 452.27 | 95.50 |
1.0×10-8 | 728.50 | 354.17 | 105.69 |
5.0×10-8 | 642.26 | 289.72 | 121.68 |
1.0×10-7 | 609.29 | 262.40 | 132.20 |
从以上图表中可以看出在一定范围内随着收敛标准降低,初值外插法的效率会提高。但即使对于高精度的数值求解,初值外插法的效率仍然可以将效率提高将近一倍。因此初值外插法对于不同精度的求解均具有非常明显的效果。
2.3 流场误差传播分析一般情况下,外插精度越高,流场迭代初值与收敛解越接近,初值外插方法计算效率就越高。然而三阶精度外插的计算效率却最低,甚至还没有不用初值外插方法的效果好。产生这种现象的原因归结于误差的2个分量。
不同的初值与真实解的误差可以按照误差各个傅里叶分量的衰减程度不同,分为高频振荡误差和低频光滑误差。
在线性空间内任何初始误差均可表示为
在迭代过程中
很明显,特征向量代表误差分量,其对应的特征值代表误差分量的放大系数。误差每一个分量vi都对应一个误差频率k,当1≤k <N/2时,vk称为低频分量,由于频率低,变化缓慢,看上去比较光滑,又称为光滑分量。当N/2≤k≤N-1时,称为高频分量,由于它们频率高,摆动快,又称为摆动分量。详细定义见文献[12]。
高频振荡误差是局部行为,来源于附近几个网格点之间的相互耦合,与边界距离较远的网格点的信息无关;而低频光滑误差是全局行为,主要来源于边界信息。通过局部松弛后误差呈现光滑性,此时误差主要来自于边界。可以设想二维N×N网格上的点松驰方法,将边界信息传播到所有点至少需O(N)次迭代,因此收敛速度极慢。低频光滑误差量阶很小,但在细网格上很难收敛[13]。
在迭代过程中低频误差的大小与迭代初值和迭代格式有关。采用初值外插方法会对迭代初值进行改变,低频误差的大小也会随之发生变化。在外插过程中,不仅仅是对精确解的插值,由于收敛标准的限制,存在一定的误差,因此外插的对象还包括残余的误差。高精度外插对插值对象的连续性要求较高,而残余误差的变化本身具有随机性与不连续性。低频误差由于在迭代过程中难以收敛在残余误差中占有较高的比例。通过高阶精度外插格式得到的初始迭代误差比低阶精度外插格式小,但是高阶精度外插格式只可以大幅衰减高频误差,低频误差反而被放大,低频误差收敛速度极慢因此高精度外插会出现收敛困难情况。
一阶外插能够保证误差变化的单调性,不会出现低频误差的放大。但是精度较低,高频误差分量衰减幅度小;高阶外插虽然可以大幅度衰减高频误差,但是会出现低频误差的放大。通过采用三步一阶三步三阶的交替外插格式,既可以大幅衰减高频误差分量,又能保证低频误差不会被放大。
3 结 论本文在非定常流场双时间步法求解的基础上,提出了将前2~3个物理时刻的流场信息外插求解出待求时刻的流场作为伪时间步迭代初场的初值外插方法,并通过绕圆柱非定常流动算例验证了本方法的有效性,研究结果显示:
1) 相对统一精度的外插格式,交替精度的外插格式具有更高的鲁棒性和计算效率,可普遍将计算效率提高一倍左右。
2) 不同的空间离散格式对初值外插方法的效果有一定的影响,相对于迎风格式,中心格式采用初值外插方法的效率会更高。
3) 初值外插方法的效率提升量会受时间步长和收敛标准影响,但在较宽的时间步长和收敛标准范围内,效率提升仍十分可观。
[1] | 燕振国. 高精度混合线性紧致格式的隐式时间推进方法研究[D]. 成都: 中国空气动力研究与发展中心, 2013 Yan Zhenguo. Investigation of Implicit Time Integration Methods with HDCS Schemes[D]. Chengdu: China Aerodynamics Research and Development Center Graduate School, 2013 (in Chinese) |
Cited By in Cnki (4) | |
[2] | Pulliam T H, Chaussee D S. A Diagonal Form of an Implicit Approximate-Factorization Algorithm[J]. Journal of Computational Physics, 1981, 39: 347-363 |
Click to display the text | |
[3] | Young D M. A Historical Review of Iterative Methods[C]//The ACM Conference on History of Scientific and Numeric Computation, 1987 |
[4] | Lomax H, Pulliam T H, Zingg D W, et al. Fundamentals of Computational Fluid Dynamics[M]. Springer, Berlin, 2001 |
[5] | Yoon S, Jameson A. Lower-Upper Symmetric Gauss-Seidel Method for the Euler and Navier-Stokes Equations[J]. AIAA Journal, 1988, 26: 1025-1026 |
Click to display the text | |
[6] | Chen R F, Wang Z J. Fast, Block Lower-Upper Symmetric Gauss-Seidel Scheme for Arbitrary Grids[J]. AIAA Journal, 2000, 38: 2238-2245 |
Click to display the text | |
[7] | Michalak C, Ollivier-Gooch C. Globalized Matrix-Explicit Newton-GMRES for the High-Order Accurate Solution of the Euler Equations[J]. Computers & Fluids, 2010, 39: 1156-1167 |
Click to display the text | |
[8] | Knoll D A, Keyes D E. Jacobian-Free Newton-Krylov Methods: a Survey of Approaches and Applications[J]. Journal of Computational Physics, 2004, 193: 357-397 |
Click to display the text | |
[9] | 赵慧勇, 乐嘉陵. 双时间步法的应用分析[J], 计算物理, 2008, 25(3): 253-258 Zhao Huiyong, Le Jialing. Application Analysis on Dual-Time Stepping[J]. Computational Physics, 2008, 25(3): 253-258 (in Chinese) |
[10] | Jameson A. Time Dependent Calculations Using Multigrid with Applications to Unsteady Flows Past Airfoils and Wings[C]. AIAA Paper 91-1596, 1991 |
[11] | 王刚. 复杂流动的网格技术及高效、高精度算法研究[D]. 西安: 西北工业大学, 2005 Wang Gang. New Type of Grid Generation Technique together with the High Efficiency and High Accuracy Scheme Researches for Complex Flow Simulation[D]. Xi'an, Northwestern Polytechnical University, 2005 (in Chinese) |
[12] | 刘超群. 多重网格法及其在计算流体力学中的应用 [M]. 北京: 清华大学出版社, 1995 Liu Chaoqun. Multigrid method with Application to Computational Fluid Dynamics[M]. Beijing, Tsinghua University Press, 1995 (in Chinese) |
[13] | 李晓梅, 莫则尧. 多重网格算法综述[J]. 中国科学基金, 1996, 10(1): 4-11 Li Xiaomei, Mo Zeyao. Viewpoints of Multigrid Algorithms[J]. Bulletin of National Science Foundation of China, 1996, 10(1): 4-11 (in Chinese) |
Cited By in Cnki (30) |