subroutine Jacobi(n)
      integer n              ! size including boundaries
      integer i, j, iters
      real grid(n,n), new(n,n), maxdiff

!HPF$ PROCESSORS pr(PR)      ! use PR processors
!HPF$ ALIGN grid(i,j) WITH new(i,j)
!HPF$ DISTRIBUTE grid(BLOCK) ONTO pr

      ! need to initialize grid and new, including boundaries
      do iters = 1, MAXITERS
        FORALL (i=2:n-1, j=2:n-1)
          new(i,j) = (grid(i-1,j) + grid(i+1,j) +
                      grid(i,j-1) + grid(i,j+1)) / 4
        grid = new           ! copies array in parallel
      end do

      maxdiff = MAXVAL(ABS(grid-new))  ! reduction
      end