资源描述：
《利用有限体积算法三阶迎风型quick离散格式和人工压缩算法求解方腔流动问题fortran语言》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。
1、语言源程序!————————————————————————————————————!利用有限体积算法三阶迎风型离散格式和!人工压缩算法求解方腔流动问题（语言版本）!————————————————————————————————————programQUICK_cavityparameter(mx=101,my=101)implicitdoubleprecision(a-h,o-z)dimensionu(mx,my+1),v(mx+1,my),p(mx+1,my+1)dimensionun(mx,m
2、y+1),vn(mx+1,my),pn(mx+1,my+1)common/ini/u,v,pc2=2.25re=100.0dt=0.0005dx=1.0/float(mx-1)dy=1.0/float(my-1)!----------------------------------------------------------------------------------------!u、v、p为t时刻值，un、vn、pn为t+1时刻值，!mx、my为最大网格数，c2为虚拟压缩系数的平方，re为雷
3、诺数。!-----------------------------------------------------------------------------------------num=0err=100.00!nun，计数器；err，判断人工压缩法求解收敛的标准callinitial!调入初始条件，以下为人工压缩算法求解dowhile(err.gt.1e-4.and.num.lt.1e6)err=0.0callquick(u,v,p,un,vn,mx,my,dx,dy,dt,re)!QUICK
4、离散格式求解动量方程，得到un、vncallcalp(p,un,vn,pn,mx,my,dx,dy,dt,c2)!求压强pncallcheck(u,v,p,un,vn,pn,mx,my,dx,dy,dt,c2,err)!校验流场信息，判断是否收敛，同时更新u、v、pwrite(*,*)'error=',errnum=num+1write(*,*)num!屏幕跟踪输出enddocalloutput(u,v,p,mx,my,dx,dy)!输出结果文件End!-F.7-!subroutineinitial!
5、初始化流场parameter(mx=101,my=101)doubleprecisionu(mx,my+1),v(mx+1,my),p(mx+1,my+1)common/ini/u,v,pdoi=1,mx+1doj=1,my+1p(i,j)=1.0enddoenddodoi=1,mxdoj=1,my+1u(i,j)=0.0if(j.eq.my+1)u(i,j)=4.0/3.0if(j.eq.my)u(i,j)=2.0/3.0enddoenddodoi=1,mx+1doj=1,myv(i,j)=0.0e
6、nddoenddoendsubroutine!!subroutinequick(u,v,p,un,vn,mx,my,dx,dy,dt,re)!以QUICK格式离散动量方程implicitdoubleprecision(a-h,o-z)dimensionu(mx,my+1),v(mx+1,my),p(mx+1,my+1),un(mx,my+1),vn(mx+1,my)doubleprecisionmiumiu=1.0/re!以下求解x方向速度un-----------------------------
7、-----------------------------------------------doi=3,mx-2doj=3,my-1fw=0.5*(u(i-1,j)+u(i,j))*dyfe=0.5*(u(i,j)+u(i+1,j))*dyfs=0.5*(v(i,j-1)+v(i+1,j-1))*dxfn=0.5*(v(i,j)+v(i+1,j))*dxdf=fe-fw+fn-fsaw=miu+0.750*alfa(fw)*fw+0.125*alfa(fe)*fe+0.375*(1.0-alfa(f
8、w))*fwaww=-0.125*alfa(fw)*fwae=miu-0.375*alfa(fe)*fe-0.750*(1.0-alfa(fe))*fe-0.125*(1.0-alfa(fw))*fwaee=0.125*(1.0-alfa(fe))*fe-F.7-as=miu+0.750*alfa(fs)*fs+0.125*alfa(fn)*fn+0.375*(1.0-alfa(fs))*fsass=-0.125*alfa(fs)*fsan=miu-0.
