How to get the derivative of an interpolated function?

General Usage
, ,
1

I have the following code:

using DataInterpolations

Cp=[0.024384377390186646,0.03811829845979822,0.05731175742000688,0.0803120051369546,0.10526323380995135,0.1312349486430945,
    0.1588590278697544,0.1882995422205645,0.21937359147848456,0.25269290703920616,0.2882162631147021,0.3240355349501659,
    0.3607290350346774,0.39295261085207084,0.4141821251393675,0.4274331262027119,0.43636422894857274,0.4414767530303278,
    0.443412806515125,0.44365385445693295,0.44317448744425714,0.44218641629727234,0.4407405317795181,0.43888039054963845,
    0.4362875461540461,0.4325702155989231,0.4278606704093836,0.4224263790651815,0.4164272616252002,0.40987027258773945,
    0.4027832291503407,0.39518785578723936,0.38719847687832065]
TSR = 2.0:0.25:10.0

cp = CubicSpline(Cp, TSR)

Now I want to get the derivative of the function cp at a given point. How can I achieve that?

EDIT:
I am using now this function, but it does not look very accurate:

function cp_der(λ)
    dλ = 0.0001
    y = cp(λ-dλ)
    y1 = cp(λ+dλ)
    (y1 - y) / 2dλ
end
1 Like
2

Perhaps you could use ForwardDiff for automatic differentiation?

julia> ForwardDiff.derivative(cp, 5)
0.14320579478665185
2 Likes
3

There doesn’t seem to be a full docs page with tutorials or API for DataInterpolations, but the source code has a derivatives.jl. If that function derivative is what it sounds like, I wonder how it performs compared to automatic differentiation.

1 Like
4 
julia> DataInterpolations.derivative(cp,5) - ForwardDiff.derivative(cp, 5)
0.0

# CubicSpline Interpolation
function derivative(A::CubicSpline{<:AbstractVector}, t::Number, iguess)
  i = searchsortedfirstcorrelated(A.t, t, iguess)
  isnothing(i) ? i = length(A.t) - 1 : i -= 1
  i == 0 ? i += 1 : nothing
  dI = -3A.z[i] * (A.t[i + 1] - t)^2 / (6A.h[i + 1]) + 3A.z[i + 1] * (t - A.t[i])^2 / (6A.h[i + 1])
  dC = A.u[i + 1] / A.h[i + 1] - A.z[i + 1] * A.h[i + 1] / 6
  dD = -(A.u[i] / A.h[i + 1] - A.z[i] * A.h[i + 1] / 6)
  dI + dC + dD, i
end
5

Any AD that uses ChainRules.jl will be the same as DataInterpolations.jl.derivative is called in custom forward and reverse diff rules. See: DataInterpolations.jl/ext/DataInterpolationsChainRulesCoreExt.jl at master · PumasAI/DataInterpolations.jl · GitHub

1 Like
6

There are probably some cases where the DataInterpolations derivative call may be more numerically stable than just the ForwardDiff version, but the ForwardDiff one is fine for many of the interpolations in there. We should probably just do a dual overload to make it call the internal derivative function, that can be done in an extension if anyone’s up for it.