用有限体积方法求解欧拉方程.docx
- 文档编号:11930977
- 上传时间:2023-04-16
- 格式:DOCX
- 页数:14
- 大小:292.84KB
用有限体积方法求解欧拉方程.docx
《用有限体积方法求解欧拉方程.docx》由会员分享,可在线阅读,更多相关《用有限体积方法求解欧拉方程.docx(14页珍藏版)》请在冰豆网上搜索。
用有限体积方法求解欧拉方程
实用标准文案
有限体积法求解二维可压缩Euler方程
——计算流体力学课程大作业
老师:
夏健、刘学强
学生:
徐锡虎
学号:
Q0901*******
日期:
2010年2月5日
文档
一、内容摘要2)(
二、流动控制方程2)(
三、有限体积法的空间离散22222222222222
(2)
四、人工耗散222222222222222222232(
五、时间离散222222222222222222242(
六、边界条件222222222222222222252(
七、计算结果222222222222222222282(
八、结论与展望222222222222222222211()参考文献222222222222222222222121(
一、内容摘要
本文通过运用JAMESON有限体积法求解了二维定常和非定常可压缩Euler方程。
程
序实现语言为C++。
其中,使用的网格是三角形非结构网格。
在时间推进上使用的是四步龙一库塔推进格式。
推进的时间步长取的是当地的时间步长。
为了消除迭代误差、round-of
等误差,本文采用了添加人工耗散项的办法。
另外,本文计算了NACA0012翼型在跨音速
下不同迎角的情况,并与fluent软件的计算结果进行了比较,来验证程序的准确性。
二、流动控制方程
守恒形式的Euler方程:
其中x和y代表笛卡儿坐标系。
W是守恒变量。
U
⑵
W
V
E
F,G表示通量
U
V
U2P
UV
F
G
⑶
UV
V2P
UH
VH
P,H和E表示密度,压强,单元总焓和单元总能量。
U,V表示笛卡儿坐标系下
的速度矢量。
这些量由理想气体的单位体积的总能量和总焓相互联系。
EP
(1)(U2V2)/2
三、有限体积法的空间离散
计算域被划分为互不重叠的单元。
在每个单元运用守恒形式的
Euler方程。
由于每个单
元相对于时间都是不变的,所以等式
(1)可以写成:
(FdyGdx)
S
d
(8)
其中和S是单元的体积和边界。
W是单元的平均值。
在对上述方程进行时间离散前,先对空间进行离散,则方程(
6)可以写为:
(9)
其中k表示第k个单元的体积,Wk是第k个单元的守恒变量。
Qk表示第k个单元
的通量。
方程(7)的右边项可以写成:
kedges
Qk(FyGx)i
i1
其中
XiXbXa,yiyya
(10)
(11)
(8)式中的求和是对第
k个单元的所有边进行的。
守恒参数的量是单元中心值,在
求通量时,第I条边的守恒参数值是用左右单元的平均来表示的:
Wi(WkWp)/2
引入变量:
(12)
ZiUiyiViXi
(13)
则第k单元的Euler方程可以写为:
Z
dWk丄kedgesZUPy
dtki1ZVPx
在本文中,采用的是JAMENSON有限体积法,为了减少存储的相关信息的量,其存
储的方式选择的是按边存储的方法。
在存储的每条边的信息中,包含了这条边的边号,左右
单元号和边的端点。
在计算通量时采用按边循环的方式:
do1=1,nedge
k=connmatrix(l,1)
a=connmatrix(l,2)
b=connmatrix(l,3)
p=connmatrix(l,4)flux=function(k,a,b,p)sum(k)=sum(k)+fluxsum(p)=sum(p)_fluxenddo
这里给出的是FOTRAN语言的形式,我编写采用的是C,具体表现在上交的程序中。
在计算时间步长、人工耗散项等也可用象这样按边循环,从此处我们可以看出求解时与
单元的形状无关。
四、人工耗散
人工粘性模型对方法的成功应用起着关键作用,人工粘性抑制解在激波附近的振荡,又
阻尼解在光滑区域内的高阶误差,对解的线性稳定和收敛于定态是很重要的。
本文在方程
(14)的右端加入了人工耗散项,如对于单元k,其表达式可以表示为:
在有限体积法中,耗散项的公式可以表示为:
kedgeskedges
(2)(4)
Dkdidi
i1i1
其中:
di⑵
ii⑵Wp
Wki
(17)
di⑷
ii⑷2Wp2Wki
其中的
I表示单元k和p的公共边,2定义为:
kedges
2Wk
WjW”
i1
(18)
上面的
j表示与
k相邻的单元。
kedges
(PpPk)i
i1
(19)
k
kedges
(PpPk)i
i1
(2)
i
k
(2)max(p
k)i
(20)
(4)
i
max(O,k⑷
i⑵)
其中的量k
(2),k(4)的范围是:
1256k⑷132,12k⑺1.0。
在计算时发现上面方法得到的人工耗散项并不太适合。
其在光滑区域耗散项太大,而
在大剃度区域又显得太小,为了弥补上面的不足,作下面的修改:
(21)
PpPk
i
PpPk
自适应系数为:
尺度系数为:
(23)
UyVxc、:
x2y2
其中的U,V表示边上的值,C表示当地声速。
五、时间离散
方程最后的稳定解是通过时间上的迭代得到的,可以写为:
为了减小计算时间,人工耗散项的计算只在第一步进行,在下面几步的迭代中保持不
变。
运用上面的方法计算,可以发现CFL数可以取到2、2,本文中使用的是2.0。
使用显示格式迭代的主要缺点是由于稳定区域的限制,所以不能使用过大的时间步长。
可以用近似的方法估算时间步长,对于任意形状的网格,可以使用下面的方法:
六、边界条件
1固面边界条件
对与无粘流动,固面边界条件无穿透条件,设其法向的速度通量为零,即乙0。
由
于压强项的影响,x—向和y—向的动量通量并不为零。
固面的压强近似的取为其相邻的单
元的单元中心压强。
广泛的数值研究证明,如果贴近壁面的单元足够小,并且人工耗散项运
用正确,用这种方法取得的压强对结果的精度不会产生太大的影响。
2远场边界条件
本文提到的远场,实际是人为的有界边界,对于流场中的扰动会传到很远的地方。
因
而对于远场边界条件,情况比较复杂,它不能直接给定具体的流场值,需要与流场内的值来
共同确定远边场的流场值。
如果边界取得过小,则通常采用环量修正。
一般情况下,我们采
用无反射边界条件。
方法,来建立无反射的远场边界条件。
一维均熵流动的Euler方程可写成:
(30)
(31)
tUUx0
2a
UtUUxx0
这里,动量方程中除去了压强项,将上式写成矩阵形式为:
这里:
(32)
(33)
A的特征值为:
因此,上式的两族特征线为:
线是常数,可以用流场内部向外插值计算:
况:
(45)
B、亚音速出流条件Un0,它有一条入流特征线,需规定一个条件
2a
2a
(46)
qn
1
qn
1
2a
2ae
(47)
qn
1
qne
1
C、超音速出流条件Un0,无入流特征线,不需在边界上规定边界条件
(50)
WWi
其中的W表示边上的守恒变量,W|表示与此边相邻元素的守恒变量值。
D、超音速入流条件Un0,它有四条入流特征线,需规定全部四个条件
WW
(51)
其中的W表示边上的守恒变量,W表示来流值。
七、计算结果
本文计算了三个算例,一个是攻角0,马赫数0.80的情况,二是攻角1.25,马赫数
0.80的情况,三是攻角2.5°,马赫数1.5的情况
翼型网格示意图
1.=0,Ma0.80
表面压强系数分布
iA
等压线
压力云图
等马赫线
马赫数云图
升力系数
CL=-3.27408E-005
阻力系数
CD=1.14998670E-002
2.1.25,Ma0.80
jn
■20-
I5r
-IJ]-
flD-
o-i«
上下表面压强系数分布图
等压线图
马赫数云图
等马赫图
压力云图
升力系数为
CL=2.741130755E-001
阻力系数为
CD=2.15905376E-002
5
a
M
5
2
DI
fl0
1I3
lllil鱼贰耳A
n?
上下表面压强系数分布
压力云图马赫数云图
升力系数为CL=0.1352582564
八、结论与展望
本文主要介绍了求解二维非结构网格欧拉方程定常解的问题。
采用的是有限体积空间离散方法,单元形状可以是任意多边形。
使用显式多步推进过程,求解非定常方程得到定常解。
一些标准的加速技术可以加快收敛速度。
Jameson格式是学习CFD编程的入门格式,其精度和收敛速度都不高,但通过实际编写程序可以学到很多方法和技巧。
因此,学习Jameson格式为以
后编写其他格式如Roe,AUSM等等以及求解三维方程和求解N-S方程都打下了良好的基础。
本次编程采用的是C语言,并且存在着一个比较大的问题,CFL数不能取到大于1.2的数值,在超音速计算时候CFL数能取的大值变的更小,估计是时间步长上的问题但是目前还没有找到问题所在。
忘老师指点。
参考文献
1、L.STOLCIS*andL.J.JOHNSTON*,Solution
oftheEulerequationson
unstructuredgrids
for
two-dimensional
compressible
flow,
Aeronautics/Aerospace
Department
VonKarman
InstituteforFluid
Dynamics
Rhode-St-Genese,Belgium.
2、J.Blazek,ComputationalFluidDynamicsPrinciplesandApplications,ELSEVIER.
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 有限 体积 方法 求解 方程
![提示](https://static.bdocx.com/images/bang_tan.gif)