Fix num_vecjac

Author: avik-pal

Project.toml

@@ -1,7 +1,7 @@
 name = "SparseDiffTools"
 uuid = "47a9eef4-7e08-11e9-0b38-333d64bd3804"
 authors = ["Pankaj Mishra <pankajmishra1511@gmail.com>", "Chris Rackauckas <contact@chrisrackauckas.com>"]
-version = "2.15.0"
+version = "2.15.1"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"

ext/SparseDiffToolsZygoteExt.jl

@@ -45,7 +45,8 @@ end
 
 ### Jac, Hes products
 
-function numback_hesvec!(dy, f::F, x, v, cache1 = similar(v), cache2 = similar(v), cache3 = similar(v)) where {F}
+function numback_hesvec!(dy, f::F, x, v, cache1 = similar(v), cache2 = similar(v),
+        cache3 = similar(v)) where {F}
     g = let f = f
         (dx, x) -> dx .= first(Zygote.gradient(f, x))
     end

src/differentiation/jaches_products.jl

@@ -32,8 +32,8 @@ function auto_jacvec(f, x, v)
     vec(partials.(vec(f(y)), 1))
 end
 
-function num_jacvec!(dy, f, x, v, cache1 = similar(v), cache2 = similar(v), cache3 = similar(v);
-        compute_f0 = true)
+function num_jacvec!(dy, f, x, v, cache1 = similar(v), cache2 = similar(v),
+        cache3 = similar(v); compute_f0 = true)
     vv = reshape(v, axes(x))
     compute_f0 && (f(cache1, x))
     T = eltype(x)
@@ -134,7 +134,8 @@ function autonum_hesvec(f, x, v)
     partials.(g(Dual{DeivVecTag}.(x, v)), 1)
 end
 
-function num_hesvecgrad!(dy, g, x, v, cache1 = similar(v), cache2 = similar(v), cache3 = similar(v))
+function num_hesvecgrad!(dy, g, x, v, cache1 = similar(v), cache2 = similar(v),
+        cache3 = similar(v))
     T = eltype(x)
     # Should it be min? max? mean?
     ϵ = sqrt(eps(real(T))) * max(one(real(T)), abs(norm(x)))

src/differentiation/vecjac_products.jl

@@ -1,5 +1,5 @@
-function num_vecjac!(du, f::F, x, v, cache1 = similar(v), cache2 = similar(v), cache3 = similar(v);
-        compute_f0 = true) where {F}
+function num_vecjac!(du, f::F, x, v, cache1 = similar(v), cache2 = similar(v),
+        cache3 = similar(x); compute_f0 = true) where {F}
     compute_f0 && (f(cache1, x))
     T = eltype(x)
     # Should it be min? max? mean?
@@ -22,10 +22,11 @@ function num_vecjac(f::F, x, v, f0 = nothing) where {F}
     # Should it be min? max? mean?
     ϵ = sqrt(eps(real(T))) * max(one(real(T)), abs(norm(x)))
     du = similar(x)
-    cache = copy(x)
+    cache = similar(x)
+    copyto!(cache, x)
     for i in 1:length(x)
         cache[i] += ϵ
-        f0 = f(x)
+        f0 = f(cache)
         cache[i] = x[i]
         du[i] = (((f0 .- _f0) ./ ϵ)' * vv)[1]
     end
@@ -93,7 +94,7 @@ function VecJac(f, u::AbstractArray, p = nothing, t = nothing; fu = nothing,
 end
 
 function _vecjac(f::F, fu, u, autodiff::AutoFiniteDiff) where {F}
-    cache = (similar(fu), similar(fu), similar(fu))
+    cache = (similar(fu), similar(fu), similar(u))
     pullback = nothing
     return AutoDiffVJP(f, u, cache, autodiff, pullback)
 end