[CUSOLVER] Interface sparse Cholesky and QR factorizations

[amontoison]

CUSOLVER has sparse Cholesky and QR factorizations but the API is not documented...
Performing these factorizations on GPU is something that many Julia users want.

[maleadt]

As seen on CI: https://buildkite.com/julialang/cuda-dot-jl/builds/4483#018b7ff0-eaf5-41cb-8e60-a5f7e112a3cc

Worker 3 failed running test libraries/cusolver/sparse_factorizations:
Some tests did not pass: 0 passed, 0 failed, 8 errored, 0 broken.
libraries/cusolver/sparse_factorizations: Error During Test at /var/lib/buildkite-agent/builds/gpuci-13/julialang/cuda-dot-jl/test/libraries/cusolver/sparse_factorizations.jl:7
  Got exception outside of a @test
  UndefVarError: `sprand` not defined

[maleadt]

Ah, I guess CI didn't fully run on here because you marked it as 'ready' only after pushing. I'll push a fix.

[shakedregev]

Very nice! One comment that I have is that the functions use British English (as opposed to American English) spelling for words such as analyze/analyse and factorize/factorise. This is sure to spark a heated debate, though I am not sure what currently the default is in Julia. Maybe it would be better to stay away from words that have these differences in function names.

[ytdHuang]

Would it be possible to support the type of CuSparseMatrixCSC ?

[amontoison]

@shakedregev
What I can do is to definite both analyse/analyze and factorise/factorize.
No need to start an heated debate if we have both versions :)

@ytdHuang
It is only possible for the Cholesky decomposition.
The CuSparseMatrixCSC should be identical to the CuSparseMatrixCSR because the matrix is symmetric in the real case. If the matrix is complex Hermitian, we can do the same trick if we conjugate the non-zeros with conj(A.nzVal).

[ytdHuang]

@amontoison
Thank you for the information.
But the problem I would like to solve is complex non-Hermitian, so I have to use QR algorithm.

[amontoison]

I suggest to try a Krylov solver with an incomplete ILU(0) preconditioner in that case.
You can find how to do that here.

You can also solve A'Ax = Ab to recover a complex Hermitian problem.
The product A' * A will return a CuSparseMatrixCSC but the matrix will be more ill-conditioned and will have more non-zeros.

[ytdHuang]

@amontoison
Thank you again for the information~
Due to some reason, the problem I'm solving doesn't satisfy A'Ax = A'b.

[amontoison]

What is the problem that you want to solve?
A'Ax = A'b is an optimality condition for least-squares problems.

[ytdHuang]

Kinda hard to explain ...
But this is our package https://github.com/NCKU-QFort/HierarchicalEOM.jl
We tried this before A'Ax = A'b, but it never gives us correct answers.

[ChrisRackauckas]

Why is this not for example overloading qr(A:: CuSparseMatrixCSC) to do the right thing and instead defining things that are not in the standard linear algebra API?