I’ve learned that.
(But the title of the topic is very general, which is tempting for people to discuss things.)
But I find that I cannot perform U \ UQ
(see the 🟠
line) due to
ERROR: MethodError: no method matching /(::LinearAlgebra.AdjointQ{Float64, LinearAlgebra.QRCompactWYQ{…}}, ::LowerTriangular{Float64, Adjoint{…}})
.
Is this normal?
Another question is, why is the R
factor not of type UpperTriangular
? (From the perspective of LinearAlgebra.jl, UpperTriangular
is not considered a sparse matrix? since I notice that lu
has the same behavior) Would it be beneficial if I do it myself (see the 🟤
line)?
using LinearAlgebra
N = 4
P = rand(-9:.17:9, N, N); P = P' * P
U = cholesky(P).U # U' * U == P, Moreover P is posdef if no ERROR occurs
A = rand(-9:.17:9, N, N) # columns in A are linearly independent, verified as follows
rank(A) == N || error("A is not full column rank")
UQ, R = qr(U * A)
R = UpperTriangular(R) # 🟤
Q = U \ Matrix(UQ) # 🟠
Q' * P * Q - I |> norm
Q * R .- A |> norm
In this example, columns of A
are represented with the basis being the columns of Q
, which are “P-conjugate”—a concept in the context of conjugate gradient method—as testified in the last two lines.