search

SciML / LinearSolve.jl

LinearSolve.jl: High-Performance Unified Interface for Linear Solvers in Julia. Easily switch between factorization and Krylov methods, add preconditioners, and all in one interface.
https://docs.sciml.ai/LinearSolve/stable/
Other
244 stars 52 forks source link

LinearSolve + Krylove + CUDA #341

Open ASLabadi opened 1 year ago

ASLabadi commented 1 year ago

Hi guys, i understood from the documentation that I can use LinearSolve with Krylov in both CPU and GPU. i am trying to solve a sparse linear system Ax=b as the following:

Pl_cpu = ilu(A_cpu) # the preconditioner prob_cpu = LinearProblem(A_cpu, b_cpu) sol = solve(prob_cpu, KrylovJL_MINRES(rtol=1e-8),Pl =Pl_cpu, maxiters =200, verbose = true, log = true);

This version works super fine.

But when I try to perform it on the GPU, I convert A, b to CuSparseMatrixCSC.

P_gpul=ilu02(A_gpu) prob_gpu=LinearProblem(A_gpu, b_gpu) sol = solve(prob_gpu, KrylovJL_MINRES(rtol=1e-8),Pl =Pl_gpu, maxiters =200, verbose = true, log = true);

I get an error:

ERROR: MethodError: no method matching ldiv!(::CuSparseMatrixCSR{Float64, Int32}, ::CuArray{Float64, 1, CUDA.Mem.DeviceBuffer})

Closest candidates are: ldiv!(::AbstractVecOrMat, ::AbstractVecOrMat, ::SciMLOperators.AbstractSciMLOperator)
@ SciMLOperators C:\Users\ahmad.julia\packages\SciMLOperators\M8FAe\src\left.jl:27
ldiv!(::Any, ::Sparspak.SpkSparseSolver.SparseSolver{IT, FT}, ::Any) where {IT, FT}
@ Sparspak C:\Users\ahmad.julia\packages\Sparspak\oqBYl\src\SparseCSCInterface\SparseCSCInterface.jl:263 ldiv!(::Any, ::ChainRulesCore.AbstractThunk) @ ChainRulesCore C:\Users\ahmad.julia\packages\ChainRulesCore\0t04l\src\tangent_types\thunks.jl:95 ...

ChrisRackauckas commented 1 year ago

Sparse GPU is something that's still a bit in progress. The current state is https://github.com/SciML/OrdinaryDiffEq.jl/issues/1566

vpuri3 commented 1 year ago
  1. does it work without preconditioner on gpu?
  2. where is the error happening? can you share a little more of the stacktrace?
ChrisRackauckas commented 1 year ago

Last I checked, without a preconditioner should work https://github.com/SciML/OrdinaryDiffEq.jl/issues/1566#issuecomment-1048737240 . It's the preconditioners that cause an issue. They need ldiv! wrapped in CUDA.jl.

amontoison commented 1 year ago

If it can help, I explain how we can do the backward and forward sweeps with IC(0) and ILU(0) preconditioners here.

Note that the sparse CSR format is more relevant if you solve square system with a Krylov method. A * v is faster with this format.