ssiev.f 3.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113
  1. *DECK SSIEV
  2. SUBROUTINE SSIEV (A, LDA, N, E, WORK, JOB, INFO)
  3. C***BEGIN PROLOGUE SSIEV
  4. C***PURPOSE Compute the eigenvalues and, optionally, the eigenvectors
  5. C of a real symmetric matrix.
  6. C***LIBRARY SLATEC
  7. C***CATEGORY D4A1
  8. C***TYPE SINGLE PRECISION (SSIEV-S, CHIEV-C)
  9. C***KEYWORDS COMPLEX HERMITIAN, EIGENVALUES, EIGENVECTORS, MATRIX,
  10. C SYMMETRIC
  11. C***AUTHOR Kahaner, D. K., (NBS)
  12. C Moler, C. B., (U. of New Mexico)
  13. C Stewart, G. W., (U. of Maryland)
  14. C***DESCRIPTION
  15. C
  16. C Abstract
  17. C SSIEV computes the eigenvalues and, optionally, the eigenvectors
  18. C of a real symmetric matrix.
  19. C
  20. C Call Sequence Parameters-
  21. C (The values of parameters marked with * (star) will be changed
  22. C by SSIEV.)
  23. C
  24. C A* REAL (LDA,N)
  25. C real symmetric input matrix.
  26. C Only the diagonal and upper triangle of A must be input,
  27. C as SSIEV copies the upper triangle to the lower.
  28. C That is, the user must define A(I,J), I=1,..N, and J=I,.
  29. C ..,N.
  30. C On return from SSIEV, if the user has set JOB
  31. C = 0 the lower triangle of A has been altered.
  32. C = nonzero the N eigenvectors of A are stored in its
  33. C first N columns. See also INFO below.
  34. C
  35. C LDA INTEGER
  36. C set by the user to
  37. C the leading dimension of the array A.
  38. C
  39. C N INTEGER
  40. C set by the user to
  41. C the order of the matrix A and
  42. C the number of elements in E.
  43. C
  44. C E* REAL (N)
  45. C on return from SSIEV, E contains the N
  46. C eigenvalues of A. See also INFO below.
  47. C
  48. C WORK* REAL (2*N)
  49. C temporary storage vector. Contents changed by SSIEV.
  50. C
  51. C JOB INTEGER
  52. C set by user on input
  53. C = 0 only calculate eigenvalues of A.
  54. C = nonzero calculate eigenvalues and eigenvectors of A.
  55. C
  56. C INFO* INTEGER
  57. C on return from SSIEV, the value of INFO is
  58. C = 0 for normal return.
  59. C = K if the eigenvalue iteration fails to converge.
  60. C eigenvalues and vectors 1 through K-1 are correct.
  61. C
  62. C
  63. C Error Messages-
  64. C No. 1 recoverable N is greater than LDA
  65. C No. 2 recoverable N is less than one
  66. C
  67. C***REFERENCES (NONE)
  68. C***ROUTINES CALLED IMTQL2, TQLRAT, TRED1, TRED2, XERMSG
  69. C***REVISION HISTORY (YYMMDD)
  70. C 800808 DATE WRITTEN
  71. C 890831 Modified array declarations. (WRB)
  72. C 890831 REVISION DATE from Version 3.2
  73. C 891214 Prologue converted to Version 4.0 format. (BAB)
  74. C 900315 CALLs to XERROR changed to CALLs to XERMSG. (THJ)
  75. C 900326 Removed duplicate information from DESCRIPTION section.
  76. C (WRB)
  77. C***END PROLOGUE SSIEV
  78. INTEGER INFO,JOB,LDA,N
  79. REAL A(LDA,*),E(*),WORK(*)
  80. C***FIRST EXECUTABLE STATEMENT SSIEV
  81. IF (N .GT. LDA) CALL XERMSG ('SLATEC', 'SSIEV', 'N .GT. LDA.',
  82. + 1, 1)
  83. IF(N .GT. LDA) RETURN
  84. IF (N .LT. 1) CALL XERMSG ('SLATEC', 'SSIEV', 'N .LT. 1', 2, 1)
  85. IF(N .LT. 1) RETURN
  86. C
  87. C CHECK N=1 CASE
  88. C
  89. E(1) = A(1,1)
  90. INFO = 0
  91. IF(N .EQ. 1) RETURN
  92. C
  93. C COPY UPPER TRIANGLE TO LOWER
  94. C
  95. DO 10 J=1,N
  96. DO 10 I=1,J
  97. A(J,I)=A(I,J)
  98. 10 CONTINUE
  99. C
  100. IF(JOB.NE.0) GO TO 20
  101. C
  102. C EIGENVALUES ONLY
  103. C
  104. CALL TRED1(LDA,N,A,E,WORK(1),WORK(N+1))
  105. CALL TQLRAT(N,E,WORK(N+1),INFO)
  106. RETURN
  107. C
  108. C EIGENVALUES AND EIGENVECTORS
  109. C
  110. 20 CALL TRED2(LDA,N,A,E,WORK,A)
  111. CALL IMTQL2(LDA,N,E,WORK,A,INFO)
  112. RETURN
  113. END