program simplegesv
      use HPF_LAPACK
      integer, parameter :: N=500, NRHS=20, NB=64, NBRHS=64, P=1, Q=3
      integer, parameter :: DP=kind(0.0D0)
      integer :: IPIV(N)
      real(DP) :: A(N, N), X(N, NRHS), B(N, NRHS)
!HPF$ PROCESSORS PROC(P,Q)
!HPF$ DISTRIBUTE A(cyclic(NB), cyclic(NB)) ONTO PROC
!HPF$ DISTRIBUTE (cyclic(NB), cyclic(NBRHS)) ONTO PROC :: B, X

!
!     Randomly generate the coefficient matrix A and the solution
!     matrix X.  Set the right hand side matrix B such that B = A * X.
!
      call random_number(A)
      call random_number(X)
      B = matmul(A, X)
!
!     Solve the linear system; the computed solution overwrites B
!
      call la_gesv(A, B, IPIV)
!
!     As a simple test, print the largest difference (in absolute value)
!     between the computed solution (B) and the generated solution (X).
!
      print*,'MAX( ABS(X~ - X) ) = ',maxval( abs(B - X) )
!
!     Shutdown the ScaLAPACK system, I'm done
!
      call SLhpf_exit()

      stop
      end