# orthit2.s                         J F Monahan, June 2001
 # orthogonal iteration for a few eigenvalues and vectors
# second problem
 n <- 4   # dimension of matrix
 p <- 4   # of cols or eigenvals to get
 A <- c(1, -3, -2, 1, -3, 10, -3, 6, -2, -3, 3, -2, 1, 6, -2, 1)
 A <- matrix(A,n,n)
 "matrix"
 A
 I <- diag(1,n)[,1:p]      # just p columns
 Z <- I
 for( i in 1:50) {
 V <- A %*% Z
 Z <- qr.qy(qr(V),I)
 R <- t(Z) %*% V 
 if( (i %% 5) == 0 ) {     # every fifth iteration
                      print("iteration"); print(i)
                      print("Z matrix"); print(Z)
                      print("R matrix"); print(R)
                     }
                }
 rm(list=ls())