Skip to content

Commit 7b67107

Browse files
eeshan9815eeshan9815
committed
Incorporate review comments
1 parent f837caf commit 7b67107

File tree

3 files changed

+112
-108
lines changed

3 files changed

+112
-108
lines changed

‎perf/slopes.jl

Lines changed: 5 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -9,12 +9,15 @@ function benchmark_results()
99
push!(f, x->(exp(x^2)))
1010
push!(f, x->(x^4 - 12x^3 + 47x^2 - 60x - 20exp(-x)))
1111
push!(f, x->(x^6 - 15x^4 + 27x^2 + 250))
12+
13+
s = interval(0.75, 1.75)
1214
df = DataFrame()
15+
1316
df[:Method] = ["Automatic Differentiation", "Slope Expansion"]
1417
for n in 1:length(f)
1518

16-
t1 = ForwardDiff.derivative(f[n], 0.75..1.75)
17-
t2 = slope(f[n], 0.75..1.75, 1.25)
19+
t1 = ForwardDiff.derivative(f[n], s)
20+
t2 = slope(f[n], s, mid(s))
1821
df[Symbol("f" * "$n")] = [t1, t2]
1922
end
2023
a = []

‎src/slopes.jl

Lines changed: 98 additions & 99 deletions
Original file line numberDiff line numberDiff line change
@@ -1,199 +1,198 @@
11
# Reference : Dietmar Ratz - An Optimized Interval Slope Arithmetic and its Application
22
using IntervalArithmetic, ForwardDiff
33
import Base: +, -, *, /, ^, sqrt, exp, log, sin, cos, tan, asin, acos, atan
4+
import IntervalArithmetic: mid, interval
45

56
function slope(f::Function, x::Interval, c::Real)
67
try
7-
f(SlopeVar(x, c)).fs
8+
f(slope_var(x, c)).s
89
catch y
910
if isa(y, MethodError)
1011
ForwardDiff.derivative(f, x)
1112
end
1213
end
1314
end
1415

15-
struct SlopeType
16-
fx::Interval
17-
fc::Interval
18-
fs::Interval
16+
struct Slope{T}
17+
x::Interval{T}
18+
c::Interval{T}
19+
s::Interval{T}
20+
Slope{T}(a, b, c) where T = new(a, b, c)
21+
Slope{T}(c) where T = Slope{T}(c, c, 0)
1922
end
2023

21-
function SlopeConst(c::Union{Real, Interval})
22-
SlopeType(c, c, 0)
24+
function slope_var(v::Real)
25+
Slope{Float64}(v, v, 1)
2326
end
2427

25-
function SlopeVar(v::Real)
26-
SlopeType(v, v, 1)
28+
function slope_var(v::Interval, c::Real)
29+
Slope{Float64}(v, c, 1)
2730
end
2831

29-
function SlopeVar(v::Interval, c::Real)
30-
SlopeType(v, c, 1)
32+
function interval(u::Slope)
33+
u.x
3134
end
3235

33-
function fxValue(u::SlopeType)
34-
u.fx
36+
function mid(u::Slope)
37+
u.c
3538
end
3639

37-
function fcValue(u::SlopeType)
38-
u.fc
40+
function slope(u::Slope)
41+
u.s
3942
end
4043

41-
function fsValue(u::SlopeType)
42-
u.fs
44+
function +(u::Slope, v::Slope)
45+
Slope{Float64}(u.x + v.x, u.c + v.c, u.s + v.s)
4346
end
4447

45-
function +(u::SlopeType, v::SlopeType)
46-
SlopeType(u.fx + v.fx, u.fc + v.fc, u.fs + v.fs)
48+
function -(u::Slope, v::Slope)
49+
Slope{Float64}(u.x - v.x, u.c - v.c, u.s - v.s)
4750
end
4851

49-
function -(u::SlopeType, v::SlopeType)
50-
SlopeType(u.fx - v.fx, u.fc - v.fc, u.fs - v.fs)
52+
function *(u::Slope, v::Slope)
53+
Slope{Float64}(u.x * v.x, u.c * v.c, u.s * v.c + u.x * v.s)
5154
end
5255

53-
function *(u::SlopeType, v::SlopeType)
54-
SlopeType(u.fx * v.fx, u.fc * v.fc, u.fs * v.fc + u.fx * v.fs)
56+
function /(u::Slope, v::Slope)
57+
Slope{Float64}(u.x / v.x, u.c / v.c, (u.s - (u.c / v.c) * v.s) / v.x)
5558
end
5659

57-
function /(u::SlopeType, v::SlopeType)
58-
SlopeType(u.fx / v.fx, u.fc / v.fc, (u.fs - (u.fc / v.fc) * v.fs) / v.fx)
60+
function +(u::Union{Interval, Real}, v::Slope)
61+
Slope{Float64}(u + v.x, u + v.c, v.s)
5962
end
6063

61-
function +(u::Union{Interval, Real}, v::SlopeType)
62-
SlopeType(u + v.fx, u + v.fc, v.fs)
64+
function -(u::Union{Interval, Real}, v::Slope)
65+
Slope{Float64}(u - v.x, u - v.c, -v.s)
6366
end
6467

65-
function -(u::Union{Interval, Real}, v::SlopeType)
66-
SlopeType(u - v.fx, u - v.fc, -v.fs)
68+
function *(u::Union{Interval, Real}, v::Slope)
69+
Slope{Float64}(u * v.x, u * v.c, u * v.s)
6770
end
6871

69-
function *(u::Union{Interval, Real}, v::SlopeType)
70-
SlopeType(u * v.fx, u * v.fc, u * v.fs)
72+
function /(u::Union{Interval, Real}, v::Slope)
73+
Slope{Float64}(u / v.x, u / v.c, -(u / v.c) * (v.s / v.x))
7174
end
7275

73-
function /(u::Union{Interval, Real}, v::SlopeType)
74-
SlopeType(u / v.fx, u / v.fc, -(u / v.fc) * (v.fs / v.fx))
75-
end
76-
77-
+(v::SlopeType, u::Union{Interval, Real}) = u + v
76+
+(v::Slope, u::Union{Interval, Real}) = u + v
7877

79-
-(v::SlopeType, u::Union{Interval, Real}) = u - v
80-
-(u::SlopeType) = u * -1
78+
-(v::Slope, u::Union{Interval, Real}) = u - v
79+
-(u::Slope) = u * -1
8180

82-
*(v::SlopeType, u::Union{Interval, Real}) = u * v
81+
*(v::Slope, u::Union{Interval, Real}) = u * v
8382

84-
/(v::SlopeType, u::Union{Interval, Real}) = u / v
83+
/(v::Slope, u::Union{Interval, Real}) = u / v
8584

86-
function sqr(u::SlopeType)
87-
SlopeType(u.fx ^ 2, u.fc ^ 2, (u.fx + u.fc) * u.fs)
85+
function sqr(u::Slope)
86+
Slope{Float64}(u.x ^ 2, u.c ^ 2, (u.x + u.c) * u.s)
8887
end
8988

90-
function ^(u::SlopeType, k::Integer)
89+
function ^(u::Slope, k::Integer)
9190
if k == 0
92-
return SlopeConst(1)
91+
return Slope{Float64}(1)
9392
elseif k == 1
9493
return u
9594
elseif k == 2
9695
return sqr(u)
9796
else
98-
hxi = interval(u.fx.lo) ^ k
99-
hxs = interval(u.fx.hi) ^ k
97+
hxi = interval(u.x.lo) ^ k
98+
hxs = interval(u.x.hi) ^ k
10099
hx = hull(hxi, hxs)
101100

102-
if (k % 2 == 0) && (0 ∈ u.fx)
101+
if (k % 2 == 0) && (0 ∈ u.x)
103102
hx = interval(0, hx.hi)
104103
end
105104

106-
hc = u.fc ^ k
105+
hc = u.c ^ k
107106

108-
i = u.fx.lo - u.fc.lo
109-
s = u.fx.hi - u.fc.hi
107+
i = u.x.lo - u.c.lo
108+
s = u.x.hi - u.c.hi
110109

111-
if ((i == 0) || (s == 0) || (k % 2 == 1 && Interval(0) ⪽ u.fx))
112-
h1 = k * (u.fx ^ (k - 1))
110+
if ((i == 0) || (s == 0) || (k % 2 == 1 && Interval(0) ⪽ u.x))
111+
h1 = k * (u.x ^ (k - 1))
113112
else
114-
if k % 2 == 0 || u.fx.lo >= 0
113+
if k % 2 == 0 || u.x.lo >= 0
115114
h1 = interval((hxi.hi - hc.lo) / i, (hxs.hi - hc.lo) / s)
116115
else
117116
h1 = interval((hxs.lo - hc.hi) / s, (hxi.lo - hc.hi) / i)
118117
end
119118
end
120-
return SlopeType(hx, hc, h1 * u.fs)
119+
return Slope{Float64}(hx, hc, h1 * u.s)
121120
end
122121
end
123122

124-
function sqrt(u::SlopeType)
125-
SlopeType(sqrt(u.fx), sqrt(u.fc), u.fs / (sqrt(u.fx) + sqrt(u.fc)))
123+
function sqrt(u::Slope)
124+
Slope{Float64}(sqrt(u.x), sqrt(u.c), u.s / (sqrt(u.x) + sqrt(u.c)))
126125
end
127126

128-
function exp(u::SlopeType)
129-
hx = exp(u.fx)
130-
hc = exp(u.fc)
127+
function exp(u::Slope)
128+
hx = exp(u.x)
129+
hc = exp(u.c)
131130

132-
i = u.fx.lo - u.fc.lo
133-
s = u.fx.hi - u.fc.hi
131+
i = u.x.lo - u.c.lo
132+
s = u.x.hi - u.c.hi
134133

135134
if (i == 0 || s == 0)
136135
h1 = hx
137136
else
138137
h1 = interval((hx.lo - hc.lo) / i, (hx.hi - hc.hi) / s)
139138
end
140139

141-
SlopeType(hx, hc, h1 * u.fs)
140+
Slope{Float64}(hx, hc, h1 * u.s)
142141
end
143142

144-
function log(u::SlopeType)
145-
hx = log(u.fx)
146-
hc = log(u.fc)
143+
function log(u::Slope)
144+
hx = log(u.x)
145+
hc = log(u.c)
147146

148-
i = u.fx.lo - u.fc.lo
149-
s = u.fx.hi - u.fc.hi
147+
i = u.x.lo - u.c.lo
148+
s = u.x.hi - u.c.hi
150149

151150
if (i == 0 || s == 0)
152-
h1 = 1 / u.fx
151+
h1 = 1 / u.x
153152
else
154153
h1 = interval((hx.hi - hc.hi) / s, (hx.lo - hc.lo) / i)
155154
end
156-
SlopeType(hx, hc, h1 * u.fs)
155+
Slope{Float64}(hx, hc, h1 * u.s)
157156
end
158157

159-
function sin(u::SlopeType) # Using derivative to upper bound the slope expansion for now
160-
hx = sin(u.fx)
161-
hc = sin(u.fc)
162-
hs = cos(u.fx)
163-
SlopeType(hx, hc, hs)
158+
function sin(u::Slope) # Using derivative to upper bound the slope expansion for now
159+
hx = sin(u.x)
160+
hc = sin(u.c)
161+
hs = cos(u.x)
162+
Slope{Float64}(hx, hc, hs)
164163
end
165164

166-
function cos(u::SlopeType) # Using derivative to upper bound the slope expansion for now
167-
hx = cos(u.fx)
168-
hc = cos(u.fc)
169-
hs = -sin(u.fx)
170-
SlopeType(hx, hc, hs)
165+
function cos(u::Slope) # Using derivative to upper bound the slope expansion for now
166+
hx = cos(u.x)
167+
hc = cos(u.c)
168+
hs = -sin(u.x)
169+
Slope{Float64}(hx, hc, hs)
171170
end
172171

173-
function tan(u::SlopeType) # Using derivative to upper bound the slope expansion for now
174-
hx = tan(u.fx)
175-
hc = tan(u.fc)
176-
hs = (sec(u.fx)) ^ 2
177-
SlopeType(hx, hc, hs)
172+
function tan(u::Slope) # Using derivative to upper bound the slope expansion for now
173+
hx = tan(u.x)
174+
hc = tan(u.c)
175+
hs = (sec(u.x)) ^ 2
176+
Slope{Float64}(hx, hc, hs)
178177
end
179178

180-
function asin(u::SlopeType)
181-
hx = asin(u.fx)
182-
hc = asin(u.fc)
183-
hs = 1 / sqrt(1 - (u.fx ^ 2))
184-
SlopeType(hx, hc, hs)
179+
function asin(u::Slope)
180+
hx = asin(u.x)
181+
hc = asin(u.c)
182+
hs = 1 / sqrt(1 - (u.x ^ 2))
183+
Slope{Float64}(hx, hc, hs)
185184
end
186185

187-
function acos(u::SlopeType)
188-
hx = acos(u.fx)
189-
hc = acos(u.fc)
190-
hs = -1 / sqrt(1 - (u.fx ^ 2))
191-
SlopeType(hx, hc, hs)
186+
function acos(u::Slope)
187+
hx = acos(u.x)
188+
hc = acos(u.c)
189+
hs = -1 / sqrt(1 - (u.x ^ 2))
190+
Slope{Float64}(hx, hc, hs)
192191
end
193192

194-
function atan(u::SlopeType)
195-
hx = atan(u.fx)
196-
hc = atan(u.fc)
197-
hs = 1 / 1 + (u.fx ^ 2)
198-
SlopeType(hx, hc, hs)
193+
function atan(u::Slope)
194+
hx = atan(u.x)
195+
hc = atan(u.c)
196+
hs = 1 / 1 + (u.x ^ 2)
197+
Slope{Float64}(hx, hc, hs)
199198
end

‎test/slopes.jl

Lines changed: 9 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -10,14 +10,16 @@ struct Slopes{T}
1010
end
1111

1212
@testset "Automatic slope expansion" begin
13+
14+
s = interval(0.75, 1.75)
1315
rts = Slopes{Float64}[]
14-
push!(rts, Slopes(x->((x + sin(x)) * exp(-x^2)), interval(0.75, 1.75), 1.25, interval(-2.8, 0.1)))
15-
push!(rts, Slopes(x->(x^4 - 10x^3 + 35x^2 - 50x + 24), interval(0.75, 1.75), 1.25, interval(-44, 38.5)))
16-
push!(rts, Slopes(x->((log(x + 1.25) - 0.84x) ^ 2), interval(0.75, 1.75), 1.25, interval(-0.16, 0.44)))
17-
push!(rts, Slopes(x->(0.02x^2 - 0.03exp(-(20(x - 0.875))^2)), interval(0.75, 1.75), 1.25, interval(0.04, 0.33)))
18-
push!(rts, Slopes(x->(exp(x^2)), interval(0.75, 1.75), 1.25, interval(6.03, 33.23)))
19-
push!(rts, Slopes(x->(x^4 - 12x^3 + 47x^2 - 60x - 20exp(-x)), interval(0.75, 1.75), 1.25, interval(-39, 65.56)))
20-
push!(rts, Slopes(x->(x^6 - 15x^4 + 27x^2 + 250), interval(0.75, 1.75), 1.25, interval(-146.9, 67.1)))
16+
push!(rts, Slopes(x->((x + sin(x)) * exp(-x^2)), s, mid(s), interval(-2.8, 0.1)))
17+
push!(rts, Slopes(x->(x^4 - 10x^3 + 35x^2 - 50x + 24), s, mid(s), interval(-44, 38.5)))
18+
push!(rts, Slopes(x->((log(x + 1.25) - 0.84x) ^ 2), s, mid(s), interval(-0.16, 0.44)))
19+
push!(rts, Slopes(x->(0.02x^2 - 0.03exp(-(20(x - 0.875))^2)), s, mid(s), interval(0.04, 0.33)))
20+
push!(rts, Slopes(x->(exp(x^2)), s, mid(s), interval(6.03, 33.23)))
21+
push!(rts, Slopes(x->(x^4 - 12x^3 + 47x^2 - 60x - 20exp(-x)), s, mid(s), interval(-39, 65.56)))
22+
push!(rts, Slopes(x->(x^6 - 15x^4 + 27x^2 + 250), s, mid(s), interval(-146.9, 67.1)))
2123

2224
for i in 1:length(rts)
2325
@test slope(rts[i].f, rts[i].x, rts[i].c) ⊆ rts[i].sol

0 commit comments

Comments
 (0)