I am trying to use the ModelingToolkit together with Unitful to model a coupled system of two ODEs with different units. The following code fails on ODEProblem with a DimensionError unless both units in u0 are the same:
using ModelingToolkit using DifferentialEquations using Plots using Unitful @parameters t k1 k2 @variables y(t) z(t) @derivatives D'~t eqs = [D(y) ~ k1*y, D(z) ~ k2*z] de = ODESystem(eqs) u0 = [y => 1.5u"V", z => 3.1u"N"] p = [k1 => 1.0u"1/s", k2 => 1.0u"1/s"] tspan = (0.0u"s",1e-5u"s") prob = ODEProblem(de,u0,tspan, p, jac=false) # ERROR: DimensionError: V and 3.1 N are not dimensionally compatible. sol = solve(prob, Euler(), dt=1e-6u"s")
I know that heterogeneous variables are an issue, as discussed here: DifferentialEquations.jl and systems with heterogenous units
The problem in that case was the adaptive stepping that can’t calculate a vector norm if the vector has heterogeneous units.
The suggested solution was to use ArrayPartition:
ODE systems with heterogenous units · Issue #147 · SciML/DifferentialEquations.jl · GitHub
It suspect that I am running into a similar problem, even though in my case the code fails at ODEProblem(), not on solve().
What is the best approach to use units in my case?