home *** CD-ROM | disk | FTP | other *** search
/ Power-Programmierung / CD1.mdf / fortran / library / ssp / mateigen / nroot.for < prev   
Text File  |  1985-11-29  |  4KB  |  144 lines
  1. C
  2. C     ..................................................................
  3. C
  4. C        SUBROUTINE NROOT
  5. C
  6. C        PURPOSE
  7. C           COMPUTE EIGENVALUES AND EIGENVECTORS OF A REAL NONSYMMETRIC
  8. C           MATRIX OF THE FORM B-INVERSE TIMES A.  THIS SUBROUTINE IS
  9. C           NORMALLY CALLED BY SUBROUTINE CANOR IN PERFORMING A
  10. C           CANONICAL CORRELATION ANALYSIS.
  11. C
  12. C        USAGE
  13. C           CALL NROOT (M,A,B,XL,X)
  14. C
  15. C        DESCRIPTION OF PARAMETERS
  16. C           M  - ORDER OF SQUARE MATRICES A, B, AND X.
  17. C           A  - INPUT MATRIX (M X M).
  18. C           B  - INPUT MATRIX (M X M).
  19. C           XL - OUTPUT VECTOR OF LENGTH M CONTAINING EIGENVALUES OF
  20. C                B-INVERSE TIMES A.
  21. C           X  - OUTPUT MATRIX (M X M) CONTAINING EIGENVECTORS COLUMN-
  22. C                WISE.
  23. C
  24. C        REMARKS
  25. C           NONE
  26. C
  27. C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
  28. C           EIGEN
  29. C
  30. C        METHOD
  31. C           REFER TO W. W. COOLEY AND P. R. LOHNES, 'MULTIVARIATE PRO-
  32. C           CEDURES FOR THE BEHAVIORAL SCIENCES', JOHN WILEY AND SONS,
  33. C           1962, CHAPTER 3.
  34. C
  35. C     ..................................................................
  36. C
  37.       SUBROUTINE NROOT (M,A,B,XL,X)
  38.       DIMENSION A(1),B(1),XL(1),X(1)
  39. C
  40. C        ...............................................................
  41. C
  42. C        IF A DOUBLE PRECISION VERSION OF THIS ROUTINE IS DESIRED, THE
  43. C        C IN COLUMN 1 SHOULD BE REMOVED FROM THE DOUBLE PRECISION
  44. C        STATEMENT WHICH FOLLOWS.
  45. C
  46. C     DOUBLE PRECISION A,B,XL,X,SUMV
  47. C
  48. C        THE C MUST ALSO BE REMOVED FROM DOUBLE PRECISION STATEMENTS
  49. C        APPEARING IN OTHER ROUTINES USED IN CONJUNCTION WITH THIS
  50. C        ROUTINE.
  51. C
  52. C        THE DOUBLE PRECISION VERSION OF THIS SUBROUTINE MUST ALSO
  53. C        CONTAIN DOUBLE PRECISION FORTRAN FUNCTIONS.  SQRT IN STATEMENTS
  54. C        110 AND 175 MUST BE CHANGED TO DSQRT.  ABS IN STATEMENT 110
  55. C        MUST BE CHANGED TO DABS.
  56. C
  57. C        ...............................................................
  58. C
  59. C     COMPUTE EIGENVALUES AND EIGENVECTORS OF B
  60. C
  61.       K=1
  62.       DO 100 J=2,M
  63.       L=M*(J-1)
  64.       DO 100 I=1,J
  65.       L=L+1
  66.       K=K+1
  67.   100 B(K)=B(L)
  68. C
  69. C        THE MATRIX B IS A REAL SYMMETRIC MATRIX.
  70. C
  71.       MV=0
  72.       CALL EIGEN (B,X,M,MV)
  73. C
  74. C     FORM RECIPROCALS OF SQUARE ROOT OF EIGENVALUES.  THE RESULTS
  75. C     ARE PREMULTIPLIED BY THE ASSOCIATED EIGENVECTORS.
  76. C
  77.       L=0
  78.       DO 110 J=1,M
  79.       L=L+J
  80.   110 XL(J)=1.0/ SQRT( ABS(B(L)))
  81.       K=0
  82.       DO 115 J=1,M
  83.       DO 115 I=1,M
  84.       K=K+1
  85.   115 B(K)=X(K)*XL(J)
  86. C
  87. C     FORM (B**(-1/2))PRIME * A * (B**(-1/2))
  88. C
  89.       DO 120 I=1,M
  90.       N2=0
  91.       DO 120 J=1,M
  92.       N1=M*(I-1)
  93.       L=M*(J-1)+I
  94.       X(L)=0.0
  95.       DO 120 K=1,M
  96.       N1=N1+1
  97.       N2=N2+1
  98.   120 X(L)=X(L)+B(N1)*A(N2)
  99.       L=0
  100.       DO 130 J=1,M
  101.       DO 130 I=1,J
  102.       N1=I-M
  103.       N2=M*(J-1)
  104.       L=L+1
  105.       A(L)=0.0
  106.       DO 130 K=1,M
  107.       N1=N1+M
  108.       N2=N2+1
  109.   130 A(L)=A(L)+X(N1)*B(N2)
  110. C
  111. C     COMPUTE EIGENVALUES AND EIGENVECTORS OF A
  112. C
  113.       CALL EIGEN (A,X,M,MV)
  114.       L=0
  115.       DO 140 I=1,M
  116.       L=L+I
  117.   140 XL(I)=A(L)
  118. C
  119. C     COMPUTE THE NORMALIZED EIGENVECTORS
  120. C
  121.       DO 150 I=1,M
  122.       N2=0
  123.       DO 150 J=1,M
  124.       N1=I-M
  125.       L=M*(J-1)+I
  126.       A(L)=0.0
  127.       DO 150 K=1,M
  128.       N1=N1+M
  129.       N2=N2+1
  130.   150 A(L)=A(L)+B(N1)*X(N2)
  131.       L=0
  132.       K=0
  133.       DO 180 J=1,M
  134.       SUMV=0.0
  135.       DO 170 I=1,M
  136.       L=L+1
  137.   170 SUMV=SUMV+A(L)*A(L)
  138.   175 SUMV= SQRT(SUMV)
  139.       DO 180 I=1,M
  140.       K=K+1
  141.   180 X(K)=A(K)/SUMV
  142.       RETURN
  143.       END
  144.