Change existing scalar rule tests to use FiniteDifferences and correct the incorrect acshc asec acsc, asecd acscd rules

Author: oxinabox

src/rulesets/Base/base.jl

@@ -4,54 +4,70 @@
 @scalar_rule(log2(x), inv(x) / log(oftype(x, 2)))
 @scalar_rule(log1p(x), inv(x + 1))
 @scalar_rule(expm1(x), exp(x))
+
 @scalar_rule(sin(x), cos(x))
 @scalar_rule(cos(x), -sin(x))
 @scalar_rule(sinpi(x), π * cospi(x))
 @scalar_rule(cospi(x), -π * sinpi(x))
 @scalar_rule(sind(x), (π / oftype(x, 180)) * cosd(x))
 @scalar_rule(cosd(x), -(π / oftype(x, 180)) * sind(x))
+
 @scalar_rule(asin(x), inv(sqrt(1 - x^2)))
 @scalar_rule(acos(x), -inv(sqrt(1 - x^2)))
 @scalar_rule(atan(x), inv(1 + x^2))
-@scalar_rule(asec(x), inv(abs(x) * sqrt(x^2 - 1)))
-@scalar_rule(acsc(x), -inv(abs(x) * sqrt(x^2 - 1)))
+@scalar_rule(asec(x::Real), inv(abs(x) * sqrt(x^2 - 1)))
+@scalar_rule(asec(x), inv(x^2 * sqrt(1 - x^-2)))
+@scalar_rule(acsc(x::Real), -inv(abs(x) * sqrt(x^2 - 1)))
+@scalar_rule(acsc(x), -inv(x^2 * sqrt(1 - x^-2)))
 @scalar_rule(acot(x), -inv(1 + x^2))
+
 @scalar_rule(asind(x), oftype(x, 180) / π / sqrt(1 - x^2))
 @scalar_rule(acosd(x), -oftype(x, 180) / π / sqrt(1 - x^2))
 @scalar_rule(atand(x), oftype(x, 180) / π / (1 + x^2))
-@scalar_rule(asecd(x), oftype(x, 180) / π / abs(x) / sqrt(x^2 - 1))
-@scalar_rule(acscd(x), -oftype(x, 180) / π / abs(x) / sqrt(x^2 - 1))
+@scalar_rule(asecd(x::Real), oftype(x, 180) / π / abs(x) / sqrt(x^2 - 1))
+@scalar_rule(asecd(x), oftype(x, 180) / π / x^2 / sqrt(1 - x^-2))
+@scalar_rule(acscd(x::Real), -oftype(x, 180) / π / abs(x) / sqrt(x^2 - 1))
+@scalar_rule(acscd(x), -oftype(x, 180) / π / x^2 / sqrt(1 - x^-2))
 @scalar_rule(acotd(x), -oftype(x, 180) / π / (1 + x^2))
+
 @scalar_rule(sinh(x), cosh(x))
 @scalar_rule(cosh(x), sinh(x))
 @scalar_rule(tanh(x), sech(x)^2)
 @scalar_rule(coth(x), -(csch(x)^2))
+
 @scalar_rule(asinh(x), inv(sqrt(x^2 + 1)))
 @scalar_rule(acosh(x), inv(sqrt(x^2 - 1)))
 @scalar_rule(atanh(x), inv(1 - x^2))
 @scalar_rule(asech(x), -inv(x * sqrt(1 - x^2)))
-@scalar_rule(acsch(x), -inv(abs(x) * sqrt(1 + x^2)))
+@scalar_rule(acsch(x::Real), -inv(abs(x) * sqrt(1 + x^2)))
+@scalar_rule(acsch(x), -inv(x^2 * sqrt(1 + x^-2)))
 @scalar_rule(acoth(x), inv(1 - x^2))
+
 @scalar_rule(deg2rad(x), π / oftype(x, 180))
 @scalar_rule(rad2deg(x), oftype(x, 180) / π)
+
 @scalar_rule(conj(x), Wirtinger(Zero(), One()))
 @scalar_rule(adjoint(x), Wirtinger(Zero(), One()))
 @scalar_rule(transpose(x), One())
+
 @scalar_rule(abs(x), sign(x))
 @scalar_rule(rem2pi(x, r::RoundingMode), (One(), DNE()))
+
 @scalar_rule(+(x), One())
 @scalar_rule(-(x), -1)
 @scalar_rule(+(x, y), (One(), One()))
 @scalar_rule(-(x, y), (One(), -1))
 @scalar_rule(/(x, y), (inv(y), -(x / y / y)))
 @scalar_rule(\(x, y), (-(y / x / x), inv(x)))
 @scalar_rule(^(x, y), (y * x^(y - 1), Ω * log(x)))
+
 @scalar_rule(inv(x), -abs2(Ω))
 @scalar_rule(sqrt(x), inv(2 * Ω))
 @scalar_rule(cbrt(x), inv(3 * Ω^2))
 @scalar_rule(exp(x), Ω)
 @scalar_rule(exp2(x), Ω * log(oftype(x, 2)))
 @scalar_rule(exp10(x), Ω * log(oftype(x, 10)))
+
 @scalar_rule(tan(x), 1 + Ω^2)
 @scalar_rule(sec(x), Ω * tan(x))
 @scalar_rule(csc(x), -Ω * cot(x))
@@ -62,9 +78,11 @@
 @scalar_rule(cotd(x), -(π / oftype(x, 180)) * (1 + Ω^2))
 @scalar_rule(sech(x), -tanh(x) * Ω)
 @scalar_rule(csch(x), -coth(x) * Ω)
+
 @scalar_rule(hypot(x, y), (x / Ω, y / Ω))
 @scalar_rule(sincos(x), @setup((sinx, cosx) = Ω), cosx, -sinx)
 @scalar_rule(atan(x, y), @setup(u = x^2 + y^2), (y / u, -x / u))
+
 @scalar_rule(max(x, y), @setup(gt = x > y), (gt, !gt))
 @scalar_rule(min(x, y), @setup(gt = x > y), (!gt, gt))
 @scalar_rule(mod(x, y), @setup((u, nan) = promote(x / y, NaN16)),

test/rulesets/Base/base.jl

@@ -1,66 +1,56 @@
-function test_scalar(f, f′, xs...)
-    for r = (rrule, frule)
-        rr = r(f, xs...)
-        @test rr !== nothing
-        fx, ∂x = rr
-        @test fx == f(xs...)
-        @test ∂x(1) ≈ f′(xs...) atol=1e-5
-    end
-end
-
 @testset "base" begin
     @testset "Trig" begin
-        @testset "Basics" for x = (Float64(π), Complex(π, π/2))
-            test_scalar(sin, cos, x)
-            test_scalar(cos, x -> -sin(x), x)
-            test_scalar(tan, x -> 1 + tan(x)^2, x)
-            test_scalar(sec, x -> sec(x) * tan(x), x)
-            test_scalar(csc, x -> -csc(x) * cot(x), x)
-            test_scalar(cot, x -> -1 - cot(x)^2, x)
-            test_scalar(sinpi, x -> π * cospi(x), x)
-            test_scalar(cospi, x -> -π * sinpi(x), x)
+        @testset "Basics" for x = (Float64(π)-0.01, Complex(π, π/2))
+            test_scalar(sin, x)
+            test_scalar(cos, x)
+            test_scalar(tan, x)
+            test_scalar(sec, x)
+            test_scalar(csc, x)
+            test_scalar(cot, x)
+            test_scalar(sinpi, x)
+            test_scalar(cospi, x)
         end
-        @testset "Hyperbolic" for x = (Float64(π), Complex(π, π/2))
-            test_scalar(sinh, cosh, x)
-            test_scalar(cosh, sinh, x)
-            test_scalar(tanh, x -> sech(x)^2, x)
-            test_scalar(sech, x -> -tanh(x) * sech(x), x)
-            test_scalar(csch, x -> -coth(x) * csch(x), x)
-            test_scalar(coth, x -> -csch(x)^2, x)
+        @testset "Hyperbolic" for x = (Float64(π)-0.01, Complex(π-0.01, π/2))
+            test_scalar(sinh, x)
+            test_scalar(cosh, x)
+            test_scalar(tanh, x)
+            test_scalar(sech, x)
+            test_scalar(csch, x)
+            test_scalar(coth, x)
         end
         @testset "Degrees" begin
             x = 45.0
-            test_scalar(sind, x -> (π / 180) * cosd(x), x)
-            test_scalar(cosd, x -> (-π / 180) * sind(x), x)
-            test_scalar(tand, x -> (π / 180) * (1 + tand(x)^2), x)
-            test_scalar(secd, x -> (π / 180) * secd(x) * tand(x), x)
-            test_scalar(cscd, x -> (-π / 180) * cscd(x) * cotd(x), x)
-            test_scalar(cotd, x -> (-π / 180) * (1 + cotd(x)^2), x)
+            test_scalar(sind, x)
+            test_scalar(cosd, x)
+            test_scalar(tand, x)
+            test_scalar(secd, x)
+            test_scalar(cscd, x)
+            test_scalar(cotd, x)
         end
-        @testset "Inverses" for x = (1.0, Complex(1.0, 0.25))
-            test_scalar(asin, x -> 1 / sqrt(1 - x^2), x)
-            test_scalar(acos, x -> -1 / sqrt(1 - x^2), x)
-            test_scalar(atan, x -> 1 / (1 + x^2), x)
-            test_scalar(asec, x -> 1 / (abs(x) * sqrt(x^2 - 1)), x)
-            test_scalar(acsc, x -> -1 / (abs(x) * sqrt(x^2 - 1)), x)
-            test_scalar(acot, x -> -1 / (1 + x^2), x)
+        @testset "Inverses" for x = (0.5, Complex(0.5, 0.25))
+            test_scalar(asin, x)
+            test_scalar(acos, x)
+            test_scalar(atan, x)
+            test_scalar(asec, 1/x)
+            test_scalar(acsc, 1/x)
+            test_scalar(acot, 1/x)
         end
-        @testset "Inverse hyperbolic" for x = (0.0, Complex(0.0, 0.25))
-            test_scalar(asinh, x -> 1 / sqrt(x^2 + 1), x)
-            test_scalar(acosh, x -> 1 / sqrt(x^2 - 1), x + 1)  # +1 accounts for domain
-            test_scalar(atanh, x -> 1 / (1 - x^2), x)
-            test_scalar(asech, x -> -1 / x / sqrt(1 - x^2), x)
-            test_scalar(acsch, x -> -1 / abs(x) / sqrt(1 + x^2), x)
-            test_scalar(acoth, x -> 1 / (1 - x^2), x + 1)
+        @testset "Inverse hyperbolic" for x = (0.5, Complex(0.5, 0.25))
+            test_scalar(asinh, x)
+            test_scalar(acosh, x + 1)  # +1 accounts for domain
+            test_scalar(atanh, x)
+            test_scalar(asech, x)
+            test_scalar(acsch, x)
+            test_scalar(acoth, x + 1)
         end
         @testset "Inverse degrees" begin
-            x = 1.0
-            test_scalar(asind, x -> 180 / π / sqrt(1 - x^2), x)
-            test_scalar(acosd, x -> -180 / π / sqrt(1 - x^2), x)
-            test_scalar(atand, x -> 180 / π / (1 + x^2), x)
-            test_scalar(asecd, x -> 180 / π / abs(x) / sqrt(x^2 - 1), x)
-            test_scalar(acscd, x -> -180 / π / abs(x) / sqrt(x^2 - 1), x)
-            test_scalar(acotd, x -> -180 / π / (1 + x^2), x)
+            x = 0.5
+            test_scalar(asind, x)
+            test_scalar(acosd, x)
+            test_scalar(atand, x)
+            test_scalar(asecd, 1/x)
+            test_scalar(acscd, 1/x)
+            test_scalar(acotd, 1/x)
         end
         @testset "Multivariate" begin
             x, y = rand(2)

test/test_util.jl

@@ -1,13 +1,43 @@
+using FiniteDifferences, Test
 using FiniteDifferences: jvp, j′vp
+using ChainRules
 
 const _fdm = central_fdm(5, 1)
 
+"""
+    test_scalar(f, x; rtol=1e-9, atol=1e-9, fdm=central_fdm(5, 1), kwargs...)
+
+Given a function `f` with scalar input an scalar output, perform finite differencing checks,
+at input point `x` to confirm that there are correct ChainRules provided.
+
+# Arguments
+- `f`: Function for which the `frule` and `rrule` should be tested.
+- `x`: input at which to evaluate `f` (should generally be set to an arbitary point in the domain).
+
+All keyword arguments except for `fdm` are passed to `isapprox`.
+"""
+function test_scalar(f, x; rtol=1e-9, atol=1e-9, fdm=_fdm, kwargs...)
+    @testset "$f at $x, $(nameof(rule))" for rule in (rrule, frule)
+        res = rule(f, x)
+        @test res !== nothing  # Check the rule was defined
+        fx, ∂x = res
+        @test fx == f(x)  # Check we still get the normal value, right
+
+        # Check that we get the derivative right:
+        @test isapprox(
+            ∂x(1), fdm(f, x);
+            rtol=rtol, atol=atol, kwargs...
+        )
+    end
+end
+
+
 """
     frule_test(f, (x, ẋ)...; rtol=1e-9, atol=1e-9, fdm=central_fdm(5, 1), kwargs...)
 
 # Arguments
 - `f`: Function for which the `frule` should be tested.
-- `x`: input at which to evaluate `f` (should generally be set randomly).
+- `x`: input at which to evaluate `f` (should generally be set to an arbitary point in the domain).
 - `ẋ`: differential w.r.t. `x` (should generally be set randomly).
 
 All keyword arguments except for `fdm` are passed to `isapprox`.
@@ -31,7 +61,7 @@ end
 # Arguments
 - `f`: Function to which rule should be applied.
 - `ȳ`: adjoint w.r.t. output of `f` (should generally be set randomly).
-- `x`: input at which to evaluate `f` (should generally be set randomly).
+- `x`: input at which to evaluate `f` (should generally be set to an arbitary point in the domain).
 - `x̄`: currently accumulated adjoint (should generally be set randomly).
 
 All keyword arguments except for `fdm` are passed to `isapprox`.