二維水箱的流場情況實驗報告

中山大學本科生實驗報告書二維水箱的流場情況院系工學院應用力學與工程系專業(yè)班級理論與應用力學10級實驗課程計算流體力學姓名程彬學號10332054實驗地點中山大學東校區(qū)實驗中心實驗時間2012年5月7日指導教師詹杰民一、問題描述：一個其右上角有一平板堵截的二維矩形水箱，堵截的平板與Ox軸成45度，水箱底面有一個入流孔和一個出流孔，當理想不可壓流體以流量1從入流孔流入，而從出流孔流出，作定常流動時，計算水箱內(nèi)的流場。設流函數(shù)為，水箱底面AB為零流線，在其上=0，在其余水箱壁上=1，即滿足：=0，在AB上=1，在其余箱壁上二、實驗方法:取正方形網(wǎng)格，網(wǎng)格格子的邊長與出入孔尺寸相等，即h=0.25。這里為簡單起見,我們?nèi)∷涑叽缂俺鋈肟孜恢?，使得邊壁正好與垂直網(wǎng)格線重合，上下壁正好與水平網(wǎng)格線重合，出入孔正好與一個網(wǎng)格格子的邊重合。由于網(wǎng)格格子是正方形的，而水箱右上角堵截平板與X軸的夾角為45°，因而堵截平板與網(wǎng)格線的交點也完全落在網(wǎng)格結點上，見下圖。所以說差分方程求解的區(qū)域邊界也就是原微分方程求解區(qū)域的邊界。網(wǎng)格線在Ox軸方向共25條(包括兩邊壁)即m=25，在Oy方向共17條即n=17。采用超松弛迭代法，因此內(nèi)點差分格式為：i=2,3,…,24j=2,3,…,9j=10時，i=2，3，…，23j=11時，i=2，3，…，22j=12時，i=2，3，…，21j=13時，i=2，3，…，20j=14時，i=2，3，…，19j=15時，i=2，3，…，18j=16時，i=2，3，…，17松弛因子ω按（2）式取最佳值其中=0.9875故ω最佳=1.7269邊界點值ψi,1=0i=8，9，…，21ψi,1=1i=1，2，…，7，22，23，24，25ψ1,j=1j=1，2，…，17ψ25,j=1j=1，2，…，9ψi,17=1i=1，2，…，17ψ18,16=ψ19,15=ψ20,14=ψ21,13=ψ22,12=ψ23,11=ψ24,10=1于是就可以通過邊界上節(jié)點的值，由內(nèi)點差分格式順序掃描，一次一次迭代，直至達到要求的精度為止，最后求得所有內(nèi)點的流函數(shù)值ψi，j。三、實驗程序！-----------------------主程序----------------！i=網(wǎng)格橫坐標，j=網(wǎng)格縱坐標，k=迭代次數(shù)變量，st=收斂的迭代次數(shù)！time=迭代次數(shù)的字符變量！w=最佳松弛因子，d=迭代精度！F=速度勢，dF=速度勢增量!flag1=循環(huán)判斷因子！---------------------定義變量----------------programwaterboximplicitnoneintegeri,j,k,stlogicialflag1character(len=2)::timereal(8)d！----------------------對變量賦初值----------------real(8)F(1:25,1:17,0:1000)，dF(1:25,1:17,0:1000)w=1.7296F(1:25,1:17,0)=1！對全體變量賦初值F(1:25,1,0:1000)=1 ！對水槽下壁賦初值F(1:17,17,0:1000)=1 ！對水槽上壁賦初值F(1,1:17,0:1000)=1 ！對水槽左側壁賦初值F(25,1:9,0:1000)=1 ！對水槽右側壁賦初值F(8:21,1,0:1000)=0 ！對AB邊賦初值doi=17,25！對水槽斜壁賦初值doj=34-i,34-iF(i,j,0:1000)=1enddoenddoflag1=.true.st=0k=0！----------------------迭代計算---------------------print*,'請輸入迭代精度：'read*,stdowhile(flag1.eqv..true.)doi=2,17doj=2,16dF(i,j,k)=-w*F(i,j,k)+.25*w*(F(i-1,j,k+1)+F(i+1,j,k)&+F(i,j-1,k+1)+F(i,j+1,k))F(i,j,k+1)=dF(i,j,k)+F(i,j,k)if(abs(dF(i,j,k))>d)then！判斷是否達到計算精度flag1=.false.endifenddoenddodoi=18,24doj=2,33-idF(i,j,k)=-w*F(i,j,k)+.25*w*(F(i-1,j,k+1)+F(i+1,j,k)&+F(i,j-1,k+1)+F(i,j+1,k))F(i,j,k+1)=dF(i,j,k)+F(i,j,k)if(abs(dF(i,j,k))>d)thenflag1=.false.endifenddoenddoflag1=.not.flag1！若每個網(wǎng)格點都達到計算精度，返回false，退出循環(huán)k=k+1st=st+1enddo！-------------------------輸出txt格式文件----------------open(10,file='test0506.txt')dok=0,stdoj=1,9write(10,100)(F(i,j,k),i=1,25)！輸出1-9行網(wǎng)格數(shù)據(jù)enddodoj=10,17write(10,100)(F(i,j,k),i=1,34-j)！輸出10-17行網(wǎng)格數(shù)據(jù)enddowrite(10,*)'','k=',k！顯示結果的迭代次數(shù)enddo100format(1x,25G10.2)！設置輸出格式：每25個數(shù)據(jù)換行！---------------------輸出dat格式文件---------------------dok=0,stwrite(time,'(i2)')kopen(unit=20,file='test0506'//trim(time)//'.dat')write(20,*)'variables=','"x"',',','"y"',',','"z"'doj=1,9doi=1,25write(20,*)i,j,F(i,j,k)！輸出1-9行網(wǎng)格數(shù)據(jù)enddoenddodoj=10,17doi=1,34-jwrite(20,*)i,j,F(i,j,k)！輸出10-17行網(wǎng)格數(shù)據(jù)enddoenddoenddo！-----------------