Support general eltypes in matrix division and SVD #1453
Conversation
Codecov Report
@@ Coverage Diff @@
## master #1453 +/- ##
==========================================
- Coverage 77.95% 77.35% -0.60%
==========================================
Files 120 120
Lines 9227 9230 +3
==========================================
- Hits 7193 7140 -53
- Misses 2034 2090 +56
Continue to review full report at Codecov.
|
Wouldn't it be surprising for |
|
I don't think so for Integers at least. Factorizing an Integer matrix is a perfectly reasonable thing to do, but requires promoting to a float eltype. Base does it all the time: julia> A = rand(1:10, 2, 3)
2×3 Matrix{Int64}:
[...]
julia> B = rand(1:10, 2, 3)
2×3 Matrix{Int64}:
[...]
julia> A \ B
3×3 Matrix{Float64}:
[...]
julia> svd(A)
SVD{Float64, Float64, Matrix{Float64}}
[...]However, the mutating methods don't promote. They won't make extra copies behind the scenes so they need a type they can actually work with: julia> svd!(A)
ERROR: MethodError: no method matching svd!(::Matrix{Int64})
[...]The more open question is for Float16 and ComplexF16. It would definitely be possible to convert the result back to Float16 for those, like, e.g., julia> A = rand(Float16, 2, 3)
2×3 Matrix{Float16}:
[...]
julia> B = rand(Float16, 2, 3)
2×3 Matrix{Float16}:
[...]
julia> A \ B
3×3 Matrix{Float16}:
[...]
julia> svd(A)
SVD{Float32, Float32, Matrix{Float32}}
[...]
julia> lu(A)
LU{Float16, Matrix{Float16}}
[...]I suppose the cleanest behavior would be to accept float-typed returns for Integer matrices, but convert back down for Float16/ComplexF16---just like |
|
On further thought, converting back to F16 after the CUSOLVER call requires an extra copy, which should be optional (unlike functions with scalar output like |
In that case, I'm happy to keep it as is, even without converting back to FP16. Thanks for the PR! |
maleadt
added
enhancement
Hi again,
Since
Base.:\, LinearAlgebra.svd, LinearAlgebra.svdvalsneed to make copies anyway, they can support general eltypes at no extra cost: just promote to a cublas/cusolver-compatible eltype when copying.For now, I only implemented this for division and SVD since CUDA.CUSOLVER already specializes the non-mutating functions for these, so it was a drop-in replacement. Ideally, it would be extended to QR, LU, and Cholesky, which currently rely on Base.LinearAlgebra's non-mutating methods for making copies as necessary---this mostly works, but eagerly promotes to Float64 rather than Float32, which is obviously not ideal. Shouldn't be too hard to add the necessary methods to the non-mutating functions, I may get around to it later.