      PROGRAM LAPLACE
c     solving the laplace pde
      REAL F(:,:), DF(:,:)
      INTEGER MAXX, MAXY
      REAL FMAX
      INTEGER ITER
c     read in sizes 
      PRINT *,'MAXX = (e.g. 64) '
      READ *,MAXX
      PRINT *,'MAXY = (e.g. 64) '
      READ *,MAXY
c     allocating memory
      ALLOCATE (F(MAXX,MAXY), DF(1:MAXX,MAXY))
c     initialization
      F = 2.
      F(:,MAXY) = 1.
      forall (i=2:maxx-1,j=2:maxy-1)
         f(i,j) =  0
      end forall
      ITER = 0
      FMAX = 1
c     iteration
      DF = 0
      DO WHILE (FMAX .gt. 0.001)
         ITER = ITER + 1
         forall (i=2:maxx-1,j=2:maxy-1)
           df(i,j) = (f(i,j+1) + f(i,j-1) + f(i-1,j) + f(i+1,j))*0.25
     $               - f(i,j)
         end forall
         forall (i=2:maxx-1,j=2:maxy-1)
            f(i,j) = f(i,j) + df(i,j)
         end forall
         DF = ABS(DF)
         FMAX = MAXVAL (DF)
         PRINT *,'iteration ',ITER,'  Max = ',FMAX
      END DO
      PRINT *, ITER, ' iterations needed'
      deallocate (df, f)
      END 
