
(FPCore (a b) :precision binary64 (* (* (/ PI 2.0) (/ 1.0 (- (* b b) (* a a)))) (- (/ 1.0 a) (/ 1.0 b))))
double code(double a, double b) {
return ((((double) M_PI) / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b));
}
public static double code(double a, double b) {
return ((Math.PI / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b));
}
def code(a, b): return ((math.pi / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b))
function code(a, b) return Float64(Float64(Float64(pi / 2.0) * Float64(1.0 / Float64(Float64(b * b) - Float64(a * a)))) * Float64(Float64(1.0 / a) - Float64(1.0 / b))) end
function tmp = code(a, b) tmp = ((pi / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b)); end
code[a_, b_] := N[(N[(N[(Pi / 2.0), $MachinePrecision] * N[(1.0 / N[(N[(b * b), $MachinePrecision] - N[(a * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * N[(N[(1.0 / a), $MachinePrecision] - N[(1.0 / b), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\left(\frac{\pi}{2} \cdot \frac{1}{b \cdot b - a \cdot a}\right) \cdot \left(\frac{1}{a} - \frac{1}{b}\right)
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 4 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (a b) :precision binary64 (* (* (/ PI 2.0) (/ 1.0 (- (* b b) (* a a)))) (- (/ 1.0 a) (/ 1.0 b))))
double code(double a, double b) {
return ((((double) M_PI) / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b));
}
public static double code(double a, double b) {
return ((Math.PI / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b));
}
def code(a, b): return ((math.pi / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b))
function code(a, b) return Float64(Float64(Float64(pi / 2.0) * Float64(1.0 / Float64(Float64(b * b) - Float64(a * a)))) * Float64(Float64(1.0 / a) - Float64(1.0 / b))) end
function tmp = code(a, b) tmp = ((pi / 2.0) * (1.0 / ((b * b) - (a * a)))) * ((1.0 / a) - (1.0 / b)); end
code[a_, b_] := N[(N[(N[(Pi / 2.0), $MachinePrecision] * N[(1.0 / N[(N[(b * b), $MachinePrecision] - N[(a * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * N[(N[(1.0 / a), $MachinePrecision] - N[(1.0 / b), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\left(\frac{\pi}{2} \cdot \frac{1}{b \cdot b - a \cdot a}\right) \cdot \left(\frac{1}{a} - \frac{1}{b}\right)
\end{array}
NOTE: a and b should be sorted in increasing order before calling this function. (FPCore (a b) :precision binary64 (if (<= b 8.8e+113) (* (/ 0.5 a) (/ PI (* b (+ a b)))) (* (/ 0.5 (* a b)) (/ PI b))))
assert(a < b);
double code(double a, double b) {
double tmp;
if (b <= 8.8e+113) {
tmp = (0.5 / a) * (((double) M_PI) / (b * (a + b)));
} else {
tmp = (0.5 / (a * b)) * (((double) M_PI) / b);
}
return tmp;
}
assert a < b;
public static double code(double a, double b) {
double tmp;
if (b <= 8.8e+113) {
tmp = (0.5 / a) * (Math.PI / (b * (a + b)));
} else {
tmp = (0.5 / (a * b)) * (Math.PI / b);
}
return tmp;
}
[a, b] = sort([a, b]) def code(a, b): tmp = 0 if b <= 8.8e+113: tmp = (0.5 / a) * (math.pi / (b * (a + b))) else: tmp = (0.5 / (a * b)) * (math.pi / b) return tmp
a, b = sort([a, b]) function code(a, b) tmp = 0.0 if (b <= 8.8e+113) tmp = Float64(Float64(0.5 / a) * Float64(pi / Float64(b * Float64(a + b)))); else tmp = Float64(Float64(0.5 / Float64(a * b)) * Float64(pi / b)); end return tmp end
a, b = num2cell(sort([a, b])){:}
function tmp_2 = code(a, b)
tmp = 0.0;
if (b <= 8.8e+113)
tmp = (0.5 / a) * (pi / (b * (a + b)));
else
tmp = (0.5 / (a * b)) * (pi / b);
end
tmp_2 = tmp;
end
NOTE: a and b should be sorted in increasing order before calling this function. code[a_, b_] := If[LessEqual[b, 8.8e+113], N[(N[(0.5 / a), $MachinePrecision] * N[(Pi / N[(b * N[(a + b), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(N[(0.5 / N[(a * b), $MachinePrecision]), $MachinePrecision] * N[(Pi / b), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
[a, b] = \mathsf{sort}([a, b])\\
\\
\begin{array}{l}
\mathbf{if}\;b \leq 8.8 \cdot 10^{+113}:\\
\;\;\;\;\frac{0.5}{a} \cdot \frac{\pi}{b \cdot \left(a + b\right)}\\
\mathbf{else}:\\
\;\;\;\;\frac{0.5}{a \cdot b} \cdot \frac{\pi}{b}\\
\end{array}
\end{array}
if b < 8.80000000000000041e113Initial program 79.6%
*-commutative79.6%
sub-neg79.6%
distribute-neg-frac79.6%
metadata-eval79.6%
*-commutative79.6%
associate-*r/79.6%
associate-*l/79.7%
*-lft-identity79.7%
difference-of-squares86.9%
associate-/r*87.4%
Simplified87.4%
*-commutative87.4%
associate-/l/87.4%
frac-add87.4%
frac-times99.1%
+-commutative99.1%
*-un-lft-identity99.1%
Applied egg-rr99.1%
associate-/l*99.1%
+-commutative99.1%
fma-define99.1%
associate-*l*99.1%
Applied egg-rr99.1%
Taylor expanded in a around 0 99.6%
frac-times99.2%
Applied egg-rr99.2%
times-frac99.6%
*-commutative99.6%
associate-/l/99.6%
*-commutative99.6%
times-frac95.1%
associate-/l*95.0%
Simplified95.0%
if 8.80000000000000041e113 < b Initial program 75.7%
*-commutative75.7%
sub-neg75.7%
distribute-neg-frac75.7%
metadata-eval75.7%
*-commutative75.7%
associate-*r/75.7%
associate-*l/75.7%
*-lft-identity75.7%
difference-of-squares88.2%
associate-/r*88.1%
Simplified88.1%
*-commutative88.1%
associate-/l/88.1%
frac-add88.1%
frac-times99.8%
+-commutative99.8%
*-un-lft-identity99.8%
Applied egg-rr99.8%
associate-/l*99.8%
+-commutative99.8%
fma-define99.8%
associate-*l*99.8%
Applied egg-rr99.8%
Taylor expanded in a around 0 99.8%
Taylor expanded in a around 0 99.8%
Final simplification95.6%
NOTE: a and b should be sorted in increasing order before calling this function. (FPCore (a b) :precision binary64 (let* ((t_0 (/ 0.5 (* a b)))) (if (<= a -2.7e+14) (* t_0 (/ PI a)) (* t_0 (/ PI b)))))
assert(a < b);
double code(double a, double b) {
double t_0 = 0.5 / (a * b);
double tmp;
if (a <= -2.7e+14) {
tmp = t_0 * (((double) M_PI) / a);
} else {
tmp = t_0 * (((double) M_PI) / b);
}
return tmp;
}
assert a < b;
public static double code(double a, double b) {
double t_0 = 0.5 / (a * b);
double tmp;
if (a <= -2.7e+14) {
tmp = t_0 * (Math.PI / a);
} else {
tmp = t_0 * (Math.PI / b);
}
return tmp;
}
[a, b] = sort([a, b]) def code(a, b): t_0 = 0.5 / (a * b) tmp = 0 if a <= -2.7e+14: tmp = t_0 * (math.pi / a) else: tmp = t_0 * (math.pi / b) return tmp
a, b = sort([a, b]) function code(a, b) t_0 = Float64(0.5 / Float64(a * b)) tmp = 0.0 if (a <= -2.7e+14) tmp = Float64(t_0 * Float64(pi / a)); else tmp = Float64(t_0 * Float64(pi / b)); end return tmp end
a, b = num2cell(sort([a, b])){:}
function tmp_2 = code(a, b)
t_0 = 0.5 / (a * b);
tmp = 0.0;
if (a <= -2.7e+14)
tmp = t_0 * (pi / a);
else
tmp = t_0 * (pi / b);
end
tmp_2 = tmp;
end
NOTE: a and b should be sorted in increasing order before calling this function.
code[a_, b_] := Block[{t$95$0 = N[(0.5 / N[(a * b), $MachinePrecision]), $MachinePrecision]}, If[LessEqual[a, -2.7e+14], N[(t$95$0 * N[(Pi / a), $MachinePrecision]), $MachinePrecision], N[(t$95$0 * N[(Pi / b), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
[a, b] = \mathsf{sort}([a, b])\\
\\
\begin{array}{l}
t_0 := \frac{0.5}{a \cdot b}\\
\mathbf{if}\;a \leq -2.7 \cdot 10^{+14}:\\
\;\;\;\;t\_0 \cdot \frac{\pi}{a}\\
\mathbf{else}:\\
\;\;\;\;t\_0 \cdot \frac{\pi}{b}\\
\end{array}
\end{array}
if a < -2.7e14Initial program 83.3%
*-commutative83.3%
sub-neg83.3%
distribute-neg-frac83.3%
metadata-eval83.3%
*-commutative83.3%
associate-*r/83.3%
associate-*l/83.4%
*-lft-identity83.4%
difference-of-squares86.7%
associate-/r*88.5%
Simplified88.5%
*-commutative88.5%
associate-/l/88.6%
frac-add88.5%
frac-times98.9%
+-commutative98.9%
*-un-lft-identity98.9%
Applied egg-rr98.9%
associate-/l*98.9%
+-commutative98.9%
fma-define98.9%
associate-*l*98.9%
Applied egg-rr98.9%
Taylor expanded in a around 0 99.7%
Taylor expanded in a around inf 88.4%
if -2.7e14 < a Initial program 77.8%
*-commutative77.8%
sub-neg77.8%
distribute-neg-frac77.8%
metadata-eval77.8%
*-commutative77.8%
associate-*r/77.8%
associate-*l/77.9%
*-lft-identity77.9%
difference-of-squares87.1%
associate-/r*87.2%
Simplified87.2%
*-commutative87.2%
associate-/l/87.2%
frac-add87.1%
frac-times99.3%
+-commutative99.3%
*-un-lft-identity99.3%
Applied egg-rr99.3%
associate-/l*99.3%
+-commutative99.3%
fma-define99.3%
associate-*l*99.3%
Applied egg-rr99.3%
Taylor expanded in a around 0 99.6%
Taylor expanded in a around 0 66.4%
Final simplification71.6%
NOTE: a and b should be sorted in increasing order before calling this function. (FPCore (a b) :precision binary64 (* (/ PI (+ a b)) (/ 0.5 (* a b))))
assert(a < b);
double code(double a, double b) {
return (((double) M_PI) / (a + b)) * (0.5 / (a * b));
}
assert a < b;
public static double code(double a, double b) {
return (Math.PI / (a + b)) * (0.5 / (a * b));
}
[a, b] = sort([a, b]) def code(a, b): return (math.pi / (a + b)) * (0.5 / (a * b))
a, b = sort([a, b]) function code(a, b) return Float64(Float64(pi / Float64(a + b)) * Float64(0.5 / Float64(a * b))) end
a, b = num2cell(sort([a, b])){:}
function tmp = code(a, b)
tmp = (pi / (a + b)) * (0.5 / (a * b));
end
NOTE: a and b should be sorted in increasing order before calling this function. code[a_, b_] := N[(N[(Pi / N[(a + b), $MachinePrecision]), $MachinePrecision] * N[(0.5 / N[(a * b), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
[a, b] = \mathsf{sort}([a, b])\\
\\
\frac{\pi}{a + b} \cdot \frac{0.5}{a \cdot b}
\end{array}
Initial program 79.1%
*-commutative79.1%
sub-neg79.1%
distribute-neg-frac79.1%
metadata-eval79.1%
*-commutative79.1%
associate-*r/79.1%
associate-*l/79.2%
*-lft-identity79.2%
difference-of-squares87.0%
associate-/r*87.5%
Simplified87.5%
*-commutative87.5%
associate-/l/87.5%
frac-add87.5%
frac-times99.2%
+-commutative99.2%
*-un-lft-identity99.2%
Applied egg-rr99.2%
associate-/l*99.2%
+-commutative99.2%
fma-define99.2%
associate-*l*99.2%
Applied egg-rr99.2%
Taylor expanded in a around 0 99.6%
Final simplification99.6%
NOTE: a and b should be sorted in increasing order before calling this function. (FPCore (a b) :precision binary64 (* (/ 0.5 (* a b)) (/ PI a)))
assert(a < b);
double code(double a, double b) {
return (0.5 / (a * b)) * (((double) M_PI) / a);
}
assert a < b;
public static double code(double a, double b) {
return (0.5 / (a * b)) * (Math.PI / a);
}
[a, b] = sort([a, b]) def code(a, b): return (0.5 / (a * b)) * (math.pi / a)
a, b = sort([a, b]) function code(a, b) return Float64(Float64(0.5 / Float64(a * b)) * Float64(pi / a)) end
a, b = num2cell(sort([a, b])){:}
function tmp = code(a, b)
tmp = (0.5 / (a * b)) * (pi / a);
end
NOTE: a and b should be sorted in increasing order before calling this function. code[a_, b_] := N[(N[(0.5 / N[(a * b), $MachinePrecision]), $MachinePrecision] * N[(Pi / a), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
[a, b] = \mathsf{sort}([a, b])\\
\\
\frac{0.5}{a \cdot b} \cdot \frac{\pi}{a}
\end{array}
Initial program 79.1%
*-commutative79.1%
sub-neg79.1%
distribute-neg-frac79.1%
metadata-eval79.1%
*-commutative79.1%
associate-*r/79.1%
associate-*l/79.2%
*-lft-identity79.2%
difference-of-squares87.0%
associate-/r*87.5%
Simplified87.5%
*-commutative87.5%
associate-/l/87.5%
frac-add87.5%
frac-times99.2%
+-commutative99.2%
*-un-lft-identity99.2%
Applied egg-rr99.2%
associate-/l*99.2%
+-commutative99.2%
fma-define99.2%
associate-*l*99.2%
Applied egg-rr99.2%
Taylor expanded in a around 0 99.6%
Taylor expanded in a around inf 64.4%
Final simplification64.4%
herbie shell --seed 2024048
(FPCore (a b)
:name "NMSE Section 6.1 mentioned, B"
:precision binary64
(* (* (/ PI 2.0) (/ 1.0 (- (* b b) (* a a)))) (- (/ 1.0 a) (/ 1.0 b))))