program stubeic c Computes a uniformly-spaced, planar grid and an initial c solution for the shock-tube problem. Assumptions: c 1) Planar domain, 2) Perfect gas. c c----------------------------------------------------------------------- Implicit none Integer nim Integer njm Parameter ( nim = 200 ) Parameter ( njm = 50 ) Integer iprob Integer i, j Integer m Integer ni, nj Real q (nim,njm,4) Real x (nim,njm) Real y (nim,njm) Real a1, p1, r1, u1 Real a4, p4, r4, u4 Real gam Real xl, xr, xd, yh Real rmach, alpha, reyn, time Real gascon Real aref, eref, pref, rref, tref c----------------------------------------------------------------------- c.....Set gas variables. gam = 1.4 gascon = 1716.0 c.....Set reference conditions (pressure = 10 psi, temperature = 520 K). pref = 10.0 * 144.0 tref = 520.0 rref = pref / ( tref * gascon ) aref = sqrt( gam * gascon * tref ) eref = rref * aref * aref c.....Set problem to 10/1 pressure ratio. p4 = 1.0 * pref r4 = 1.0 * rref u4 = 0.0 p1 = 0.1 * pref r1 = 0.125 * rref u1 = 0.0 xl = 0.0 xr = 1.0 xd = 0.5 yh = 0.1 c.....Compute acoustic velocities. a1 = sqrt( gam * p1 / r1 ) a4 = sqrt( gam * p4 / r4 ) c.....Generate an evenly spaced, planar grid. write(*,*) ' ' write(*,*) 'Enter ni and nj' read (*,*) ni, nj do i = 1, ni do j = 1, nj x(i,j) = xl + ( xr - xl ) * ( i - 1.0 ) / ( ni - 1.0 ) y(i,j) = yh * ( j - 1.0 ) / ( nj - 1.0 ) enddo enddo open ( unit=7, file='stube.x.fmt' ) write ( 7, * ) ni, nj write ( 7, * ) (( x(i,j), i=1,ni), j=1,nj ), & (( y(i,j), i=1,ni), j=1,nj ) c.....Generate initial solution. do i = 1, ni do j = 1, nj q(i,j,2) = 0.0 q(i,j,3) = 0.0 if ( x(i,j) .lt. xd ) then q(i,j,1) = r4 / rref q(i,j,4) = ( p4 / ( gam - 1.0 ) ) / eref else q(i,j,1) = r1 / rref q(i,j,4) = ( p1 / ( gam - 1.0 ) ) / eref endif enddo enddo c.....Output the solution into file st0q.fmt (2d plot3d format) open ( unit=8, file='stubeic.q.fmt' ) rmach = 0.0 alpha = 0.0 reyn = 1.0 time = 0.0 write ( 8, * ) ni, nj write ( 8, * ) rmach, alpha, reyn, time write ( 8, * ) ((( q(i,j,m), i=1,ni), j=1,nj), m=1,4 ) c----------------------------------------------------------------------- stop end