pro setup_slab
  fin='glass_10x10x10'
  readnew,fin,h,'HEAD'
  readnew,fin,xin,'POS'
  xin[*,*] /= h.BOXSIZE
  Nin = TOTAL(h.npart,/int)

  Ndup = 20L
  fout='glass_200x10x10_dvy'

  Np = Ndup * Nin
  x=fltarr(3,Np)

  i=0L
  FOR ix=0,Ndup-1 DO BEGIN
     FOR j=0L,Nin-1 DO BEGIN
        x[0,i] = xin[0,j] + ix
        x[1,i] = xin[1,j]
        x[2,i] = xin[2,j]
        i++
     END
  END

  v = fltarr(3,Np)
  u=fltarr(Np)
  m=fltarr(Np)
  id=lindgen(Np)

  h.BOXSIZE=1
  h.MASSARR[*]=0

  gamma=5./3.
  rho=0.125
  P = 1.0

  u[*] = P/(gamma-1)/rho
  m[*] = (rho *  h.BOXSIZE^3)/Np

   h.npart[*]=0
   h.npart[0]=Np

   h.PARTTOTAL[*]=0
   h.PARTTOTAL[0]=Np

;Add wave with an amplitude and wavelength
   du = 0.01 * P/(gamma-1)/rho
   dv = 0.1
   l=2.0

   FOR j=0L,Np-1 DO BEGIN
      if(x[0,j] LT l) THEN begin
         v[0,j] = dv * sin(!Pi * x[0,j]/l)
      end
      if(x[0,j] GT Ndup-l) THEN begin
         v[0,j] = -1 * dv * sin(!Pi * (Ndup-x[0,j])/l)
      end
   end

   write_head,fout,h
   add_block,fout,x,'POS '
   add_block,fout,v,'VEL '
   add_block,fout,id,'ID  '
   add_block,fout,m,'MASS'
   add_block,fout,u,'U   '

end