I am trying to estimate the posterior distributions of a couple of parameters of a system of ODEs. The probabilistic model is data ~ MvNormal(y_sim, σ² * I), where data contains sensor measurements of a quantity of interest from several experiments:
exp_data = map((date) -> get_exp_data(date), dates)
CO2_exp = vcat(exp_data...)
and y_sim is computed from get_sim_data which solves the ODE system for inputs from the corresponding experiments and a given 𝜃.
The Turing model is given here:
@everywhere @model function outlet_CO2_fit(data)
#𝑝(𝜃): Priors of parameters 𝜃
γ ~ filldist(truncated(Normal(-5.5, 0.3); lower = -7.0, upper = -4.0),1)
α ~ filldist(truncated(Normal(1.6, 0.05); lower = 1.5, upper = 2.0),1)
σ² ~ filldist(truncated(Normal(0.01, 0.003), 0, 0.1), 1)
𝜃 = [α[1], 10^γ[1]]
outlet_mol_frac_sim = pmap((date) -> get_sim_data(𝜃, date), dates)
y_sim = 1e3.*vcat(outlet_mol_frac_sim...)
data = 1e3.*data
#𝑝(𝒟│𝜃): Likelihood pdf of data 𝒟 given parameters 𝜃
data ~ MvNormal(y_sim, σ²[1] * I)
return nothing
end
model = outlet_CO2_fit(CO2_exp)
advi = ADVI(1, 1000)
opt = Variational.TruncatedADAGrad(0.01, 1.0, 10)
q = vi(model, advi; optimizer = opt)
Now when I run this variational inference model fitting program in a serial manner, the optimization runs without any problem. However, when I run the code with julia -p 4 to distribute the ODE solves over multiple cores, I get this error:
First call to automatic differentiation for the Jacobian
failed. This means that the userffunction is not compatible
with automatic differentiation.
MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{DiffEqBase.OrdinaryDiffEqTag,
tacktrace:
[1] jacobian!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/derivative_wrappers.jl:230
[2] calc_J!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/derivative_utils.jl:169 [inlined]
[3] calc_W!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/derivative_utils.jl:716
[4] update_W!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/derivative_utils.jl:824 [inlined]
[5] update_W!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/derivative_utils.jl:823 [inlined]
[6] nlsolve!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/nlsolve/nlsolve.jl:25
[7] perform_step!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/perform_step/sdirk_perform_step.jl:481
[8] perform_step!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/perform_step/sdirk_perform_step.jl:458 [inlined]
[9] solve!
@ ~/.julia/packages/OrdinaryDiffEq/tAXVW/src/solve.jl:514
[10] macro expansion
@ ./timing.jl:380 [inlined]
[11] #evolve_system#469
@ ~/ToolsCA.jl/src/evolve.jl:51
[12] get_sim_data
@ ~/ToolsCA.jl/validation/param_fitting/mass_conv_fit/dist_fit_frost_2.jl:54
[13] #9
@ ~/ToolsCA.jl/validation/param_fitting/mass_conv_fit/dist_fit_frost_2.jl:80
[14] #110
@ /opt/julia-1.7.3/share/julia/stdlib/v1.7/Distributed/src/process_messages.jl:278
[15] run_work_thunk
@ /opt/julia-1.7.3/share/julia/stdlib/v1.7/Distributed/src/process_messages.jl:63
[16] macro expansion
@ /opt/julia-1.7.3/share/julia/stdlib/v1.7/Distributed/src/process_messages.jl:278 [inlined]
[17] #109
@ ./task.jl:429
Stacktrace:
[1] (::Base.var"#898#900")(x::Task)
@ Base ./asyncmap.jl:177
[2] foreach(f::Base.var"#898#900", itr::Vector{Any})
@ Base ./abstractarray.jl:2712
[3] maptwice(wrapped_f::Function, chnl::Channel{Any}, worker_tasks::Vector{Any}, c::Vector{String})
@ Base ./asyncmap.jl:177
[4] wrap_n_exec_twice
@ ./asyncmap.jl:153 [inlined]
[5] #async_usemap#883
@ ./asyncmap.jl:103 [inlined]
[6] #asyncmap#882
@ ./asyncmap.jl:81 [inlined]
[7] pmap(f::Function, p::WorkerPool, c::Vector{String}; distributed::Bool, batch_size::Int64, on_error::Nothing, retry_delays::Vector{Any}, retry_check::Nothing)
@ Distributed /opt/julia-1.7.3/share/julia/stdlib/v1.7/Distributed/src/pmap.jl:126
[8] pmap(f::Function, p::WorkerPool, c::Vector{String})
@ Distributed /opt/julia-1.7.3/share/julia/stdlib/v1.7/Distributed/src/pmap.jl:101
[9] pmap(f::Function, c::Vector{String}; kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
@ Distributed /opt/julia-1.7.3/share/julia/stdlib/v1.7/Distributed/src/pmap.jl:156
[10] pmap(f::Function, c::Vector{String})
@ Distributed /opt/julia-1.7.3/share/julia/stdlib/v1.7/Distributed/src/pmap.jl:156
[11] outlet_CO2_fit(model::DynamicPPL.Model