SUBROUTINE ABQMAIN
C
C     PROGRAM GRIDFILE_3D
C
C This program reads an ASCII formatted gridded wave data file and
C writes the corresponding unformatted binary gridded wave file
C required by ABAQUS.
C
C This program is valid for 3D gridded wave data only.
C
C Description of variables:
C
C First Record:
C
C  NCOMP  = number of wave components
C  DTG    = time increment at which the data are given
C  NWGX   = number of grid points in the grid's x-direction
C  NWGY   = number of grid points in the grid's y-direction 
C           (if NWGY=1, ABAQUS assumes data are constant in the 
C           local y-direction)
C  NWGZ   = number of grid points in the grid's z-direction 
C           (if NWGZ=0 or 1, the analysis is two-dimensional and the 
C           y-direction is vertical)
C  IPDYN  = flag indicated whether or not dynamic pressure 
C           information is stored (IPDYN=1) or not stored (IPDYN=0) 
C           on the gridded wave file.
C
C Second Record:  This record is not used by ABAQUS.  It is provided  
C                 for information in user subroutine UEL via the  
C                 interface GETWAVE. The meaning of arrays AMP and WXL 
C                 is left to the user, but PHI is converted to radians.
C
C  AMP    = array containing "wave amplitudes" of each wave component
C  WXL    = array containing "wavelengths" of these components
C  PHI    = array containing "phase angles" (in degrees) of these 
C           components
C
C Third record:
C
C  WGX    = array containing local x-coordinates of the grid points
C  WGY    = array containing local y-coordinates of the grid points
C  WGZ    = array containing local z-coordinates of the grid points
C
C Remaining records (one for each time at which the wave data is 
C given):
C
C  WGVX   = array containing local x-components of the wave particle 
C           velocity 
C  WGVY   = array containing local x-components of the wave particle  
C           velocity 
C  WGVZ   = array containing local x-components of the wave particle 
C           velocity 
C  WGAX   = array containing local x-components of the wave particle
C           acceleration
C  WGAY   = array containing local x-components of the wave particle
C           acceleration
C  WGAZ   = array containing local x-components of the wave particle
C           acceleration
C  WZCRST = array containing wave surface elevation
C  NCRST  = array containing the index for the vertical grid level 
C           just above the instantaneous water surface
C
C
      PARAMETER(MCOMP=10,MDIMX=10,MDIMY=100,MDIMZ=100,MX=10,MY=10)
C
      DIMENSION AMP(MCOMP),WXL(MCOMP),PHI(MCOMP),
     $     WGX(MDIMX),WGY(MDIMY),WGZ(MDIMZ),
     $     WGVX(MDIMX,MDIMY,MDIMZ),WGVY(MDIMX,MDIMY,MDIMZ),
     $     WGVZ(MDIMX,MDIMY,MDIMZ),WGAX(MDIMX,MDIMY,MDIMZ),
     $     WGAY(MDIMX,MDIMY,MDIMZ),WGAZ(MDIMX,MDIMY,MDIMZ),
     $     P(MDIMX,MDIMY,MDIMZ),DPDZ(MDIMX,MDIMY,MDIMZ),
     $     WZCRST(MX,MY),NCRST(MX,MY)
C
      OPEN(UNIT=9, FILE='gridwave_3d.inp',FORM='FORMATTED',
     $     STATUS='OLD')
      OPEN(UNIT=10,FILE='gridwave_3d.binary',FORM='UNFORMATTED',
     $     STATUS='NEW')

      READ(9,*) NCOMP,DTG,NWGX,NWGY,NWGZ,IPDYN
      WRITE(10) NCOMP,DTG,NWGX,NWGY,NWGZ,IPDYN

      READ(9,*) (AMP(K1),WXL(K1),PHI(K1),K1=1,NCOMP)
      WRITE(10) (AMP(K1),WXL(K1),PHI(K1),K1=1,NCOMP)
      
      READ(9,*) (WGX(K1),K1=1,NWGX), (WGY(K1),K1=1,NWGY),
     $     (WGZ(K1),K1=1,NWGZ)
      WRITE(10) (WGX(K1),K1=1,NWGX), (WGY(K1),K1=1,NWGY),
     $     (WGZ(K1),K1=1,NWGZ)
C
C NWHILE is a dummy flag for DO WHILE loop. 
C
      NWHILE = 1
      NCOUNT = 0
      IF(IPDYN.EQ.0) THEN
         DO WHILE (NWHILE.NE.0)
            READ(9,*,END=100) (((WGVX(K1,K2,K3), WGVY(K1,K2,K3),
     $           WGVZ(K1,K2,K3),WGAX(K1,K2,K3),WGAY(K1,K2,K3),
     $           WGAZ(K1,K2,K3),K3=1,NWGZ),WZCRST(K1,K2),
     $           NCRST(K1,K2),K1=1,NWGX),K2=1,NWGY)
            NCOUNT=NCOUNT+1
            WRITE(10) (((WGVX(K1,K2,K3), WGVY(K1,K2,K3),
     $           WGVZ(K1,K2,K3),WGAX(K1,K2,K3),WGAY(K1,K2,K3),
     $           WGAZ(K1,K2,K3),K3=1,NWGZ),WZCRST(K1,K2),
     $           NCRST(K1,K2),K1=1,NWGX),K2=1,NWGY)
         ENDDO
      ELSE
         DO WHILE (NWHILE.NE.0)
            READ(9,*,END=100) (((WGVX(K1,K2,K3), WGVY(K1,K2,K3),
     $           WGVZ(K1,K2,K3),WGAX(K1,K2,K3),WGAY(K1,K2,K3),
     $           WGAZ(K1,K2,K3),P(K1,K2,K3),DPDZ(K1,K2,K3),
     $           K3=1,NWGZ),WZCRST(K1,K2),NCRST(K1,K2),K1=1,NWGX),
     $           K2=1,NWGY)
            NCOUNT=NCOUNT+1
            WRITE(10) (((WGVX(K1,K2,K3), WGVY(K1,K2,K3),
     $           WGVZ(K1,K2,K3),WGAX(K1,K2,K3),WGAY(K1,K2,K3),
     $           WGAZ(K1,K2,K3),P(K1,K2,K3),DPDZ(K1,K2,K3),
     $           K3=1,NWGZ),WZCRST(K1,K2),NCRST(K1,K2),
     $           K1=1,NWGX),K2=1,NWGY)
         ENDDO
      ENDIF
 100  CONTINUE
      WRITE(*,'(I6," COMPLETE TIME POINTS READ")') NCOUNT

      END