From 03ee831ab5b3e2e6802ab230d9910f98ffca59ef Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Thu, 17 Sep 2026 23:41:24 +0200 Subject: [PATCH 1/4] Divide instead of multiplying by cached reciprocals (#239) The reverse rules for `mean` and elementwise `/` and `\` multiplied the adjoint by a precomputed reciprocal. That reciprocal was formed in the type of the denominator, so an integer or `Float64` denominator capped the gradient at `Float64` accuracy regardless of the deriv type, and `inv(d)` overflowed to `Inf` for subnormal `d`. Both rules now divide. `mean` divides by `length(input)` directly, and the `/` and `\` instructions no longer cache the numerator partial at all: `broadcast_increment_div_deriv!` divides by `value(b)` during propagation, mirroring how `broadcast_mul` already reads `value(b)` in its reverse pass. This drops one cache array per instruction. Scalar accumulations in `propagation.jl` now go through `sum`, whose pairwise reduction is both more accurate and several times faster than the sequential loops they replace. Co-Authored-By: Claude Opus 5 (1M context) --- Project.toml | 2 +- ext/StatisticsExt.jl | 2 +- src/derivatives/elementwise.jl | 37 ++++-------- src/derivatives/propagation.jl | 103 ++++++++++++++++----------------- test/PrecisionTests.jl | 47 +++++++++++++++ test/runtests.jl | 6 ++ 6 files changed, 117 insertions(+), 80 deletions(-) create mode 100644 test/PrecisionTests.jl diff --git a/Project.toml b/Project.toml index 596d5dd1..a24ae60d 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "ReverseDiff" uuid = "37e2e3b7-166d-5795-8a7a-e32c996b4267" -version = "1.18.0-dev" +version = "1.18.0" [deps] ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" diff --git a/ext/StatisticsExt.jl b/ext/StatisticsExt.jl index 727b9a6c..d856448d 100644 --- a/ext/StatisticsExt.jl +++ b/ext/StatisticsExt.jl @@ -14,7 +14,7 @@ end @noinline function ReverseDiff.special_reverse_exec!(instruction::SpecialInstruction{typeof(mean)}) input = instruction.input output = instruction.output - istracked(input) && increment_deriv!(input, inv(length(input)) * deriv(output)) + istracked(input) && increment_deriv!(input, deriv(output) / length(input)) unseed!(output) return nothing end diff --git a/src/derivatives/elementwise.jl b/src/derivatives/elementwise.jl index 7e588c02..7c18c711 100644 --- a/src/derivatives/elementwise.jl +++ b/src/derivatives/elementwise.jl @@ -539,26 +539,18 @@ end # ./ # #----# -numer_partials(d::Real) = Ref(inv(d)) -numer_partials(d::AbstractArray) = broadcast(inv, d) -numer_partials!(out::Ref, d) = (out[] = inv(d); nothing) -numer_partials!(out::AbstractArray, d) = (broadcast!(inv, out, d); nothing) - denom_partials_kernel(n::Real, d::Real) = -(n / (d * d)) denom_partials(n::Real, d::Real) = Ref(denom_partials_kernel(n, d)) denom_partials(n, d) = broadcast(denom_partials_kernel, n, d) denom_partials!(out::Ref, n, d) = (out[] = denom_partials_kernel(n, d); nothing) denom_partials!(out::AbstractArray, n, d) = (broadcast!(denom_partials_kernel, out, n, d); nothing) -rdiv_cache(x, y) = (numer_partials(value(y)), denom_partials(value(x), value(y))) - function broadcast_rdiv(x, y, ::Type{D}) where D tp = tape(x, y) out = track(value(x) ./ value(y), D, tp) - n_partials, d_partials = rdiv_cache(x, y) - cache = (n_partials, d_partials, - index_bound(x, out), index_bound(y, out), - index_bound(n_partials, out), index_bound(d_partials, out)) + d_partials = denom_partials(value(x), value(y)) + cache = (d_partials, index_bound(x, out), index_bound(y, out), + index_bound(d_partials, out)) record!(tp, SpecialInstruction, (broadcast, /), (x, y), out, cache) return out end @@ -566,12 +558,11 @@ end @noinline function special_forward_exec!(instruction::SpecialInstruction{Tuple{typeof(broadcast),typeof(/)}}) a, b = instruction.input a_value, b_value = value(a), value(b) - n_partials, d_partials = instruction.cache + d_partials = first(instruction.cache) output = instruction.output pull_value!(a) pull_value!(b) broadcast!(/, value(output), a_value, b_value) - istracked(a) && numer_partials!(n_partials, b_value) istracked(b) && denom_partials!(d_partials, a_value, b_value) return nothing end @@ -580,9 +571,8 @@ end a, b = instruction.input output = instruction.output output_deriv = deriv(output) - n_partials, d_partials, a_bound, b_bound, - n_partials_bound, d_partials_bound = instruction.cache - istracked(a) && broadcast_increment_deriv!(a, output_deriv, n_partials, a_bound, n_partials_bound) + d_partials, a_bound, b_bound, d_partials_bound = instruction.cache + istracked(a) && broadcast_increment_div_deriv!(a, output_deriv, value(b), a_bound, b_bound) istracked(b) && broadcast_increment_deriv!(b, output_deriv, d_partials, b_bound, d_partials_bound) unseed!(output) return nothing @@ -594,10 +584,9 @@ end function broadcast_ldiv(x, y, ::Type{D}) where D tp = tape(x, y) out = track(value(x) .\ value(y), D, tp) - n_partials, d_partials = rdiv_cache(y, x) - cache = (n_partials, d_partials, - index_bound(x, out), index_bound(y, out), - index_bound(n_partials, out), index_bound(d_partials, out)) + d_partials = denom_partials(value(y), value(x)) + cache = (d_partials, index_bound(x, out), index_bound(y, out), + index_bound(d_partials, out)) record!(tp, SpecialInstruction, (broadcast, \), (x, y), out, cache) return out end @@ -605,12 +594,11 @@ end @noinline function special_forward_exec!(instruction::SpecialInstruction{Tuple{typeof(broadcast),typeof(\)}}) a, b = instruction.input a_value, b_value = value(a), value(b) - n_partials, d_partials = instruction.cache + d_partials = first(instruction.cache) output = instruction.output pull_value!(a) pull_value!(b) broadcast!(\, value(output), a_value, b_value) - istracked(b) && numer_partials!(n_partials, a_value) istracked(a) && denom_partials!(d_partials, b_value, a_value) return nothing end @@ -619,10 +607,9 @@ end a, b = instruction.input output = instruction.output output_deriv = deriv(output) - n_partials, d_partials, a_bound, b_bound, - n_partials_bound, d_partials_bound = instruction.cache + d_partials, a_bound, b_bound, d_partials_bound = instruction.cache istracked(a) && broadcast_increment_deriv!(a, output_deriv, d_partials, a_bound, d_partials_bound) - istracked(b) && broadcast_increment_deriv!(b, output_deriv, n_partials, b_bound, n_partials_bound) + istracked(b) && broadcast_increment_div_deriv!(b, output_deriv, value(a), b_bound, a_bound) unseed!(output) return nothing end diff --git a/src/derivatives/propagation.jl b/src/derivatives/propagation.jl index c680a6a7..6867f4f9 100644 --- a/src/derivatives/propagation.jl +++ b/src/derivatives/propagation.jl @@ -7,10 +7,8 @@ form of normal scalars/arrays, or `Dual` numbers) needs to accounted for as well Often, a function's input and output values are not similarly shaped. To account for these cases, `broadcast` and `reduce` versions of the propagation functions have been implemented. -These special cases should be decently fast, but some might be prone to numerical error if -the accumulated derivative becomes too large compared to the individual terms being added to -it. This can be overcome by using the divide-and-conquer strategy from Base.mapreducedim, -but that strategy is less cache efficient and more complicated to implement. +Accumulations into a scalar derivative go through `sum`, whose pairwise reduction is both more +accurate and faster than a sequential loop. A lot of the code here is pretty repetitive, and only covers the use cases that actually arise from the derivative definitions implemented elsewhere in ReverseDiff. At some point, @@ -81,7 +79,7 @@ end function diffresult_increment_deriv!(input::AbstractArray, x::AbstractArray, results, p::Int) - for i in eachindex(x) + for i in eachindex(x, results) increment_deriv!(input, x[i] * getpartial(results[i], p), i) end return nothing @@ -89,20 +87,20 @@ end function diffresult_increment_deriv!(input::AbstractArray, x::AbstractArray, results, p::Int, bound::CartesianIndex) - for i in CartesianIndices(size(x)) - increment_deriv!(input, x[i] * getpartial(results[i], p), min(bound, i)) + axes(x) == axes(results) || + throw(DimensionMismatch("`x` and `results` must have the same indices")) + for (xi, r, i) in zip(x, results, CartesianIndices(size(x))) + increment_deriv!(input, xi * getpartial(r, p), min(bound, i)) end return nothing end function diffresult_increment_deriv!(input::TrackedReal, x::AbstractArray, results, p::Int, ::Nothing) + inds = eachindex(x, results) + isempty(inds) && return nothing pull_deriv!(input) - input_deriv = input.deriv - for i in eachindex(results) - input_deriv += x[i] * getpartial(results[i], p) - end - input.deriv = input_deriv + input.deriv += sum(i -> x[i] * getpartial(results[i], p), inds) push_deriv!(input) return nothing end @@ -116,45 +114,50 @@ end function broadcast_increment_deriv!(input::AbstractArray, x::AbstractArray, bound::CartesianIndex) - for i in CartesianIndices(size(x)) - increment_deriv!(input, x[i], min(bound, i)) + for (xi, i) in zip(x, CartesianIndices(size(x))) + increment_deriv!(input, xi, min(bound, i)) end return nothing end function broadcast_increment_deriv!(input::TrackedReal, x::AbstractArray, ::Nothing) + isempty(x) && return nothing pull_deriv!(input) - input_deriv = input.deriv - for i in x - input_deriv += i - end - input.deriv = input_deriv + input.deriv += sum(x) push_deriv!(input) return nothing end +# with partials # +#---------------# + +@inline broadcast_increment_deriv!(input, x, partials, input_bound, partials_bound) = + _broadcast_increment_deriv!(*, input, x, partials, input_bound, partials_bound) + +@inline broadcast_increment_div_deriv!(input, x, partials, input_bound, partials_bound) = + _broadcast_increment_deriv!(/, input, x, partials, input_bound, partials_bound) + # with partial array # #--------------------# -function broadcast_increment_deriv!(input::AbstractArray, x::AbstractArray, - partials::AbstractArray, - input_bound::CartesianIndex, - partials_bound::CartesianIndex) - for i in CartesianIndices(size(x)) - current_deriv = x[i] * partials[min(partials_bound, i)] +function _broadcast_increment_deriv!(op::F, input::AbstractArray, x::AbstractArray, + partials::AbstractArray, + input_bound::CartesianIndex, + partials_bound::CartesianIndex) where {F} + for (xi, i) in zip(x, CartesianIndices(size(x))) + current_deriv = op(xi, partials[min(partials_bound, i)]) increment_deriv!(input, current_deriv, min(input_bound, i)) end return nothing end -function broadcast_increment_deriv!(input::TrackedReal, x::AbstractArray, - partials::AbstractArray, ::Nothing, ::CartesianIndex) +function _broadcast_increment_deriv!(op::F, input::TrackedReal, x::AbstractArray, + partials::AbstractArray, ::Nothing, + ::CartesianIndex) where {F} + inds = eachindex(x, partials) + isempty(inds) && return nothing pull_deriv!(input) - input_deriv = input.deriv - for i in eachindex(x) - input_deriv += x[i] * partials[i] - end - input.deriv = input_deriv + input.deriv += sum(i -> op(x[i], partials[i]), inds) push_deriv!(input) return nothing end @@ -162,27 +165,24 @@ end # with partial scalar # #---------------------# -@inline function broadcast_increment_deriv!(input, x, partial::Ref, i, j) - return broadcast_increment_deriv!(input, x, partial[], i, j) +@inline function _broadcast_increment_deriv!(op::F, input, x, partial::Ref, i, j) where {F} + return _broadcast_increment_deriv!(op, input, x, partial[], i, j) end -function broadcast_increment_deriv!(input::AbstractArray, x::AbstractArray, - partial::Real, input_bound::CartesianIndex, - ::Nothing) - for i in CartesianIndices(size(x)) - increment_deriv!(input, x[i] * partial, min(input_bound, i)) +function _broadcast_increment_deriv!(op::F, input::AbstractArray, x::AbstractArray, + partial::Real, input_bound::CartesianIndex, + ::Nothing) where {F} + for (xi, i) in zip(x, CartesianIndices(size(x))) + increment_deriv!(input, op(xi, partial), min(input_bound, i)) end return nothing end -function broadcast_increment_deriv!(input::TrackedReal, x::AbstractArray, - partial::Real, ::Nothing, ::Nothing) +function _broadcast_increment_deriv!(op::F, input::TrackedReal, x::AbstractArray, + partial::Real, ::Nothing, ::Nothing) where {F} + isempty(x) && return nothing pull_deriv!(input) - input_deriv = input.deriv - for i in eachindex(x) - input_deriv += x[i] * partial - end - input.deriv = input_deriv + input.deriv += sum(xi -> op(xi, partial), x) push_deriv!(input) return nothing end @@ -193,19 +193,16 @@ end function broadcast_decrement_deriv!(input::AbstractArray, x::AbstractArray, bound::CartesianIndex) - for i in CartesianIndices(size(x)) - decrement_deriv!(input, x[i], min(bound, i)) + for (xi, i) in zip(x, CartesianIndices(size(x))) + decrement_deriv!(input, xi, min(bound, i)) end return nothing end function broadcast_decrement_deriv!(input::TrackedReal, x::AbstractArray, ::Nothing) + isempty(x) && return nothing pull_deriv!(input) - input_deriv = input.deriv - for i in x - input_deriv -= i - end - input.deriv = input_deriv + input.deriv -= sum(x) push_deriv!(input) return nothing end diff --git a/test/PrecisionTests.jl b/test/PrecisionTests.jl new file mode 100644 index 00000000..08a88a73 --- /dev/null +++ b/test/PrecisionTests.jl @@ -0,0 +1,47 @@ +module PrecisionTests + +using ReverseDiff, Test, LinearAlgebra, Statistics + +include("utils.jl") + +# gradients must be accurate to the working precision (#239) +setprecision(BigFloat, 256) do + x = big.(rand(5)) ./ 3 + n = length(x) + T = eltype(x) + rtol = 8 * eps(T) + + for (f, ∇f) in ((mean, y -> fill(one(T) / n, n)), + (y -> mean(y .^ 3), y -> 3 .* y .^ 2 ./ n), + (y -> sum(y ./ 3), y -> fill(one(T) / 3, n)), + (y -> sum(y ./ 3.0), y -> fill(one(T) / 3, n)), + (y -> sum(3 .\ y), y -> fill(one(T) / 3, n)), + (y -> sum(3.0 .\ y), y -> fill(one(T) / 3, n)), + (sum, y -> fill(one(T), n)), + (prod, y -> prod(y) ./ y), + (norm, y -> y ./ norm(y)), + (var, y -> 2 .* (y .- mean(y)) ./ (n - 1)), + (std, y -> (y .- mean(y)) ./ ((n - 1) * std(y))), + (y -> dot(y, y), y -> 2 .* y), + (y -> sum(3 ./ y), y -> -3 ./ y .^ 2), + (y -> sum(y .\ 3), y -> -3 ./ y .^ 2), + (y -> sum(y ./ (y .+ 1)), y -> 1 ./ (y .+ 1) .^ 2), + (y -> sum(exp.(y)), y -> exp.(y))) + test_println("BigFloat gradients", f) + @test ReverseDiff.gradient(f, x) ≈ ∇f(x) rtol=rtol + end + + # fused broadcasts lose precision in `ForwardDiff._div_partials` (upstream) + let f = y -> sum(y .^ 2 ./ 7) + test_println("BigFloat gradients", f) + @test_broken ReverseDiff.gradient(f, x) ≈ 2 .* x ./ 7 rtol=rtol + end +end + +# a subnormal denominator must not overflow the gradient +let f = y -> 1e-10 * sum(y ./ 1e-310) + test_println("Float64 gradients", f) + @test ReverseDiff.gradient(f, [1.0, 2.0]) == fill(1e-10 / 1e-310, 2) +end + +end # module diff --git a/test/runtests.jl b/test/runtests.jl index 5dc605c5..5da72799 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -30,6 +30,12 @@ test_println(kind, f, pad=" ") = println(pad, "testing $(kind): `$(f)`...") println("done (took $t seconds).") end + @testset "PrecisionTests" begin + println("running PrecisionTests...") + t = @elapsed include("PrecisionTests.jl") + println("done (took $t seconds).") + end + @testset "ScalarTests" begin println("running ScalarTests...") t = @elapsed include(joinpath(TESTDIR, "derivatives/ScalarTests.jl")) From da6a171c1c9d07b742a6cfaea9df29459a3e4959 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Fri, 18 Sep 2026 00:25:41 +0200 Subject: [PATCH 2/4] Divide twice instead of squaring the denominator `denom_partials_kernel` computed `-n / (d * d)`, but `d * d` overflows to `Inf` (and underflows to zero) while `-n / d^2` is still finite and normal, so the denominator adjoint was wrong at the extremes of the range: ReverseDiff.gradient(v -> sum([1e-200] ./ v), [1e-200]) # [-Inf] ReverseDiff.gradient(v -> sum([1e200] ./ v), [1e200]) # [0.0] Dividing twice never forms the square and is the same two roundings, so accuracy in the normal range is unchanged: over 200k random `(n, d)` pairs the mean relative error against a `BigFloat` reference is 4.007e-17 before and 4.017e-17 after. The kernel backs both `./` and `.\`, so one change covers the numerator adjoint's counterpart in each. `PrecisionTests` also moves after `LinAlgTests`, which asserts that `StatisticsExt` is still unloaded when it starts; `PrecisionTests` loads Statistics and so has to run later. Co-Authored-By: Claude Opus 5 (1M context) --- src/derivatives/elementwise.jl | 2 +- test/PrecisionTests.jl | 9 +++++++++ test/runtests.jl | 12 ++++++------ 3 files changed, 16 insertions(+), 7 deletions(-) diff --git a/src/derivatives/elementwise.jl b/src/derivatives/elementwise.jl index 7c18c711..d35a7a7e 100644 --- a/src/derivatives/elementwise.jl +++ b/src/derivatives/elementwise.jl @@ -539,7 +539,7 @@ end # ./ # #----# -denom_partials_kernel(n::Real, d::Real) = -(n / (d * d)) +denom_partials_kernel(n::Real, d::Real) = -(n / d) / d denom_partials(n::Real, d::Real) = Ref(denom_partials_kernel(n, d)) denom_partials(n, d) = broadcast(denom_partials_kernel, n, d) denom_partials!(out::Ref, n, d) = (out[] = denom_partials_kernel(n, d); nothing) diff --git a/test/PrecisionTests.jl b/test/PrecisionTests.jl index 08a88a73..29d777aa 100644 --- a/test/PrecisionTests.jl +++ b/test/PrecisionTests.jl @@ -44,4 +44,13 @@ let f = y -> 1e-10 * sum(y ./ 1e-310) @test ReverseDiff.gradient(f, [1.0, 2.0]) == fill(1e-10 / 1e-310, 2) end +# the denominator adjoint is `-n / d^2`, which must not overflow while it is still finite +for (n, d) in ((1e-200, 1e-200), (1e200, 1e200)) + ∂ = Float64(-big(n) / big(d)^2) + for f in (y -> sum([n] ./ y), y -> sum(y .\ [n])) + test_println("Float64 gradients", f) + @test ReverseDiff.gradient(f, [d]) ≈ [∂] rtol=8*eps(Float64) + end +end + end # module diff --git a/test/runtests.jl b/test/runtests.jl index 5da72799..b53f0e80 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -30,12 +30,6 @@ test_println(kind, f, pad=" ") = println(pad, "testing $(kind): `$(f)`...") println("done (took $t seconds).") end - @testset "PrecisionTests" begin - println("running PrecisionTests...") - t = @elapsed include("PrecisionTests.jl") - println("done (took $t seconds).") - end - @testset "ScalarTests" begin println("running ScalarTests...") t = @elapsed include(joinpath(TESTDIR, "derivatives/ScalarTests.jl")) @@ -48,6 +42,12 @@ test_println(kind, f, pad=" ") = println(pad, "testing $(kind): `$(f)`...") println("done (took $t seconds).") end + @testset "PrecisionTests" begin + println("running PrecisionTests...") + t = @elapsed include(joinpath(TESTDIR, "PrecisionTests.jl")) + println("done (took $t seconds).") + end + @testset "ElementWiseTests" begin println("running ElementwiseTests...") t = @elapsed include(joinpath(TESTDIR, "derivatives/ElementwiseTests.jl")) From f39ee64bb5fd15f0fc4939d29e8eef8e7310ff83 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Fri, 18 Sep 2026 00:42:35 +0200 Subject: [PATCH 3/4] Update version from 1.18.0 to 1.18.1 --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index a24ae60d..5bdf1447 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "ReverseDiff" uuid = "37e2e3b7-166d-5795-8a7a-e32c996b4267" -version = "1.18.0" +version = "1.18.1" [deps] ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" From 3f50fcb7398110cde5427006d744836cdfa3f047 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Fri, 18 Sep 2026 01:51:36 +0200 Subject: [PATCH 4/4] Remove the dead scalar-`Ref` partial machinery MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `numer_partials(d::Real) = Ref(inv(d))` was the last producer of a `Ref` partial, and dividing by `value(b)` directly retired it. Nothing constructs a `Ref` any more: the remaining `Ref` methods of `denom_partials`, `base_partials` and `exp_partials` need both operand values to be `<:Real`, and none of the nine operand combinations generated at `elementwise.jl:419` (or the matching `materialize` overloads in `broadcast.jl`) is scalar/scalar — at least one side is always array-valued. So the `partial::Ref` forwarder is unreachable, as is `_broadcast_increment_deriv!` for a `TrackedReal` input with a scalar partial: when the input is a tracked scalar the other operand is an array, so its partial is an `AbstractArray`. Codecov flagged exactly these two methods on #298. Dropping them is fail-loud rather than silently wrong: no less specific method survives for either signature, so a mistaken deletion is a `MethodError`, never a wrong gradient. `ElementwiseTests` covers all nine combinations. Co-Authored-By: Claude Opus 5 (1M context) --- src/derivatives/elementwise.jl | 6 ------ src/derivatives/propagation.jl | 13 ------------- 2 files changed, 19 deletions(-) diff --git a/src/derivatives/elementwise.jl b/src/derivatives/elementwise.jl index d35a7a7e..130c7aa3 100644 --- a/src/derivatives/elementwise.jl +++ b/src/derivatives/elementwise.jl @@ -540,9 +540,7 @@ end #----# denom_partials_kernel(n::Real, d::Real) = -(n / d) / d -denom_partials(n::Real, d::Real) = Ref(denom_partials_kernel(n, d)) denom_partials(n, d) = broadcast(denom_partials_kernel, n, d) -denom_partials!(out::Ref, n, d) = (out[] = denom_partials_kernel(n, d); nothing) denom_partials!(out::AbstractArray, n, d) = (broadcast!(denom_partials_kernel, out, n, d); nothing) function broadcast_rdiv(x, y, ::Type{D}) where D @@ -618,15 +616,11 @@ end #----# base_partials_kernel(b::Real, e::Real) = e * b^(e - 1) -base_partials(b::Real, e::Real) = Ref(base_partials_kernel(b, e)) base_partials(b, e) = broadcast(base_partials_kernel, b, e) -base_partials!(out::Ref, b, e) = (out[] = base_partials_kernel(b, e); nothing) base_partials!(out::AbstractArray, b, e) = (broadcast!(base_partials_kernel, out, b, e); nothing) exp_partials_kernel(b::Real, e::Real) = log(b) * b^e -exp_partials(b::Real, e::Real) = Ref(exp_partials_kernel(b, e)) exp_partials(b, e) = broadcast(exp_partials_kernel, b, e) -exp_partials!(out::Ref, b, e) = (out[] = exp_partials_kernel(b, e); nothing) exp_partials!(out::AbstractArray, b, e) = (broadcast!(exp_partials_kernel, out, b, e); nothing) function pow_cache(x, y) diff --git a/src/derivatives/propagation.jl b/src/derivatives/propagation.jl index 6867f4f9..f7dd2e55 100644 --- a/src/derivatives/propagation.jl +++ b/src/derivatives/propagation.jl @@ -165,10 +165,6 @@ end # with partial scalar # #---------------------# -@inline function _broadcast_increment_deriv!(op::F, input, x, partial::Ref, i, j) where {F} - return _broadcast_increment_deriv!(op, input, x, partial[], i, j) -end - function _broadcast_increment_deriv!(op::F, input::AbstractArray, x::AbstractArray, partial::Real, input_bound::CartesianIndex, ::Nothing) where {F} @@ -178,15 +174,6 @@ function _broadcast_increment_deriv!(op::F, input::AbstractArray, x::AbstractArr return nothing end -function _broadcast_increment_deriv!(op::F, input::TrackedReal, x::AbstractArray, - partial::Real, ::Nothing, ::Nothing) where {F} - isempty(x) && return nothing - pull_deriv!(input) - input.deriv += sum(xi -> op(xi, partial), x) - push_deriv!(input) - return nothing -end - ############################## # broadcast_decrement_deriv! # ##############################