search

SciML / DiffEqGPU.jl

GPU-acceleration routines for DifferentialEquations.jl and the broader SciML scientific machine learning ecosystem
https://docs.sciml.ai/DiffEqGPU/stable/
MIT License
272 stars 27 forks source link

EnsembleGPUArray() incompatibility with SecondOrderODEProblem() #271

Open shahabahreini opened 1 year ago

shahabahreini commented 1 year ago

Hi, I attempted to solve a second-order differential equation using a Julia library. However, I encountered an issue with the EnsembleGPUArray function, as the generated ensemble is not working as expected. Below, you can find the code that I used:

using DifferentialEquations, DataFrames, CSV, DelimitedFiles, RecursiveArrayTools, DiffEqGPU, CUDA

function lorenz_ddu(ddu, du, u, p, t)
    ddu[1] = p[1]*(du[2]-u[1])
    ddu[2] = u[1]*(p[2]-u[3]) - u[2]
    ddu[3] = u[1]*u[2] - p[3]*u[3]
end

function ensembleOutputExtracter(simulationResult, numTrajectories)
    for j = 1:min(numTrajectories, length(simulationResult))
        solnew = vcat(j, simulationResult[j].t', simulationResult[:, :, j])
        open("shahab.txt", "a+") do file
            writedlm(file, solnew', ',')
        end
    end
end

CUDA.set_runtime_version!("local")

u0 = [1.0, 1.0, 1.0]
du0 = [1.1, 1.2, 1.3]
tspan = (0.0, 100.0)
p = [10.0, 28.0, 8/3]
prob = SecondOrderODEProblem(lorenz_ddu, du0, u0, tspan, p)
prob_func = (prob, i, repeat) -> remake(prob, u0=rand(3).*u0, p=rand(3).*p)
monteprob = EnsembleProblem(prob, prob_func=prob_func, safetycopy=false)
sol = solve(monteprob, Tsit5(), EnsembleGPUArray(CUDA.CUDABackend(), 0.0), save_everystep=false, trajectories=10, saveat=1.0f0)
ensembleOutputExtracter(sol, 10)

Am I missing something? Thanks,

shahabahreini commented 1 year ago

Here is the error:

ERROR: LoadError: type DynamicalODEFunction has no field f
Stacktrace:
 [1] getproperty
   @ ./Base.jl:38 [inlined]
 [2] generate_problem(prob::ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, DynamicalODEFunction{true, SciMLBase.FullSpecialize, ODEFunction{true, SciMLBase.FullSpecialize, typeof(lorenz_ddu), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#283#285", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SecondOrderODEProblem{true}}, u0::CuArray{Float64, 2, CUDA.Mem.DeviceBuffer}, p::CuArray{Float64, 2, CUDA.Mem.DeviceBuffer}, jac_prototype::Nothing, colorvec::Nothing)
   @ DiffEqGPU ~/.julia/packages/DiffEqGPU/4SgCk/src/DiffEqGPU.jl:937
 [3] batch_solve_up(ensembleprob::EnsembleProblem{ODEProblem{ArrayPartition{Float64, Tuple{Vector{Float64}, Vector{Float64}}}, Tuple{Float64, Float64}, true, Vector{Float64}, DynamicalODEFunction{true, SciMLBase.FullSpecialize, ODEFunction{true, SciMLBase.FullSpecialize, typeof(lorenz_ddu), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#283#285", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SecondOrderODEProblem{true}}, var"#3#4", typeof(SciMLBase.DEFAULT_OUTPUT_FUNC), typeof(SciMLBase.DEFAULT_REDUCTION), Nothing}, probs::Vector{ODEProblem{Vector{Float64}, Tuple{Float64, Float64}, true, Vector{Float64}, DynamicalODEFunction{true, SciMLBase.FullSpecialize, ODEFunction{true, SciMLBase.FullSpecialize, typeof(lorenz_ddu), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#283#285", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SecondOrderODEProblem{true}}}, alg::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, ensemblealg::EnsembleGPUArray{CUDABackend}, I::UnitRange{Int64}, u0::Matrix{Float64}, p::Matrix{Float64}; kwargs::Base.Pairs{Symbol, Any, NTuple{4, Symbol}, NamedTuple{(:adaptive, :unstable_check, :save_everystep, :saveat), Tuple{Bool, DiffEqGPU.var"#19#25", Bool, Float32}}})
   @ DiffEqGPU ~/.julia/packages/DiffEqGPU/4SgCk/src/DiffEqGPU.jl:829
 [4] batch_solve(ensembleprob::EnsembleProblem{ODEProblem{ArrayPartition{Float64, Tuple{Vector{Float64}, Vector{Float64}}}, Tuple{Float64, Float64}, true, Vector{Float64}, DynamicalODEFunction{true, SciMLBase.FullSpecialize, ODEFunction{true, SciMLBase.FullSpecialize, typeof(lorenz_ddu), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#283#285", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SecondOrderODEProblem{true}}, var"#3#4", typeof(SciMLBase.DEFAULT_OUTPUT_FUNC), typeof(SciMLBase.DEFAULT_REDUCTION), Nothing}, alg::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, ensemblealg::EnsembleGPUArray{CUDABackend}, I::UnitRange{Int64}, adaptive::Bool; kwargs::Base.Pairs{Symbol, Any, Tuple{Symbol, Symbol, Symbol}, NamedTuple{(:unstable_check, :save_everystep, :saveat), Tuple{DiffEqGPU.var"#19#25", Bool, Float32}}})
   @ DiffEqGPU ~/.julia/packages/DiffEqGPU/4SgCk/src/DiffEqGPU.jl:725
 [5] macro expansion
   @ ./timing.jl:382 [inlined]
 [6] __solve(ensembleprob::EnsembleProblem{ODEProblem{ArrayPartition{Float64, Tuple{Vector{Float64}, Vector{Float64}}}, Tuple{Float64, Float64}, true, Vector{Float64}, DynamicalODEFunction{true, SciMLBase.FullSpecialize, ODEFunction{true, SciMLBase.FullSpecialize, typeof(lorenz_ddu), LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, ODEFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#283#285", LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing, Nothing}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}, SecondOrderODEProblem{true}}, var"#3#4", typeof(SciMLBase.DEFAULT_OUTPUT_FUNC), typeof(SciMLBase.DEFAULT_REDUCTION), Nothing}, alg::Tsit5{typeof(OrdinaryDiffEq.trivial_limiter!), typeof(OrdinaryDiffEq.trivial_limiter!), Static.False}, ensemblealg::EnsembleGPUArray{CUDABackend}; trajectories::Int64, batch_size::Int64, unstable_check::Function, adaptive::Bool, kwargs::Base.Pairs{Symbol, Real, Tuple{Symbol, Symbol}, NamedTuple{(:save_everystep, :saveat), Tuple{Bool, Float32}}})
   @ DiffEqGPU ~/.julia/packages/DiffEqGPU/4SgCk/src/DiffEqGPU.jl:552
 [7] #solve#31
   @ ~/.julia/packages/DiffEqBase/ihYDa/src/solve.jl:956 [inlined]
 [8] top-level scope
   @ ~/julia_pratice/GPU5.jl:27
utkarsh530 commented 1 year ago

No, SecondOrderODEProblem is not supported yet.