Skip to content

gradient! writes to the wrong entries for structured inputs #838

Description

@devmotion

extract_gradient! and extract_gradient_chunk! take their target positions from structural_eachindex(result), i.e. from the result container, while the seeds are laid out according to structural_eachindex(x). When the two disagree, the derivatives land in the wrong entries.

julia> using ForwardDiff, LinearAlgebra, DiffResults

julia> x = UpperTriangular([1.0 2.0 3.0; 0.0 4.0 5.0; 0.0 0.0 6.0]);

julia> f(z) = sum(abs2, z) / 2;   # df/dz[i,j] == z[i,j] on the structural entries

julia> ForwardDiff.gradient(f, x)   # correct -- the result is allocated as `similar(x)`
3×3 UpperTriangular{Float64, Matrix{Float64}}:
 1.0  2.0  3.0
     4.0  5.0
         6.0

julia> out = fill(NaN, 3, 3); ForwardDiff.gradient!(out, f, x); out
3×3 Matrix{Float64}:
 1.0  3.0  NaN
 2.0  5.0  NaN
 4.0  6.0  NaN

The six derivatives are written to linear positions 1:6 of the dense result instead of to the upper triangle. No error, no warning, and it happens at every chunk size including vector mode.

With a DiffResult it errors instead, because DiffResults.gradient! does a linear copyto! into the UpperTriangular buffer that GradientResult allocated:

julia> ForwardDiff.gradient!(DiffResults.GradientResult(x), f, x)
ERROR: ArgumentError: cannot set index (2, 1) in the lower triangular part of an UpperTriangular matrix to a nonzero value (2.0)

Same for LowerTriangular and Diagonal. Introduced in #739.

The fix is presumably to derive the positions from x rather than from resultgradient! has x right there. That is what I ended up doing for the Hessian in #837, where the same delegation would otherwise have mis-scattered the gradient stored in a HessianResult (whose gradient buffer is dense even when x is not).

ForwardDiff v1.4.5, Julia 1.12.6.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions