Something off with broadcast of ArrayPartition with NTuple

[tomaklutfu]

I tried something with OrdinaryDiffEq.jl bare with me.

julia> f = (u, p, t) -> .- u
#38 (generic function with 1 method)

julia> prb = ODEProblem(f, (1.0:5.0...,), (0.0,2.0))
ODEProblem with uType NTuple{5,Float64} and tType Float64. In-place: false
timespan: (0.0, 2.0)
u0: (1.0, 2.0, 3.0, 4.0, 5.0)

julia> sprb = solve(prb, Tsit5(), abstol=1e-17, reltol=1e-17, save_on=false)
retcode: Success
Interpolation: specialized 4th order "free" interpolation
t: 2-element Array{Float64,1}:
 0.0
 2.0
u: 2-element Array{RecursiveArrayTools.ArrayPartition{Float64,NTuple{5,Float64}},1}:
 1.02.03.04.05.0                                                                             
 0.135335283236613170.270670566473226350.40600584970983820.54134113294645270.6766764161830637

julia> sprb.u[2] .- ((1:5).*exp(-2)) #Check with the real result
(4.718447854656915e-16, 0.13533528323661365, 0.2706705664732255, 0.40600584970984, 0.541341132946451)

It seems like second argument of the broadcast uses only the first element.

julia> sprb.u[2] .- 1
(-0.8646647167633869, -0.7293294335267737, -0.5939941502901618, -0.4586588670535473, -0.32332358381693627)

julia> sprb.u[2] .- (1:5)
(-0.8646647167633869, -0.7293294335267737, -0.5939941502901618, -0.4586588670535473, -0.32332358381693627)

julia> sprb.u[2] .- collect(1:5)
(-0.8646647167633869, -0.7293294335267737, -0.5939941502901618, -0.4586588670535473, -0.32332358381693627)

Except, a value with the same type

julia> sprb.u[2] .- sprb.u[2]
(0.0, 0.0, 0.0, 0.0, 0.0)

julia> sprb.u[2] .- sprb.u[1]
(-0.8646647167633869, -1.7293294335267737, -2.593994150290162, -3.4586588670535474, -4.323323583816936)

I do not know where to blame (here or Base) but normal NTuple works as expected

julia> (1:5...,) .- (1:5)
5-element Array{Int64,1}:
 0
 0
 0
 0
 0

[ChrisRackauckas]

Interesting. I can't say I know why this would occur.

[tomaklutfu]

Let me chase down a bit more then. Now I see there are some discrepancies for broadcast here and Base

julia> Base.Broadcast.combine_styles(((NTuple{5, Float64})))
Base.Broadcast.Style{Tuple}()

julia> RecursiveArrayTools.combine_styles(((NTuple{5, Float64}).parameters...,))
Base.Broadcast.DefaultArrayStyle{0}()

May be this is a peculiar and desired situation of ArrayPartition. I think, OrdinaryDiffEq should have just returned Array{NTuple{5,Float64}, 1} similar to normal array initial conditions (u0::Array{Float64, 1} gives sol.u::Array{Array{Float64, 1}, 1}). If that's the case I should open an issue to OrdinaryDiffEq.

[tomaklutfu]

I meant to ask first. Should I open issue to OrdinaryDiffEq.jl?

[ChrisRackauckas]

This is a DifferentialEquations.jl issue, that we shouldn't be changing from tuple in DiffEqBase. Please open that issue. I hope to fix it soon.