
(FPCore (x) :precision binary64 (- (sqrt (+ 1.0 x)) (sqrt (- 1.0 x))))
double code(double x) {
return sqrt((1.0 + x)) - sqrt((1.0 - x));
}
real(8) function code(x)
real(8), intent (in) :: x
code = sqrt((1.0d0 + x)) - sqrt((1.0d0 - x))
end function
public static double code(double x) {
return Math.sqrt((1.0 + x)) - Math.sqrt((1.0 - x));
}
def code(x): return math.sqrt((1.0 + x)) - math.sqrt((1.0 - x))
function code(x) return Float64(sqrt(Float64(1.0 + x)) - sqrt(Float64(1.0 - x))) end
function tmp = code(x) tmp = sqrt((1.0 + x)) - sqrt((1.0 - x)); end
code[x_] := N[(N[Sqrt[N[(1.0 + x), $MachinePrecision]], $MachinePrecision] - N[Sqrt[N[(1.0 - x), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\sqrt{1 + x} - \sqrt{1 - x}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 4 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (x) :precision binary64 (- (sqrt (+ 1.0 x)) (sqrt (- 1.0 x))))
double code(double x) {
return sqrt((1.0 + x)) - sqrt((1.0 - x));
}
real(8) function code(x)
real(8), intent (in) :: x
code = sqrt((1.0d0 + x)) - sqrt((1.0d0 - x))
end function
public static double code(double x) {
return Math.sqrt((1.0 + x)) - Math.sqrt((1.0 - x));
}
def code(x): return math.sqrt((1.0 + x)) - math.sqrt((1.0 - x))
function code(x) return Float64(sqrt(Float64(1.0 + x)) - sqrt(Float64(1.0 - x))) end
function tmp = code(x) tmp = sqrt((1.0 + x)) - sqrt((1.0 - x)); end
code[x_] := N[(N[Sqrt[N[(1.0 + x), $MachinePrecision]], $MachinePrecision] - N[Sqrt[N[(1.0 - x), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\sqrt{1 + x} - \sqrt{1 - x}
\end{array}
(FPCore (x) :precision binary64 (+ x (+ (* 0.0546875 (pow x 5.0)) (* 0.125 (pow x 3.0)))))
double code(double x) {
return x + ((0.0546875 * pow(x, 5.0)) + (0.125 * pow(x, 3.0)));
}
real(8) function code(x)
real(8), intent (in) :: x
code = x + ((0.0546875d0 * (x ** 5.0d0)) + (0.125d0 * (x ** 3.0d0)))
end function
public static double code(double x) {
return x + ((0.0546875 * Math.pow(x, 5.0)) + (0.125 * Math.pow(x, 3.0)));
}
def code(x): return x + ((0.0546875 * math.pow(x, 5.0)) + (0.125 * math.pow(x, 3.0)))
function code(x) return Float64(x + Float64(Float64(0.0546875 * (x ^ 5.0)) + Float64(0.125 * (x ^ 3.0)))) end
function tmp = code(x) tmp = x + ((0.0546875 * (x ^ 5.0)) + (0.125 * (x ^ 3.0))); end
code[x_] := N[(x + N[(N[(0.0546875 * N[Power[x, 5.0], $MachinePrecision]), $MachinePrecision] + N[(0.125 * N[Power[x, 3.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
x + \left(0.0546875 \cdot {x}^{5} + 0.125 \cdot {x}^{3}\right)
\end{array}
Initial program 8.6%
Taylor expanded in x around 0 99.9%
Final simplification99.9%
(FPCore (x) :precision binary64 (/ (* x 2.0) (fma -0.25 (pow x 2.0) 2.0)))
double code(double x) {
return (x * 2.0) / fma(-0.25, pow(x, 2.0), 2.0);
}
function code(x) return Float64(Float64(x * 2.0) / fma(-0.25, (x ^ 2.0), 2.0)) end
code[x_] := N[(N[(x * 2.0), $MachinePrecision] / N[(-0.25 * N[Power[x, 2.0], $MachinePrecision] + 2.0), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\frac{x \cdot 2}{\mathsf{fma}\left(-0.25, {x}^{2}, 2\right)}
\end{array}
Initial program 8.6%
flip--8.6%
div-inv8.6%
add-sqr-sqrt8.6%
add-sqr-sqrt8.7%
associate--r-21.2%
add-exp-log21.2%
log1p-udef21.2%
expm1-udef100.0%
expm1-log1p-u100.0%
Applied egg-rr100.0%
associate-*r/100.0%
*-rgt-identity100.0%
count-2100.0%
associate-/l*99.7%
remove-double-neg99.7%
sub-neg99.7%
sub-neg99.7%
remove-double-neg99.7%
+-commutative99.7%
Simplified99.7%
Taylor expanded in x around 0 99.5%
add-log-exp8.5%
*-un-lft-identity8.5%
log-prod8.5%
metadata-eval8.5%
div-inv8.5%
add-log-exp99.5%
clear-num99.8%
+-commutative99.8%
fma-def99.8%
Applied egg-rr99.8%
+-lft-identity99.8%
*-commutative99.8%
associate-*l/99.8%
Simplified99.8%
Final simplification99.8%
(FPCore (x) :precision binary64 (+ x (* 0.125 (pow x 3.0))))
double code(double x) {
return x + (0.125 * pow(x, 3.0));
}
real(8) function code(x)
real(8), intent (in) :: x
code = x + (0.125d0 * (x ** 3.0d0))
end function
public static double code(double x) {
return x + (0.125 * Math.pow(x, 3.0));
}
def code(x): return x + (0.125 * math.pow(x, 3.0))
function code(x) return Float64(x + Float64(0.125 * (x ^ 3.0))) end
function tmp = code(x) tmp = x + (0.125 * (x ^ 3.0)); end
code[x_] := N[(x + N[(0.125 * N[Power[x, 3.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
x + 0.125 \cdot {x}^{3}
\end{array}
Initial program 8.6%
Taylor expanded in x around 0 99.8%
Final simplification99.8%
(FPCore (x) :precision binary64 x)
double code(double x) {
return x;
}
real(8) function code(x)
real(8), intent (in) :: x
code = x
end function
public static double code(double x) {
return x;
}
def code(x): return x
function code(x) return x end
function tmp = code(x) tmp = x; end
code[x_] := x
\begin{array}{l}
\\
x
\end{array}
Initial program 8.6%
Taylor expanded in x around 0 98.9%
Final simplification98.9%
(FPCore (x) :precision binary64 (/ (* 2.0 x) (+ (sqrt (+ 1.0 x)) (sqrt (- 1.0 x)))))
double code(double x) {
return (2.0 * x) / (sqrt((1.0 + x)) + sqrt((1.0 - x)));
}
real(8) function code(x)
real(8), intent (in) :: x
code = (2.0d0 * x) / (sqrt((1.0d0 + x)) + sqrt((1.0d0 - x)))
end function
public static double code(double x) {
return (2.0 * x) / (Math.sqrt((1.0 + x)) + Math.sqrt((1.0 - x)));
}
def code(x): return (2.0 * x) / (math.sqrt((1.0 + x)) + math.sqrt((1.0 - x)))
function code(x) return Float64(Float64(2.0 * x) / Float64(sqrt(Float64(1.0 + x)) + sqrt(Float64(1.0 - x)))) end
function tmp = code(x) tmp = (2.0 * x) / (sqrt((1.0 + x)) + sqrt((1.0 - x))); end
code[x_] := N[(N[(2.0 * x), $MachinePrecision] / N[(N[Sqrt[N[(1.0 + x), $MachinePrecision]], $MachinePrecision] + N[Sqrt[N[(1.0 - x), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\frac{2 \cdot x}{\sqrt{1 + x} + \sqrt{1 - x}}
\end{array}
herbie shell --seed 2024020
(FPCore (x)
:name "bug333 (missed optimization)"
:precision binary64
:pre (and (<= -1.0 x) (<= x 1.0))
:herbie-target
(/ (* 2.0 x) (+ (sqrt (+ 1.0 x)) (sqrt (- 1.0 x))))
(- (sqrt (+ 1.0 x)) (sqrt (- 1.0 x))))