使用 fortran 包装的右侧时，Julia 的 DifferentialEquations 包失败

Question

我正在尝试用朱莉娅的方法求解常微分方程组DifferentialEquations。我的 ODE 的右侧是Fortran 90包装的。这是我的朱莉娅代码：

using DifferentialEquations
function rhs(dNdt,N,p,t)
    ccall((:__atmos_MOD_rhs, "./EvolveAtmFort.so"), Cvoid,(Ref{Float64}, Ref{Float64},Ref{Float64}),t,N,dNdt)
end

N0 = [0.0,298.9,0.0562,22.9,0.0166,35.96,0.0,0.0,0.0,0.0]*6.022e23
tspan = [0.0,1.0e6*365.0*24.0*60.0*60.0]
prob = ODEProblem(rhs,N0,tspan)
sol = solve(prob,Rodas5());

这会产生以下长错误，该错误与计算右侧的导数/雅可比有关。下面，我只列出了一些看起来很重要的部分Stacktrace。

MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.TimeGradientWrapper{ODEFunction{true,typeof(rhs),LinearAlgebra.UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Float64,1},DiffEqBase.NullParameters},Float64},Float64,1})
Closest candidates are:
  Float64(::Real, !Matched::RoundingMode) where T<:AbstractFloat at rounding.jl:200
  Float64(::T) where T<:Number at boot.jl:715
  Float64(!Matched::Int8) at float.jl:60
  ...

Stacktrace:
[1] convert(::Type{Float64},...
[2] Base.RefValue{Float64}...
[3] convert(::Type{Ref{Float64}}, ...
[4] cconvert(::Type{T} where T, 
[5] rhs(::Array{ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.TimeGradientWrapper{...
[6] (::ODEFunction{true,typeof(rhs),LinearAlgebra...
[7] (::DiffEqBase.TimeGradientWrapper{ODEFunction{true,typeof(rhs),LinearAlgebra...
[8] derivative!(::Array{Float64,1},...
[9] calc_tderivative!(::OrdinaryDiffEq...
[10] calc_rosenbrock_differentiation! at...
[etc...]

当我使用无雅可比方法时，例如Tsit5()，集成工作得很好。只有需要雅可比计算的方法才会失败。我究竟做错了什么？如何调整 Fortran 包装器以便可以使用隐式方法？谢谢！

Answer 1

这是由于某些隐式求解器中的隐式自动微分造成的。您需要将其关闭，即Rodas5(autodiff=false)。