
subroutine setup(nx, ny, nz, maxitr, dosmooth, width, ampl, &
     vx, vy, vz, epsilon, bbox, lb, cfl, dt, h, time, interval)
  implicit none

  integer nx, ny, nz, maxitr, dosmooth, interval

  real*8 width, ampl, vx,vy,vz, epsilon, bbox(6), lb(3), cfl, dt, h, time

  open(unit=1, file='wave3d.in')

  read(1,*) bbox
  read(1,*) width
  read(1,*) ampl
  read(1,*) vx, vy, vz
  read(1,*) cfl
  read(1,*) nx, ny, nz
  read(1,*) maxitr
  read(1,*) epsilon
  read(1,*) dosmooth
  read(1,*) interval
  
  close (unit=1)
  
  lb(:) = 0.0d0

  h = (bbox(4) - bbox(1)) / float(nx-1)

  dt = h * cfl

  time = 0.0d0

  return
end subroutine setup


subroutine carefulswap(nx,ny,nz, a, b)
  implicit none

  integer nx,ny,nz

  real*8, dimension(nx,ny,nz) :: a,b

  ! Locals
  
  integer i,j,k

  real*8 r

  do i=1, nx
     do j=1, ny
        do k=1, nz

           r = a(i,j,k)
           a(i,j,k) = b(i,j,k)
           b(i,j,k) = r
   
        end do
     end do
  end do
  
end
