Thank you for these answers.
@ctkelley @cgeoga Having already tried to implement the Gram-Schmidt algorithm in the past, I fear I won’t be able to match the stability and the performance of the library’s QR method. Thus my question was toward using the existing function of
LinearAlgebra to achieve my goal, if possible.
@stevengj Here is a usecase
I want to create an orthogonal matrix distributed on different processors, but this matrix can be very large, so it shouldn’t be generated all at once and then split the columns.
I generate the different partitions on the different processors, and then I orthogonalize the partitions against one another.
for p in workers()
# Compute the local orthonormal partition
@spawnat p A[:l] = Matrix(LinearAlgebra.qr(A[:l]).Q)
# On processors p+1 .. P, orthogonalize against the basis of p
for p2 in workers()[p:length(workers())]
@spawnat p2 begin
A_on_p = @fetchfrom p A[:l]
A[:l] = (I - A_on_p * A_on_p') * A[:l]
U = drand((m,m), workers(), [1,4])
See the line 4, instead I would like to perform the inplace orthonormalization of
A[:l] (the local part of the input matrix).
If I use
@spawnat p LinearAlgebra.qr!(A[:l]) on line 4 instead, I don’t understand how to use
A[:l] afterwards to orthogonalize the other partitions, because it contains only the reflectors. I cannot do
A_on_p * A_on_p' for example.