subroutine init
C
C Sod shock tube problem (a whimpy test) in 1, 2, or 3 dimensions
C 24jan92 blondin
C=======================================================================
C
C GLOBALS
C
      include 'global.h'
      include 'zone.h'
C
C LOCALS
C
      integer i, j, k
C
C Set up geometry and boundary conditions of grid
C
C Boundary condition flags : nleft, nright
C   = 0  :  reflecting boudary condition
C   = 1  :  inflow/outflow boundary condition
C   = 2  :  fixed inflow boundary condition
C   = 5  :  ``special'' boundary (cook up your own)
C Geometry flag : ngeom                         |  Cartesian:
C   = 0  :  planar                              |    gx = 0, gy = 0, gz = 0
C   = 1  :  cylindrical radial                  |  Cylindrical:
C   = 2  :  spherical   radial             3D= {     gx = 0, gy = 1, gz = 3
C   = 3  :  cylindrical angle                   |
C   = 4  :  spherical polar angle (theta)       |  Spherical:
C   = 5  :  spherical azimu angle (phi)         |    gx = 2, gy = 4, gz = 5
C
      ndim   = 1
      ngeomx = 0
      ngeomy = 0
      ngeomz = 0
      nleftx = 0
      nrightx= 0
      nlefty = 0
      nrighty= 0
      nleftz = 0
      nrightz= 0
C
C=======================================================================
C
C constants that define problem
C
      data pright / 0.1 /
      data dright / 0.125 /
      data pleft  / 1.0 /
      data dleft  / 1.0 /
      cmin = 0.0
      cmax = 1.0
      time = 0.0
      dt   = 1.0
C
C=======================================================================
C
C initial state  - The left hand state of the standard Sod shock
C tube problem is set up in one corner of the grid.
C
      write (8,*) 
      write (8,*) 'Oblique Sod shock tube in ',ndim,' dimensions.'
      write (8,*) 
      write (8,*) 'Adiabatic index, gamma = ', gam
      write (8,*) 'Pressure ratio is ', pright
      write (8,*) 'Density ratio is ', dright
C
C initialize grid:
C
      i0 = imax / 2  
      j0 = jmax / 2
      k0 = kmax / 2
c
      do 1 k = 1, kmax
      do 1 j = 1, jmax
      do 1 i = 1, imax
         zrho(i,j,k) = dright
         zp  (i,j,k) = pright
         zux (i,j,k) = 0.
         zuy (i,j,k) = 0.
         zuz (i,j,k) = 0.
 1    continue
C
      do 2 i = 1, i0
	 js = max(1, j0 - i + 1)
         do 3 j = 1, js
            ks = max(1, k0 - i - j + 2)
            do 4 k = 1, ks
               zrho(i,j,k) = dleft
               zp  (i,j,k) = pleft
 4          continue
 3       continue
 2    continue
C
C########################################################################
C
        return
        end