      PROGRAM Shear
      include 'undefined.h'
      INTEGER tmax,tprint,tsave,tmovie
      INTEGER t,i,tgraph
      DOUBLE PRECISION tphys
      INCLUDE 'size.h'
      INCLUDE 'global.f'
      INCLUDE 'pointers.f'
      write(6,*) 
     c "$Id: Shear.f,v 1.2 2001/07/09 15:02:58 zaleski Exp $"
      call initlog()
      CALL machineid()
      open(unit=11)
      write(11,*) 
     c "$Id: Shear.f,v 1.2 2001/07/09 15:02:58 zaleski Exp $"
 
      call read_stdin()
      call initialize(tphys,tmax,tprint)
      call closelog()

      i = ifetch("tsave",tsave)
      i = ifetch("tmovie",tmovie)
      call flush(6)
      do t=t0,tmax
         if(mod(t,tsave).eq.0) then
            call save(u,v,z(icc(ng)),z(iu(ng)),
     $                tau,nx,ny)
	    write (*,*) ' --- Writing dump.dat at t = ',t,' ---'
            call flush(6)
         endif
         call timestep(tphys,t,tprint,tmovie)
      enddo
c
c  Close movie file
c
      close(12)
c
c  Close spurious currents file
c
      close(67) 
      STOP
      end

      SUBROUTINE initialize(tphys,tmax,tprint)
      include 'undefined.h'
      INTEGER tmax,tprint,restart,i,t0
      DOUBLE PRECISION tphys
      INCLUDE 'size.h'
      INCLUDE 'global.f'
      INCLUDE 'pointers.f'
      call readglob(tmax,tprint)
           write(6,*) "tauuuu",tau	
      restart=0
      i=ifetch("restart",restart)
      write(6,*) "restart=",restart
      call initpvel(nx,ny)
c
      call printhist(tprint)
      call initpressure(nx,ny)

      if(restart.eq.0) then
	t0=0;
c
c Include here relevant initialisation module
c
c It should be called ini"something". 
c

         call inivort(u,v,z(icc(ng)),z(iu(ng)),timescalesi,tau,du,nx,ny)
c$$$        write(6,*) "x symmetry check for disk initialisation c"
c$$$        call  chkx(cc,nx,ny,0,1)
c$$$        write(6,*) "y symmetry check for disk initialisation c"
c$$$        call  chky(cc,nx,ny,0,1)
         write(6,*) v(12,70)
      else
         call iniget(u,v,z(icc(ng)),z(iu(ng)),nx,ny,t0)
      endif

      call makebccf(z(icc(ng)),nx,ny)
      call makebcf(u,v,nx,ny)
c
c Used in rhouv version
c
c      call acfromcc(z(icc(ng)),z(ia(ng)),z(ic(ng)),
c     %              z(iei(ng)),z(ires(ng)),nx,ny)
c      call rhouv(rhou,rhov,u,v,z(ia(ng)),z(ic(ng)),rhogas,nx,ny)
c      call makebcf(rhou,rhov,nx,ny)

      return
      end
