Conversion issue when combining StaticArrays, ForwardDif, and IntervalArithmetic

[maltezfaria]

In trying to combine ForwardDiff with StaticArrays and IntervalArithmetic I ran into the following issue:

using StaticArrays, ForwardDiff, IntervalArithmetic

x  = SVector(1.0, 2.0)
xI = SVector(interval(0.0, 1.0), interval(0.0, 1.0))
g  = (x) -> sum(insert(x, 1, 1.0))
∇g = x -> ForwardDiff.gradient(g, x)
∇g(x) # OK
∇g(xI) # error

I did not look too deep, but it seems that somewhere in the conversion pipeline StaticArrays calls a conversion method, which in turn calls the constructor Interval{Float64}(x). I saw that this constructor was intentionally commented out

IntervalArithmetic.jl/src/intervals/construction.jl

Line 371 in 429ad16

# Interval{T}(x) where {T<:NumTypes} = convert(Interval{T}, x)

Uncommenting the line above does "solve" the issue, but the IntervalArithmetic tests fail because it is explicitly tested that Interval{T}(x) should throw

IntervalArithmetic.jl/test/interval_tests/construction.jl

Line 47 in 429ad16

@test_throws MethodError Interval{Float64}(1)

So I see that this is intentional. I was wondering if there is a workaround?

Thanks.

[OlivierHnt]

Thx for reporting this issue. If one defines g = (x) -> sum(insert(x, 1, interval(1))), then ∇g(xI) does work.
We removed Interval and Interval{T} as constructors because of implicit conversion. But now that we have the "NG flag" we could/should bring these back I suppose 🤔.
This would help with #668.

Your example also pointed out that ExactReal and ForwardDiff.Dual do not play well together. The following works

using StaticArrays, ForwardDiff, IntervalArithmetic

julia> Base.convert(::Type{ForwardDiff.Dual{T,V,N}}, x::ExactReal) where {T,V,N} = ForwardDiff.Dual{T}(convert(V, x), zero(ForwardDiff.Partials{N,V})) # missing conversion!

julia> x  = SVector(1.0, 2.0)
2-element SVector{2, Float64} with indices SOneTo(2):
 1.0
 2.0

julia> xI = SVector(interval(0.0, 1.0), interval(0.0, 1.0))
2-element SVector{2, Interval{Float64}} with indices SOneTo(2):
 [0.0, 1.0]_com
 [0.0, 1.0]_com

julia> @exact g  = (x) -> sum(insert(x, 1, 1.0)); # wraps `1.0` into `ExactReal(1.0)`

julia> ∇g = x -> ForwardDiff.gradient(g, x);

julia> ∇g(x) # OK (relies on `Float64(::ExactReal)`)
2-element SVector{2, Float64} with indices SOneTo(2):
 1.0
 1.0

julia> ∇g(xI) # also OK (relies on `Interval{Float64}(::ExactReal)`)
2-element SVector{2, Interval{Float64}} with indices SOneTo(2):
 [1.0, 1.0]_com
 [1.0, 1.0]_com

[maltezfaria]

Thanks for looking into this. In my application I can't really use g = (x) -> sum(insert(x,1,interval(1))) since that breaks ∇g(x)...

I don't know the internals of IntervalArithmetic, but if bringing back implicit conversations is not a major problem, it would certainly fix my issue (currently I am just committing type piracy and defining Interval{T}(x)).

[OlivierHnt]

I think we can bring back such implicit conversions, it will just be marked with a "not guaranteed" flag. E.g.

julia> insert(xI, 1, 1.0)
3-element SVector{3, Interval{Float64}} with indices SOneTo(3):
 [1.0, 1.0]_com_NG
 [0.0, 1.0]_com
 [0.0, 1.0]_com

Note that for functions expected to work for both floats and intervals, the current way is to use ExactReal and @exact (to not have the "NG flag").
The proposed PR gives:

julia> g  = (x) -> sum(insert(x, 1, ExactReal(1.0)))

julia> g(x)
4.0

julia> g(xI)
[1.0, 3.0]_com

julia> ∇g = x -> ForwardDiff.gradient(g, x);

julia> ∇g(x)
2-element SVector{2, Float64} with indices SOneTo(2):
 1.0
 1.0

julia> ∇g(xI)
2-element SVector{2, Interval{Float64}} with indices SOneTo(2):
 [1.0, 1.0]_com
 [1.0, 1.0]_com