ChrisRackauckas
commented
6 months ago We have a derivative overload, but it must only capture and handle the chunk size 1 case. @sathvikbhagavan could you take a look?
Closed Leebre closed 6 months ago
ChrisRackauckas
commented
6 months ago We have a derivative overload, but it must only capture and handle the chunk size 1 case. @sathvikbhagavan could you take a look?
sathvikbhagavan
commented
6 months ago @Leebre, I am not sure why you would do ForwardDiff.gradient? This is a not a vector valued function.
Doing:
julia> grad = ForwardDiff.derivative.((interp,), myvec)
20-element Vector{Float64}:
0.34909318900754527
0.5810773018949766
0.5810773018949766
0.34909318900754527
0.34909318900754527
-0.44018781608409474
0.5810773018949766
-0.28332889046138166
⋮
0.34909318900754527
0.34909318900754527
-0.44018781608409474
0.34909318900754527
0.34909318900754527
0.5810773018949766
0.34909318900754527works. You can also use derivative method in DataInterpolations as well:
julia> DataInterpolations.derivative.((interp,), myvec)
20-element Vector{Float64}:
0.34909318900754527
0.5810773018949766
0.5810773018949766
0.34909318900754527
0.34909318900754527
-0.44018781608409474
0.5810773018949766
-0.28332889046138166
⋮
0.34909318900754527
0.34909318900754527
-0.44018781608409474
0.34909318900754527
0.34909318900754527
0.5810773018949766
0.34909318900754527
Leebre
commented
6 months ago @sathvikbhagavan well ... it can be vector-valued though, in the sense that the function can accept a vector argument and output a vector.
To provide more context for my use case: I am not calling ForwardDiff.gradient directly. I am solving a PDE over a discretized mesh domain using OrdinaryDiffEq, and I am using interpolated functions for several model parameters, which are built into my f-function (derivative function). Those parameters are interpolated based on data. Part of the OrdinaryDiffEq solver process requires computation of the jacobian of the f-function, which (if using AD) requires a vector of dual numbers to be passed through the f-function. So, for that to work, the interpolation functions also need to be able to pass through this vector of dual numbers. I used ForwardDiff.gradient as a simple example in the MWE, since it seems to be doing essentially the same thing and seems to reproduce the same error.
I suppose I could add a for loop in my f-function, to pass each entry of the vector individually to the interpolation function. However, I would be concerned that adding for loops to my inner f-function may have a performance impact on my code.
The interpolation functions in Interpolations.jl do allow for a vector of dual numbers to be passed through. I noticed this issue after switching from Interpolations.jl to DataInterpolations.
Perhaps @ChrisRackauckas could provide some additional comments on whether this functionality should be present?
ChrisRackauckas
commented
6 months ago It needs to be a Jacobian otherwise you get a ForwardDiff error:
using DataInterpolations
using ForwardDiff
# create an interpolated function from some discrete data:
u = rand(5)
t = 0:4
interp = LinearInterpolation(u, t,extrapolate=true)
# ForwardDiff derivative works with a single scalar argument:
grad1 = ForwardDiff.derivative(interp,2.4)
# However, the AD seems to fail, if I try to compute a gradient from a vector input (which my more complex f-function is doing):
myvec = rand(20).*4.0
interp(myvec) # (works)
grad = ForwardDiff.jacobian(interp,myvec) # works!but this needs #230
Leebre
commented
6 months ago @ChrisRackauckas thanks for addressing this. I will give it a try when I get a chance over the next few days!
Describe the bug 🐞
ForwardDiff.gradient seems to fail (AD error) when trying to differentiate an interpolated function using a vector input.
Expected behavior
ForwardDiff should work and output a valid gradient vector.
Minimal Reproducible Example 👇
MWE to duplicate the AD error I'm seeing with DataInterpolations.jl
Error & Stacktrace ⚠️
Environment (please complete the following information):
using Pkg; Pkg.status()using Pkg; Pkg.status(; mode = PKGMODE_MANIFEST)versioninfo()Additional context
The issue was discovered when trying to use an interpolated function inside an f-function being used with OrdinaryDiffEQ (KenCarp4 solver).