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

Fix default rows_index and cols_index in Jacobian #302

Open
wants to merge 5 commits into
base: master
Choose a base branch
from
Open

Conversation

gdalle
Copy link
Member

@gdalle gdalle commented Oct 7, 2024

Partial fix for SciML/NonlinearSolve.jl#473

julia> x = rand(2, 3);

julia> cartind = vec(CartesianIndices(x))
6-element reshape(::CartesianIndices{2, Tuple{Base.OneTo{Int64}, Base.OneTo{Int64}}}, 6) with eltype CartesianIndex{2}:
 CartesianIndex(1, 1)
 CartesianIndex(2, 1)
 CartesianIndex(1, 2)
 CartesianIndex(2, 2)
 CartesianIndex(1, 3)
 CartesianIndex(2, 3)

julia> rows_index = Base.Iterators.map(first  Tuple, cartind);

julia> cols_index = Base.Iterators.map(last  Tuple, cartind);

julia> collect(rows_index), collect(cols_index)
([1, 2, 1, 2, 1, 2], [1, 1, 2, 2, 3, 3])

@gdalle gdalle changed the title Set default rows_index and cols_index to nothing Fix default rows_index and cols_index in Jacobian Oct 7, 2024
@ChrisRackauckas
Copy link
Member

This isn't quite correct either? It's correct for Matrix, but other structured sparse matrices like banded matrices do still have a non-trivial sparsity pattern, but don't use structurednz. However, in that case these are just sentinals, so I guess it's fine to have this, especially since this is legacy code at this point that should get archived in the very near future.

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

It's correct for Matrix, but other structured sparse matrices like banded matrices do still have a non-trivial sparsity pattern, but don't use structurednz.

The previous version was incorrect too, it basically assumed the matrix was diagonal.
At least with this version, rows_index and cols_index have the same length, which means we'll get a cleaner error with other structured matrices (trying to access a nonexistent element).
To ensure that, I would also like to remove the @inbounds in this package and in the FiniteDiff extension.

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

For testing this it would be great to have code coverage (#303), any idea why we don't?

@ChrisRackauckas
Copy link
Member

No idea what happened to code coverage. That's news to me that it's not working. Though it has always been a bit wonky with test groups.

@ChrisRackauckas
Copy link
Member

The previous version was incorrect too, it basically assumed the matrix was diagonal

Agreed, that's why I still approved it. It's a sentinel in the other cases, so it should be fine? Wanted to make sure that was documented in case anyone looked back at this PR as to why it's not the right value.

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

Alternately we could put nothing as you had suggested earlier. Without code coverage, it's hard to know whether this branch is even hit

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

Apparently we're hitting the rate limit of codecov but that's weird, this is JuliaDiff not SciML. Anyone else pushing updates all the time?
EDIT: no, we need to set an upload token and update to codecov v4.

@ChrisRackauckas
Copy link
Member

Codecov's rate limit on forks is effectively zero. What's weird is that this isn't from a fork though...

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

Fixed in #305

@ChrisRackauckas
Copy link
Member

Looks like this breaks GPU

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

We could use map instead of Iterators.map to allocate the vector of indices, although I'm unsure whether CartesianIndices(::CuArray) lives on the GPU.
Of course it would be inefficient but anything going through that branch at the moment is incorrect anyway.

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

@gdalle
Copy link
Member Author

gdalle commented Oct 7, 2024

What do you want to do here @ChrisRackauckas?
We somehow have to generate two sequences that look like [1, 2, 3, 1, 2, 3] and [1, 1, 2, 2, 3, 3]. If we use lazy iterables, GPU is not happy, but if we allocate, well we allocate.

@gdalle
Copy link
Member Author

gdalle commented Oct 16, 2024

Gentle bump

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

Successfully merging this pull request may close these issues.

2 participants