QR-Iteration Code

New to Julia
1

Hello everyone,

I have to program the qr iteration for the university and compare it with the commands eigvals and eigvecs. Unfortunately I can’t get the code to output the correct and same vectors as eigvecs. eigvecs gives me a complex vector and the Qr algorithm only a real number. Can someone maybe help me with this? Or are there certain matrix defaults why the qr code does not work?

Thanks in advance :slight_smile:

My code:
using LinearAlgebra
using BenchmarkTools

function qr_algorithm(matrix)
A = copy(matrix)
eigenvalues =
eigenvectors = Matrix{Float64}(I, size(A, 1), size(A, 2))

for k = 1:1000
    Q, R = qr(A)
    A = R * Q
    
    eigenvalues = diag(A)
    
    if norm(diag(A, 1), 2) < 1e-6
        break
    end
    
    eigenvectors *= Q
end

eigenvalues, eigenvectors

end

matrix1 = [7 9 3; 5 2 1; 8 4 6]

@time (qr_eigenvalues1, qr_eigenvectors1) = qr_algorithm(matrix1)
@time eigvals(matrix1)
@time eigvecs(matrix1)

matrix2 = [5 4 6 9 2 4 5 6;
6 8 0 5 3 4 3 2;
4 3 5 7 4 8 9 6;
2 4 5 7 8 7 3 2;
9 8 6 4 3 4 5 7;
3 4 5 8 2 6 3 5;
3 4 6 8 5 3 5 6;
5 6 8 4 2 5 7 5]

@time (qr_eigenvalues2, qr_eigenvectors2) = qr_algorithm(matrix2)
@time eigvals(matrix2)
@time eigvecs(matrix2)

println(“Qr Algorithmus Eigenvalues:”)
display(qr_eigenvalues1)
println(“\n”)

println(“QR Algorithm Eigenvalues:”)
display(qr_eigenvalues2)
println(“\n”)

println(“eigvals Eigenvalues:”)
display(eigvals(matrix1))
println(“\n”)

println(“eigvals Eigenvalues:”)
display(eigvals(matrix2))
println(“\n”)

2

here with fixed code formatting:

 using LinearAlgebra
 using BenchmarkTools
 
 function qr_algorithm(matrix)
 A = copy(matrix)
 eigenvalues = []
 eigenvectors = Matrix{Float64}(I, size(A, 1), size(A, 2))
 
 for k = 1:1000
     Q, R = qr(A)
     A = R * Q
 
     eigenvalues = diag(A)
 
     if norm(diag(A, 1), 2) < 1e-6
         break
     end
 
     eigenvectors *= Q
 end
 
 eigenvalues, eigenvectors
 
 end
 
 matrix1 = [7 9 3; 5 2 1; 8 4 6]
 
 @time (qr_eigenvalues1, qr_eigenvectors1) = qr_algorithm(matrix1)
 @time eigvals(matrix1)
 @time eigvecs(matrix1)
 
 matrix2 = [5 4 6 9 2 4 5 6;
 6 8 0 5 3 4 3 2;
 4 3 5 7 4 8 9 6;
 2 4 5 7 8 7 3 2;
 9 8 6 4 3 4 5 7;
 3 4 5 8 2 6 3 5;
 3 4 6 8 5 3 5 6;
 5 6 8 4 2 5 7 5]
 
 @time (qr_eigenvalues2, qr_eigenvectors2) = qr_algorithm(matrix2)
 @time eigvals(matrix2)
 @time eigvecs(matrix2)
 
 println("Qr Algorithmus Eigenvalues:")
 display(qr_eigenvalues1)
 println("\n")
 
 println("QR Algorithm Eigenvalues:")
 display(qr_eigenvalues2)
 println("\n")
 
 println("eigvals Eigenvalues:")
 display(eigvals(matrix1))
 println("\n")
 
 println("eigvals Eigenvalues:")
 display(eigvals(matrix2))
 println("\n")
1 Like
3

If this is a homework assignment, I’m not sure what sort of help is appropriate. An explicit unshifted QR iteration like this is not what is normally considered a practical QR iteration. But it doing this sort thing is instructional when learning about the properties of the QR iteration.

I have assigned this sort of thing before and would not want my students getting overly specific help on the internet. But even for homework , it seems reasonable to give you two general things to consider that I think you’ll need to look at if you want to understand what’s going on. Both of them relate to the fact that you need to look more closely at what the QR iteration is supposed to do to the matrix A:

  1. In the convergence criteria you are using, you are looking only at the first superdiagonal. This is sufficient if the matrix is symmetric tridiagonal, but not in general. As it is, your termination condition doesn’t do anything and both examples run for 1000 iterations.

  2. You might want to note that you appear to be getting the real eigenvalues, just not the conjugate pairs of eigenvalues. You should look up what the QR iteration does to a real, non-symmetric matrix. Specifically, it does not in general have all the eigenvalues on the diagonal of A. A good search term is the “real Schur decomposition.” Not all books in numerical linear algebra cover this. The popular book by Trefethen and Bau focuses on the symmetric/Hermitian case.

3 Likes