Method error composing quadqk and forwarddiff

Specific Domains Numerics
,
1

Greetings! Listed below is a small piece of code that composes QuadGK and ForwardDiff and that fails. Including the resulting method error message for the function extract_derivative(). I have two questions regarding this code.

1/ Is there a way to see that the quadqk() call results in the function u(x) = x - 0.5 other then by sampling the output u(x)?

2/ How do I make the ForwardDiff.derivative() call work (i.e. avoid the method error)?

Thx!

using QuadGK
using ForwardDiff

# define two input integrand
integrand(x,y) = x - y

# compute integral by quadrature over second input - results in u(x) = x - 0.5
u(x) = quadgk(y -> integrand(x,y), 0, 1)

# compute derivative - should result in du/dx = 1
# toy model ForwardDiff.derivative(x -> x^2, 1) works just fine
# alternative attempts ForwardDiff.derivative(x -> u(x), 1) or
# ForwardDiff.derivative(x -> u(x), 1) result in same method error
ForwardDiff.derivative(u, 1)

ERROR: MethodError: no method matching extract_derivative(::Type{ForwardDiff.Tag{var"#9#10", Int64}}, ::typeof(u))

Closest candidates are:
extract_derivative(::Type{T}, ::ForwardDiff.Dual) where T at ~/.julia/packages/ForwardDiff/PcZ48/src/derivative.jl:84
extract_derivative(::Type{T}, ::Real) where T at ~/.julia/packages/ForwardDiff/PcZ48/src/derivative.jl:82
extract_derivative(::Type{T}, ::AbstractArray) where T at ~/.julia/packages/ForwardDiff/PcZ48/src/derivative.jl:85

…

Stacktrace:

[1] derivative(f::var"#9#10", x::Int64)

@ ForwardDiff ~/.julia/packages/ForwardDiff/PcZ48/src/derivative.jl:14

[2] top-level scope

@ REPL[27]:1

2

I’m not sure that is the only issue but IIRC quadgk returns a tuple of outputs: the integral and the error estimate. You are only interested in the first one

3

Let me check …

4

Replacing

ForwardDiff.derivative(x -> u(x) ,5.)

by

ForwardDiff.derivative(x -> u(x)[1] ,5.)

seems to do the trick (for now - requires more elaborate check). Thx!

Q: Is there a way to verify the output of the quadqk() call that avoids sampling?

5

what do you mean by verify?

6

By “verify” I mean checking that u(x) = x - 0.5 after the quadgk() call (in the toy example provided above). Maybe something like @code_lowered u(1), but more informative. Not sure that I am making sense here (?).

code_lowered u(1)
CodeInfo(
1 ─ %1 = Main.:(var"#17#18")
│   %2 = Core.typeof(x)
│   %3 = Core.apply_type(%1, %2)
│        #17 = %new(%3, x)
│   %5 = #17
│   %6 = Main.quadgk(%5, 0, 1)
└──      return %6
)
7

not completely ^^ for which x do you want to verify that?

8

I will make some homework and get back here.

Hurray for your input for now!

9

do you mean getting the value along with the derivative? there’s no built in way to do that with ForwardDiff but other packages can offer that option