C------------------------------------------------------------------------------
C       Example program to show the use of the "wgsmp" routine.
C------------------------------------------------------------------------------
C.. This program can be obtained from:
C
C   http://www.cs.umn.edu/~agupta/wsmp
C
C   (C) IBM Corporation, 2000.
C
c       Acceptance and use of this program constitutes the user's understanding
c       that he/she will have no recourse to IBM for any actual or consequential
c       damages, including, but not limited to, lost profits or savings, arising
c       out of the use or inability to use this program.
C------------------------------------------------------------------------------
        program wgsmp_ex1
        implicit none

        integer n, ldb, nrhs
        integer iparm(64)
        integer ia(10)
        integer ja(32)
        double precision dparm(64)
        double precision avals(32)
        double precision b(9)
        double precision darray          ! just a placeholder in this program

        integer i
        double precision waltime, wsmprtc

C.. Fill all arrays containing matrix data.

        data n /9/, ldb /9/, nrhs /1/

        data ia /1,5,9,12,14,17,22,24,29,33/

        data ja
     1        /1,          3,                      7,    8,
     2               2,    3,          5,                      9,
     3         1,          3,                      7,            
     4                           4,                      8,
     5                                 5,    6,                9,
     6               2,          4,          6,    7,    8,      
     7                           4,                7,      
     8         1,          3,          5,          7,    8,      
     9               2,    3,                            8,    9/
        data avals
     1      /14.d0,      -5.d0,                  -1.d0,-6.d0,
     2             14.d0,-1.d0,      -3.d0,                  -1.d0,
     3       -1.d0,      16.d0,                  -2.d0,
     4                         14.d0,                  -3.d0,
     5                               14.d0,-1.d0,            -1.d0,
     6             -2.d0,      -1.d0,      16.d0,-2.d0,-4.d0,
     7                         -1.d0,            16.d0,
     8       -3.d0,      -4.d0,      -3.d0,      -4.d0,71.d0,
     9             -1.d0,-2.d0,                        -4.d0,16.d0/

C.. Executable part of the program starts here.

      waltime = wsmprtc()

C.. Fill 'iparm' and 'dparm' arrays with default values.

C.. As an alternative to this step, the values in 'iparm' and 'dparm' can be
c   filled with values suitable for the application either manually or by
c   using a "data" statement according to the description in the User's guide.

      iparm(1) = 0
      iparm(2) = 0
      iparm(3) = 0
      call wgsmp (n, ia, ja, avals, b, ldb, nrhs, darray, iparm, dparm)      

      print *,'Initialization complete in time - ',wsmprtc()-waltime
      if (iparm(64) .ne. 0) then
        print *,'The following ERROR was detected: ',iparm(64)
        stop
      end if

C.. Analysis.

      waltime = wsmprtc()
      iparm(2) = 1
      iparm(3) = 1
      call wgsmp (n, ia, ja, avals, b, ldb, nrhs, darray, iparm, dparm)      

      print *,'Analysis complete in time - ',wsmprtc()-waltime

      if (iparm(64) .ne. 0) then
        print *,'The following ERROR was detected: ',iparm(64)
        stop
      end if

      print *,'Number of nonzeros in LU factors = ',iparm(24)
      print *,'Number of FLOPS in factorization = ',dparm(24)

C.. Factorization.

      waltime = wsmprtc()
      iparm(2) = 2
      iparm(3) = 2
      call wgsmp (n, ia, ja, avals, b, ldb, nrhs, darray, iparm, dparm)

      waltime = wsmprtc() - waltime
      print *,'Factorization complete in time - ',waltime
      if (iparm(64) .ne. 0) then
        print *,'The following ERROR was detected: ',iparm(64)
        stop
      end if
      print *,'Factorization MegaFlops = ',(dparm(23)*1.d-6)/waltime

C.. Back substitution.

      do i = 1, n
        b(i) = 1.d0
      end do

      waltime = wsmprtc()
      iparm(2) = 3
      iparm(3) = 3
      call wgsmp (n, ia, ja, avals, b, ldb, nrhs, darray, iparm, dparm)

      print *,'Back substitution complete in time - ',wsmprtc()-waltime
      if (iparm(64) .ne. 0) then
        print *,'The following ERROR was detected: ',iparm(64)
        stop
      end if

C.. Iterative refinement.

      waltime = wsmprtc()
      iparm(2) = 4
      iparm(3) = 4
      call wgsmp (n, ia, ja, avals, b, ldb, nrhs, darray, iparm, dparm)

      print *,'Iterative refinement done in time - ',wsmprtc()-waltime
      if (iparm(64) .ne. 0) then
        print *,'The following ERROR was detected: ',iparm(64)
        stop
      end if

      print *,'Maximum relative error = ',dparm(7)
      print *,'The solution of the system is as follows:'
      do i = 1, n
        print *,i,' : ',b(i)
      end do

      stop
      end