Skip to content

Commit

Permalink
Add safe guard for domain violation.
Browse files Browse the repository at this point in the history
  • Loading branch information
@mewilhel
mewilhel committed Apr 9, 2020
1 parent 20b9ee0 commit 51c1d0f
Show file tree
Hide file tree
Showing 4 changed files with 60 additions and 29 deletions.
2 changes: 2 additions & 0 deletions src/McCormick.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -77,6 +77,8 @@ const MC_ENV_TOL = 1E-10
const MC_DIFF_MU = 1
const MC_MV_TOL = 1E-8
const MC_DEGEN_TOL = 1E-14
const MC_DOMAIN_TOL = 1E-12
const MC_DOMAIN_CATCH = true

const IntervalConstr = interval
const Half64 = Float64(0.5)
Expand Down
23 changes: 17 additions & 6 deletions src/forward_operators/division.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,8 @@ function div_MV(x::MC{N,Diff}, y::MC{N,Diff}, z::Interval{Float64}) where N
cc = -cc
cc_grad = -cc_grad
else
return MC{N,Diff}(union(x.Intv/Interval{Float64}(y.Intv.lo, -MC_DOMAIN_TOL),
x.Intv/Interval{Float64}(MC_DOMAIN_TOL, y.Intv.hi)))
error("Division (x/y) is unbounded on intervals y containing 0.")
end
return MC{N,Diff}(cv, cc, z, cv_grad, cc_grad, x.cnst && y.cnst)
Expand Down Expand Up @@ -87,16 +89,25 @@ function div_kernel(x::MC{N,Diff}, y::MC{N,Diff}, z::Interval{Float64}) where {N
elseif ~(degen1||degen2)
zMC = div_MV(x, y, z)
else
q = inv(y)
zMC = mult_kernel(x, q, z)
if (y.Intv.lo <= 0.0 <= y.Intv.hi)
q = inv(y)
zMC = mult_kernel(x, q, x.Intv/q.Intv)
else
q = inv(y)
zMC = mult_kernel(x, q, z)
end
end
return zMC
end

function /(x::MC{N,T}, y::MC{N,T}) where {N,T<:RelaxTag}
@assert ~(y.Intv.lo <= 0.0 <= y.Intv.hi) "Domain Error: When computing the relaxations of x/y the
interval bounds of y contained zero. As such the x/y is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this."
if ~(y.Intv.lo <= 0.0 <= y.Intv.hi)
if ~MC_DOMAIN_CATCH
error("Domain Error: When computing the relaxations of x/y the
interval bounds of y contained zero. As such the x/y is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this.")
end
end
return div_kernel(x, y, x.Intv/y.Intv)
end
48 changes: 33 additions & 15 deletions src/forward_operators/power.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -243,7 +243,7 @@ function neg_powneg_even(x::MC{N,Diff}, c::Z, y::Interval{Float64}) where {N, Z<
return MC{N,Diff}(cv, cc, y, cv_grad, cc_grad, x.cnst)
end

function pow_kernel(x::MC, c::Z, y::Interval{Float64}) where {Z<:Integer}
function pow_kernel(x::MC{N,T}, c::Z, y::Interval{Float64}) where {Z<:Integer,N,T<:RelaxTag}
if (c == 0)
z = one(x)
elseif (c == 1)
Expand All @@ -266,20 +266,30 @@ function pow_kernel(x::MC, c::Z, y::Interval{Float64}) where {Z<:Integer}
elseif (x.Intv.lo > 0.0)
z = npp_or_pow4(x, c, y)
else
error("Domain Error: When computing the relaxations of y^c (c < 0) the
interval bounds of y contained zero. As such the y^c is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this.")
if ~MC_DOMAIN_CATCH
error("Domain Error: When computing the relaxations of y^c (c < 0) the
interval bounds of y contained zero. As such the y^c is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this.")
else
z = MC{N,T}(union(Interval{Float64}(x.Intv.lo, -MC_DOMAIN_TOL)^c,
Interval{Float64}(MC_DOMAIN_TOL, x.Intv.hi)^c))
end
end
end
return z
return z
end
function pow(x::MC, c::Z) where {Z<:Integer}
function pow(x::MC{N,T}, c::Z) where {Z<:Integer,N,T<:RelaxTag}
if (x.Intv.lo <= 0.0 <= x.Intv.hi) && (c < 0)
error("Domain Error: When computing the relaxations of y^c (c < 0) the
interval bounds of y contained zero. As such the y^c is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this.")
if ~MC_DOMAIN_CATCH
error("Domain Error: When computing the relaxations of y^c (c < 0) the
interval bounds of y contained zero. As such the y^c is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this.")
else
return MC{N,T}(union(Interval{Float64}(x.Intv.lo, -MC_DOMAIN_TOL)^c,
Interval{Float64}(MC_DOMAIN_TOL, x.Intv.hi)^c))
end
end
return pow_kernel(x, c, x.Intv^c)
end
Expand Down Expand Up @@ -366,15 +376,23 @@ function inv1(x::MC{N,NS}, y::Interval{Float64}) where N
end
function inv_kernel(x::MC{N,T}, y::Interval{Float64}) where {N,T<:RelaxTag}
if (x.Intv.lo <= 0.0 <= x.Intv.hi)
error("Domain Error: When computing the relaxations of inv(y) the
interval bounds of y contained zero. As such the y is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this.")
if ~MC_DOMAIN_CATCH
error("Domain Error: When computing the relaxations of inv(x) the
interval bounds of x contained zero. As such the y is unbounded
and the relaxations do not exist. This may occur due to the expansiveness
in long calculations. Reformulating the function may remedy this. Also,
consider setting the MC_DOMAIN_CATCH flag to true. This will treat
inv(x) as union(inv(Interval(x.Intv.lo,-MC_DOMAIN_TOL),
inv(Interval(MC_DOMAIN_TOL, x.Intv.hi)")
end
end
if (x.Intv.hi < 0.0)
x = neg_powneg_odd(x, -1, y)
elseif (x.Intv.lo > 0.0)
x = inv1(x, y)
elseif MC_DOMAIN_CATCH
x = MC{N,T}(union(inv(Interval{Float64}(x.Intv.lo,-MC_DOMAIN_TOL)),
inv(Interval{Float64}(MC_DOMAIN_TOL, x.Intv.hi))))
end
return x
end
Expand Down
16 changes: 8 additions & 8 deletions test/forward_mccormick.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -700,26 +700,26 @@ end
out1 = McCormick.mul1_u1pos_u2mix(x1a, x2b, z1, false)
@test out1.cv == 1.0
@test out1.cc == 1.5
@test out1.Intv.lo == 0.5
@test out1.Intv.hi == 2.0
@test isapprox(out1.Intv.lo, 0.5, atol =1E-6)
@test isapprox(out1.Intv.hi, 2.0, atol =1E-6)

out2 = McCormick.mul1_u1pos_u2mix(x1a, x2c, z2, false)
@test out2.cv == 4.0
@test out2.cc == 6.0
@test out2.Intv.lo == 2.0
@test out2.Intv.hi == 8.0
@test isapprox(out2.Intv.lo, 2.0, atol =1E-6)
@test isapprox(out2.Intv.hi, 8.0, atol =1E-6)

out3 = McCormick.mul2_u1pos_u2mix(x1a, x2b, z1, false)
@test out3.cv == 1.0
@test out3.cc == 1.5
@test out3.Intv.lo == 0.5
@test out3.Intv.hi == 2.0
@test isapprox(out3.Intv.lo, 0.5, atol =1E-6)
@test isapprox(out3.Intv.hi, 2.0, atol =1E-6)

out4 = McCormick.mul2_u1pos_u2mix(x1a, x2c, z2, false)
@test out4.cv == 4.0
@test out4.cc == 6.0
@test out4.Intv.lo == 2.0
@test out4.Intv.hi == 8.0
@test isapprox(out4.Intv.lo, 2.0, atol =1E-6)
@test isapprox(out4.Intv.hi, 8.0, atol =1E-6)

x1a = MC{2,NS}(1.1, 2.3, Interval(0.1,3.3))
x2a = MC{2,NS}(2.1, 3.3, Interval(1.1,4.3))
Expand Down

0 comments on commit 51c1d0f

Please sign in to comment.