
(FPCore (re im) :precision binary64 (* 0.5 (sqrt (* 2.0 (+ (sqrt (+ (* re re) (* im im))) re)))))
double code(double re, double im) {
return 0.5 * sqrt((2.0 * (sqrt(((re * re) + (im * im))) + re)));
}
real(8) function code(re, im)
real(8), intent (in) :: re
real(8), intent (in) :: im
code = 0.5d0 * sqrt((2.0d0 * (sqrt(((re * re) + (im * im))) + re)))
end function
public static double code(double re, double im) {
return 0.5 * Math.sqrt((2.0 * (Math.sqrt(((re * re) + (im * im))) + re)));
}
def code(re, im): return 0.5 * math.sqrt((2.0 * (math.sqrt(((re * re) + (im * im))) + re)))
function code(re, im) return Float64(0.5 * sqrt(Float64(2.0 * Float64(sqrt(Float64(Float64(re * re) + Float64(im * im))) + re)))) end
function tmp = code(re, im) tmp = 0.5 * sqrt((2.0 * (sqrt(((re * re) + (im * im))) + re))); end
code[re_, im_] := N[(0.5 * N[Sqrt[N[(2.0 * N[(N[Sqrt[N[(N[(re * re), $MachinePrecision] + N[(im * im), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] + re), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
0.5 \cdot \sqrt{2 \cdot \left(\sqrt{re \cdot re + im \cdot im} + re\right)}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 5 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (re im) :precision binary64 (* 0.5 (sqrt (* 2.0 (+ (sqrt (+ (* re re) (* im im))) re)))))
double code(double re, double im) {
return 0.5 * sqrt((2.0 * (sqrt(((re * re) + (im * im))) + re)));
}
real(8) function code(re, im)
real(8), intent (in) :: re
real(8), intent (in) :: im
code = 0.5d0 * sqrt((2.0d0 * (sqrt(((re * re) + (im * im))) + re)))
end function
public static double code(double re, double im) {
return 0.5 * Math.sqrt((2.0 * (Math.sqrt(((re * re) + (im * im))) + re)));
}
def code(re, im): return 0.5 * math.sqrt((2.0 * (math.sqrt(((re * re) + (im * im))) + re)))
function code(re, im) return Float64(0.5 * sqrt(Float64(2.0 * Float64(sqrt(Float64(Float64(re * re) + Float64(im * im))) + re)))) end
function tmp = code(re, im) tmp = 0.5 * sqrt((2.0 * (sqrt(((re * re) + (im * im))) + re))); end
code[re_, im_] := N[(0.5 * N[Sqrt[N[(2.0 * N[(N[Sqrt[N[(N[(re * re), $MachinePrecision] + N[(im * im), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] + re), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
0.5 \cdot \sqrt{2 \cdot \left(\sqrt{re \cdot re + im \cdot im} + re\right)}
\end{array}
im_m = (fabs.f64 im) (FPCore (re im_m) :precision binary64 (if (<= re -1.9e-7) (* 0.5 (* im_m (* (sqrt 2.0) (sqrt (/ -0.5 re))))) (sqrt (* 0.5 (+ re (hypot re im_m))))))
im_m = fabs(im);
double code(double re, double im_m) {
double tmp;
if (re <= -1.9e-7) {
tmp = 0.5 * (im_m * (sqrt(2.0) * sqrt((-0.5 / re))));
} else {
tmp = sqrt((0.5 * (re + hypot(re, im_m))));
}
return tmp;
}
im_m = Math.abs(im);
public static double code(double re, double im_m) {
double tmp;
if (re <= -1.9e-7) {
tmp = 0.5 * (im_m * (Math.sqrt(2.0) * Math.sqrt((-0.5 / re))));
} else {
tmp = Math.sqrt((0.5 * (re + Math.hypot(re, im_m))));
}
return tmp;
}
im_m = math.fabs(im) def code(re, im_m): tmp = 0 if re <= -1.9e-7: tmp = 0.5 * (im_m * (math.sqrt(2.0) * math.sqrt((-0.5 / re)))) else: tmp = math.sqrt((0.5 * (re + math.hypot(re, im_m)))) return tmp
im_m = abs(im) function code(re, im_m) tmp = 0.0 if (re <= -1.9e-7) tmp = Float64(0.5 * Float64(im_m * Float64(sqrt(2.0) * sqrt(Float64(-0.5 / re))))); else tmp = sqrt(Float64(0.5 * Float64(re + hypot(re, im_m)))); end return tmp end
im_m = abs(im); function tmp_2 = code(re, im_m) tmp = 0.0; if (re <= -1.9e-7) tmp = 0.5 * (im_m * (sqrt(2.0) * sqrt((-0.5 / re)))); else tmp = sqrt((0.5 * (re + hypot(re, im_m)))); end tmp_2 = tmp; end
im_m = N[Abs[im], $MachinePrecision] code[re_, im$95$m_] := If[LessEqual[re, -1.9e-7], N[(0.5 * N[(im$95$m * N[(N[Sqrt[2.0], $MachinePrecision] * N[Sqrt[N[(-0.5 / re), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[Sqrt[N[(0.5 * N[(re + N[Sqrt[re ^ 2 + im$95$m ^ 2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]]
\begin{array}{l}
im_m = \left|im\right|
\\
\begin{array}{l}
\mathbf{if}\;re \leq -1.9 \cdot 10^{-7}:\\
\;\;\;\;0.5 \cdot \left(im\_m \cdot \left(\sqrt{2} \cdot \sqrt{\frac{-0.5}{re}}\right)\right)\\
\mathbf{else}:\\
\;\;\;\;\sqrt{0.5 \cdot \left(re + \mathsf{hypot}\left(re, im\_m\right)\right)}\\
\end{array}
\end{array}
if re < -1.90000000000000007e-7Initial program 13.3%
sqr-neg13.3%
+-commutative13.3%
sqr-neg13.3%
+-commutative13.3%
distribute-rgt-in13.3%
cancel-sign-sub13.3%
distribute-rgt-out--13.3%
sub-neg13.3%
remove-double-neg13.3%
+-commutative13.3%
hypot-define33.0%
Simplified33.0%
*-commutative33.0%
hypot-define13.3%
+-commutative13.3%
*-commutative13.3%
add-cbrt-cube13.2%
pow1/312.6%
Applied egg-rr23.3%
unpow1/324.3%
Simplified24.3%
Taylor expanded in re around -inf 38.2%
*-commutative38.2%
associate-*l/38.2%
Simplified38.2%
pow1/336.9%
pow-pow53.6%
metadata-eval53.6%
pow1/253.6%
*-commutative53.6%
sqrt-prod53.5%
associate-/l*53.5%
sqrt-prod62.8%
sqrt-pow139.9%
metadata-eval39.9%
pow139.9%
Applied egg-rr39.9%
associate-*l*39.9%
*-commutative39.9%
Simplified39.9%
if -1.90000000000000007e-7 < re Initial program 45.4%
sqr-neg45.4%
+-commutative45.4%
sqr-neg45.4%
+-commutative45.4%
distribute-rgt-in45.4%
cancel-sign-sub45.4%
distribute-rgt-out--45.4%
sub-neg45.4%
remove-double-neg45.4%
+-commutative45.4%
hypot-define94.5%
Simplified94.5%
*-commutative94.5%
hypot-define45.4%
+-commutative45.4%
*-commutative45.4%
add-sqr-sqrt45.1%
sqrt-unprod45.4%
*-commutative45.4%
*-commutative45.4%
swap-sqr45.4%
Applied egg-rr94.5%
*-commutative94.5%
associate-*r*94.5%
metadata-eval94.5%
Simplified94.5%
im_m = (fabs.f64 im) (FPCore (re im_m) :precision binary64 (if (<= re -1.35e+53) (* 0.5 (sqrt (/ (pow im_m 2.0) (- re)))) (sqrt (* 0.5 (+ re (hypot re im_m))))))
im_m = fabs(im);
double code(double re, double im_m) {
double tmp;
if (re <= -1.35e+53) {
tmp = 0.5 * sqrt((pow(im_m, 2.0) / -re));
} else {
tmp = sqrt((0.5 * (re + hypot(re, im_m))));
}
return tmp;
}
im_m = Math.abs(im);
public static double code(double re, double im_m) {
double tmp;
if (re <= -1.35e+53) {
tmp = 0.5 * Math.sqrt((Math.pow(im_m, 2.0) / -re));
} else {
tmp = Math.sqrt((0.5 * (re + Math.hypot(re, im_m))));
}
return tmp;
}
im_m = math.fabs(im) def code(re, im_m): tmp = 0 if re <= -1.35e+53: tmp = 0.5 * math.sqrt((math.pow(im_m, 2.0) / -re)) else: tmp = math.sqrt((0.5 * (re + math.hypot(re, im_m)))) return tmp
im_m = abs(im) function code(re, im_m) tmp = 0.0 if (re <= -1.35e+53) tmp = Float64(0.5 * sqrt(Float64((im_m ^ 2.0) / Float64(-re)))); else tmp = sqrt(Float64(0.5 * Float64(re + hypot(re, im_m)))); end return tmp end
im_m = abs(im); function tmp_2 = code(re, im_m) tmp = 0.0; if (re <= -1.35e+53) tmp = 0.5 * sqrt(((im_m ^ 2.0) / -re)); else tmp = sqrt((0.5 * (re + hypot(re, im_m)))); end tmp_2 = tmp; end
im_m = N[Abs[im], $MachinePrecision] code[re_, im$95$m_] := If[LessEqual[re, -1.35e+53], N[(0.5 * N[Sqrt[N[(N[Power[im$95$m, 2.0], $MachinePrecision] / (-re)), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[Sqrt[N[(0.5 * N[(re + N[Sqrt[re ^ 2 + im$95$m ^ 2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]]
\begin{array}{l}
im_m = \left|im\right|
\\
\begin{array}{l}
\mathbf{if}\;re \leq -1.35 \cdot 10^{+53}:\\
\;\;\;\;0.5 \cdot \sqrt{\frac{{im\_m}^{2}}{-re}}\\
\mathbf{else}:\\
\;\;\;\;\sqrt{0.5 \cdot \left(re + \mathsf{hypot}\left(re, im\_m\right)\right)}\\
\end{array}
\end{array}
if re < -1.3500000000000001e53Initial program 8.8%
sqr-neg8.8%
+-commutative8.8%
sqr-neg8.8%
+-commutative8.8%
distribute-rgt-in8.8%
cancel-sign-sub8.8%
distribute-rgt-out--8.8%
sub-neg8.8%
remove-double-neg8.8%
+-commutative8.8%
hypot-define30.1%
Simplified30.1%
Taylor expanded in re around -inf 60.7%
mul-1-neg60.7%
distribute-neg-frac260.7%
Simplified60.7%
if -1.3500000000000001e53 < re Initial program 44.3%
sqr-neg44.3%
+-commutative44.3%
sqr-neg44.3%
+-commutative44.3%
distribute-rgt-in44.3%
cancel-sign-sub44.3%
distribute-rgt-out--44.3%
sub-neg44.3%
remove-double-neg44.3%
+-commutative44.3%
hypot-define91.0%
Simplified91.0%
*-commutative91.0%
hypot-define44.3%
+-commutative44.3%
*-commutative44.3%
add-sqr-sqrt44.1%
sqrt-unprod44.3%
*-commutative44.3%
*-commutative44.3%
swap-sqr44.3%
Applied egg-rr91.0%
*-commutative91.0%
associate-*r*91.0%
metadata-eval91.0%
Simplified91.0%
im_m = (fabs.f64 im) (FPCore (re im_m) :precision binary64 (sqrt (* 0.5 (+ re (hypot re im_m)))))
im_m = fabs(im);
double code(double re, double im_m) {
return sqrt((0.5 * (re + hypot(re, im_m))));
}
im_m = Math.abs(im);
public static double code(double re, double im_m) {
return Math.sqrt((0.5 * (re + Math.hypot(re, im_m))));
}
im_m = math.fabs(im) def code(re, im_m): return math.sqrt((0.5 * (re + math.hypot(re, im_m))))
im_m = abs(im) function code(re, im_m) return sqrt(Float64(0.5 * Float64(re + hypot(re, im_m)))) end
im_m = abs(im); function tmp = code(re, im_m) tmp = sqrt((0.5 * (re + hypot(re, im_m)))); end
im_m = N[Abs[im], $MachinePrecision] code[re_, im$95$m_] := N[Sqrt[N[(0.5 * N[(re + N[Sqrt[re ^ 2 + im$95$m ^ 2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
im_m = \left|im\right|
\\
\sqrt{0.5 \cdot \left(re + \mathsf{hypot}\left(re, im\_m\right)\right)}
\end{array}
Initial program 37.1%
sqr-neg37.1%
+-commutative37.1%
sqr-neg37.1%
+-commutative37.1%
distribute-rgt-in37.1%
cancel-sign-sub37.1%
distribute-rgt-out--37.1%
sub-neg37.1%
remove-double-neg37.1%
+-commutative37.1%
hypot-define78.7%
Simplified78.7%
*-commutative78.7%
hypot-define37.1%
+-commutative37.1%
*-commutative37.1%
add-sqr-sqrt36.9%
sqrt-unprod37.1%
*-commutative37.1%
*-commutative37.1%
swap-sqr37.1%
Applied egg-rr78.7%
*-commutative78.7%
associate-*r*78.7%
metadata-eval78.7%
Simplified78.7%
im_m = (fabs.f64 im) (FPCore (re im_m) :precision binary64 (if (<= re 1.55e-15) (* 0.5 (sqrt (* im_m 2.0))) (sqrt re)))
im_m = fabs(im);
double code(double re, double im_m) {
double tmp;
if (re <= 1.55e-15) {
tmp = 0.5 * sqrt((im_m * 2.0));
} else {
tmp = sqrt(re);
}
return tmp;
}
im_m = abs(im)
real(8) function code(re, im_m)
real(8), intent (in) :: re
real(8), intent (in) :: im_m
real(8) :: tmp
if (re <= 1.55d-15) then
tmp = 0.5d0 * sqrt((im_m * 2.0d0))
else
tmp = sqrt(re)
end if
code = tmp
end function
im_m = Math.abs(im);
public static double code(double re, double im_m) {
double tmp;
if (re <= 1.55e-15) {
tmp = 0.5 * Math.sqrt((im_m * 2.0));
} else {
tmp = Math.sqrt(re);
}
return tmp;
}
im_m = math.fabs(im) def code(re, im_m): tmp = 0 if re <= 1.55e-15: tmp = 0.5 * math.sqrt((im_m * 2.0)) else: tmp = math.sqrt(re) return tmp
im_m = abs(im) function code(re, im_m) tmp = 0.0 if (re <= 1.55e-15) tmp = Float64(0.5 * sqrt(Float64(im_m * 2.0))); else tmp = sqrt(re); end return tmp end
im_m = abs(im); function tmp_2 = code(re, im_m) tmp = 0.0; if (re <= 1.55e-15) tmp = 0.5 * sqrt((im_m * 2.0)); else tmp = sqrt(re); end tmp_2 = tmp; end
im_m = N[Abs[im], $MachinePrecision] code[re_, im$95$m_] := If[LessEqual[re, 1.55e-15], N[(0.5 * N[Sqrt[N[(im$95$m * 2.0), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[Sqrt[re], $MachinePrecision]]
\begin{array}{l}
im_m = \left|im\right|
\\
\begin{array}{l}
\mathbf{if}\;re \leq 1.55 \cdot 10^{-15}:\\
\;\;\;\;0.5 \cdot \sqrt{im\_m \cdot 2}\\
\mathbf{else}:\\
\;\;\;\;\sqrt{re}\\
\end{array}
\end{array}
if re < 1.5499999999999999e-15Initial program 38.3%
sqr-neg38.3%
+-commutative38.3%
sqr-neg38.3%
+-commutative38.3%
distribute-rgt-in38.3%
cancel-sign-sub38.3%
distribute-rgt-out--38.3%
sub-neg38.3%
remove-double-neg38.3%
+-commutative38.3%
hypot-define70.9%
Simplified70.9%
Taylor expanded in re around 0 27.3%
*-commutative27.3%
Simplified27.3%
if 1.5499999999999999e-15 < re Initial program 33.8%
sqr-neg33.8%
+-commutative33.8%
sqr-neg33.8%
+-commutative33.8%
distribute-rgt-in33.8%
cancel-sign-sub33.8%
distribute-rgt-out--33.8%
sub-neg33.8%
remove-double-neg33.8%
+-commutative33.8%
hypot-define100.0%
Simplified100.0%
Taylor expanded in re around inf 75.1%
*-commutative75.1%
unpow275.1%
rem-square-sqrt76.6%
associate-*r*76.6%
metadata-eval76.6%
Simplified76.6%
Final simplification40.4%
im_m = (fabs.f64 im) (FPCore (re im_m) :precision binary64 (sqrt re))
im_m = fabs(im);
double code(double re, double im_m) {
return sqrt(re);
}
im_m = abs(im)
real(8) function code(re, im_m)
real(8), intent (in) :: re
real(8), intent (in) :: im_m
code = sqrt(re)
end function
im_m = Math.abs(im);
public static double code(double re, double im_m) {
return Math.sqrt(re);
}
im_m = math.fabs(im) def code(re, im_m): return math.sqrt(re)
im_m = abs(im) function code(re, im_m) return sqrt(re) end
im_m = abs(im); function tmp = code(re, im_m) tmp = sqrt(re); end
im_m = N[Abs[im], $MachinePrecision] code[re_, im$95$m_] := N[Sqrt[re], $MachinePrecision]
\begin{array}{l}
im_m = \left|im\right|
\\
\sqrt{re}
\end{array}
Initial program 37.1%
sqr-neg37.1%
+-commutative37.1%
sqr-neg37.1%
+-commutative37.1%
distribute-rgt-in37.1%
cancel-sign-sub37.1%
distribute-rgt-out--37.1%
sub-neg37.1%
remove-double-neg37.1%
+-commutative37.1%
hypot-define78.7%
Simplified78.7%
Taylor expanded in re around inf 27.6%
*-commutative27.6%
unpow227.6%
rem-square-sqrt28.1%
associate-*r*28.1%
metadata-eval28.1%
Simplified28.1%
Final simplification28.1%
(FPCore (re im)
:precision binary64
(let* ((t_0 (sqrt (+ (* re re) (* im im)))))
(if (< re 0.0)
(* 0.5 (* (sqrt 2.0) (sqrt (/ (* im im) (- t_0 re)))))
(* 0.5 (sqrt (* 2.0 (+ t_0 re)))))))
double code(double re, double im) {
double t_0 = sqrt(((re * re) + (im * im)));
double tmp;
if (re < 0.0) {
tmp = 0.5 * (sqrt(2.0) * sqrt(((im * im) / (t_0 - re))));
} else {
tmp = 0.5 * sqrt((2.0 * (t_0 + re)));
}
return tmp;
}
real(8) function code(re, im)
real(8), intent (in) :: re
real(8), intent (in) :: im
real(8) :: t_0
real(8) :: tmp
t_0 = sqrt(((re * re) + (im * im)))
if (re < 0.0d0) then
tmp = 0.5d0 * (sqrt(2.0d0) * sqrt(((im * im) / (t_0 - re))))
else
tmp = 0.5d0 * sqrt((2.0d0 * (t_0 + re)))
end if
code = tmp
end function
public static double code(double re, double im) {
double t_0 = Math.sqrt(((re * re) + (im * im)));
double tmp;
if (re < 0.0) {
tmp = 0.5 * (Math.sqrt(2.0) * Math.sqrt(((im * im) / (t_0 - re))));
} else {
tmp = 0.5 * Math.sqrt((2.0 * (t_0 + re)));
}
return tmp;
}
def code(re, im): t_0 = math.sqrt(((re * re) + (im * im))) tmp = 0 if re < 0.0: tmp = 0.5 * (math.sqrt(2.0) * math.sqrt(((im * im) / (t_0 - re)))) else: tmp = 0.5 * math.sqrt((2.0 * (t_0 + re))) return tmp
function code(re, im) t_0 = sqrt(Float64(Float64(re * re) + Float64(im * im))) tmp = 0.0 if (re < 0.0) tmp = Float64(0.5 * Float64(sqrt(2.0) * sqrt(Float64(Float64(im * im) / Float64(t_0 - re))))); else tmp = Float64(0.5 * sqrt(Float64(2.0 * Float64(t_0 + re)))); end return tmp end
function tmp_2 = code(re, im) t_0 = sqrt(((re * re) + (im * im))); tmp = 0.0; if (re < 0.0) tmp = 0.5 * (sqrt(2.0) * sqrt(((im * im) / (t_0 - re)))); else tmp = 0.5 * sqrt((2.0 * (t_0 + re))); end tmp_2 = tmp; end
code[re_, im_] := Block[{t$95$0 = N[Sqrt[N[(N[(re * re), $MachinePrecision] + N[(im * im), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]}, If[Less[re, 0.0], N[(0.5 * N[(N[Sqrt[2.0], $MachinePrecision] * N[Sqrt[N[(N[(im * im), $MachinePrecision] / N[(t$95$0 - re), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(0.5 * N[Sqrt[N[(2.0 * N[(t$95$0 + re), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := \sqrt{re \cdot re + im \cdot im}\\
\mathbf{if}\;re < 0:\\
\;\;\;\;0.5 \cdot \left(\sqrt{2} \cdot \sqrt{\frac{im \cdot im}{t\_0 - re}}\right)\\
\mathbf{else}:\\
\;\;\;\;0.5 \cdot \sqrt{2 \cdot \left(t\_0 + re\right)}\\
\end{array}
\end{array}
herbie shell --seed 2024107
(FPCore (re im)
:name "math.sqrt on complex, real part"
:precision binary64
:alt
(if (< re 0.0) (* 0.5 (* (sqrt 2.0) (sqrt (/ (* im im) (- (sqrt (+ (* re re) (* im im))) re))))) (* 0.5 (sqrt (* 2.0 (+ (sqrt (+ (* re re) (* im im))) re)))))
(* 0.5 (sqrt (* 2.0 (+ (sqrt (+ (* re re) (* im im))) re)))))