飞行器流动仿真讲稿第7章-拟一维喷管流动的数值模拟.ppt

程序编写举例 U1m(n)=U1(n) U2m(n)=U2(n) U3m(n)=U3(n) pm(n)=p(n) C dU1dtm=-(U2m(i)-U2m(i-1))/deltx cJ2=g7*pm(i)*(a(i)-a(i-1))/deltx dU2dtm=-(g5*(U2m(i)**2/U1m(i)-U2m(i-1)**2/U1m(i-1))+ + g6*(U3m(i)-U3m(i-1)))/deltx+cJ2 dU3dtm=-(g0*(U2m(i)*U3m(i)/U1m(i)-U2m(i-1)*U3m(i-1)/U1m(i-1)) + -g4*(U2m(i)**3/U1m(i)**2-U2m(i-1)**3/U1m(i-1)**2))/deltx C dU1dtav=(dU1dt+dU1dtm)/2.0 ！计算时间导数的算术平均值 dU2dtav=(dU2dt+dU2dtm)/2.0 dU3dtav=(dU3dt+dU3dtm)/2.0 程序编写举例 if (icase.eq.2) then U1(i)=U1(i)+dU1dtav*dt ！计算校正步变量值 U2(i)=U2(i)+dU2dtav*dt U3(i)=U3(i)+dU3dtav*dt else dU1m=U1m(i+1)-2.0*U1m(i)+U1m(i-1) dU2m=U2m(i+1)-2.0*U2m(i)+U2m(i-1) dU3m=U3m(i+1)-2.0*U3m(i)+U3m(i-1) dpm=abs(pm(i+1)-2.0*pm(i)+pm(i-1))/(pm(i+1)+2.0*pm(i)+pm(i-1)) s1m=Cx*dpm*dU1m s2m=Cx*dpm*dU2m s3m=Cx*dpm*dU3m C U1(i)=U1(i)+dU1dtav*dt+s1m ！计算校正步变量值（含人工黏性） U2(i)=U2(i)+dU2dtav*dt+s2m U3(i)=U3(i)+dU3dtav*dt+s3m end if 程序编写举例 r(i)=U1(i)/a(i) u(i)=U2(i)/U1(i) t(i)=g1*U3(i)/U1(i)-g4*u(i)**2 p(i)=r(i)*t(i) am(i)=u(i)/sqrt(t(i)) fm(i)=r(i)*a(i)*u(i) if (i.eq.it) write(12,22) j,r(i),t(i),p(i),am(i) 400 continue C进口边界条件i=1 r(1)=1.0 ！物理边界条件 t(1)=1.0 U1(1)=r(1)*a(1) U2(1)=2.0*U2(2)-U2(3) ！一阶外推 u(1)=U2(1)/U1(1) U3(1)=U1(1)*(g2*t(1)+g3*u(1)**2) p(1)=r(1)*t(1) am(1)=u(1)/sqrt(t(1)) fm(1)=r(1)*a(1)*