Fix NG flag

src/contractors.jl

@@ -25,9 +25,9 @@ end
 
 function contract(::Type{Krawczyk}, f, derivative, X::Interval)
     m = interval(mid(X))
-    Y = 1 / derivative(m)
+    Y = inv(derivative(m))
 
-    return m - Y*f(m) + (1 - Y*derivative(X)) * (X - m)
+    return m - Y*f(m) + (interval(1) - Y*derivative(X)) * (X - m)
 end
 
 function contract(::Type{Krawczyk}, f, derivative, X::AbstractVector)
@@ -40,10 +40,12 @@ function contract(::Type{Krawczyk}, f, derivative, X::AbstractVector)
     dm = derivative(m)
     det(dm) == 0 && return interval(-Inf, Inf, trv)
 
-    Y = inv(dm)
+    Y = interval.(inv(dm))
     mm = interval.(m)
 
-    return mm - Y*f(mm) + (I - Y*J) * (X - mm)
+    Id = interval.(Matrix(I, length(m), length(m)))
+
+    return mm - Y*f(mm) + (Id - Y*J) * (X - mm)
 end
 
 function contract(root_problem::RootProblem{C}, X) where C

src/linear_eq.jl

@@ -4,7 +4,7 @@ Preconditions the matrix A and b with the inverse of mid(A)
 function preconditioner(A::AbstractMatrix, b::AbstractArray)
     Aᶜ = mid.(A)
     det(Aᶜ) == 0 && return A, b
-    B = inv(Aᶜ)
+    B = interval.(inv(Aᶜ))
     return B*A, B*b
 end
 
@@ -85,7 +85,7 @@ function gauss_elimination_interval(A::AbstractMatrix, b::AbstractArray ; precon
     n = size(A, 1)
 
     p = similar(b)
-    p .= 0
+    p .= interval(0)
 
     if any(a -> in_interval(0, a), diag(A))
         p .= interval(-Inf, Inf)

test/roots.jl

@@ -171,30 +171,38 @@ end
 end
 
 @testset "@exact" begin
-    @exact f(x) = x^2 - 2.25
-    @exact g(x) = sin(1.57x)
+    @testset "1D" begin
+        @exact f(x) = x^2 - 2.25
+        @exact g(x) = sin(1.57x)
 
-    f0(x) = x^2 - 2.25
-    g0(x) = sin(1.57x)
+        f0(x) = x^2 - 2.25
+        g0(x) = sin(1.57x)
 
-    X = interval(-2.2, 2.2)
+        X = interval(-2.2, 2.2)
 
-    for (func, func0, nsol) in [(f, f0, 2), (g, g0, 3)]
-        for method in newtonlike_methods
-            rts = roots(func, X ; contractor = method)
-            @test length(rts) == nsol
+        @testset for method in newtonlike_methods
+            for (func, func0, nsol) in [(f, f0, 2), (g, g0, 3)]
+                rts = roots(func, X ; contractor = method)
+                @test length(rts) == nsol
 
-            rts0 = roots(func0, X ; contractor = method)
-            for (rt, rt0) in zip(rts, rts0)
-                @test isequal_interval(root_region(rt), root_region(rt0))
-            end
+                rts0 = roots(func0, X ; contractor = method)
+                for (rt, rt0) in zip(rts, rts0)
+                    @test isequal_interval(root_region(rt), root_region(rt0))
+                end
 
-            # Guarantee is currently broken with Krawczyk
-            if method == Krawczyk
-                @test_broken all(isguaranteed.(root_region.(rts)))
-            else
                 @test all(isguaranteed.(root_region.(rts)))
             end
         end
     end
+
+    @testset "2D" begin
+        @exact f(xx) = [xx[1]^2 - 1, xx[2]^2 - 5]
+        X = [interval(-2, 2), interval(-10, 10)]
+        @testset for method in newtonlike_methods
+            rts = roots(f, X ; contractor = method)
+            @test all(rts) do rt
+                all(isguaranteed.(root_region(rt)))
+            end
+        end
+    end
 end