JuliaDiff / JuliaDiff/ForwardDiff.jl

`gradient!` leaves part of the result unwritten when `f` does not depend on its argument

Open
#847 0 comments 0 reactions 0 assignees View on GitHub

Nobody has claimed this yet.

Dominant language
Julia
Stars
1k
Forks
160
PR merge metrics
No merged PRs in 30d

Description

When f does not depend on its argument but carries a perturbation from an enclosing differentiation, ForwardDiff.gradient! writes only part of the result and leaves the rest exactly as it found it. The gradient is zero in this situation, so every entry should be written.

using ForwardDiff, DiffResults
x = [1.0, 2.0, 3.0]

ForwardDiff.derivative(1.0) do a
    out = fill(a * 111.0, 3)                    # a buffer with recognisable junk
    ForwardDiff.gradient!(out, z -> a * 2.0, x)
    @show ForwardDiff.value.(out)               # [0.0, 111.0, 111.0], should be [0.0, 0.0, 0.0]
    return zero(a)
end

The allocating form has the same defect against the buffer it just allocated, so it returns uninitialized memory. That is often zero, which hides the bug, but not always:

julia> ForwardDiff.derivative(1.0) do a
           @show ForwardDiff.gradient(z -> a * 2.0, x)
           return zero(a)
       end
gradient = Dual{…}[Dual(0.0,0.0), Dual(2.3611913415e-314,2.361521164e-314), Dual(5.0e-324,2.260386021e-314)]

The DiffResult form errors instead:

ForwardDiff.derivative(1.0) do a
    r = DiffResults.GradientResult(fill(a * 0.0, 3))
    ForwardDiff.gradient!(r, z -> a * 2.0, x)
    return zero(a)
end
# ERROR: MethodError: Cannot `convert` an object of type Dual{…} to an object of type Float64

The StaticArray path is correct and returns [0.0, 0.0, 0.0], because its @generated extract_gradient is built from length(x) and uses the indexed partials(T, y, i).

Cause

src/gradient.jl#L65-L72:

extract_gradient!(::Type{T}, result::AbstractArray, y::Real) where {T} = fill!(result, zero(y))
function extract_gradient!(::Type{T}, result::AbstractArray, dual::Dual) where {T}
    idxs = structural_eachindex(result)
    for (i, idx) in zip(1:npartials(dual), idxs)
        result[idx] = partials(T, dual, i)
    end
    return result
end

A Dual{S,V,N} with S ≺ T carries no T-perturbation, so its gradient is zero everywhere and the ::Real method above is what it needs — but it is a Dual, so it takes the second method, where npartials reports N for the wrong layer. Whenever that N is smaller than length(x) the tail of the result is never written.

L59-L63 shares the dispatch and then passes partials(T, dual) to DiffResults.gradient!; for this case that is zero(dual), a scalar Dual rather than a Partials, which is the same shape problem as in #846.

Suggested fix

Dispatch on Dual{T} rather than on Dual, so a value carrying no T-perturbation falls through to the ::Real method, which already fills the whole result with zeros.

Version: ForwardDiff master (v1.4.5), Julia 1.12.7, StaticArrays 1.9.19.

Contributor guide

No contributing guide indexed for this repository

First steps

  1. Read the whole issue, then the project's contributing guide.
  2. Comment on the issue to say you are picking it up — it saves two people doing the same work.
  3. Fork the repository and make your change on a branch.
  4. Open a pull request that references the issue number.

Research direction

Start in src/gradient.jl around extract_gradient! and inspect how nested Dual values are dispatched for array and DiffResults gradient extraction. Reproduce the provided nested-differentiation examples, then verify that gradients write zero to every entry and that the DiffResult form completes without a conversion error. Add or run focused regression tests for these cases.

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
78/100

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.