I you don’t care too much about performance, and easier way to do this would be to just not have a strictly typed array.
rv0 = Vector{Number}()
r0 = [1131.340, -2282.343, 6672.423]u"km"
v0 = [-5.64305, 4.30333, 2.42879]u"km/s"
append!(rv0,r0); append!(rv0,v0)
But I noticed that
Δt = 86400.0*365u"s"
mu = 398600.4418u"km^3/s^2"
prob = ODEProblem((t, y, dy) -> newton!(t, y, dy, mu), rv0, (0.0u"s", Δt))
solve(prob, Vern8())[end]
this errors.
ERROR: DimensionError: 1//1000000 and 1.13134 km are not dimensionally compatible.
in .+(::Rational{Int64}, ::Unitful.Quantity{Float64,Unitful.Dimensions{(Unitful.Dimension{:Length}(1//1),)},Unitful.FreeUnits{(Unitful.Unit{:Meter,Unitful.Dimensions{(Unitful.Dimension{:Length}(1//1),)}}(3,1//1),),Unitful.Dimensions{(Unitful.Dimension{:Length}(1//1),)}}}) at .\operators.jl:152
in collect(::Base.Generator{Array{Unitful.Quantity{Float64,D,U},1},Base.##143#145{Rational{Int64}}}) at .\array.jl:307
in .+(::Rational{Int64}, ::Array{Unitful.Quantity{Float64,D,U},1}) at .\arraymath.jl:90
in +(::Rational{Int64}, ::Array{Unitful.Quantity{Float64,D,U},1}) at .\arraymath.jl:130
in ode_determine_initdt(::Array{Number,1}, ::Unitful.Quantity{Float64,Unitful.Dimensions{(Unitful.Dimension{:Time}(1//1),)},Unitful.FreeUnits{(Unitful.Unit{:Second,Unitful.Dimensions{(Unitful.Dimension{:Time}(1//1),)}}(0,1//1),),Unitful.Dimensions{(Unitful.Dimension{:Time}(1//1),)}}}, ::Float64, ::Unitful.Quantity{Float64,Unitful.Dimensions{(Unitful.Dimension{:Time}(1//1),)},Unitful.FreeUnits{(Unitful.Unit{:Second,Unitful.Dimensions{(Unitful.Dimension{:Time}(1//1),)}}(0,1//1),),Unitful.Dimensions{(Unitful.Dimension{:Time}(1//1),)}}}, ::Rational{Int64}, ::Float64, ::OrdinaryDiffEq.#ODE_DEFAULT_NORM, ::DiffEqBase.ODEProblem{Array{Number,1},Quantity{Float64, Dimensions:{�}, Units:{s}},true,##3#4}, ::Int64) at C:\Users\Chris\.julia\v0.5\OrdinaryDiffEq\src\initdt.jl:5
This makes me think that the previous solution won’t work either. This isn’t going to work because it doesn’t know what to type some of the internal constants, and it usually types them using eltype(u0). In this case, eltype(u0) is not concrete, so it dies. To handle this case, the solvers will need to generalize how it does the typing.
Would you mind opening an issue? It’ll take some work to make heterogenous arrays like this work. And I think that with a proper broadcast (making them an iterator in the right), they can even work well. But it’ll take some work.