
(FPCore (p x) :precision binary64 (sqrt (* 0.5 (+ 1.0 (/ x (sqrt (+ (* (* 4.0 p) p) (* x x))))))))
double code(double p, double x) {
return sqrt((0.5 * (1.0 + (x / sqrt((((4.0 * p) * p) + (x * x)))))));
}
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
code = sqrt((0.5d0 * (1.0d0 + (x / sqrt((((4.0d0 * p) * p) + (x * x)))))))
end function
public static double code(double p, double x) {
return Math.sqrt((0.5 * (1.0 + (x / Math.sqrt((((4.0 * p) * p) + (x * x)))))));
}
def code(p, x): return math.sqrt((0.5 * (1.0 + (x / math.sqrt((((4.0 * p) * p) + (x * x)))))))
function code(p, x) return sqrt(Float64(0.5 * Float64(1.0 + Float64(x / sqrt(Float64(Float64(Float64(4.0 * p) * p) + Float64(x * x))))))) end
function tmp = code(p, x) tmp = sqrt((0.5 * (1.0 + (x / sqrt((((4.0 * p) * p) + (x * x))))))); end
code[p_, x_] := N[Sqrt[N[(0.5 * N[(1.0 + N[(x / N[Sqrt[N[(N[(N[(4.0 * p), $MachinePrecision] * p), $MachinePrecision] + N[(x * x), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\sqrt{0.5 \cdot \left(1 + \frac{x}{\sqrt{\left(4 \cdot p\right) \cdot p + x \cdot x}}\right)}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 6 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (p x) :precision binary64 (sqrt (* 0.5 (+ 1.0 (/ x (sqrt (+ (* (* 4.0 p) p) (* x x))))))))
double code(double p, double x) {
return sqrt((0.5 * (1.0 + (x / sqrt((((4.0 * p) * p) + (x * x)))))));
}
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
code = sqrt((0.5d0 * (1.0d0 + (x / sqrt((((4.0d0 * p) * p) + (x * x)))))))
end function
public static double code(double p, double x) {
return Math.sqrt((0.5 * (1.0 + (x / Math.sqrt((((4.0 * p) * p) + (x * x)))))));
}
def code(p, x): return math.sqrt((0.5 * (1.0 + (x / math.sqrt((((4.0 * p) * p) + (x * x)))))))
function code(p, x) return sqrt(Float64(0.5 * Float64(1.0 + Float64(x / sqrt(Float64(Float64(Float64(4.0 * p) * p) + Float64(x * x))))))) end
function tmp = code(p, x) tmp = sqrt((0.5 * (1.0 + (x / sqrt((((4.0 * p) * p) + (x * x))))))); end
code[p_, x_] := N[Sqrt[N[(0.5 * N[(1.0 + N[(x / N[Sqrt[N[(N[(N[(4.0 * p), $MachinePrecision] * p), $MachinePrecision] + N[(x * x), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\sqrt{0.5 \cdot \left(1 + \frac{x}{\sqrt{\left(4 \cdot p\right) \cdot p + x \cdot x}}\right)}
\end{array}
(FPCore (p x) :precision binary64 (sqrt (* 0.5 (exp (log1p (/ x (hypot x (* p 2.0))))))))
double code(double p, double x) {
return sqrt((0.5 * exp(log1p((x / hypot(x, (p * 2.0)))))));
}
public static double code(double p, double x) {
return Math.sqrt((0.5 * Math.exp(Math.log1p((x / Math.hypot(x, (p * 2.0)))))));
}
def code(p, x): return math.sqrt((0.5 * math.exp(math.log1p((x / math.hypot(x, (p * 2.0)))))))
function code(p, x) return sqrt(Float64(0.5 * exp(log1p(Float64(x / hypot(x, Float64(p * 2.0))))))) end
code[p_, x_] := N[Sqrt[N[(0.5 * N[Exp[N[Log[1 + N[(x / N[Sqrt[x ^ 2 + N[(p * 2.0), $MachinePrecision] ^ 2], $MachinePrecision]), $MachinePrecision]], $MachinePrecision]], $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\sqrt{0.5 \cdot e^{\mathsf{log1p}\left(\frac{x}{\mathsf{hypot}\left(x, p \cdot 2\right)}\right)}}
\end{array}
Initial program 98.4%
add-exp-log98.4%
log1p-define98.4%
div-inv95.1%
+-commutative95.1%
associate-*r*95.1%
fma-undefine95.2%
div-inv98.4%
fma-undefine98.4%
associate-*r*98.4%
add-sqr-sqrt98.4%
hypot-define98.5%
associate-*r*98.5%
*-commutative98.5%
sqrt-prod98.5%
sqrt-prod48.7%
add-sqr-sqrt98.5%
metadata-eval98.5%
Applied egg-rr98.5%
Final simplification98.5%
(FPCore (p x) :precision binary64 (pow (pow (* 0.5 (+ (/ x (hypot x (* p 2.0))) 1.0)) 1.5) 0.3333333333333333))
double code(double p, double x) {
return pow(pow((0.5 * ((x / hypot(x, (p * 2.0))) + 1.0)), 1.5), 0.3333333333333333);
}
public static double code(double p, double x) {
return Math.pow(Math.pow((0.5 * ((x / Math.hypot(x, (p * 2.0))) + 1.0)), 1.5), 0.3333333333333333);
}
def code(p, x): return math.pow(math.pow((0.5 * ((x / math.hypot(x, (p * 2.0))) + 1.0)), 1.5), 0.3333333333333333)
function code(p, x) return (Float64(0.5 * Float64(Float64(x / hypot(x, Float64(p * 2.0))) + 1.0)) ^ 1.5) ^ 0.3333333333333333 end
function tmp = code(p, x) tmp = ((0.5 * ((x / hypot(x, (p * 2.0))) + 1.0)) ^ 1.5) ^ 0.3333333333333333; end
code[p_, x_] := N[Power[N[Power[N[(0.5 * N[(N[(x / N[Sqrt[x ^ 2 + N[(p * 2.0), $MachinePrecision] ^ 2], $MachinePrecision]), $MachinePrecision] + 1.0), $MachinePrecision]), $MachinePrecision], 1.5], $MachinePrecision], 0.3333333333333333], $MachinePrecision]
\begin{array}{l}
\\
{\left({\left(0.5 \cdot \left(\frac{x}{\mathsf{hypot}\left(x, p \cdot 2\right)} + 1\right)\right)}^{1.5}\right)}^{0.3333333333333333}
\end{array}
Initial program 98.4%
add-exp-log98.4%
log1p-define98.4%
div-inv95.1%
+-commutative95.1%
associate-*r*95.1%
fma-undefine95.2%
div-inv98.4%
fma-undefine98.4%
associate-*r*98.4%
add-sqr-sqrt98.4%
hypot-define98.5%
associate-*r*98.5%
*-commutative98.5%
sqrt-prod98.5%
sqrt-prod48.7%
add-sqr-sqrt98.5%
metadata-eval98.5%
Applied egg-rr98.5%
add-cbrt-cube98.4%
pow1/398.5%
Applied egg-rr98.5%
Final simplification98.5%
(FPCore (p x) :precision binary64 (sqrt (+ 0.5 (/ (* 0.5 x) (hypot x (* p 2.0))))))
double code(double p, double x) {
return sqrt((0.5 + ((0.5 * x) / hypot(x, (p * 2.0)))));
}
public static double code(double p, double x) {
return Math.sqrt((0.5 + ((0.5 * x) / Math.hypot(x, (p * 2.0)))));
}
def code(p, x): return math.sqrt((0.5 + ((0.5 * x) / math.hypot(x, (p * 2.0)))))
function code(p, x) return sqrt(Float64(0.5 + Float64(Float64(0.5 * x) / hypot(x, Float64(p * 2.0))))) end
function tmp = code(p, x) tmp = sqrt((0.5 + ((0.5 * x) / hypot(x, (p * 2.0))))); end
code[p_, x_] := N[Sqrt[N[(0.5 + N[(N[(0.5 * x), $MachinePrecision] / N[Sqrt[x ^ 2 + N[(p * 2.0), $MachinePrecision] ^ 2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\sqrt{0.5 + \frac{0.5 \cdot x}{\mathsf{hypot}\left(x, p \cdot 2\right)}}
\end{array}
Initial program 98.4%
+-commutative98.4%
sqr-neg98.4%
associate-*l*98.4%
sqr-neg98.4%
fma-define98.4%
sqr-neg98.4%
fma-define98.4%
associate-*l*98.4%
+-commutative98.4%
Simplified98.4%
*-commutative98.4%
associate-*l/98.4%
fma-undefine98.4%
associate-*r*98.4%
add-sqr-sqrt98.4%
hypot-define98.5%
associate-*r*98.5%
*-commutative98.5%
sqrt-prod98.5%
sqrt-prod48.7%
add-sqr-sqrt98.5%
metadata-eval98.5%
Applied egg-rr98.5%
Final simplification98.5%
(FPCore (p x) :precision binary64 (if (<= p 7e-7) (+ (- 1.0 (/ p x)) -1.0) (sqrt 0.5)))
double code(double p, double x) {
double tmp;
if (p <= 7e-7) {
tmp = (1.0 - (p / x)) + -1.0;
} else {
tmp = sqrt(0.5);
}
return tmp;
}
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
real(8) :: tmp
if (p <= 7d-7) then
tmp = (1.0d0 - (p / x)) + (-1.0d0)
else
tmp = sqrt(0.5d0)
end if
code = tmp
end function
public static double code(double p, double x) {
double tmp;
if (p <= 7e-7) {
tmp = (1.0 - (p / x)) + -1.0;
} else {
tmp = Math.sqrt(0.5);
}
return tmp;
}
def code(p, x): tmp = 0 if p <= 7e-7: tmp = (1.0 - (p / x)) + -1.0 else: tmp = math.sqrt(0.5) return tmp
function code(p, x) tmp = 0.0 if (p <= 7e-7) tmp = Float64(Float64(1.0 - Float64(p / x)) + -1.0); else tmp = sqrt(0.5); end return tmp end
function tmp_2 = code(p, x) tmp = 0.0; if (p <= 7e-7) tmp = (1.0 - (p / x)) + -1.0; else tmp = sqrt(0.5); end tmp_2 = tmp; end
code[p_, x_] := If[LessEqual[p, 7e-7], N[(N[(1.0 - N[(p / x), $MachinePrecision]), $MachinePrecision] + -1.0), $MachinePrecision], N[Sqrt[0.5], $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;p \leq 7 \cdot 10^{-7}:\\
\;\;\;\;\left(1 - \frac{p}{x}\right) + -1\\
\mathbf{else}:\\
\;\;\;\;\sqrt{0.5}\\
\end{array}
\end{array}
if p < 6.99999999999999968e-7Initial program 97.9%
+-commutative97.9%
sqr-neg97.9%
associate-*l*97.9%
sqr-neg97.9%
fma-define97.9%
sqr-neg97.9%
fma-define97.9%
associate-*l*97.9%
+-commutative97.9%
Simplified97.9%
*-commutative97.9%
associate-*l/97.9%
fma-undefine97.9%
associate-*r*97.9%
add-sqr-sqrt97.9%
hypot-define98.0%
associate-*r*98.0%
*-commutative98.0%
sqrt-prod98.0%
sqrt-prod29.9%
add-sqr-sqrt98.0%
metadata-eval98.0%
Applied egg-rr98.0%
expm1-log1p-u97.4%
expm1-undefine97.4%
+-commutative97.4%
associate-/l*93.4%
fma-define66.3%
Applied egg-rr66.3%
Taylor expanded in x around -inf 34.3%
mul-1-neg34.3%
unsub-neg34.3%
Simplified34.3%
if 6.99999999999999968e-7 < p Initial program 99.8%
Taylor expanded in x around 0 93.1%
Final simplification50.2%
(FPCore (p x) :precision binary64 (+ (- 1.0 (/ p x)) -1.0))
double code(double p, double x) {
return (1.0 - (p / x)) + -1.0;
}
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
code = (1.0d0 - (p / x)) + (-1.0d0)
end function
public static double code(double p, double x) {
return (1.0 - (p / x)) + -1.0;
}
def code(p, x): return (1.0 - (p / x)) + -1.0
function code(p, x) return Float64(Float64(1.0 - Float64(p / x)) + -1.0) end
function tmp = code(p, x) tmp = (1.0 - (p / x)) + -1.0; end
code[p_, x_] := N[(N[(1.0 - N[(p / x), $MachinePrecision]), $MachinePrecision] + -1.0), $MachinePrecision]
\begin{array}{l}
\\
\left(1 - \frac{p}{x}\right) + -1
\end{array}
Initial program 98.4%
+-commutative98.4%
sqr-neg98.4%
associate-*l*98.4%
sqr-neg98.4%
fma-define98.4%
sqr-neg98.4%
fma-define98.4%
associate-*l*98.4%
+-commutative98.4%
Simplified98.4%
*-commutative98.4%
associate-*l/98.4%
fma-undefine98.4%
associate-*r*98.4%
add-sqr-sqrt98.4%
hypot-define98.5%
associate-*r*98.5%
*-commutative98.5%
sqrt-prod98.5%
sqrt-prod48.7%
add-sqr-sqrt98.5%
metadata-eval98.5%
Applied egg-rr98.5%
expm1-log1p-u97.7%
expm1-undefine97.7%
+-commutative97.7%
associate-/l*94.4%
fma-define74.2%
Applied egg-rr74.2%
Taylor expanded in x around -inf 26.8%
mul-1-neg26.8%
unsub-neg26.8%
Simplified26.8%
Final simplification26.8%
(FPCore (p x) :precision binary64 (/ p (- x)))
double code(double p, double x) {
return p / -x;
}
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
code = p / -x
end function
public static double code(double p, double x) {
return p / -x;
}
def code(p, x): return p / -x
function code(p, x) return Float64(p / Float64(-x)) end
function tmp = code(p, x) tmp = p / -x; end
code[p_, x_] := N[(p / (-x)), $MachinePrecision]
\begin{array}{l}
\\
\frac{p}{-x}
\end{array}
Initial program 98.4%
+-commutative98.4%
sqr-neg98.4%
associate-*l*98.4%
sqr-neg98.4%
fma-define98.4%
sqr-neg98.4%
fma-define98.4%
associate-*l*98.4%
+-commutative98.4%
Simplified98.4%
Taylor expanded in x around -inf 6.5%
mul-1-neg6.5%
distribute-neg-frac26.5%
Simplified6.5%
Final simplification6.5%
(FPCore (p x) :precision binary64 (sqrt (+ 0.5 (/ (copysign 0.5 x) (hypot 1.0 (/ (* 2.0 p) x))))))
double code(double p, double x) {
return sqrt((0.5 + (copysign(0.5, x) / hypot(1.0, ((2.0 * p) / x)))));
}
public static double code(double p, double x) {
return Math.sqrt((0.5 + (Math.copySign(0.5, x) / Math.hypot(1.0, ((2.0 * p) / x)))));
}
def code(p, x): return math.sqrt((0.5 + (math.copysign(0.5, x) / math.hypot(1.0, ((2.0 * p) / x)))))
function code(p, x) return sqrt(Float64(0.5 + Float64(copysign(0.5, x) / hypot(1.0, Float64(Float64(2.0 * p) / x))))) end
function tmp = code(p, x) tmp = sqrt((0.5 + ((sign(x) * abs(0.5)) / hypot(1.0, ((2.0 * p) / x))))); end
code[p_, x_] := N[Sqrt[N[(0.5 + N[(N[With[{TMP1 = Abs[0.5], TMP2 = Sign[x]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision] / N[Sqrt[1.0 ^ 2 + N[(N[(2.0 * p), $MachinePrecision] / x), $MachinePrecision] ^ 2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\sqrt{0.5 + \frac{\mathsf{copysign}\left(0.5, x\right)}{\mathsf{hypot}\left(1, \frac{2 \cdot p}{x}\right)}}
\end{array}
herbie shell --seed 2024045
(FPCore (p x)
:name "Given's Rotation SVD example"
:precision binary64
:pre (and (< 1e-150 (fabs x)) (< (fabs x) 1e+150))
:herbie-target
(sqrt (+ 0.5 (/ (copysign 0.5 x) (hypot 1.0 (/ (* 2.0 p) x)))))
(sqrt (* 0.5 (+ 1.0 (/ x (sqrt (+ (* (* 4.0 p) p) (* x x))))))))