ForwardDiff gradient iterated quadqk

Now the problem is here. ForwardDiff propagates a special type of number called Dual through your code to compute derivatives (see Limitations of ForwardDiff · ForwardDiff). It requires your code to be open to such generality, but here you create a vector bfluxvec of Float64. As a result, ForwardDiff cannot fill it with Dual numbers.
The fix here is to make the eltype of bfluxvec depend on the input:

bfluxvec = zeros(eltype(m), N)
1 Like