How to promote a factorization?

Specific Domains Numerics
1

I’m trying to promote an LU factorization in one precision to a higher
precision. The permutation vector does not map correctly and I can’t
figure out why.

Any ideas?

MWE:

julia> A=rand(Float32,3,3);

julia> AF=lu!(A);

julia> AF64=LU(Float64.(A), AF.p, AF.info);

julia> norm(AF64.L - AF.L)
0.00000e+00

julia> norm(AF64.U - AF.U)
0.00000e+00

julia> # so the factors look good, but

julia> AF.p
3-element Vector{Int64}:
 2
 1
 3

julia> AF64.p
3-element Vector{Int64}:
 1
 2
 3
2

Too me these look like permutation vectors, for which Int should be perfect on 64-bit systems…

3

You want:

AF64=LU(Float64.(A), AF.ipiv, AF.info)

(The p property is not actually stored, but rather is computed by LinearAlgebra.ipiv2perm.)

3 Likes
4

Sorry, I’m too stupid (I didn’t even realize the intent of original question). So what are you trying to do here: solving in Float32 and pretending to have it solved in Float64?

5

I’m computing a factorization in low precision (Float16 is the target) and using that factorization as a preconditioner for GMRES in a higher precision (Float32 is the target). This only works right if I promote the factorization before I use it in the preconditioner.

I could define the preconditioiner via brute force using L, U, and P, but it is far cleaner if I can simply promote the factorization. @stevengj fixed me up, as he as done in the past.

1 Like
6

That’s it! Thanks.