C$TEST PRSF
C TO RUN AS A MAIN PROGRAM REMOVE NEXT LINE
      SUBROUTINE PRSF
C***********************************************************************
C
C  EXAMPLE OF USE OF THE PORT PROGRAM SPFLE
C
C***********************************************************************
       INTEGER N1, N2, N, I, ITIME, ITIME1, ITIME2, ILAPSZ
       INTEGER I1MACH, IREAD, IWRITE, ISTAK(25000), ISIZE, ISIZE2
       REAL RSTAK(25000), B(500), A1, A2
       EXTERNAL QUEA
       COMMON /QUE/ A1, A2, N1, N2, N
       COMMON /CSTAK/ ISTAK
       EQUIVALENCE (ISTAK(1),RSTAK(1))
       IREAD = I1MACH(1)
       IWRITE = I1MACH(2)
       CALL ISTKIN(25000,2)
  10   READ(IREAD, 11)N1, N2, A1, A2
  11   FORMAT(2I3,2E15.5)
       IF (N1 .EQ. 0) STOP
       WRITE(IWRITE,12)N1, N2, A1, A2
  12   FORMAT(/4H N1=,I3,4H N2=,I3,4H A1=,E15.5,4H A2=,E15.5)
       N = (N1+1)*(N2+1)
C
C GENERATE THE RIGHT HAND SIDE
C
       DO 20 I=1, N
          B(I) = 0.0
  20   CONTINUE
       B(N) = 1.0
C
C SOLVE THE SYSTEM WITH PIVOTING FOR SPARSITY AND TIME IT
C
       ITIME=ILAPSZ(0)
       CALL SPFLE(N, .TRUE., QUEA, ISIZE, B, 500, 1)
       ITIME1=ILAPSZ(0)-ITIME
C
C FIND THE PROBABILITIES
C
       IT = (N2+1)*N1+1
       WRITE(IWRITE,21)B(N),SASUM(N1+1,B(N2+1),N2+1),SASUM(N2,B(IT),1)
  21   FORMAT(6H L1 = ,E15.5,6H L2 = ,E15.5,7H I12 = ,E15.5)
       WRITE(IWRITE,22)
  22   FORMAT(22H PIVOTING FOR SPARSITY)
       WRITE(IWRITE,23)ISIZE
  23   FORMAT(34H SPACE NEEDED FOR DECOMPOSITION - ,I8)
       WRITE(IWRITE,24)ITIME1
  24   FORMAT(15H TIME NEEDED - ,I8)
C
C FOR COMPARISON, REDO PROBLEM WITHOUT REQUESTING PIVOTING
C
      DO 30 I= 1,N
         B(I)=0.0
  30  CONTINUE
      B(N)=1.0
      ITIM=ILAPSZ(0)
      CALL SPFLE(N,.FALSE.,QUEA, ISIZE2, B, 500, 1)
      ITIME2=ILAPSZ(0)-ITIME
      WRITE(IWRITE, 31)
  31  FORMAT(25H NO PIVOTING FOR SPARSITY)
      WRITE(IWRITE,23)ISIZE2
      WRITE(IWRITE, 24)ITIME2
       GO TO 10
       END
       SUBROUTINE QUEA(I, ROW, JCOL, NUM)
       INTEGER I, NUM, JCOL(100), N, N1, N2, II, JJ, J
       REAL ROW(100)
       COMMON /QUE/ A1,A2, N1,N2, N
       IF (I.NE.N) GO TO 20
C TREAT LAST ROW AS SPECIAL CASE
       DO 10 J=1, N
          JCOL(J) = J
          ROW(J) = 1.0
  10   CONTINUE
       NUM = N
       RETURN
  20   N2P1=N2+1
C DETERMINE WHICH MAJOR BLOCK
       II=(I-1)/N2P1
C DETERMINE POSITION IN BLOCK
       JJ = MOD(I-1, N2P1)
C FILL IN DIAGONAL
       JCOL(1) = I
       ROW(1) = -A1-A2-FLOAT(II+JJ)
       IF (JJ.EQ.N2) ROW(1)=ROW(1)+A2
       NUM = 1
       IF (II .EQ. 0) GO TO 30
C THIS IS NOT THE FIRST BLOCK
       JCOL(2) = I-N2P1
       ROW(2) = A1
       NUM = 2
  30   IF (JJ.EQ.0) GO TO 40
C THIS IS NOT FIRST ROW IN THE BLOCK
       NUM = NUM+1
       JCOL(NUM) = I-1
       ROW(NUM) = A2
  40   IF (JJ.EQ. N2) GO TO 50
C THIS IS NOT LAST ROW IN THE BLOCK
       NUM=NUM+1
       JCOL(NUM)= I+1
       ROW(NUM) = JJ+1
  50   IF (II .EQ. N1) RETURN
C THIS IS NOT THE LAST BLOCK
       NUM = NUM +1
       JCOL(NUM) = I+N2P1
       ROW(NUM)= II+1
       RETURN
       END
      INTEGER FUNCTION ILAPSZ(N)
      INTEGER N
      ILAPSZ = 0
      RETURN
      END
C
C DATA FOR THE EXAMPLE IN THE PORT SHEET...  (REMOVE THE C
C IN COLUMN 1 BEFORE FEEDING THIS DATA TO THE PROGRAM ABOVE.)
C$DATA
C10 10   9.             9.
C20 20  19.            19.
C 0  0   0.             0.