Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

sparse Jacobians fail with Julia 1.7 beta3 #315

Closed
sjdaines opened this issue Jul 14, 2021 · 4 comments
Closed

sparse Jacobians fail with Julia 1.7 beta3 #315

sjdaines opened this issue Jul 14, 2021 · 4 comments

Comments

@sjdaines
Copy link
Contributor

Visible as a Sundials.jl package test failure:

julia> versioninfo()
Julia Version 1.7.0-beta3.0
Commit e76c9dad42 (2021-07-07 08:12 UTC)
Platform Info:
  OS: Linux (x86_64-pc-linux-gnu)
  CPU: Intel(R) Core(TM) i7-8650U CPU @ 1.90GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-12.0.0 (ORCJIT, skylake)

(sundialstest) pkg> st
      Status `~/software/julia/sundialstest/Project.toml`
  [c3572dad] Sundials v4.5.1

julia> Pkg.test("Sundials")


Test Summary:    | Pass  Error  Total
Common Interface |   45      1     46
  CVODE          |   25            25
  ARKODE         |    4             4
  IDA            |    7             7
  Jacobians      |    1      1      2
  Callbacks      |    6             6
  Iterator       |              No tests
  Errors         |    1             1
  Mass Matrix    |    1             1

Looks like this might just be that the initial sparse Jacobian construction fails the new SparseMatrixCSC consistency checks introduced in Julia 1.7 ?

Jacobians: Error During Test at /home/sd336/.julia/packages/Sundials/cdlY7/test/runtests.jl:63
  Got exception outside of a @test
  LoadError: ArgumentError: Illegal buffers for SparseMatrixCSC construction 2 [1, 1, 1] [0, 0, 0, 0] [0.0, 0.0, 0.0, 0.0]
  Stacktrace:
    [1] SparseMatrixCSC
      @ /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/SparseArrays/src/sparsematrix.jl:29 [inlined]
    [2] SparseMatrixCSC(m::Int64, n::Int64, colptr::Vector{Int64}, rowval::Vector{Int64}, nzval::Vector{Float64})
      @ SparseArrays /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.7/SparseArrays/src/sparsematrix.jl:44
    [3] convert(#unused#::Type{SparseMatrixCSC}, J::Ptr{Sundials._generic_SUNMatrix})
      @ Sundials ~/.julia/packages/Sundials/cdlY7/src/types_and_consts_additions.jl:54
    [4] cvodejac(t::Float64, u::Ptr{Sundials._generic_N_Vector}, du::Ptr{Sundials._generic_N_Vector}, 
_J::Ptr{Sundials._generic_SUNMatrix}, funjac::Sundials.FunJac{ODEFunction{true, typeof(Lotka), 
LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, typeof(Lotka_jac), Nothing, Nothing, SparseMatrixCSC{Float64, Int64}, 
SparseMatrixCSC{Float64, Int64}, Nothing, Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), 
Nothing}, Nothing, typeof(Lotka_jac), SciMLBase.NullParameters, Nothing, SparseMatrixCSC{Float64, Int64}, Vector{Float64}, 
Nothing, Nothing, Nothing}, tmp1::Ptr{Sundials._generic_N_Vector}, tmp2::Ptr{Sundials._generic_N_Vector}, 
tmp3::Ptr{Sundials._generic_N_Vector})
      @ Sundials ~/.julia/packages/Sundials/cdlY7/src/common_interface/function_types.jl:68
    [5] CVode
      @ ~/.julia/packages/Sundials/cdlY7/src/API/cvodes.jl:166 [inlined]
    [6] CVode(cvode_mem::Sundials.Handle{Sundials.CVODEMem}, tout::Float64, yout::Vector{Float64}, tret::Vector{Float64}, itask::Int32)
@ChrisRackauckas
Copy link
Member

JuliaLang/julia#40523 I guess might be that? It might be related to the 1-based indexing correction we do?

https://github.com/SciML/Sundials.jl/blob/master/src/types_and_consts_additions.jl#L48-L53
https://github.com/SciML/Sundials.jl/blob/master/src/common_interface/function_types.jl#L67-L78

that's all of the code involved in the process in case you have any ideas.

@sjdaines
Copy link
Contributor Author

OK, that was it.

Proximate cause was that the sparse SUNMatrix isn't fully initialised (?).

More fundamental problem is that
Base.convert(::Type{SparseArrays.SparseMatrixCSC}, J::SUNMatrix)
was inherently unsafe as an in-place 1-based indexing correction creates an invalid SUNMatrix
(the Sundials.jl cvodejac etc did overwrite again with zero-based indices, but this presumably could cause invalid memory access in the Sundials library if there were code paths that didn't subsequently
undo the correction, and that could happen in external code as this was in Base).

I've created a PR that replaces the problematic Base.convert with a Base.copyto!
#316
Tests now pass on Julia 1.6 and 1.7.

ChrisRackauckas pushed a commit that referenced this issue Jul 14, 2021
…ix) with copyto!

Fix for #315

Proximate cause was that the SUNMatrix sparse Jacobian was not fully initialised,
hence failed the SparseMatrixCSC consistency checks introduces with
Julia 1.7

More fundamentally, given that Sundials uses zero-based index arrays,
convert(::Type{SparseArrays.SparseMatrixCSC}, J::SUNMatrix)
would either need to allocate, or is inherently unsafe as an in-place
version that applied a 1-based index conversion would create an
invalid SUNMatrix.
Also the convert code here looks like it modified colptr but
not rowptr, so couldn't have worked anyway... (cvodejac and similar
used to overwrite this, which is why those calls worked).
@ChrisRackauckas
Copy link
Member

Thanks! Should all be good on latest patches now.

@sjdaines
Copy link
Contributor Author

Confirm looks good with v4.5.2. Thank you!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants