From da5557576255261d71d65dba6281a5f681b14f5d Mon Sep 17 00:00:00 2001 From: Andreas Noack Date: Sat, 25 May 2019 21:40:27 +0200 Subject: [PATCH] Introduce absstep to specify the stepsize to use close to zero and rename epsilon_factor to relstep. --- src/derivatives.jl | 34 ++++++++++++++++++++++------------ src/finitediff.jl | 14 +++++++------- src/gradients.jl | 36 ++++++++++++++++++++---------------- src/jacobians.jl | 17 ++++++++++------- test/finitedifftests.jl | 30 +++++++++++++++--------------- 5 files changed, 74 insertions(+), 57 deletions(-) diff --git a/src/derivatives.jl b/src/derivatives.jl index 6aac033..a91d661 100644 --- a/src/derivatives.jl +++ b/src/derivatives.jl @@ -1,10 +1,16 @@ #= Single-point derivatives of scalar->scalar maps. =# -function finite_difference_derivative(f, x::T, fdtype::Type{T1}=Val{:central}, - returntype::Type{T2}=eltype(x), f_x::Union{Nothing,T}=nothing) where {T<:Number,T1,T2} +function finite_difference_derivative( + f, + x::T, + fdtype::Type{T1}=Val{:central}, + returntype::Type{T2}=eltype(x), + f_x::Union{Nothing,T}=nothing; + relstep=default_relstep(fdtype, T), + absstep=relstep) where {T<:Number,T1,T2} - epsilon = compute_epsilon(fdtype, x) + epsilon = compute_epsilon(fdtype, x, relstep, absstep) if fdtype==Val{:forward} return (f(x+epsilon) - f(x)) / epsilon elseif fdtype==Val{:central} @@ -69,10 +75,11 @@ function finite_difference_derivative( returntype :: Type{T2} = eltype(x), # return type of f fx :: Union{Nothing,AbstractArray{<:Number}} = nothing, epsilon :: Union{Nothing,AbstractArray{<:Real}} = nothing; - epsilon_factor = compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2} df = fill(zero(returntype), size(x)) - finite_difference_derivative!(df, f, x, fdtype, returntype, fx, epsilon; epsilon_factor=epsilon_factor) + finite_difference_derivative!(df, f, x, fdtype, returntype, fx, epsilon; relstep=relstep, absstep=absstep) end function finite_difference_derivative!( @@ -83,10 +90,11 @@ function finite_difference_derivative!( returntype :: Type{T2} = eltype(x), fx :: Union{Nothing,AbstractArray{<:Number}} = nothing, epsilon :: Union{Nothing,AbstractArray{<:Real}} = nothing; - epsilon_factor = compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2} cache = DerivativeCache(x, fx, epsilon, fdtype, returntype) - finite_difference_derivative!(df, f, x, cache; epsilon_factor=epsilon_factor) + finite_difference_derivative!(df, f, x, cache; relstep=relstep, absstep=absstep) end function finite_difference_derivative!( @@ -94,11 +102,12 @@ function finite_difference_derivative!( f, x::AbstractArray{<:Number}, cache::DerivativeCache{T1,T2,fdtype,returntype}; - epsilon_factor = compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,fdtype,returntype} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,fdtype,returntype} fx, epsilon = cache.fx, cache.epsilon if typeof(epsilon) != Nothing - @. epsilon = compute_epsilon(fdtype, x, epsilon_factor) + @. epsilon = compute_epsilon(fdtype, x, relstep, absstep) end if fdtype == Val{:forward} if typeof(fx) == Nothing @@ -127,12 +136,13 @@ function finite_difference_derivative!( f, x::StridedArray, cache::DerivativeCache{T1,T2,fdtype,returntype}; - epsilon_factor = compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,fdtype,returntype} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,fdtype,returntype} if fdtype == Val{:forward} fx = cache.fx @inbounds for i ∈ eachindex(x) - epsilon = compute_epsilon(Val{:forward}, x[i], epsilon_factor) + epsilon = compute_epsilon(Val{:forward}, x[i], relstep, absstep) x_plus = x[i] + epsilon if typeof(fx) == Nothing df[i] = (f(x_plus) - f(x[i])) / epsilon @@ -142,7 +152,7 @@ function finite_difference_derivative!( end elseif fdtype == Val{:central} @inbounds for i ∈ eachindex(x) - epsilon = compute_epsilon(Val{:central}, x[i], epsilon_factor) + epsilon = compute_epsilon(Val{:central}, x[i], relstep, absstep) epsilon_double_inv = one(typeof(epsilon)) / (2*epsilon) x_plus, x_minus = x[i]+epsilon, x[i]-epsilon df[i] = (f(x_plus) - f(x_minus)) * epsilon_double_inv diff --git a/src/finitediff.jl b/src/finitediff.jl index cdd07b5..b170324 100644 --- a/src/finitediff.jl +++ b/src/finitediff.jl @@ -6,19 +6,19 @@ Very heavily inspired by Calculus.jl, but with an emphasis on performance and Di Compute the finite difference interval epsilon. Reference: Numerical Recipes, chapter 5.7. =# -@inline function compute_epsilon(::Type{Val{:forward}}, x::T, eps_sqrt=sqrt(eps(real(T)))) where T<:Number - eps_sqrt * max(one(real(T)), abs(x)) +@inline function compute_epsilon(::Type{Val{:forward}}, x::T, relstep::Real, absstep::Real) where T<:Number + return relstep*abs(x) + absstep end -@inline function compute_epsilon(::Type{Val{:central}}, x::T, eps_cbrt=cbrt(eps(real(T)))) where T<:Number - eps_cbrt * max(one(real(T)), abs(x)) +@inline function compute_epsilon(::Type{Val{:central}}, x::T, relstep::Real, absstep::Real) where T<:Number + return relstep*abs(x) + absstep end -@inline function compute_epsilon(::Type{Val{:complex}}, x::T, ::Union{Nothing,T}=nothing) where T<:Real - eps(T) +@inline function compute_epsilon(::Type{Val{:complex}}, x::T, ::Union{Nothing,T}=nothing, ::Union{Nothing,T}=nothing) where T<:Real + return eps(T) end -@inline function compute_epsilon_factor(fdtype::DataType, ::Type{T}) where T<:Number +@inline function default_relstep(fdtype::DataType, ::Type{T}) where T<:Number if fdtype==Val{:forward} return sqrt(eps(real(T))) elseif fdtype==Val{:central} diff --git a/src/gradients.jl b/src/gradients.jl index 5b888ca..205cc04 100644 --- a/src/gradients.jl +++ b/src/gradients.jl @@ -114,7 +114,8 @@ function finite_difference_gradient( fx::Union{Nothing,AbstractArray{<:Number}}=nothing, c1::Union{Nothing,AbstractArray{<:Number}}=nothing, c2::Union{Nothing,AbstractArray{<:Number}}=nothing; - epsilon_factor=compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,T3} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3} if typeof(x) <: AbstractArray df = fill(zero(returntype), size(x)) @@ -133,7 +134,7 @@ function finite_difference_gradient( end end cache = GradientCache(df, x, fdtype, returntype, inplace) - finite_difference_gradient!(df, f, x, cache, epsilon_factor=epsilon_factor) + finite_difference_gradient!(df, f, x, cache, relstep=relstep, absstep=absstep) end function finite_difference_gradient!( @@ -146,24 +147,26 @@ function finite_difference_gradient!( fx::Union{Nothing,AbstractArray{<:Number}}=nothing, c1::Union{Nothing,AbstractArray{<:Number}}=nothing, c2::Union{Nothing,AbstractArray{<:Number}}=nothing; - epsilon_factor=compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,T3} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3} cache = GradientCache(df, x, fdtype, returntype, inplace) - finite_difference_gradient!(df, f, x, cache, epsilon_factor=epsilon_factor) + finite_difference_gradient!(df, f, x, cache, relstep=relstep, absstep=absstep) end function finite_difference_gradient( f, x, cache::GradientCache{T1,T2,T3,fdtype,returntype,inplace}; - epsilon_factor=compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,T3,fdtype,returntype,inplace} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3,fdtype,returntype,inplace} if typeof(x) <: AbstractArray df = fill(zero(returntype), size(x)) else df = zero(cache.c1) end - finite_difference_gradient!(df, f, x, cache, epsilon_factor=epsilon_factor) + finite_difference_gradient!(df, f, x, cache, relstep=relstep, absstep=absstep) df end @@ -174,13 +177,14 @@ function finite_difference_gradient!( f, x::AbstractArray{<:Number}, cache::GradientCache{T1,T2,T3,fdtype,returntype,inplace}; - epsilon_factor=compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,T3,fdtype,returntype,inplace} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3,fdtype,returntype,inplace} # NOTE: in this case epsilon is a vector, we need two arrays for epsilon and x1 # c1 denotes x1, c2 is epsilon fx, c1, c2 = cache.fx, cache.c1, cache.c2 if fdtype != Val{:complex} - @. c2 = compute_epsilon(fdtype, x, epsilon_factor) + @. c2 = compute_epsilon(fdtype, x, relstep, absstep) copyto!(c1,x) end if fdtype == Val{:forward} @@ -246,8 +250,8 @@ function finite_difference_gradient!( f, x::StridedVector{<:Number}, cache::GradientCache{T1,T2,T3,fdtype,returntype,inplace}; - epsilon_factor=compute_epsilon_factor(fdtype, eltype(x)), - ) where {T1,T2,T3,fdtype,returntype,inplace} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3,fdtype,returntype,inplace} # c1 is x1 if we need a complex copy of x, otherwise Nothing # c2 is Nothing @@ -259,7 +263,7 @@ function finite_difference_gradient!( end if fdtype == Val{:forward} for i ∈ eachindex(x) - epsilon = compute_epsilon(fdtype, x[i], epsilon_factor) + epsilon = compute_epsilon(fdtype, x[i], relstep, absstep) x_old = x[i] if typeof(fx) != Nothing x[i] += epsilon @@ -296,7 +300,7 @@ function finite_difference_gradient!( end elseif fdtype == Val{:central} @inbounds for i ∈ eachindex(x) - epsilon = compute_epsilon(fdtype, x[i], epsilon_factor) + epsilon = compute_epsilon(fdtype, x[i], relstep, absstep) x_old = x[i] x[i] += epsilon dfi = f(x) @@ -344,15 +348,15 @@ function finite_difference_gradient!( f, x::Number, cache::GradientCache{T1,T2,T3,fdtype,returntype,inplace}; - epsilon_factor=compute_epsilon_factor(fdtype, eltype(x)) - ) where {T1,T2,T3,fdtype,returntype,inplace} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3,fdtype,returntype,inplace} # NOTE: in this case epsilon is a scalar, we need two arrays for fx1 and fx2 # c1 denotes fx1, c2 is fx2, sizes guaranteed by the cache constructor fx, c1, c2 = cache.fx, cache.c1, cache.c2 if fdtype == Val{:forward} - epsilon = compute_epsilon(Val{:forward}, x, epsilon_factor) + epsilon = compute_epsilon(Val{:forward}, x, relstep, absstep) if inplace == Val{true} f(c1, x+epsilon) else @@ -369,7 +373,7 @@ function finite_difference_gradient!( @. df = (c1 - c2) / epsilon end elseif fdtype == Val{:central} - epsilon = compute_epsilon(Val{:central}, x, epsilon_factor) + epsilon = compute_epsilon(Val{:central}, x, relstep, absstep) if inplace == Val{true} f(c1, x+epsilon) f(c2, x-epsilon) diff --git a/src/jacobians.jl b/src/jacobians.jl index 2ebfe88..cd83bbf 100644 --- a/src/jacobians.jl +++ b/src/jacobians.jl @@ -90,10 +90,11 @@ function finite_difference_jacobian(f, x::AbstractArray{<:Number}, returntype :: Type{T2}=eltype(x), inplace :: Type{Val{T3}}=Val{true}, f_in :: Union{T2,Nothing}=nothing; - epsilon_factor = compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,T3} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3} cache = JacobianCache(x, fdtype, returntype, inplace) - finite_difference_jacobian(f, x, cache, f_in; epsilon_factor=epsilon_factor) + finite_difference_jacobian(f, x, cache, f_in; relstep=relstep, absstep=absstep) end function finite_difference_jacobian( @@ -101,10 +102,11 @@ function finite_difference_jacobian( x, cache::JacobianCache{T1,T2,T3,fdtype,returntype,inplace}, f_in=nothing; - epsilon_factor = compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,T3,fdtype,returntype,inplace} + relstep=default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3,fdtype,returntype,inplace} J = fill(zero(eltype(x)), length(x), length(x)) - finite_difference_jacobian!(J, f, x, cache, f_in; epsilon_factor=epsilon_factor) + finite_difference_jacobian!(J, f, x, cache, f_in; relstep=relstep, absstep=absstep) J end @@ -114,7 +116,8 @@ function finite_difference_jacobian!( x::AbstractArray{<:Number}, cache::JacobianCache{T1,T2,T3,fdtype,returntype,inplace}, f_in::Union{T2,Nothing}=nothing; - epsilon_factor = compute_epsilon_factor(fdtype, eltype(x))) where {T1,T2,T3,fdtype,returntype,inplace} + relstep = default_relstep(fdtype, eltype(x)), + absstep=relstep) where {T1,T2,T3,fdtype,returntype,inplace} m, n = size(J) x1, fx, fx1 = cache.x1, cache.fx, cache.fx1 @@ -123,7 +126,7 @@ function finite_difference_jacobian!( if fdtype == Val{:forward} vfx1 = vec(fx1) @inbounds for i ∈ 1:n - epsilon = compute_epsilon(Val{:forward}, x[i], epsilon_factor) + epsilon = compute_epsilon(Val{:forward}, x[i], relstep, absstep) x1_save = x1[i] x1[i] += epsilon if inplace == Val{true} @@ -148,7 +151,7 @@ function finite_difference_jacobian!( elseif fdtype == Val{:central} vfx1 = vec(fx1) @inbounds for i ∈ 1:n - epsilon = compute_epsilon(Val{:central}, x[i], epsilon_factor) + epsilon = compute_epsilon(Val{:central}, x[i], relstep, absstep) x1_save = x1[i] x_save = x[i] x1[i] += epsilon diff --git a/test/finitedifftests.jl b/test/finitedifftests.jl index 0dcd321..1b56068 100644 --- a/test/finitedifftests.jl +++ b/test/finitedifftests.jl @@ -96,31 +96,31 @@ end @time @testset "Derivative StridedArray f : R -> C tests" begin @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}), df_ref) < 1e-3 - @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, epsilon_factor=sqrt(eps())), df_ref) < 1e-3 + @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, relstep=sqrt(eps())), df_ref) < 1e-3 @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:central}, Complex{eltype(x)}), df_ref) < 1e-6 @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, y), df_ref) < 1e-3 - @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon_factor=sqrt(eps())), df_ref) < 1e-3 + @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, y, relstep=sqrt(eps())), df_ref) < 1e-3 @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:central}, Complex{eltype(x)}, y), df_ref) < 1e-6 @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon), df_ref) < 1e-3 - @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon, epsilon_factor=sqrt(eps())), df_ref) < 1e-3 + @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon, relstep=sqrt(eps())), df_ref) < 1e-3 @test err_func(DiffEqDiffTools.finite_difference_derivative(f, x, Val{:central}, Complex{eltype(x)}, y, epsilon), df_ref) < 1e-6 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}), df_ref) < 1e-3 - @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, epsilon_factor=sqrt(eps())), df_ref) < 1e-3 + @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, relstep=sqrt(eps())), df_ref) < 1e-3 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:central}, Complex{eltype(x)}), df_ref) < 1e-6 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, y), df_ref) < 1e-3 - @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon_factor=sqrt(eps())), df_ref) < 1e-3 + @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, y, relstep=sqrt(eps())), df_ref) < 1e-3 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:central}, Complex{eltype(x)}, y), df_ref) < 1e-6 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon), df_ref) < 1e-3 - @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon, epsilon_factor=sqrt(eps())), df_ref) < 1e-3 + @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:forward}, Complex{eltype(x)}, y, epsilon, relstep=sqrt(eps())), df_ref) < 1e-3 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, Val{:central}, Complex{eltype(x)}, y, epsilon), df_ref) < 1e-6 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, forward_cache), df_ref) < 1e-3 - @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, forward_cache, epsilon_factor=sqrt(eps())), df_ref) < 1e-3 + @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, forward_cache, relstep=sqrt(eps())), df_ref) < 1e-3 @test err_func(DiffEqDiffTools.finite_difference_derivative!(df, f, x, central_cache), df_ref) < 1e-6 end @@ -176,17 +176,17 @@ complex_cache = DiffEqDiffTools.GradientCache(df,x,Val{:complex}) @time @testset "Gradient of f:vector->scalar real-valued tests" begin @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:forward}), df_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:forward}, epsilon_factor=sqrt(eps())), df_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:forward}, relstep=sqrt(eps())), df_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:central}), df_ref) < 1e-8 @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:complex}), df_ref) < 1e-15 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:forward}), df_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:forward}, epsilon_factor=sqrt(eps())), df_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:forward}, relstep=sqrt(eps())), df_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:central}), df_ref) < 1e-8 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:complex}), df_ref) < 1e-15 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, forward_cache), df_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, forward_cache, epsilon_factor=sqrt(eps())), df_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, forward_cache, relstep=sqrt(eps())), df_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, central_cache), df_ref) < 1e-8 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, complex_cache), df_ref) < 1e-15 end @@ -201,15 +201,15 @@ central_cache = DiffEqDiffTools.GradientCache(df,x,Val{:central}) @time @testset "Gradient of f : C^N -> C tests" begin @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:forward}), df_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:forward}, epsilon_factor=sqrt(eps())), df_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:forward}, relstep=sqrt(eps())), df_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_gradient(f, x, Val{:central}), df_ref) < 1e-8 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:forward}), df_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:forward}, epsilon_factor=sqrt(eps())), df_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:forward}, relstep=sqrt(eps())), df_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, Val{:central}), df_ref) < 1e-8 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, forward_cache), df_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, forward_cache, epsilon_factor=sqrt(eps())), df_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, forward_cache, relstep=sqrt(eps())), df_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_gradient!(df, f, x, central_cache), df_ref) < 1e-8 end @@ -316,7 +316,7 @@ f_in = copy(y) @time @testset "Jacobian StridedArray real-valued tests" begin @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache), J_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache, epsilon_factor=sqrt(eps())), J_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache, relstep=sqrt(eps())), J_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache, f_in), J_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, central_cache), J_ref) < 1e-8 @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x), J_ref) < 1e-8 @@ -341,7 +341,7 @@ f_in = copy(y) @time @testset "Jacobian StridedArray f : C^N -> C^N tests" begin @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache), J_ref) < 1e-4 - @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache, epsilon_factor=sqrt(eps())), J_ref) < 1e-4 + @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache, relstep=sqrt(eps())), J_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, forward_cache, f_in), J_ref) < 1e-4 @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x, central_cache), J_ref) < 1e-8 @test err_func(DiffEqDiffTools.finite_difference_jacobian(f, x), J_ref) < 1e-8