From 16e82e12310a5112a93e84399fd84bffcc91bf5d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Mon, 12 Oct 2026 15:42:21 +0200 Subject: [PATCH 1/5] Make `x^y` with `Dual` arguments with zero partials consistent with `Real` arguments A `Dual` with zero partials represents the same number as a `Real`, so the `Dual^Dual` branch of `^` and `NaNMath.pow` now uses the same check as the `Dual^Real` branch (`iszero(vy) || isconstant(x)`) and no longer returns `NaN` partials for, e.g., `Dual(0.0, 0.0)^Dual(0.5, 1.0)` and `Dual(0.0, 1.0)^Dual(0.0, 0.0)` (the latter also in NaN-safe mode). The check in the `Dual^Real` branch is rewritten equivalently as `iszero(y) || isconstant(x)` for consistency. Co-Authored-By: Claude Opus 5.5 --- src/dual.jl | 8 ++++++-- test/DerivativeTest.jl | 5 +++++ test/DualTest.jl | 30 ++++++++++++++++++++++++++++++ 3 files changed, 41 insertions(+), 2 deletions(-) diff --git a/src/dual.jl b/src/dual.jl index 6d13dec3..24891b8d 100644 --- a/src/dual.jl +++ b/src/dual.jl @@ -553,7 +553,11 @@ for (f, log) in ((:(Base.:^), :(Base.log)), (:(NaNMath.pow), :(NaNMath.log))) begin vx, vy = value(x), value(y) expv = ($f)(vx, vy) - powval = vy * ($f)(vx, vy - 1) + if iszero(vy) || isconstant(x) + powval = zero(expv) + else + powval = vy * ($f)(vx, vy - 1) + end if isconstant(y) logval = one(expv) elseif iszero(vx) && vy > 0 @@ -567,7 +571,7 @@ for (f, log) in ((:(Base.:^), :(Base.log)), (:(NaNMath.pow), :(NaNMath.log))) begin v = value(x) expv = ($f)(v, y) - if y == zero(y) || iszero(partials(x)) + if iszero(y) || isconstant(x) new_partials = zero(partials(x)) else new_partials = partials(x) * y * ($f)(v, y - 1) diff --git a/test/DerivativeTest.jl b/test/DerivativeTest.jl index d66a7cc7..22b2c521 100644 --- a/test/DerivativeTest.jl +++ b/test/DerivativeTest.jl @@ -93,6 +93,11 @@ end @test (x -> ForwardDiff.derivative(y -> x^y, 0.0))(0.0) === -Inf @test (x -> ForwardDiff.derivative(y -> x^y, 0.5))(0.0) === 0.0 @test (x -> ForwardDiff.derivative(y -> x^y, 1.5))(0.0) === 0.0 + + # base is a `Dual` with zero partials + for f in (^, NaNMath.pow), y in (0.5, 1.5) + @test ForwardDiff.derivative(y -> f(0 * y, y), y) === 0.0 + end end @testset "exponentiation with NaNMath" begin diff --git a/test/DualTest.jl b/test/DualTest.jl index 67d9c9f7..657bd099 100644 --- a/test/DualTest.jl +++ b/test/DualTest.jl @@ -679,6 +679,36 @@ end y = Dual{TestTag}(1.0, 0.0, 1.0); x = Dual{OuterTestTag}(0*y, 0*y); @test iszero(ForwardDiff.partials(ForwardDiff.partials(x^y)[1])) + + # base with zero partials and exponent with the same tag + for f in (^, NaNMath.pow), vy in (0.5, 1.5) + x = Dual{TestTag}(0.0, 0.0, 0.0) + y = Dual{TestTag}(vy, 1.0, 0.0) + @test f(x, y) === Dual{TestTag}(0.0, 0.0, 0.0) + + # nested + x = Dual{OuterTestTag}(Dual{TestTag}(0.0, 0.0), Dual{TestTag}(0.0, 0.0)) + y = Dual{OuterTestTag}(Dual{TestTag}(vy, 1.0), Dual{TestTag}(1.0, 0.0)) + z = f(x, y) + @test iszero(ForwardDiff.value(z)) + @test iszero(ForwardDiff.partials(z)) + + # base is constant with respect to the outer tag but not the inner tag + x = Dual{OuterTestTag}(Dual{TestTag}(0.25, 1.0), Dual{TestTag}(0.0, 0.0)) + y = Dual{OuterTestTag}(Dual{TestTag}(vy, 0.5), Dual{TestTag}(1.0, 0.0)) + z = f(x, y) + @test dual_isapprox(z, exp(y * log(x))) + @test !iszero(ForwardDiff.partials(ForwardDiff.value(z))) + end +end + +@testset "Exponentiation: a `Dual` with zero partials behaves like a `Real`" begin + for f in (^, NaNMath.pow), vx in (0.0, 0.5, 2.0), vy in (-0.5, 0.0, 0.5, 1.5) + x = Dual{TestTag}(vx, 1.0, 2.0) + y = Dual{TestTag}(vy, 1.0, 2.0) + @test isequal(f(Dual{TestTag}(vx, 0.0, 0.0), y), f(vx, y)) + @test isequal(f(x, Dual{TestTag}(vy, 0.0, 0.0)), f(x, vy)) + end end @testset "Type min/max" begin From d05f10a18237d9c24e557aa232c2c9a78f64ba92 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Mon, 12 Oct 2026 16:15:04 +0200 Subject: [PATCH 2/5] Update gradient tests of `x^y` at `(0, 0)` MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit With the fix, `∂x^y/∂x = 0` at `(0, 0)`: in NaN-safe mode the gradient is now `[0.0, -Inf]`, and in the default mode `[NaN, -Inf]` (the remaining `NaN` is caused by the zero partial of `y` in the first direction multiplied by `log(0) = -Inf`, which only NaN-safe mode handles). Co-Authored-By: Claude Opus 5.5 --- test/GradientTest.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/GradientTest.jl b/test/GradientTest.jl index bf121239..7403b22b 100644 --- a/test/GradientTest.jl +++ b/test/GradientTest.jl @@ -151,11 +151,11 @@ end @testset "exponential function at base zero" begin if ForwardDiff.NANSAFE_MODE_ENABLED @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, -0.5]), [-Inf, -Inf]) - @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, 0.0]), [NaN, -Inf]) + @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, 0.0]), [0.0, -Inf]) @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, 0.5]), [Inf, 0.0]) else @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, -0.5]), [NaN, NaN]) - @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, 0.0]), [NaN, NaN]) + @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, 0.0]), [NaN, -Inf]) @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, 0.5]), [Inf, NaN]) end @test isequal(ForwardDiff.gradient(t -> t[1]^t[2], [0.0, 1.5]), [0.0, 0.0]) From 409bbf3c04f4f76c8dfc47256b47818f1ed29f86 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Mon, 12 Oct 2026 16:25:57 +0200 Subject: [PATCH 3/5] Bump version from 1.4.6 to 1.4.7 --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 801792c6..0bc90893 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "ForwardDiff" uuid = "f6369f11-7733-5829-9624-2563aa707210" -version = "1.4.6" +version = "1.4.7" [deps] CommonSubexpressions = "bbf7d656-a473-5ed7-a52c-81e309532950" From 7e0e0436487bd6a16978f3c33e3af9eeaaa08e45 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Wed, 14 Oct 2026 15:34:13 +0200 Subject: [PATCH 4/5] Make `x^y` with zero-partials `Dual`s match `Real` arguments exactly If one argument of `Dual^Dual` has zero partials, defer to the `Dual^Real` or `Real^Dual` method instead of combining the partials with `_mul_partials`, which turned `-0.0` partials into `0.0` with NaN-safe mode disabled. In the `Dual^Real` method, compute the derivative before scaling the partials: `partials(x) * y * f(v, y - 1)` scaled them twice and rounded differently than the `Dual^Dual` method. Co-Authored-By: Claude Opus 5.5 --- src/dual.jl | 19 ++++++------------- test/DualTest.jl | 10 +++++----- 2 files changed, 11 insertions(+), 18 deletions(-) diff --git a/src/dual.jl b/src/dual.jl index 24891b8d..278dddec 100644 --- a/src/dual.jl +++ b/src/dual.jl @@ -551,20 +551,13 @@ for (f, log) in ((:(Base.:^), :(Base.log)), (:(NaNMath.pow), :(NaNMath.log))) @define_binary_dual_op( $f, begin + # a `Dual` with zero partials behaves like a `Real` + isconstant(y) && return ($f)(x, value(y)) + isconstant(x) && return ($f)(value(x), y) vx, vy = value(x), value(y) expv = ($f)(vx, vy) - if iszero(vy) || isconstant(x) - powval = zero(expv) - else - powval = vy * ($f)(vx, vy - 1) - end - if isconstant(y) - logval = one(expv) - elseif iszero(vx) && vy > 0 - logval = zero(vx) - else - logval = expv * ($log)(vx) - end + powval = iszero(vy) ? zero(expv) : vy * ($f)(vx, vy - 1) + logval = (iszero(vx) && vy > 0) ? zero(vx) : expv * ($log)(vx) new_partials = _mul_partials(partials(x), partials(y), powval, logval) return Dual{Txy}(expv, new_partials) end, @@ -574,7 +567,7 @@ for (f, log) in ((:(Base.:^), :(Base.log)), (:(NaNMath.pow), :(NaNMath.log))) if iszero(y) || isconstant(x) new_partials = zero(partials(x)) else - new_partials = partials(x) * y * ($f)(v, y - 1) + new_partials = (y * ($f)(v, y - 1)) * partials(x) end return Dual{Tx}(expv, new_partials) end, diff --git a/test/DualTest.jl b/test/DualTest.jl index 657bd099..9f22afa0 100644 --- a/test/DualTest.jl +++ b/test/DualTest.jl @@ -703,11 +703,11 @@ end end @testset "Exponentiation: a `Dual` with zero partials behaves like a `Real`" begin - for f in (^, NaNMath.pow), vx in (0.0, 0.5, 2.0), vy in (-0.5, 0.0, 0.5, 1.5) - x = Dual{TestTag}(vx, 1.0, 2.0) - y = Dual{TestTag}(vy, 1.0, 2.0) - @test isequal(f(Dual{TestTag}(vx, 0.0, 0.0), y), f(vx, y)) - @test isequal(f(x, Dual{TestTag}(vy, 0.0, 0.0)), f(x, vy)) + for f in (^, NaNMath.pow), vx in (0.0, 0.5, 1.3, 2.0), vy in (-0.5, 0.0, 0.5, 1.5, 2.7) + x = Dual{TestTag}(vx, 0.7, 0.3, -0.0) + y = Dual{TestTag}(vy, 0.7, 0.3, -0.0) + @test isequal(f(Dual{TestTag}(vx, 0.0, 0.0, 0.0), y), f(vx, y)) + @test isequal(f(x, Dual{TestTag}(vy, 0.0, 0.0, 0.0)), f(x, vy)) end end From 571ca6a37ea004fb39f25af667a59527d094ee07 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?David=20M=C3=BCller-Widmann?= Date: Wed, 14 Oct 2026 23:45:15 +0200 Subject: [PATCH 5/5] Fix `x^y` with constant same-tag nested `Dual`s Calling `^` again with `value(y)` (or `value(x)`) in the `Dual^Dual` branch is wrong if the value is itself a `Dual` with the same tag: it dispatches back to the `Dual^Dual` branch and mixes up the nesting levels. Call the `Dual^Real` and `Real^Dual` implementations directly instead. Also make the zero partials in the `Dual^Real` branch have the same type as the non-zero partials, so that a constant `BigInt` exponent does not make `Dual^Dual` type-unstable. Co-Authored-By: Claude Opus 5.5 --- src/dual.jl | 40 ++++++++++++++++++++++------------------ 1 file changed, 22 insertions(+), 18 deletions(-) diff --git a/src/dual.jl b/src/dual.jl index 278dddec..b1bf1cae 100644 --- a/src/dual.jl +++ b/src/dual.jl @@ -546,14 +546,32 @@ end # exponentiation # #----------------# +@inline function _pow_dual_base(f::F, x::Dual{T}, y) where {F,T} + v = value(x) + expv = f(v, y) + if iszero(y) || isconstant(x) + new_partials = zero((y * expv) * partials(x)) + else + new_partials = (y * f(v, y - 1)) * partials(x) + end + return Dual{T}(expv, new_partials) +end + +@inline function _pow_dual_exponent(f::F, log::L, x, y::Dual{T}) where {F,L,T} + v = value(y) + expv = f(x, v) + deriv = (iszero(x) && v > 0) ? zero(expv) : expv*log(oftype(expv, x)) + return Dual{T}(expv, deriv * partials(y)) +end + for (f, log) in ((:(Base.:^), :(Base.log)), (:(NaNMath.pow), :(NaNMath.log))) @eval begin @define_binary_dual_op( $f, begin # a `Dual` with zero partials behaves like a `Real` - isconstant(y) && return ($f)(x, value(y)) - isconstant(x) && return ($f)(value(x), y) + isconstant(y) && return _pow_dual_base($f, x, value(y)) + isconstant(x) && return _pow_dual_exponent($f, $log, value(x), y) vx, vy = value(x), value(y) expv = ($f)(vx, vy) powval = iszero(vy) ? zero(expv) : vy * ($f)(vx, vy - 1) @@ -561,22 +579,8 @@ for (f, log) in ((:(Base.:^), :(Base.log)), (:(NaNMath.pow), :(NaNMath.log))) new_partials = _mul_partials(partials(x), partials(y), powval, logval) return Dual{Txy}(expv, new_partials) end, - begin - v = value(x) - expv = ($f)(v, y) - if iszero(y) || isconstant(x) - new_partials = zero(partials(x)) - else - new_partials = (y * ($f)(v, y - 1)) * partials(x) - end - return Dual{Tx}(expv, new_partials) - end, - begin - v = value(y) - expv = ($f)(x, v) - deriv = (iszero(x) && v > 0) ? zero(expv) : expv*($log)(oftype(expv, x)) - return Dual{Ty}(expv, deriv * partials(y)) - end + _pow_dual_base($f, x, y), + _pow_dual_exponent($f, $log, x, y) ) end end