Sundials.jl icon indicating copy to clipboard operation
Sundials.jl copied to clipboard

sundials doesn't GPU

Open vpuri3 opened this issue 2 years ago • 4 comments

using CUDA, Sundials

N  = 4
f  = (u, p, t) -> zero(u)
u0 = CUDA.rand(N) 
ts = (0f0, 1f0)

prob = ODEProblem(f, u0, ts)
alg  = CVODE_BDF()
sol  = solve(prob, alg)
julia> include("examples/ad/sundials_gpu.jl")                                     [42/1828]
ERROR: LoadError: MethodError: no method matching Ptr{Sundials._generic_N_Vector}(::Vector{
Float32})                                                                                  
Closest candidates are:                                                                    
  Ptr{Sundials._generic_N_Vector}(::Sundials.NVector) at ~/.julia/packages/Sundials/k9hc3/s
rc/nvector_wrapper.jl:29                                                                   
  Ptr{T}() where T at boot.jl:781                                                          
  Ptr{T}(::Union{Int64, UInt64, Ptr}) where T at boot.jl:780                               
Stacktrace:                                                                                
  [1] convert(#unused#::Type{Ptr{Sundials._generic_N_Vector}}, v::Vector{Float32})         
    @ Sundials ~/.julia/packages/Sundials/k9hc3/src/nvector_wrapper.jl:70                  
  [2] __init(prob::ODEProblem{CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}, Tuple{Float32, Fl
oat32}, false, SciMLBase.NullParameters, ODEFunction{false, var"#196#197", UniformScaling{B
ool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Noth
ing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Pairs{Symbol, Uni
on{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::CVODE_BDF{:New
ton, :Dense, Nothing, Nothing}, timeseries::Vector{Any}, ts::Vector{Any}, ks::Vector{Any}; 
verbose::Bool, callback::Nothing, abstol::Float64, reltol::Float64, saveat::Vector{Float64}
, tstops::Vector{Float64}, maxiters::Int64, dt::Nothing, dtmin::Float64, dtmax::Float64, ti
meseries_errors::Bool, dense_errors::Bool, save_everystep::Bool, save_idxs::Nothing, save_o
n::Bool, save_start::Bool, save_end::Bool, dense::Bool, progress::Bool, progress_name::Stri
ng, progress_message::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), save_timeseries::Nothing
, advance_to_tstop::Bool, stop_at_next_tstop::Bool, userdata::Nothing, alias_u0::Bool, kwar
gs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})  
    @ Sundials ~/.julia/packages/Sundials/k9hc3/src/common_interface/solve.jl:232          
  [3] __init(prob::ODEProblem{CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}, Tuple{Float32, Fl
oat32}, false, SciMLBase.NullParameters, ODEFunction{false, var"#196#197", UniformScaling{B
ool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Noth
ing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Pairs{Symbol, Uni
on{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::CVODE_BDF{:New
ton, :Dense, Nothing, Nothing}, timeseries::Vector{Any}, ts::Vector{Any}, ks::Vector{Any})
    @ Sundials ~/.julia/packages/Sundials/k9hc3/src/common_interface/solve.jl:87           
  [4] __solve(prob::ODEProblem{CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}, Tuple{Float32, F
loat32}, false, SciMLBase.NullParameters, ODEFunction{false, var"#196#197", UniformScaling{
Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Not
hing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Pairs{Symbol, Un
ion{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::CVODE_BDF{:Ne
wton, :Dense, Nothing, Nothing}, timeseries::Vector{Any}, ts::Vector{Any}, ks::Vector{Any},
 recompile::Type{Val{true}}; calculate_error::Bool, kwargs::Base.Pairs{Symbol, Union{}, Tup
le{}, NamedTuple{(), Tuple{}}})                                                            
    @ Sundials ~/.julia/packages/Sundials/k9hc3/src/common_interface/solve.jl:14           
  [5] __solve(prob::ODEProblem{CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}, Tuple{Float32, F
loat32}, false, SciMLBase.NullParameters, ODEFunction{false, var"#196#197", UniformScaling{
Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Not
hing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Pairs{Symbol, Un
ion{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, alg::CVODE_BDF{:Ne
wton, :Dense, Nothing, Nothing}, timeseries::Vector{Any}, ts::Vector{Any}, ks::Vector{Any},
 recompile::Type{Val{true}}) (repeats 5 times)                      
    @ Sundials ~/.julia/packages/Sundials/k9hc3/src/common_interface/solve.jl:3            
  [6] solve_call(_prob::ODEProblem{CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}, Tuple{Float3
2, Float32}, false, SciMLBase.NullParameters, ODEFunction{false, var"#196#197", UniformScal
ing{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing,
 Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Pairs{Symbol
, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, args::CVODE_BD
F{:Newton, :Dense, Nothing, Nothing}; merge_callbacks::Bool, kwargshandle::KeywordArgError,
 kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/xOkrD/src/solve.jl:439
  [7] solve_call(_prob::ODEProblem{CuArray{Float32, 1, CUDA.Mem.DeviceBuffer}, Tuple{Float3
2, Float32}, false, SciMLBase.NullParameters, ODEFunction{false, var"#196#197", UniformScal
ing{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing,
 Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Base.Pairs{Symbol
, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SciMLBase.StandardODEProblem}, args::CVODE_BD
F{:Newton, :Dense, Nothing, Nothing})
    @ DiffEqBase ~/.julia/packages/DiffEqBase/xOkrD/src/solve.jl:409

vpuri3 avatar Aug 09 '22 20:08 vpuri3

That is expected. Sundials is a different repository from OrdinaryDiffEq.jl. It's a completely different set of solvers and in fact it's a set of wrappers to the C/Fortran SUNDIALS library. As such, it won't work with Julia GPU kernels. In theory we could wrapper the CUDA N_Vector and pass pointers or something, but at this time it's not in scope for the people developing the library.

ChrisRackauckas avatar Aug 09 '22 20:08 ChrisRackauckas

There is significant improvement in building binaries that can use the CUDA stack, which should pave the path for this if someone wants to work on it.

ViralBShah avatar Jul 15 '23 00:07 ViralBShah

@abillscmu if you are still interested.

@ViralBShah is support for distributed computing with Sundials within scope?

vpuri3 avatar Jul 19 '23 19:07 vpuri3

I believe there is now a recipe for building MPI libraries in BinaryBuilder - so it is certainly possible. Should be a separate JLL whenever we decide to do it.

ViralBShah avatar Jul 19 '23 21:07 ViralBShah