Merge pull request #23 from mewilhel/master

Update to v0.6

Project.toml

@@ -1,7 +1,7 @@
 name = "ReverseMcCormick"
 uuid = "f9d89e60-3bef-4186-aa8f-bdf0d7614dd6"
 authors = ["Matthew Wilhelm <[email protected]>"]
-version = "0.5.3"
+version = "0.6.0"
 
 [deps]
 DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
@@ -11,7 +11,7 @@ McCormick = "53c679d3-6890-5091-8386-c291e8c8aaa1"
 [compat]
 DocStringExtensions = "~0.8"
 IntervalContractors = "~0.4.2"
-McCormick = "0.5.2"
+McCormick = "0.6.0"
 julia = "~1.2, ~1.3, ~1.4"
 
 [extras]

src/ReverseMcCormick.jl

@@ -6,7 +6,7 @@ using DocStringExtensions, McCormick
 
 import IntervalContractors: sin_rev, cos_rev, tan_rev
 
-export plus_rev, mul_rev, min_rev, max_rev, minus_rev, div_rev, exp_rev,
+export plus_rev, mult_rev, min_rev, max_rev, minus_rev, div_rev, exp_rev,
        exp2_rev, exp10_rev, expm1_rev, log_rev, log2_rev, log10_rev,
        log1p_rev, sin_rev, cos_rev, tan_rev, asin_rev, acos_rev, atan_rev,
        sec_rev, csc_rev, cot_rev, asec_rev, acsc_rev, acot_rev,

src/reverse_operators/arithmetic.jl

@@ -66,10 +66,9 @@ Creates reverse McCormick contractor for `a` = `b` +`c`
 """
 function plus_rev(a::MC, b::MC, c::MC)
     isempty(a) && (return a, a, a)
-    bold = b
-    b = b ∩ (a - c)
-    c = c ∩ (a - bold)
-    a, b, c
+    b_new = b ∩ (a - c)
+    c_new = c ∩ (a - b)
+    a, b_new, c_new
 end
 function plus_rev(a::MC{N,T}, b::MC{N,T}, c::C) where {N, T<:RelaxTag, C<:NumberNotRelax}
     isempty(a) && (return a, a, a)
@@ -96,21 +95,20 @@ Creates reverse McCormick contractor for `a` = `b`- `c`
 """
 function minus_rev(a::MC, b::MC, c::MC)
     isempty(a) && (return a, a, a)
-    bold = b
-    b = b ∩ (a + c)
-    c = c ∩ (bold - a)
-    a, b, c
+    b_new = b ∩ (a + c)
+    c_new = c ∩ (b - a)
+    a, b_new, c_new
 end
 function minus_rev(a::MC{N,T}, b::MC{N,T}, c::C) where {N, T<:RelaxTag, C<:NumberNotRelax}
     isempty(a) && (return a, a, a)
-    b = b ∩ (a + c)
-    a, b, MC{N,T}(c)
+    b_new = b ∩ (a + c)
+    a, b_new, MC{N,T}(c)
 end
 
 function minus_rev(a::MC{N,T}, b::C, c::MC{N,T}) where {N, T<:RelaxTag, C<:NumberNotRelax}
     isempty(a) && (return a, a, a)
-    c = c ∩ (b - a)
-    a, MC{N,T}(b), c
+    c_new = c ∩ (b - a)
+    a, MC{N,T}(b), c_new
 end
 
 
@@ -121,40 +119,37 @@ Creates reverse McCormick contractor for `a` = `-b``
 """
 function minus_rev(a::MC, b::MC)
     isempty(a) && (return a, a)
-    b = b ∩ -a
-    return a, b
+    b_new = b ∩ (-a)
+    return a, b_new
 end
 
 """
 $(SIGNATURES)
 
 Creates reverse McCormick contractor for `a` = `b`*`c`
 """
-function mul_rev(a::MC, b::MC, c::MC)
+function mult_rev(a::MC, b::MC, c::MC)
     isempty(a) && (return a, a, a)
-    bflag = 0.0 ∉ b
-    if bflag
-        temp1 = a / b
-        ((0.0 ∉ a) || bflag) && (c = c ∩ temp1)
-    end
-    cflag = 0.0 ∉ c
-    if cflag
-        temp2 = a / c
-        ((0.0 ∉ a) || cflag) && (b = b ∩ temp2)
-    end
-    a, b, c
+    c_new = (0.0 ∉ b) ? (c ∩ (a/b)) : c
+    b_new = (0.0 ∉ c) ? (b ∩ (a/c)) : b
+    a, b_new, c_new
 end
-function mul_rev(a::MC{N,T}, b::MC{N,T}, c::C) where {N, T<:RelaxTag, C<:NumberNotRelax}
+function mult_rev(a::MC{N,T}, b::MC{N,T}, c::C) where {N, T<:RelaxTag, C<:NumberNotRelax}
+    #println("ran me 1, a = $a, b = $b, c = $c")
+    #println("a*inv(c) = $(a*inv(c))")
+    isempty(a) && (return a, a, a)
     if !iszero(c)
-        b = b ∩ a*inv(c)
+        #println("b ∩ a*inv(c) = $(b ∩ a*inv(c))")
+        b = b ∩ (a*inv(c))
     elseif 0.0 ∉ a
         return empty(a), empty(a), empty(a)
     end
     a, b, MC{N,T}(c)
 end
-function mul_rev(a::MC{N,T}, c::C, b::MC{N,T}) where {N, T<:RelaxTag, C<:NumberNotRelax}
+function mult_rev(a::MC{N,T}, c::C, b::MC{N,T}) where {N, T<:RelaxTag, C<:NumberNotRelax}
+    #println("ran me 2")
     if !iszero(c)
-        b = b ∩ a*inv(c)
+        b = b ∩ (a*inv(c))
     elseif 0.0 ∉ a
         return empty(a), empty(a), empty(a)
     end
@@ -232,14 +227,44 @@ function power_rev(a::MC{N,T}, b::C, c::MC{N,T}) where {N, T<:RelaxTag, C<:Numbe
     end
     a, MC{N,T}(b), c
 end
+
+function int_power_rev(a::MC{N,T}, b::MC{N,T}, n::Int, c::C) where {N, T<:RelaxTag, C<:NumberNotRelax}
+    if n == 2
+        root = sqrt(a)
+        if lo(b) > 0.0
+            b = b ∩ root
+        else
+            b = b ∩ -root
+        end
+    elseif iseven(n)
+        root = a^(1/n)
+        if lo(b) > 0.0
+            b = b ∩ root
+        else
+            b = b ∩ -root
+        end
+    elseif isodd(n)
+        if lo(b) > 0.0
+            b = b ∩ (a ∩ Interval(0, Inf))^(1/n)
+        else
+            b = b ∩ -(((-a) ∩ Interval(0, Inf))^(1/n))
+        end
+    end
+    a, b, MC{N,T}(c)
+end
+
 function power_rev(a::MC{N,T}, b::MC{N,T}, c::C) where {N, T<:RelaxTag, C<:NumberNotRelax}
     isempty(a) && (return a, a, a)
     (isempty(b) || 0.0 > lo(b)) && (return a, a, a)
     if isone(-c)
         a, b = inv_rev(a,b)
         return a, b, MC{N,T}(c)
     elseif !iszero(c)
-        b = b ∩ a^(1/c)
+        if isinteger(c)
+            return int_power_rev(a, b, Int(c), c)
+        elseif lo(b) > 0.0
+            b = b ∩ a^(one(C)/c)
+        end
     elseif iszero(c) && 1.0 ∉ a
         return empty(a), empty(a), empty(a)
     end
@@ -254,7 +279,7 @@ Creates reverse McCormick contractor for `a` = `sqrt(b)`. That is
 """
 function sqrt_rev(a::MC, b::MC)
     isempty(a) && (return a, a)
-    b = b ∩ a^2
+    b = b ∩ (a^2)
     a, b
 end
 sqr_rev(f, x) = power_rev(f, x, 2)

test/reverse_mccormick.jl

@@ -83,7 +83,7 @@ end
     b = MC{1,NS}(Interval{Float64}(-10.0,-1.0))
     c = MC{1,NS}(2.0, Interval{Float64}(1.1,4.5), 1)
 
-    aout1, bout1, cout1 = mul_rev(a,b,c)
+    aout1, bout1, cout1 = mult_rev(a,b,c)
 
     @test bout1.Intv.lo == Inf
     @test bout1.Intv.hi == -Inf
@@ -94,7 +94,7 @@ end
     cout1 = MC{1,NS}(Interval{Float64}(-10.0,-1.0))
     aout1 = bout1*cout1
 
-    aout1_a, bout1_a, cout1_a = mul_rev(aout1, bout1, cout1)
+    aout1_a, bout1_a, cout1_a = mult_rev(aout1, bout1, cout1)
 
     MC_1_is_equal(aout1_a, aout1, 0.00001)
     MC_1_is_equal(bout1_a, bout1, 0.00001)
@@ -104,7 +104,7 @@ end
     cout2 = MC{1,NS}(Interval{Float64}(-10.0,-1.0))
     aout2 = 0.3*bout1*cout1+1.0
 
-    aout2_a, bout2_a, cout2_a = mul_rev(aout2, bout2, cout2)
+    aout2_a, bout2_a, cout2_a = mult_rev(aout2, bout2, cout2)
 
     MC_1_is_equal(aout2_a, aout2, 0.00001)
     MC_1_is_equal(bout2_a, bout2, 0.00001)
@@ -120,18 +120,18 @@ end
     b = 0.5
     c = MC{1,NS}(Interval{Float64}(1.1,4.5))
 
-    aout2, bout2, cout2 = mul_rev(a, b, c)
-    @test cout2.cv == 2.2
-    @test cout2.cc == 6.0
+    aout2, bout2, cout2 = mult_rev(a, b, c)
+    @test cout2.cv == 1.1
+    @test cout2.cc == 4.5
 
-    aout2, bout2, cout2 = mul_rev(a, c, b)
-    @test bout2.cv == 2.2
-    @test bout2.cc == 6.0
+    aout2, bout2, cout2 = mult_rev(a, c, b)
+    @test bout2.cv == 1.1
+    @test bout2.cc == 4.5
 
-    aout2, bout2, cout2 = mul_rev(a, 0.0, c)
+    aout2, bout2, cout2 = mult_rev(a, 0.0, c)
     @test isempty(cout2)
 
-    aout2, bout2, cout2 = mul_rev(a, c, 0.0)
+    aout2, bout2, cout2 = mult_rev(a, c, 0.0)
     @test isempty(bout2)
 end
 
@@ -801,7 +801,7 @@ end
             bout, aout = f(b, a)
             @test isempty(aout)
         end
-        for f in (plus_rev, minus_rev, mul_rev, div_rev)
+        for f in (plus_rev, minus_rev, mult_rev, div_rev)
             bout, cout, aout = f(b, c, a)
             ~isempty(aout) && (println("f = $(f)"))
             ~isempty(cout) && (println("f = $(f)"))