JuliaSmoothOptimizers / JuliaSmoothOptimizers/ADNLPModels.jl
SparseReverseADHessian records its ReverseDiff tape on uninitialised memory
- Dominant language
- Julia
- Stars
- 44
- Forks
- 19
- PR merge metrics
- No merged PRs in 30d
Description
# `SparseReverseADHessian` records its ReverseDiff tape on uninitialised memory
## Summary
`SparseReverseADHessian` (and the `:optimized` predefined backend that uses it)
allocates its tape input as `Vector{ForwardDiff.Dual{…}}(undef, nvar)` and then
runs the model function on it while building the ReverseDiff tape. Uninitialised
memory can legally hold non-finite `Float64`s (e.g. when the allocator hands back
a just-freed `-Inf` bounds array). If the objective/constraints contain a
domain-restricted primitive — `cos`/`sin` (via `Base.sincos`), `sqrt`, `log`, … —
tape construction throws instead of returning a backend:
```
DomainError with -Inf:
sincos(x) is only defined for finite x.
```
The model is never solved; `ADNLPModel!` / the backend constructor itself fails.
## Where
`src/sparse_hessian.jl`, both `SparseReverseADHessian` methods
(`v0.8.13`, lines ~193 and ~207):
```julia
z = Vector{ForwardDiff.Dual{tagf, T, 1}}(undef, nvar) # not initialised
f_tape = ReverseDiff.GradientTape(f, z) # executes f(z)
...
zψ = Vector{ForwardDiff.Dual{tagψ, T, 1}}(undef, nvar) # not initialised
yψ = fill!(similar(zψ, ncon), zero(T))
cfgψ = ReverseDiff.compile(ReverseDiff.GradientTape(ψ, (zψ, yψ), ψ_tape)) # executes ψ(zψ, yψ)
```
`SparseADHessian` (forward mode) is **not** affected: `ForwardDiff.GradientConfig`
only sizes buffers, it never runs the function.
## Reproducer
`ADNLPModels 0.8.13`, `ForwardDiff 1.4.5`, `ReverseDiff 1.17.0`, Julia 1.12.
```julia
using ADNLPModels, ForwardDiff, ReverseDiff, LinearAlgebra, SparseArrays
# min sum(x) s.t. cos(x[1]) - 1 = 0 , x ∈ Rⁿ (x unbounded: no finite fallback)
f(x) = sum(x)
c!(cx, x) = (cx[1] = cos(x[1]) - 1; cx)
n, m = 4, 1
Hpat = sparse(1:n, 1:n, trues(n), n, n)
# --- 1. deterministic: the exact internal snippet, with zψ holding a value that
# `undef` legally may hold (here: -Inf, as a recycled `-Inf` bounds array).
ψ = (x, u) -> (tmp = similar(x, length(u)); c!(tmp, x); dot(tmp, u))
tagψ = ForwardDiff.Tag{typeof(ψ), Float64}
zψ = Vector{ForwardDiff.Dual{tagψ, Float64, 1}}(undef, n)
fill!(zψ, ForwardDiff.Dual{tagψ}(-Inf, 0.0))
yψ = fill!(similar(zψ, m), 0.0)
ReverseDiff.GradientTape(ψ, (zψ, yψ), ReverseDiff.GradientConfig((zψ, yψ)))
# -> ERROR: DomainError with -Inf: sincos(x) is only defined for finite x.
# --- 2. through the public constructor, no touching of internals: just make the
# small-object pool contain the -Inf pattern first (as a freed bounds
# array would), then call it. Fires ~80-90% of the time.
NEGINFDUAL = ForwardDiff.Dual{tagψ}(-Inf, -Inf)
poison() = (for _ in 1:5000; b = fill(NEGINFDUAL, n); Base.donotdelete(b); end; GC.gc(); GC.gc())
fails = 0
for _ in 1:100
poison()
try ADNLPModels.SparseReverseADHessian(n, f, m, c!, Hpat) catch; global fails += 1 end
end
@show fails # e.g. 81 / 100
# --- 3. contrast: forward-mode Hessian on the same model, same poisoning: 0 / 100.
fwd = 0
for _ in 1:100
poison()
try ADNLPModels.SparseADHessian(n, f, m, c!, Hpat) catch; global fwd += 1 end
end
@show fwd # 0 / 100
```
## Real-world hit
CTDirect.jl builds NLP variable bounds as `-Inf .* ones(nvar)`; those arrays are
freed and their storage is recycled into the `undef` `zψ` above. A free-final-time
optimal-control problem with a rotation-matrix dynamics (`cos(θ)`, `sin(θ)`, `θ`
unbounded) then fails to build on CI — deterministically on Julia 1.12, flakily on
1.10 — depending only on which dependency versions were resolved that run.
(control-toolbox/CTDirect.jl CI, `test/ci/test_all_ocp.jl` `:moonlander`.)
## Suggested fix
Initialise the tape inputs before recording — from `x0` (already a kwarg,
defaulting to `rand(nvar)`), or at least `fill!` the primal parts with a finite
value:
```julia
z = Vector{ForwardDiff.Dual{tagf, T, 1}}(undef, nvar)
@inbounds for i in 1:nvar
z[i] = ForwardDiff.Dual{tagf}(x0[i], zero(T))
end
```
and likewise for `zψ` / `lz`. Recording on a well-defined point also makes tape
compilation deterministic.
Contributor guide
Research direction
Start in src/sparse_hessian.jl at both SparseReverseADHessian methods, then run the provided reproducer with the domain-restricted cos model. Ensure the tape inputs z, zψ, and lz are initialized from the available starting values before tape construction; done means the constructor succeeds reliably with recycled -Inf memory and the existing forward-mode behavior remains unaffected.
Written by the indexing model from the issue text.
Assessment
- Tech stack
- julia
- Domain
- backend-api-design
- Issue type
- Bug
- Difficulty
- 3/5
- Estimated time
- 1-2 days
- Activity status
- Active
- Clarity
- Clearly specified
- Newbie friendliness
- 84/100