This matrix-renormalizer comes from WHAT IF. Feel free to use it any way
you want:
SUBROUTINE GVSREN (RMAT)
C+++
C-------------------------------------------------------------------------------
C----
----
C---- RMAT IS A REAL MATRIX DIMENSIONED
(3,3). ----
C---- GVSEIG can be any eigenvalue calculator. Ask me if you want the
WHAT IF one.
C----
----
C-------------------------------------------------------------------------------
C===
IMPLICIT NONE
INTEGER I, J, K, L
REAL A(6), RMAT(3,3), S(3,3), T(3,3), X(3,3), Y
C----
C---- FORM THE PRODUCT OF RMAT * RMAT(TRANSPOSE)
C----
L=0
DO 30 I=1,3
DO 20 J=1,I
L=L+1
A(L)=0.0
DO 10 K=1,3
A(L)=A(L)+RMAT(I,K)*RMAT(J,K)
10 CONTINUE
20 CONTINUE
30 CONTINUE
C----
C---- CALCULATE THE EIGENVALUES AND THE EIGENVECTORS OF A
C----
CALL GVSEIG (A,X,3,0)
C----
C---- SMALL PRECAUTION AGAINST FUTURE OVERFLOWS
C----
A(1)=1.0/SQRT(AMAX1(A(1),1.E-12))
A(2)=1.0/SQRT(AMAX1(A(3),1.E-12))
A(3)=1.0/SQRT(AMAX1(A(6),1.E-12))
C----
C---- FORM THE PRODUCT OF X * A
C----
DO 50 I=1,3
DO 40 J=1,3
S(I,J)=X(I,J)*A(J)
40 CONTINUE
50 CONTINUE
C----
C---- TRANSPOSE X
C----
DO 70 I=2,3
DO 60 J=1,I-1
Y=X(I,J)
X(I,J)=X(J,I)
X(J,I)=Y
60 CONTINUE
70 CONTINUE
C----
C---- CALCULATE THE RENORMALIZED RMAT AS S * X * RMAT
C----
CALL GVS3X3 (T,X,RMAT)
CALL GVS3X3 (RMAT,S,T)
RETURN
END
