二维水箱流畅问题1

《二维水箱流畅问题1》由会员分享，可在线阅读，更多相关《二维水箱流畅问题1（10页珍藏版）》请在金锄头文库上搜索。

1、二二维维水箱水箱的流场情况的流场情况 一一、实验目的实验目的 一个其右上角有一平板堵截的二维矩形水箱， 堵截的平板与 Ox 轴成 45 度， 水箱底面有 一个入流孔和一个出流孔，当理想不可压流体以流量 1 从入流孔流入，而从出流孔流出，作 定常流动时，计算水箱内的流场。 设流函数为，水箱底面 AB 为零流线，在其上=0，在其余水箱壁上=1，即满足： 02222 yx=0， 在 AB 上 =1， 在其余箱壁上 Figure 1.1 水箱示意图 二、二、实验实验要求要求 1. 使用不同方法对实验问题进行求解； 2. 对不同方法求出的数值解进行比较； 三三、实验实验步骤步骤 （一）以 0.250.2

2、5 为尺寸的正方形网格进行网格划分，划分出 2517 个网格，再分开计 算区域（红色）与非计算区域（蓝色） ； Figure 3.1 水箱网格划分图 1 （二）使用 SOR 方法（超松弛迭代法）进行求解 1.差分方程： ,+1= 4(1,+1+ ,1+1+ +1,+ ,+1) + (1 ),其中 k 为迭代次数，为超松弛因子。一般超松弛因子为 1 106)，松弛因子(=1.7000)迭代次数比最佳松弛因子(最佳= 1.7277)迭代次数少；而在高精度时（err 106） ，最佳松弛因子(最佳= 1.7277)迭代次数比松弛因子(=1.7000)迭代次数少。因此，在进行高精度计算时，最佳松弛因子

3、还是正确的。 下图为=1.7277，err=(106)时的二维水箱流场图： Figure 4.1 二维水箱流场图（SOR 方法数值解） （二） ADI 方法（交替方向隐格式） 利用上述 fortran 程序的计算结果， 对不同的时间步长 （） 和允许误差时， 得出计算结果的迭代步数，进行比较；然后利用 MATLAB 绘图，更形象地体现 计算结果。不同的和允许误差对迭代次数的对比数据如下： 3 允许误差允许误差(err) 迭代次数迭代次数 0.01 0.000001(106) 876 0.01 0.00001(105) 635 0.01 0.0001(104) 393 0.01 0.001(10

4、3) 155 0.001 0.0001(104) 1538 0.1 0.0001(104) 64 1 0.0001(104) 26 2 0.0001(104) 39 5 0.0001(104) 74 Table 4.3 不同时间步长、允许误差的迭代次数 从上表可以看出，当=0.01 时，允许误差越小，迭代次数越多，收敛速度 越慢；当允许误差 err=0.0001 时，从 0.001 增大到 5 时，迭代次数先减少后 增多，因此会有一个最佳时间步长能使收敛速度更快。 另外单独取允许误差 err=106和 err=107，观察在取不同值时，迭代步 数的变化情况。 允许误差允许误差(err) 迭代次

5、数迭代次数 0.9 0.000001(106) 48 0.8 0.000001(106) 44 0.7 0.000001(106) 43 0.6 0.000001(106) 44 0.7 0.0000001(107) 52 0.6 0.0000001(107) 52 0.5 0.0000001(107) 52 0.4 0.0000001(107) 51 0.3 0.0000001(107) 50 0.2 0.0000001(107) 72 Table 4.3 不同时间步长在两种允许误差下的迭代次数 从上表可以看出，在不同的精度下，最佳时间步长取值不一样。由于比较难确定最佳时间步长， 因此难以控

6、制差分方程的收敛速度。 所以， SOR 方法在解决 Laplace方程上比 ADI 方法要好。 下图为 =0.1，err=104时的二维水箱流场图： Figure 4.2 二维水箱流场图（ADI 方法数值解） 4 五、五、 实验结论实验结论 求解二维水箱的 Laplace 方程时，SOR 方法（超松弛迭代法）比ADI 方法（交替方向隐格 式）收敛得更快，计算效果更好，能更快地得出计算结果。在 SOR 方法中，最佳松弛因子对 于高精度计算会有更快的收敛速度；而在 ADI 方法中，最佳时间步长是与计算精度有关，不 同的计算精度，对应的最佳时间步长也不一样。 附件附件： （1） SOR.f90 文件

7、 !使用说明：更改参数 omega(松弛因子)和 err(允许误差) !运行成功后在程序目录下可以找到输出文件”output.txt” program main implicit none REAL,PARAMETER:XMIN=0.0 REAL,PARAMETER:XMAX=6.0 REAL,PARAMETER:YMIN=0.0 REAL,PARAMETER:YMAX=4.0 REAL,PARAMETER:DX=0.25 REAL,PARAMETER:DY=0.25 INTEGER,PARAMETER:XN=(XMAX-XMIN)/DX+1 INTEGER,PARAMETER:YN=(YMA

8、X-YMIN)/DY+1 REAL:U(XN,YN),U0(XN,YN) INTEGER:I,J,N REAL,PARAMETER:err=0.000001 REAL,PARAMETER:omega=1.7269 !计算域矩阵 DO i=1,XN DO J=1,YN IF(XMIN+DX*real(I)-DX)+(YMIN+DY*REAL(J)-DY)=5.25) THEN 5 U(I,J)=1 ELSE U(I,J)=0 END IF else !其他初值全部为1 U(I,J)=1 end if END DO end do N=0 loop1: do while(.true.) N=N+1

9、U0=U do i=2,xn-1 loop3: do j=2,yn-1 if(CA(I,J+1)=0) exit loop3 U(i,j)=(1.0-omega)*U0(i,j)+(omega/4.0)*(U(i-1,j)+U0(i+1,j)+U(i,j-1)+U0(i,j+1) end do loop3 end do loop2:do i=1,xn do j=1,yn if(abs(U(i,j)-U0(i,j)err) then exit loop2 end if end do end do loop2 end do loop1 write(*,*) “一共迭代了“,N,“次“ open(u

10、nit=8,file=“output.txt“) DO I=1,XN DO J=1,YN write(8,*) XMIN+DX*real(I)-DX, YMIN+DY*real(J)-DY,U(i,j) END DO END DO End program （2） ADI.f90 文件 !使用说明：更改参数 dt(时间步长)和 err(允许误差) !运行成功后在程序目录下可以找到输出文件”output.txt” program main 6 implicit none real:a,b,c REAL(kind=8),PARAMETER: DT=0.01 REAL(kind=8),PARAMETE

11、R:XMIN=0.0 REAL(kind=8),PARAMETER:XMAX=6.0 REAL(kind=8),PARAMETER:YMIN=0.0 REAL(kind=8),PARAMETER:YMAX=4.0 REAL(kind=8),PARAMETER:DX=0.25 REAL(kind=8),PARAMETER:DY=0.25 REAL(kind=8),PARAMETER:S=DT/(2*DX*2) REAL(kind=8),PARAMETER:P=DT/(2*DY*2) INTEGER,PARAMETER:XN=(XMAX-XMIN)/DX+1 INTEGER,PARAMETER:Y

12、N=(YMAX-YMIN)/DY+1 INTEGER:FLAG integer:CA(XN,YN) real(kind=8):U(XN,YN),U0(XN,YN),U1(XN,YN),W(XN,YN),G(XN,YN),F(XN,YN) integer:I,J,N real(kind=8),parameter:err=0.0001 !计算域矩阵 DO i=1,XN DO J=1,YN IF(XMIN+DX*real(I)-DX)+(YMIN+DY*REAL(J)-DY)=8.0 .or. i=1 .or. i=xn .or. j=1 .or. j=yn) then CA(i,j)=0 else

13、 CA(i,j)=1 end if END DO END DO !初值条件 U=1 !下边界条件 DO I=1,XN if(XMIN+DX*real(I)-DX1.5 .and. XMIN+DX*real(I)-DXerr)then goto 100 end if end do 9 end do write(*,“(A6,F8.6,A12,I6,A2)“) “精度为“,err,“,一共迭代了“,n,“次“ open(unit=8,file=“output.txt“) DO I=1,XN DO J=1,YN write(8,*) XMIN+DX*real(I)-DX, YMIN+DY*real(J)-DY,U(i,j) END DO END DO End program （3） Draw.m 文件 clear all;clc; num=importdata(output.txt); x=num(:,1); y=num(:,2); z=num(:,3); for i=1:25 for j=1:17 z1(i,j)=z(i-1)*17+j); end end contourf(y,x,z1)