搜索FLUENT不收敛案例及解决方法：压力出口导致的不收敛【转发】
2017-06-25 by:CAE仿真在线 来源:互联网
今天分享一个FLUENT的不收敛案例及其解决方法。计算的对象是一个新型的涡扇发动机加力燃烧室(图1)。在这种新型加力燃烧室中,火焰稳定器被整合到整流支板上,因此整流支板和整流锥都需要冷却。在整流支板和整流锥上开了很多小孔,冷却气从这些孔渗出,形成冷却气膜。
图1 加力燃烧室
这个算例模拟的是实验的工况。实验中没有在加力燃烧室内燃烧,而只是在“主流入口”处引入高温气体,在“冷却气入口”处引入冷却气,以检验气膜冷却的效能。
整流支板共有15块,为了减小计算量,只计算其中的一块(360°/15=24°)。主流入口和冷却气入口都采用“mass-flow-inlet”条件,其中主流入口的流量是0.8kg/s,总温是1241.3K,冷却气入口的流量是0.024kg/s,总温是490.3K。出口采用“pressure-outlet”条件,反压是101325Pa(绝对压力)。由于形状比较复杂,特别是其中有很多小孔,所以采用非结构网格,四面体单元。流体的状态方程采用理想气体(ideal-gas)模型,湍流模型采用Realizable k-ε模型。采用基于压力的求解器。
采用定常算法计算不收敛(图2;这里我们使用FLUENT默认的收敛条件,即能量方程的残差降低到1e-6以下,其余方程降低到1e-3以下)。考虑到可能是分离流诱发的非定常效应导致不收敛(见公众号先前的文章“为何我这个流动总是算不收敛?我要砸电脑!”),我们尝试使用非定常算法。但是不幸的是非定常算法仍然不能在每个时间步内收敛。非定常计算的典型残差曲线如图3所示。
图2 定常计算的残差曲线
图3 非定常计算仍然不收敛。此图是时间步长设为3×10-6秒时的结果。图中显示了5个时间步的残差曲线。
为了弄清不收敛的原因,我们用MATLAB编写两个小程序。第一个程序用来产生一个命令文件j1.jou,让FLUENT迭代30次,并把每次迭代后的流场都保存到文件里面:
clear fid=fopen('j1.jou','wt'); % 打开命令文件 n=30; fprintf(fid,'solve update-physical-time\n'); % 下一个时间步 for i=1:n fprintf(fid,'solve iterate 1\n'); % 迭代一次 fprintf(fid,'file interpolate write-data "d:\\case1-%d" yes yes \n',i); % 将计算结果保存到文件 end fclose(fid); % 关闭命令文件在FLUENT中打开非定常计算的case和data,在菜单栏选择[File]->[Read]->[Journal…],选取命令文件j1.jou,FLUENT便会更新到下一个时间步并迭代30次,并在D盘根目录下产生30个计算结果文件(图4)。
图4 计算结果文件
第二个程序分析这些计算结果,找出压力波动最剧烈的点:
clear n=30; m=375045; % 网格数量 n_points=15; fid1=fopen('j2.jou','wt'); % 打开命令文件 L=1; for i=27:n filename=sprintf('d:\\case1-%d.ip',i); fid2=fopen(filename); % 打开计算结果文件 % ASCII码中,40代表左圆括号 % FLUENT的计算结果文件中,每一组数据都是以左圆括号开头的 % 因此,可以通过查找左圆括号的方法找到每一组数据的起点 % 需要了解更多关于FLUENT计算结果文件格式的内容,请 % 阅读User's Guide中“Format of the Interpolation File”这一节 % until函数需要自己编写,见注1 until(fid2,40); arrx=fread(fid2,[m,1],'double'); % 读取x坐标 until(fid2,40); arry=fread(fid2,[m,1],'double'); % 读取y坐标 until(fid2,40); arrz=fread(fid2,[m,1],'double'); % 读取z坐标 until(fid2,40); arrp=fread(fid2,[m,1],'double'); % 读取压力场 fclose(fid2); if i>27 % 只关注最后三次迭代 % 比较相邻两次迭代的结果,找出压力波动最大的15个点 [sorted,ix]=sort(abs(arrp_old-arrp),'descend'); x_maxchange=arrx(ix(1:n_points)); y_maxchange=arry(ix(1:n_points)); z_maxchange=arrz(ix(1:n_points)); for j=1:n_points % 让FLUENT标出这些点 fprintf(fid1,'surface point-surf point_%d %e %e %e\n',L,x_maxchange(j),y_maxchange(j),z_maxchange(j)); L=L+1; end end arrp_old=arrp; end fclose(fid1); % 关闭命令文件运行这个程序后将生成一个FLUENT命令文件j2.jou,在FLUENT中执行它,便将最后三次迭代中压力波动最剧烈的一些点标记了出来(每次迭代标记15个,共45个点)。
在FLUENT菜单栏选择[Display]->[Mesh…],将这些点显示出来,可以发现,压力波动最剧烈的点都位于出口截面上(图5)。因此推测可能是出口边界设置不当导致不收敛。
图5 压力波动最剧烈的点
尝试对出口的边界条件进行修改,发现当使用“无反射”选项(Non-Reflecting Boundary,图6)的时候,不收敛的问题就得以解决(图7)。
图6 “无反射”选项。
由于参考压力设为101325Pa,所以表压(Gauge Pressure)是0。
图7 非定常计算的残差曲线。时间步长设为3×10-6秒。出口边界使用“无反射”选项。图中显示了约25个时间步的残差曲线。可以看出,在每一个时间步内,残差都能降低到默认的收敛标准以下。
究其原因,在这个算例中出口边界已经达到了壅塞状态,这可以从马赫数云图上看出来(图8)。从图中可以看出,出口附近有马赫数超过1的局部超音速区域。马赫波在出口边界的反射导致了出口截面的参数不断振荡,不能收敛。这种反射是边界条件的数学处理造成的——因为我们强制地让出口截面的压力等于指定的数值,而这是不符合物理事实的。“无反射”的具体处理方法涉及特征线理论,这里不予以叙述,有兴趣的读者可以看计算流体力学原理方面的资料(例如[1];以及FLUENT的User's Guide中的“General Non-Reflecting Boundary Conditions”这一节)。
图8 马赫数云图
如果计算域的出口没有局部超音速区域,就不必将出口边界设置为无反射的了。计算实践表明,这时不将出口设成无反射的也是可以收敛的。
作者非常感谢北京航空航天大学能源与动力工程学院的研究生王昌胜。他提供了本文算例的case文件。此外,北航航空科学与工程学院的研究生李亮阅读了本文的初稿并提出了宝贵的修改意见。
未经许可,不得转载
长按二维码关注流体那些事儿
参考文献
[1] 吴子牛. 计算流体力学基本原理. 北京: 科学出版社, 2001
注1:函数M文件until.m的代码如下,这个函数的作用是在文件中从当前位置向后查找指定的字符:
function until(fid,ch) while true a=fread(fid,1,'char*1'); if a==ch break end end
转自微信公众号:流体那些事儿 叶汉玉
相关标签搜索：FLUENT不收敛案例及解决方法：压力出口导致的不收敛【转发】 Fluent培训 Fluent流体培训 Fluent软件培训 fluent技术教程 fluent在线视频教程 fluent资料下载 fluent分析理论 fluent化学反应 fluent软件下载 UDF编程代做 Fluent、CFX流体分析 HFSS电磁分析
