home *** CD-ROM | disk | FTP | other *** search
/ Simtel MSDOS 1992 September / Simtel20_Sept92.cdr / msdos / fortran / linpklib.arc / STRSL.FOR < prev   
Text File  |  1984-01-12  |  4KB  |  135 lines
  1.       SUBROUTINE STRSL(T,LDT,N,B,JOB,INFO)
  2.       INTEGER LDT,N,JOB,INFO
  3.       REAL T(LDT,1),B(1)
  4. C
  5. C
  6. C     STRSL SOLVES SYSTEMS OF THE FORM
  7. C
  8. C                   T * X = B
  9. C     OR
  10. C                   TRANS(T) * X = B
  11. C
  12. C     WHERE T IS A TRIANGULAR MATRIX OF ORDER N. HERE TRANS(T)
  13. C     DENOTES THE TRANSPOSE OF THE MATRIX T.
  14. C
  15. C     ON ENTRY
  16. C
  17. C         T         REAL(LDT,N)
  18. C                   T CONTAINS THE MATRIX OF THE SYSTEM. THE ZERO
  19. C                   ELEMENTS OF THE MATRIX ARE NOT REFERENCED, AND
  20. C                   THE CORRESPONDING ELEMENTS OF THE ARRAY CAN BE
  21. C                   USED TO STORE OTHER INFORMATION.
  22. C
  23. C         LDT       INTEGER
  24. C                   LDT IS THE LEADING DIMENSION OF THE ARRAY T.
  25. C
  26. C         N         INTEGER
  27. C                   N IS THE ORDER OF THE SYSTEM.
  28. C
  29. C         B         REAL(N).
  30. C                   B CONTAINS THE RIGHT HAND SIDE OF THE SYSTEM.
  31. C
  32. C         JOB       INTEGER
  33. C                   JOB SPECIFIES WHAT KIND OF SYSTEM IS TO BE SOLVED.
  34. C                   IF JOB IS
  35. C
  36. C                        00   SOLVE T*X=B, T LOWER TRIANGULAR,
  37. C                        01   SOLVE T*X=B, T UPPER TRIANGULAR,
  38. C                        10   SOLVE TRANS(T)*X=B, T LOWER TRIANGULAR,
  39. C                        11   SOLVE TRANS(T)*X=B, T UPPER TRIANGULAR.
  40. C
  41. C     ON RETURN
  42. C
  43. C         B         B CONTAINS THE SOLUTION, IF INFO .EQ. 0.
  44. C                   OTHERWISE B IS UNALTERED.
  45. C
  46. C         INFO      INTEGER
  47. C                   INFO CONTAINS ZERO IF THE SYSTEM IS NONSINGULAR.
  48. C                   OTHERWISE INFO CONTAINS THE INDEX OF
  49. C                   THE FIRST ZERO DIAGONAL ELEMENT OF T.
  50. C
  51. C     LINPACK. THIS VERSION DATED 08/14/78 .
  52. C     G. W. STEWART, UNIVERSITY OF MARYLAND, ARGONNE NATIONAL LAB.
  53. C
  54. C     SUBROUTINES AND FUNCTIONS
  55. C
  56. C     BLAS SAXPY,SDOT
  57. C     FORTRAN MOD
  58. C
  59. C     INTERNAL VARIABLES
  60. C
  61.       REAL SDOT,TEMP
  62.       INTEGER CASE,J,JJ
  63. C
  64. C     BEGIN BLOCK PERMITTING ...EXITS TO 150
  65. C
  66. C        CHECK FOR ZERO DIAGONAL ELEMENTS.
  67. C
  68.          DO 10 INFO = 1, N
  69. C     ......EXIT
  70.             IF (T(INFO,INFO) .EQ. 0.0E0) GO TO 150
  71.    10    CONTINUE
  72.          INFO = 0
  73. C
  74. C        DETERMINE THE TASK AND GO TO IT.
  75. C
  76.          CASE = 1
  77.          IF (MOD(JOB,10) .NE. 0) CASE = 2
  78.          IF (MOD(JOB,100)/10 .NE. 0) CASE = CASE + 2
  79.          GO TO (20,50,80,110), CASE
  80. C
  81. C        SOLVE T*X=B FOR T LOWER TRIANGULAR
  82. C
  83.    20    CONTINUE
  84.             B(1) = B(1)/T(1,1)
  85.             IF (N .LT. 2) GO TO 40
  86.             DO 30 J = 2, N
  87.                TEMP = -B(J-1)
  88.                CALL SAXPY(N-J+1,TEMP,T(J,J-1),1,B(J),1)
  89.                B(J) = B(J)/T(J,J)
  90.    30       CONTINUE
  91.    40       CONTINUE
  92.          GO TO 140
  93. C
  94. C        SOLVE T*X=B FOR T UPPER TRIANGULAR.
  95. C
  96.    50    CONTINUE
  97.             B(N) = B(N)/T(N,N)
  98.             IF (N .LT. 2) GO TO 70
  99.             DO 60 JJ = 2, N
  100.                J = N - JJ + 1
  101.                TEMP = -B(J+1)
  102.                CALL SAXPY(J,TEMP,T(1,J+1),1,B(1),1)
  103.                B(J) = B(J)/T(J,J)
  104.    60       CONTINUE
  105.    70       CONTINUE
  106.          GO TO 140
  107. C
  108. C        SOLVE TRANS(T)*X=B FOR T LOWER TRIANGULAR.
  109. C
  110.    80    CONTINUE
  111.             B(N) = B(N)/T(N,N)
  112.             IF (N .LT. 2) GO TO 100
  113.             DO 90 JJ = 2, N
  114.                J = N - JJ + 1
  115.                B(J) = B(J) - SDOT(JJ-1,T(J+1,J),1,B(J+1),1)
  116.                B(J) = B(J)/T(J,J)
  117.    90       CONTINUE
  118.   100       CONTINUE
  119.          GO TO 140
  120. C
  121. C        SOLVE TRANS(T)*X=B FOR T UPPER TRIANGULAR.
  122. C
  123.   110    CONTINUE
  124.             B(1) = B(1)/T(1,1)
  125.             IF (N .LT. 2) GO TO 130
  126.             DO 120 J = 2, N
  127.                B(J) = B(J) - SDOT(J-1,T(1,J),1,B(1),1)
  128.                B(J) = B(J)/T(J,J)
  129.   120       CONTINUE
  130.   130       CONTINUE
  131.   140    CONTINUE
  132.   150 CONTINUE
  133.       RETURN
  134.       END
  135.