Matrix multiplication with mixed Interval SMatrix and SMatrix broken

[schillic]

This came up while working on JuliaIntervals/IntervalRootFinding.jl#192.

julia> A = SMatrix{2,2}([-12..12 -12..12; -2.. -2 1..1])
2×2 SMatrix{2, 2, Interval{Float64}, 4} with indices SOneTo(2)×SOneTo(2):
 [-12.0, 12.0]  [-12.0, 12.0]
  [-2.0, -2.0]    [1.0, 1.0]

julia> B = SMatrix{2,2}([1.0 1.0; 1.0 1.0])
2×2 SMatrix{2, 2, Float64, 4} with indices SOneTo(2)×SOneTo(2):
 1.0  1.0
 1.0  1.0

A*A still works, but A*B and B*A do not work anymore.

v0.20.9

julia> A*B
2×2 SMatrix{2, 2, Interval{Float64}, 4} with indices SOneTo(2)×SOneTo(2):
 [-24, 24]  [-24, 24]
  [-1, -1]   [-1, -1]

v0.21.0

julia> A*B
ERROR: promotion of types Interval{Float64}, Float64 and Interval{Float64} failed to change any arguments
Stacktrace:
  [1] error(::String, ::String, ::String)
    @ Base ./error.jl:44
  [2] sametype_error(input::Tuple{Interval{Float64}, Float64, Interval{Float64}})
    @ Base ./promotion.jl:405
  [3] not_sametype(x::Tuple{Interval{Float64}, Float64, Interval{Float64}}, y::Tuple{Interval{Float64}, Float64, Interval{Float64}})
    @ Base ./promotion.jl:399
  [4] promote
    @ ./promotion.jl:388 [inlined]
  [5] muladd(x::Interval{Float64}, y::Float64, z::Interval{Float64})
    @ Base ./promotion.jl:447
  [6] macro expansion
    @ ~/.julia/packages/StaticArrays/9KYrc/src/matrix_multiply.jl:221 [inlined]
  [7] mul_loop
    @ ~/.julia/packages/StaticArrays/9KYrc/src/matrix_multiply.jl:188 [inlined]
  [8] macro expansion
    @ ~/.julia/packages/StaticArrays/9KYrc/src/matrix_multiply.jl:139 [inlined]
  [9] _mul
    @ ~/.julia/packages/StaticArrays/9KYrc/src/matrix_multiply.jl:130 [inlined]
 [10] *(A::SMatrix{2, 2, Interval{Float64}, 4}, B::SMatrix{2, 2, Float64, 4})
    @ StaticArrays ~/.julia/packages/StaticArrays/9KYrc/src/matrix_multiply.jl:11
 [11] top-level scope
    @ REPL[67]:1

[OlivierHnt]

First let me say that for the 0.21 release and onward, we are being more strict to not conflate Number and Interval (note that Interval should not be a subtype of Real in the future). The rational is that this mixing sometimes leads to code, not geared towards the specifics of interval arithmetic, to still run and silently return a non rigorous result. So you should no longer expect things like Float64(::Interval), or convert(Interval, ::Number), or promote(::Number, ::Interval) to work.

So this is why you are hitting this issue you are reporting: the generic fallback method of muladd in Base tries to promote its arguments.

That being said, since we do allow *(::Interval, ::Number) and +(::Interval, ::Number) I think it is reasonable to expect muladd(::Union{Interval,Number}, ::Union{Interval,Number}, ::Union{Interval,Number}) to work.

cc-ed @Kolaru, @lbenet, @lucaferranti

[schillic]

It seems to me that IntervalRootFinding (see the link in the OP) was relying on this to work. I do not know who feels responsible for that package, but since it is part of the same organization, maybe some developers here have to decide whether you want to support v0.21 there. Note that, without support in IntervalRootFinding, there will be no support in the TaylorModels package, and that then cascades.

I just wanted to let you know that this is a consequence of the recent changes and thought that maybe this was unintended.

[OlivierHnt]

Yes thank you very much for pointing this out! 🙂
This is definitely undesirable and ought to be resolved.

In IntervalArithmetic.jl, just defining muladd for a mix of Number and Interval arguments ought to be enough to close this issue. Then, the two other problems mentioned in the OP should be addressed directly in IntervalRootFinding.jl perhaps?