Skip to content

Commit

Permalink
Merge pull request #47 from RXGottlieb/master
Browse files Browse the repository at this point in the history 
Added xabs(x) and NS odd power relaxations
  • Loading branch information
@mewilhel
mewilhel authored May 12, 2021
2 parents f569525 + 81cd02c commit 55569f3
Show file tree
Hide file tree
Showing 4 changed files with 120 additions and 3 deletions.
3 changes: 2 additions & 1 deletion src/McCormick.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -71,7 +71,8 @@ export MC, cc, cv, Intv, lo, hi, cc_grad, cv_grad, cnst, +, -, *, /, convert,
relu, param_relu, leaky_relu, maxsig, maxtanh, softplus, pentanh,
sigmoid, bisigmoid, softsign, gelu, elu, selu, swish1,
positive, negative, lower_bnd, upper_bnd, bnd, xlogx,
<, <=, ==, fma, cbrt, abs2, sinpi, cospi, arh, xexpax, trilinear
<, <=, ==, fma, cbrt, abs2, sinpi, cospi, arh, xexpax, trilinear,
xabsx

# Export kernel operators
export plus_kernel, minus_kernel, mult_kernel, div_kernel, max_kernel,
Expand Down
72 changes: 72 additions & 0 deletions src/forward_operators/other.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -565,3 +565,75 @@ function xexpax(x::MC{N,T}, a::Float64) where {N, T <: RelaxTag}
yMC, tp1, tp2 = xexpax_kernel(x, a, intv_xexpax, Inf, Inf)
return yMC
end

"""
xabsx
The function `xabsx` is defined as `xabsx(x) = x*abs(x)`.
"""
xabsx(x::Float64) = x*abs(x)
xabsx_deriv(x::Float64) = 2*abs(x)
xabsx_deriv2(x::Float64) = 2*abs(x)/x #BUT DOESN'T EXIST AT 0
@inline function cc_xabsx(x::Float64, xL::Float64, xU::Float64)
(xU <= 0.0) && (return xabsx(x), xabsx_deriv(x))
(0.0 <= xL) && (return dline_seg(xabsx, xabsx_deriv, x, xL, xU))

root_f(x, xL, xU) = (xabsx(xU) - xabsx(x)) - 2*(xU - x)*abs(x)
root_df(x, xL, xU) = xabsx_deriv(x) - xU*xabsx_deriv2(x)
inflection, flag = newton(xL, xL, 0.0, root_f, root_df, xL, xU)
flag && (inflection = golden_section(xL, xU, root_f, xL, xU))
if x <= inflection
return xabsx(x), xabsx_deriv(x)
else
return dline_seg(xabsx, xabsx_deriv, x, inflection, xU)
end
end
@inline function cv_xabsx(x::Float64, xL::Float64, xU::Float64)
(xU <= 0.0) && (return dline_seg(xabsx, xabsx_deriv, x, xL, xU))
(0.0 <= xL) && (return xabsx(x), xabsx_deriv(x))

root_f(x, xL, xU) = (xabsx(x) - xabsx(xL)) - 2*(x - xL)*abs(x)
root_df(x, xL, xU) = -xabsx_deriv(x) + xL*xabsx_deriv2(x)
inflection, flag = newton(xU, 0.0, xU, root_f, root_df, xL, xU)
flag && (inflection = golden_section(xL, xU, root_f, xL, xU))
if x <= inflection
return dline_seg(xabsx, xabsx_deriv, x, xL, inflection)
else
return xabsx(x), xabsx_deriv(x)
end
end
function xabsx(x::Interval{Float64})
xabsx_xL = Interval(x.lo)*abs(Interval(x.lo))
xabsx_xU = Interval(x.hi)*abs(Interval(x.hi))
Interval{Float64}(xabsx_xL.lo, xabsx_xU.hi)
end
@inline function xabsx_kernel(x::MC{N, T}, y::Interval{Float64}) where {N,T<:Union{NS,MV}}
xL = x.Intv.lo
xU = x.Intv.hi
eps_max = xU
eps_min = xL
midcc, cc_id = mid3(x.cc, x.cv, eps_max)
midcv, cv_id = mid3(x.cc, x.cv, eps_min)
cc, dcc = cc_xabsx(midcc, xL, xU)
cv, dcv = cv_xabsx(midcv, xL, xU)
cc_grad = mid_grad(x.cc_grad, x.cv_grad, cc_id)*dcc
cv_grad = mid_grad(x.cc_grad, x.cv_grad, cv_id)*dcv
cv, cc, cv_grad, cc_grad = cut(y.lo, y.hi, cv, cc, cv_grad, cc_grad)
return MC{N,T}(cv, cc, y, cv_grad, cc_grad, x.cnst)
end
@inline function xabsx_kernel(x::MC{N, T}, y::Interval{Float64}) where {N,T<:Diff}
xL = x.Intv.lo
xU = x.Intv.hi
eps_max = xU
eps_min = xL
midcc, cc_id = mid3(x.cc, x.cv, eps_max)
midcv, cv_id = mid3(x.cc, x.cv, eps_min)
cc, dcc = cc_xabsx(midcc, xL, xU)
cv, dcv = cv_xabsx(midcv, xL, xU)
cc_grad = mid_grad(x.cc_grad, x.cv_grad, cc_id)*dcc
cv_grad = mid_grad(x.cc_grad, x.cv_grad, cv_id)*dcv
return MC{N,T}(cv, cc, y, cv_grad, cc_grad, x.cnst)
end
function xabsx(x::MC{N,T}) where {N, T <: RelaxTag}
xabsx_kernel(x, xabsx(x.Intv))
end
32 changes: 30 additions & 2 deletions src/forward_operators/power.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -161,8 +161,8 @@ end
@inline function pos_odd(x::MC{N,T}, c::Z, y::Interval{Float64}) where {N, Z<:Integer, T<:Union{NS,MV}}
midcc, cc_id = mid3(x.cc, x.cv, x.Intv.hi)
midcv, cv_id = mid3(x.cc, x.cv, x.Intv.lo)
cc, dcc = cc_powodd(midcc, x.Intv.lo, x.Intv.hi, c)
cv, dcv = cv_powodd(midcv, x.Intv.lo, x.Intv.hi, c)
cc, dcc = cc_powodd_ns(midcc, x.Intv.lo, x.Intv.hi, c)
cv, dcv = cv_powodd_ns(midcv, x.Intv.lo, x.Intv.hi, c)
cc_grad = mid_grad(x.cc_grad, x.cv_grad, cc_id)*dcc
cv_grad = mid_grad(x.cc_grad, x.cv_grad, cv_id)*dcv
cv, cc, cv_grad, cc_grad = cut(y.lo, y.hi, cv, cc, cv_grad, cc_grad)
Expand All @@ -181,6 +181,34 @@ end
cc_grad = min(0.0, gdcc1)*x.cv_grad + max(0.0, gdcc2)*x.cc_grad
return MC{N,Diff}(cv, cc, y, cv_grad, cc_grad, x.cnst)
end
@inline function cv_powodd_ns(x::Float64, xL::Float64, xU::Float64, n::Z) where {Z <: Integer}
(xU <= 0.0) && (return dline_seg(^, pow_deriv, x, xL, xU, n))
(0.0 <= xL) && (return x^n, n*x^(n - 1))

root_f(x::Float64, xL::Float64, xU::Float64) = (x^n - xL^n) - (x-xL)*(n)*(x^(n-1));
root_df(x::Float64, xL::Float64, xU::Float64) = n*x^(n-1) + xL*n*(n-1)*x - n*x^(n-1) - n*(n-1)*x^(n-1);
inflection, flag = newton(xU, 0.0, xU, root_f, root_df, xL, xU)
flag && (inflection = golden_section(xL, xU, root_f, xL, xU))
if x <= inflection
return dline_seg(^, pow_deriv, x, xL, inflection, n)
else
return x^n, n*x^(n-1)
end
end
@inline function cc_powodd_ns(x::Float64, xL::Float64, xU::Float64, n::Z) where {Z <: Integer}
(xU <= 0.0) && (return x^n, n*x^(n - 1))
(0.0 <= xL) && (return dline_seg(^, pow_deriv, x, xL, xU, n))

root_f(x::Float64, xL::Float64, xU::Float64) = (xU^n-x^n) - (xU-x)*(n)*(x^(n-1));
root_df(x::Float64, xL::Float64, xU::Float64) = -n*x^(n-1) - xU*n*(n-1)*x + n*x^(n-1) + n*(n-1)*x^(n-1);
inflection, flag = newton(xL, xL, 0.0, root_f, root_df, xL, xU)
flag && (inflection = golden_section(xL, xU, root_f, xL, xU))
if x <= inflection
return x^n, n*x^(n-1)
else
return dline_seg(^, pow_deriv, x, inflection, xU, n)
end
end

# neg_powneg_odd computes the McComrick relaxation of x^c where x < 0.0 and c is odd
@inline function cv_neg_powneg_odd(x::Float64, xL::Float64, xU::Float64, c::Z) where Z<:Integer
Expand Down
16 changes: 16 additions & 0 deletions test/forward_mccormick.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -554,6 +554,22 @@ end
@test check_vs_ref1(atanh, x_atanh_z2, yref_atanh_d1_z2, mctol)
@test check_vs_ref1(atanh, x_atanh_z1_ns, yref_atanh_ns, mctol)

##### x Times Abs(x) #####
x_xabsx_p = MC{2,Diff}(0.3, 0.3, Interval{Float64}(0.1,0.7), seed_gradient(1,Val(2)), seed_gradient(1,Val(2)), false)
x_xabsx_n = MC{2,Diff}(-0.3, -0.3, Interval{Float64}(-0.7,-0.1), seed_gradient(1,Val(2)), seed_gradient(1,Val(2)), false)
x_xabsx_z1 = MC{2,Diff}(0.2, 0.2, Interval{Float64}(-0.7,0.5), seed_gradient(1,Val(2)), seed_gradient(1,Val(2)), false)
x_xabsx_z2 = MC{2,Diff}(0.4, 0.4, Interval{Float64}(-0.7,0.5), seed_gradient(1,Val(2)), seed_gradient(1,Val(2)), false)

yref_xabsx_d1_p = MC{2, Diff}(0.09, 0.16999999999999998, Interval{Float64}(0.01, 0.49), @SVector[0.6, 0.0], @SVector[0.7999999999999999, 0.0], false)
yref_xabsx_d1_n = MC{2, Diff}(-0.16999999999999998, -0.09, Interval{Float64}(-0.49, -0.01), @SVector[0.7999999999999999, 0.0], @SVector[0.6, 0.0], false)
yref_xabsx_d1_z1 = MC{2, Diff}(0.03190908859009977, 0.1257359312880715, Interval{Float64}(-0.49, 0.25), @SVector[0.579898987322333, 0.0], @SVector[0.4142135623730951, 0.0], false)
yref_xabsx_d1_z2 = MC{2, Diff}(0.16000000000000003, 0.20857864376269047, Interval{Float64}(-0.49, 0.25), @SVector[0.8, 0.0], @SVector[0.4142135623730951, 0.0], false)

@test check_vs_ref1(xabsx, x_xabsx_p, yref_xabsx_d1_p, mctol)
@test check_vs_ref1(xabsx, x_xabsx_n, yref_xabsx_d1_n, mctol)
@test check_vs_ref1(xabsx, x_xabsx_z1, yref_xabsx_d1_z1, mctol)
@test check_vs_ref1(xabsx, x_xabsx_z2, yref_xabsx_d1_z2, mctol)

# CONVERSION
X = MC{2,NS}(4.5, 4.5, Interval{Float64}(-3.0,8.0), seed_gradient(1,Val(2)), seed_gradient(1,Val(2)), false)
X1 = convert(MC{2,NS}, 1)
Expand Down

2 comments on commit 55569f3

@mewilhel
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@JuliaRegistrator register

@JuliaRegistrator
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Registration pull request created: JuliaRegistries/General/36612

After the above pull request is merged, it is recommended that a tag is created on this repository for the registered package version.

This will be done automatically if the Julia TagBot GitHub Action is installed, or can be done manually through the github interface, or via:

git tag -a v0.11.0 -m "<description of version>" 55569f379614a5be3cc1b2f2f43147af8a56980a
git push origin v0.11.0

Please sign in to comment.