c     These are test routines to make sure
c     my poisson solver was working correctly
c     
c     They have no other purpose beyond this.

      subroutine testPsiFromW(nX, nY)

      implicit none

      integer nX, nY, nMax

      parameter(nMax=1000)

      real*8 psi, w, x, y, wTmp
      real*8 fEps 
      parameter (fEps=1.d-5)

      dimension psi(0:nmax+1,0:nmax+1)
      dimension w(0:nmax+1,0:nmax+1)
      dimension wTmp(0:nmax+1,0:nmax+1)
      dimension x(0:nx+1), y(0:ny+1)

      call step(x, 0.d0, 1.d0, nx, .TRUE.)
      call step(y, 0.d0, 1.d0, ny, .FALSE.)

      call setTestField1(w(0,0), nx, ny, nMax,
     *     x(0),y(0), 1.d-1)

      call printResult(x(0),y,w(0,0),
     *     nx,ny,nmax,0.d0,0,
     *     .TRUE. , .TRUE.)

      call set_const_2d(psi(0,0), 
     *     nx, ny, nMax, 0.d0)

      call mul_and_copy(w, 
     *     psi, 2.d0*dacos(-1.d0), 
     *     nX, nY,
     *     nMax, nMax,
     *     .TRUE.)



      call psiFromW(psi(0,0), w,wTmp,
     *     nx, ny, nMax,
     *     x(1)-x(0), y(1)-y(0), fEps)

      call printResult(x(0), y(0), psi(0,0),
     *     nx, ny, nmax, 0.d0, 1,
     *     .TRUE. , .TRUE.)

      pause "Is all in order?"


      call wFromPsi(psi(0,0), w(0,0), 
     *     nx, ny, nMax,
     *     x(1)-x(0), y(1)-y(0), fEps)


      call printResult(x(0),y(0),w(0,0),
     *     nx,ny,nmax,0.d0,2,
     *     .TRUE. , .TRUE.)


      return
      end



      subroutine setTestField1(fieldToSet, nX, nY, nMax,
     *     x, y, fHeight)

      implicit none
      
      integer i,j, nX, nY, nMax
      real*8 fieldToSet, x, y, fHeight
      real*8 fXBounds, fYBounds
      real*8 fPi

      dimension fieldToSet(0:nMax+1, 0:nMax+1)
      dimension x(0:nX+1), y(0:nY+1)

      fXBounds = x(nX)-x(1)
      fYBounds = y(nY+1)-y(0)

      fPi = dacos(-1.d0)

      do j=1,nY
         do i=1,nX
            fieldToSet(i,j) = fHeight*
     *           (dsin(2.d0*fPi*x(i)) +
     +           dcos(2.d0*fPi*y(j)))
         enddo
      enddo

      return
      end