Hi,

I would like to complete the example PDE in the SciML News (also in the ModelingToolKit.jl documentation), section “ModelingToolkit Automatic Parallelism”, on automatically computed multithreaded Jacobian functions via `ModelingToolkit.build_function`

with the objective to speed-up ODE simulations in DifferentialEquations.jl. I follow the steps and obtain the function `multithreadedjac`

. I then go on with

```
#=
previous code, see provided links
...
multithreadedjac = eval(ModelingToolkit.build_function(vec(jac),u,multithread=true)[2])
=#
MyA = zeros(N,N)
AMx = zeros(N,N)
DA = zeros(N,N)
# Define the initial condition as normal arrays
u0 = rand(N,N,3)
#solve
using DifferentialEquations
prob_jac = ODEProblem(ODEFunction(f!,jac = multithreadedjac),u0,(0.0,10.0))
solve(prob_jac)
```

leading to a

`ERROR: MethodError: no method matching (::var"#31#32")(::Array{Float64,2}, ::Array{Float64,3}, ::DiffEqBase.NullParameters, ::Float64)...`

Any assistance on how to successfully translate back the multithreaded Jacobian for applications as DiffEq ODE is much appreciated.

I am using Julia v1.4, DifferentialEquations.jl v6.14.0 and ModelingToolKit.jl v3.10.2