Make x^y with Dual arguments with zero partials consistent with Real arguments

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"

src/dual.jl

@@ -546,40 +546,41 @@ 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 _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 = vy * ($f)(vx, vy - 1)
-                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,
-            begin
-                v = value(x)
-                expv = ($f)(v, y)
-                if y == zero(y) || iszero(partials(x))
-                    new_partials = zero(partials(x))
-                else
-                    new_partials = partials(x) * y * ($f)(v, y - 1)
-                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

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

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, 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
 
 @testset "Type min/max" begin

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])