
(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 5 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}
NOTE: p should be positive before calling this function (FPCore (p x) :precision binary64 (if (<= (/ x (sqrt (+ (* p (* 4.0 p)) (* x x)))) -0.6) (fma (pow (/ p x) 3.0) 1.5 (/ (- p) x)) (sqrt (* (+ 1.0 (/ x (hypot x (* p 2.0)))) 0.5))))
p = abs(p);
double code(double p, double x) {
double tmp;
if ((x / sqrt(((p * (4.0 * p)) + (x * x)))) <= -0.6) {
tmp = fma(pow((p / x), 3.0), 1.5, (-p / x));
} else {
tmp = sqrt(((1.0 + (x / hypot(x, (p * 2.0)))) * 0.5));
}
return tmp;
}
p = abs(p) function code(p, x) tmp = 0.0 if (Float64(x / sqrt(Float64(Float64(p * Float64(4.0 * p)) + Float64(x * x)))) <= -0.6) tmp = fma((Float64(p / x) ^ 3.0), 1.5, Float64(Float64(-p) / x)); else tmp = sqrt(Float64(Float64(1.0 + Float64(x / hypot(x, Float64(p * 2.0)))) * 0.5)); end return tmp end
NOTE: p should be positive before calling this function code[p_, x_] := If[LessEqual[N[(x / N[Sqrt[N[(N[(p * N[(4.0 * p), $MachinePrecision]), $MachinePrecision] + N[(x * x), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], -0.6], N[(N[Power[N[(p / x), $MachinePrecision], 3.0], $MachinePrecision] * 1.5 + N[((-p) / x), $MachinePrecision]), $MachinePrecision], N[Sqrt[N[(N[(1.0 + N[(x / N[Sqrt[x ^ 2 + N[(p * 2.0), $MachinePrecision] ^ 2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * 0.5), $MachinePrecision]], $MachinePrecision]]
\begin{array}{l}
p = |p|\\
\\
\begin{array}{l}
\mathbf{if}\;\frac{x}{\sqrt{p \cdot \left(4 \cdot p\right) + x \cdot x}} \leq -0.6:\\
\;\;\;\;\mathsf{fma}\left({\left(\frac{p}{x}\right)}^{3}, 1.5, \frac{-p}{x}\right)\\
\mathbf{else}:\\
\;\;\;\;\sqrt{\left(1 + \frac{x}{\mathsf{hypot}\left(x, p \cdot 2\right)}\right) \cdot 0.5}\\
\end{array}
\end{array}
if (/.f64 x (sqrt.f64 (+.f64 (*.f64 (*.f64 4 p) p) (*.f64 x x)))) < -0.599999999999999978Initial program 18.3%
distribute-lft-in18.3%
metadata-eval18.3%
fma-def18.3%
Simplified18.3%
*-commutative18.3%
fma-def18.3%
distribute-rgt1-in18.3%
+-commutative18.3%
Applied egg-rr18.3%
Taylor expanded in x around -inf 41.6%
+-commutative41.6%
fma-def41.6%
Simplified57.5%
Taylor expanded in x around -inf 54.9%
+-commutative54.9%
*-commutative54.9%
fma-def54.9%
cube-div59.5%
associate-*r/59.5%
neg-mul-159.5%
Simplified59.5%
if -0.599999999999999978 < (/.f64 x (sqrt.f64 (+.f64 (*.f64 (*.f64 4 p) p) (*.f64 x x)))) Initial program 100.0%
distribute-lft-in100.0%
metadata-eval100.0%
fma-def100.0%
Simplified100.0%
*-commutative100.0%
fma-def100.0%
distribute-rgt1-in100.0%
+-commutative100.0%
Applied egg-rr100.0%
Final simplification89.7%
NOTE: p should be positive before calling this function
(FPCore (p x)
:precision binary64
(let* ((t_0 (/ (- p) x)))
(if (<= p 7.6e-200)
t_0
(if (<= p 4.3e-185)
1.0
(if (<= p 2.5e-154) t_0 (if (<= p 1.06e-54) 1.0 (sqrt 0.5)))))))p = abs(p);
double code(double p, double x) {
double t_0 = -p / x;
double tmp;
if (p <= 7.6e-200) {
tmp = t_0;
} else if (p <= 4.3e-185) {
tmp = 1.0;
} else if (p <= 2.5e-154) {
tmp = t_0;
} else if (p <= 1.06e-54) {
tmp = 1.0;
} else {
tmp = sqrt(0.5);
}
return tmp;
}
NOTE: p should be positive before calling this function
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
real(8) :: t_0
real(8) :: tmp
t_0 = -p / x
if (p <= 7.6d-200) then
tmp = t_0
else if (p <= 4.3d-185) then
tmp = 1.0d0
else if (p <= 2.5d-154) then
tmp = t_0
else if (p <= 1.06d-54) then
tmp = 1.0d0
else
tmp = sqrt(0.5d0)
end if
code = tmp
end function
p = Math.abs(p);
public static double code(double p, double x) {
double t_0 = -p / x;
double tmp;
if (p <= 7.6e-200) {
tmp = t_0;
} else if (p <= 4.3e-185) {
tmp = 1.0;
} else if (p <= 2.5e-154) {
tmp = t_0;
} else if (p <= 1.06e-54) {
tmp = 1.0;
} else {
tmp = Math.sqrt(0.5);
}
return tmp;
}
p = abs(p) def code(p, x): t_0 = -p / x tmp = 0 if p <= 7.6e-200: tmp = t_0 elif p <= 4.3e-185: tmp = 1.0 elif p <= 2.5e-154: tmp = t_0 elif p <= 1.06e-54: tmp = 1.0 else: tmp = math.sqrt(0.5) return tmp
p = abs(p) function code(p, x) t_0 = Float64(Float64(-p) / x) tmp = 0.0 if (p <= 7.6e-200) tmp = t_0; elseif (p <= 4.3e-185) tmp = 1.0; elseif (p <= 2.5e-154) tmp = t_0; elseif (p <= 1.06e-54) tmp = 1.0; else tmp = sqrt(0.5); end return tmp end
p = abs(p) function tmp_2 = code(p, x) t_0 = -p / x; tmp = 0.0; if (p <= 7.6e-200) tmp = t_0; elseif (p <= 4.3e-185) tmp = 1.0; elseif (p <= 2.5e-154) tmp = t_0; elseif (p <= 1.06e-54) tmp = 1.0; else tmp = sqrt(0.5); end tmp_2 = tmp; end
NOTE: p should be positive before calling this function
code[p_, x_] := Block[{t$95$0 = N[((-p) / x), $MachinePrecision]}, If[LessEqual[p, 7.6e-200], t$95$0, If[LessEqual[p, 4.3e-185], 1.0, If[LessEqual[p, 2.5e-154], t$95$0, If[LessEqual[p, 1.06e-54], 1.0, N[Sqrt[0.5], $MachinePrecision]]]]]]
\begin{array}{l}
p = |p|\\
\\
\begin{array}{l}
t_0 := \frac{-p}{x}\\
\mathbf{if}\;p \leq 7.6 \cdot 10^{-200}:\\
\;\;\;\;t_0\\
\mathbf{elif}\;p \leq 4.3 \cdot 10^{-185}:\\
\;\;\;\;1\\
\mathbf{elif}\;p \leq 2.5 \cdot 10^{-154}:\\
\;\;\;\;t_0\\
\mathbf{elif}\;p \leq 1.06 \cdot 10^{-54}:\\
\;\;\;\;1\\
\mathbf{else}:\\
\;\;\;\;\sqrt{0.5}\\
\end{array}
\end{array}
if p < 7.6000000000000001e-200 or 4.3000000000000001e-185 < p < 2.5000000000000001e-154Initial program 74.6%
distribute-lft-in74.6%
metadata-eval74.6%
fma-def74.6%
Simplified74.6%
Taylor expanded in x around -inf 18.7%
associate-*r/18.7%
mul-1-neg18.7%
Simplified18.7%
if 7.6000000000000001e-200 < p < 4.3000000000000001e-185 or 2.5000000000000001e-154 < p < 1.0600000000000001e-54Initial program 79.3%
distribute-lft-in79.3%
metadata-eval79.3%
fma-def79.3%
add-log-exp79.2%
*-un-lft-identity79.2%
log-prod79.2%
metadata-eval79.2%
add-log-exp79.3%
+-commutative79.3%
Applied egg-rr79.3%
+-lft-identity79.3%
fma-udef79.3%
associate-*r/79.3%
*-commutative79.3%
associate-*r/79.3%
metadata-eval79.3%
associate-*l/79.3%
fma-udef78.5%
associate-*l/78.5%
metadata-eval78.5%
Simplified78.5%
Taylor expanded in x around inf 65.5%
if 1.0600000000000001e-54 < p Initial program 89.9%
Taylor expanded in x around 0 84.2%
Final simplification40.8%
NOTE: p should be positive before calling this function (FPCore (p x) :precision binary64 (if (<= x -3.2e-207) (/ (- p) x) 1.0))
p = abs(p);
double code(double p, double x) {
double tmp;
if (x <= -3.2e-207) {
tmp = -p / x;
} else {
tmp = 1.0;
}
return tmp;
}
NOTE: p should be positive before calling this function
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
real(8) :: tmp
if (x <= (-3.2d-207)) then
tmp = -p / x
else
tmp = 1.0d0
end if
code = tmp
end function
p = Math.abs(p);
public static double code(double p, double x) {
double tmp;
if (x <= -3.2e-207) {
tmp = -p / x;
} else {
tmp = 1.0;
}
return tmp;
}
p = abs(p) def code(p, x): tmp = 0 if x <= -3.2e-207: tmp = -p / x else: tmp = 1.0 return tmp
p = abs(p) function code(p, x) tmp = 0.0 if (x <= -3.2e-207) tmp = Float64(Float64(-p) / x); else tmp = 1.0; end return tmp end
p = abs(p) function tmp_2 = code(p, x) tmp = 0.0; if (x <= -3.2e-207) tmp = -p / x; else tmp = 1.0; end tmp_2 = tmp; end
NOTE: p should be positive before calling this function code[p_, x_] := If[LessEqual[x, -3.2e-207], N[((-p) / x), $MachinePrecision], 1.0]
\begin{array}{l}
p = |p|\\
\\
\begin{array}{l}
\mathbf{if}\;x \leq -3.2 \cdot 10^{-207}:\\
\;\;\;\;\frac{-p}{x}\\
\mathbf{else}:\\
\;\;\;\;1\\
\end{array}
\end{array}
if x < -3.2000000000000003e-207Initial program 55.4%
distribute-lft-in55.4%
metadata-eval55.4%
fma-def55.4%
Simplified55.4%
Taylor expanded in x around -inf 33.7%
associate-*r/33.7%
mul-1-neg33.7%
Simplified33.7%
if -3.2000000000000003e-207 < x Initial program 100.0%
distribute-lft-in100.0%
metadata-eval100.0%
fma-def100.0%
add-log-exp100.0%
*-un-lft-identity100.0%
log-prod100.0%
metadata-eval100.0%
add-log-exp100.0%
+-commutative100.0%
Applied egg-rr100.0%
+-lft-identity100.0%
fma-udef100.0%
associate-*r/100.0%
*-commutative100.0%
associate-*r/100.0%
metadata-eval100.0%
associate-*l/100.0%
fma-udef100.0%
associate-*l/100.0%
metadata-eval100.0%
Simplified100.0%
Taylor expanded in x around inf 57.4%
Final simplification46.4%
NOTE: p should be positive before calling this function (FPCore (p x) :precision binary64 (if (<= x -3.6e+105) (/ p x) 1.0))
p = abs(p);
double code(double p, double x) {
double tmp;
if (x <= -3.6e+105) {
tmp = p / x;
} else {
tmp = 1.0;
}
return tmp;
}
NOTE: p should be positive before calling this function
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
real(8) :: tmp
if (x <= (-3.6d+105)) then
tmp = p / x
else
tmp = 1.0d0
end if
code = tmp
end function
p = Math.abs(p);
public static double code(double p, double x) {
double tmp;
if (x <= -3.6e+105) {
tmp = p / x;
} else {
tmp = 1.0;
}
return tmp;
}
p = abs(p) def code(p, x): tmp = 0 if x <= -3.6e+105: tmp = p / x else: tmp = 1.0 return tmp
p = abs(p) function code(p, x) tmp = 0.0 if (x <= -3.6e+105) tmp = Float64(p / x); else tmp = 1.0; end return tmp end
p = abs(p) function tmp_2 = code(p, x) tmp = 0.0; if (x <= -3.6e+105) tmp = p / x; else tmp = 1.0; end tmp_2 = tmp; end
NOTE: p should be positive before calling this function code[p_, x_] := If[LessEqual[x, -3.6e+105], N[(p / x), $MachinePrecision], 1.0]
\begin{array}{l}
p = |p|\\
\\
\begin{array}{l}
\mathbf{if}\;x \leq -3.6 \cdot 10^{+105}:\\
\;\;\;\;\frac{p}{x}\\
\mathbf{else}:\\
\;\;\;\;1\\
\end{array}
\end{array}
if x < -3.5999999999999999e105Initial program 52.7%
distribute-lft-in52.7%
metadata-eval52.7%
fma-def52.7%
Simplified52.7%
*-commutative52.7%
fma-def52.7%
distribute-rgt1-in52.7%
+-commutative52.7%
Applied egg-rr52.7%
Taylor expanded in x around -inf 48.6%
unpow248.6%
unpow248.6%
associate-*r/48.9%
Simplified48.9%
Taylor expanded in p around 0 57.3%
if -3.5999999999999999e105 < x Initial program 82.2%
distribute-lft-in82.2%
metadata-eval82.2%
fma-def82.3%
add-log-exp82.2%
*-un-lft-identity82.2%
log-prod82.2%
metadata-eval82.2%
add-log-exp82.3%
+-commutative82.3%
Applied egg-rr82.3%
+-lft-identity82.3%
fma-udef82.3%
associate-*r/82.3%
*-commutative82.3%
associate-*r/82.3%
metadata-eval82.3%
associate-*l/82.3%
fma-udef81.3%
associate-*l/81.3%
metadata-eval81.3%
Simplified81.3%
Taylor expanded in x around inf 39.4%
Final simplification41.2%
NOTE: p should be positive before calling this function (FPCore (p x) :precision binary64 1.0)
p = abs(p);
double code(double p, double x) {
return 1.0;
}
NOTE: p should be positive before calling this function
real(8) function code(p, x)
real(8), intent (in) :: p
real(8), intent (in) :: x
code = 1.0d0
end function
p = Math.abs(p);
public static double code(double p, double x) {
return 1.0;
}
p = abs(p) def code(p, x): return 1.0
p = abs(p) function code(p, x) return 1.0 end
p = abs(p) function tmp = code(p, x) tmp = 1.0; end
NOTE: p should be positive before calling this function code[p_, x_] := 1.0
\begin{array}{l}
p = |p|\\
\\
1
\end{array}
Initial program 79.2%
distribute-lft-in79.2%
metadata-eval79.2%
fma-def79.2%
add-log-exp79.2%
*-un-lft-identity79.2%
log-prod79.2%
metadata-eval79.2%
add-log-exp79.2%
+-commutative79.2%
Applied egg-rr79.2%
+-lft-identity79.2%
fma-udef79.2%
associate-*r/79.2%
*-commutative79.2%
associate-*r/79.2%
metadata-eval79.2%
associate-*l/79.2%
fma-udef75.5%
associate-*l/75.5%
metadata-eval75.5%
Simplified75.5%
Taylor expanded in x around inf 36.2%
Final simplification36.2%
(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 2023192
(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))))))))