home *** CD-ROM | disk | FTP | other *** search
/ Simtel MSDOS 1992 September / Simtel20_Sept92.cdr / msdos / fortran / linpklib.arc / SSISL.FOR < prev    next >
Text File  |  1984-01-07  |  4KB  |  172 lines
  1.       SUBROUTINE SSISL(A,LDA,N,KPVT,B)
  2.       INTEGER LDA,N,KPVT(1)
  3.       REAL A(LDA,1),B(1)
  4. C
  5. C     SSISL SOLVES THE REAL SYMMETRIC SYSTEM
  6. C     A * X = B
  7. C     USING THE FACTORS COMPUTED BY SSIFA.
  8. C
  9. C     ON ENTRY
  10. C
  11. C        A       REAL(LDA,N)
  12. C                THE OUTPUT FROM SSIFA.
  13. C
  14. C        LDA     INTEGER
  15. C                THE LEADING DIMENSION OF THE ARRAY  A .
  16. C
  17. C        N       INTEGER
  18. C                THE ORDER OF THE MATRIX  A .
  19. C
  20. C        KPVT    INTEGER(N)
  21. C                THE PIVOT VECTOR FROM SSIFA.
  22. C
  23. C        B       REAL(N)
  24. C                THE RIGHT HAND SIDE VECTOR.
  25. C
  26. C     ON RETURN
  27. C
  28. C        B       THE SOLUTION VECTOR  X .
  29. C
  30. C     ERROR CONDITION
  31. C
  32. C        A DIVISION BY ZERO MAY OCCUR IF  SSICO  HAS SET RCOND .EQ. 0.0
  33. C        OR  SSIFA  HAS SET INFO .NE. 0  .
  34. C
  35. C     TO COMPUTE  INVERSE(A) * C  WHERE  C  IS A MATRIX
  36. C     WITH  P  COLUMNS
  37. C           CALL SSIFA(A,LDA,N,KPVT,INFO)
  38. C           IF (INFO .NE. 0) GO TO ...
  39. C           DO 10 J = 1, P
  40. C              CALL SSISL(A,LDA,N,KPVT,C(1,J))
  41. C        10 CONTINUE
  42. C
  43. C     LINPACK. THIS VERSION DATED 08/14/78 .
  44. C     JAMES BUNCH, UNIV. CALIF. SAN DIEGO, ARGONNE NAT. LAB.
  45. C
  46. C     SUBROUTINES AND FUNCTIONS
  47. C
  48. C     BLAS SAXPY,SDOT
  49. C     FORTRAN IABS
  50. C
  51. C     INTERNAL VARIABLES.
  52. C
  53.       REAL AK,AKM1,BK,BKM1,SDOT,DENOM,TEMP
  54.       INTEGER K,KP
  55. C
  56. C     LOOP BACKWARD APPLYING THE TRANSFORMATIONS AND
  57. C     D INVERSE TO B.
  58. C
  59.       K = N
  60.    10 IF (K .EQ. 0) GO TO 80
  61.          IF (KPVT(K) .LT. 0) GO TO 40
  62. C
  63. C           1 X 1 PIVOT BLOCK.
  64. C
  65.             IF (K .EQ. 1) GO TO 30
  66.                KP = KPVT(K)
  67.                IF (KP .EQ. K) GO TO 20
  68. C
  69. C                 INTERCHANGE.
  70. C
  71.                   TEMP = B(K)
  72.                   B(K) = B(KP)
  73.                   B(KP) = TEMP
  74.    20          CONTINUE
  75. C
  76. C              APPLY THE TRANSFORMATION.
  77. C
  78.                CALL SAXPY(K-1,B(K),A(1,K),1,B(1),1)
  79.    30       CONTINUE
  80. C
  81. C           APPLY D INVERSE.
  82. C
  83.             B(K) = B(K)/A(K,K)
  84.             K = K - 1
  85.          GO TO 70
  86.    40    CONTINUE
  87. C
  88. C           2 X 2 PIVOT BLOCK.
  89. C
  90.             IF (K .EQ. 2) GO TO 60
  91.                KP = IABS(KPVT(K))
  92.                IF (KP .EQ. K - 1) GO TO 50
  93. C
  94. C                 INTERCHANGE.
  95. C
  96.                   TEMP = B(K-1)
  97.                   B(K-1) = B(KP)
  98.                   B(KP) = TEMP
  99.    50          CONTINUE
  100. C
  101. C              APPLY THE TRANSFORMATION.
  102. C
  103.                CALL SAXPY(K-2,B(K),A(1,K),1,B(1),1)
  104.                CALL SAXPY(K-2,B(K-1),A(1,K-1),1,B(1),1)
  105.    60       CONTINUE
  106. C
  107. C           APPLY D INVERSE.
  108. C
  109.             AK = A(K,K)/A(K-1,K)
  110.             AKM1 = A(K-1,K-1)/A(K-1,K)
  111.             BK = B(K)/A(K-1,K)
  112.             BKM1 = B(K-1)/A(K-1,K)
  113.             DENOM = AK*AKM1 - 1.0E0
  114.             B(K) = (AKM1*BK - BKM1)/DENOM
  115.             B(K-1) = (AK*BKM1 - BK)/DENOM
  116.             K = K - 2
  117.    70    CONTINUE
  118.       GO TO 10
  119.    80 CONTINUE
  120. C
  121. C     LOOP FORWARD APPLYING THE TRANSFORMATIONS.
  122. C
  123.       K = 1
  124.    90 IF (K .GT. N) GO TO 160
  125.          IF (KPVT(K) .LT. 0) GO TO 120
  126. C
  127. C           1 X 1 PIVOT BLOCK.
  128. C
  129.             IF (K .EQ. 1) GO TO 110
  130. C
  131. C              APPLY THE TRANSFORMATION.
  132. C
  133.                B(K) = B(K) + SDOT(K-1,A(1,K),1,B(1),1)
  134.                KP = KPVT(K)
  135.                IF (KP .EQ. K) GO TO 100
  136. C
  137. C                 INTERCHANGE.
  138. C
  139.                   TEMP = B(K)
  140.                   B(K) = B(KP)
  141.                   B(KP) = TEMP
  142.   100          CONTINUE
  143.   110       CONTINUE
  144.             K = K + 1
  145.          GO TO 150
  146.   120    CONTINUE
  147. C
  148. C           2 X 2 PIVOT BLOCK.
  149. C
  150.             IF (K .EQ. 1) GO TO 140
  151. C
  152. C              APPLY THE TRANSFORMATION.
  153. C
  154.                B(K) = B(K) + SDOT(K-1,A(1,K),1,B(1),1)
  155.                B(K+1) = B(K+1) + SDOT(K-1,A(1,K+1),1,B(1),1)
  156.                KP = IABS(KPVT(K))
  157.                IF (KP .EQ. K) GO TO 130
  158. C
  159. C                 INTERCHANGE.
  160. C
  161.                   TEMP = B(K)
  162.                   B(K) = B(KP)
  163.                   B(KP) = TEMP
  164.   130          CONTINUE
  165.   140       CONTINUE
  166.             K = K + 2
  167.   150    CONTINUE
  168.       GO TO 90
  169.   160 CONTINUE
  170.       RETURN
  171.       END
  172.