
(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 4 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}
(FPCore (re im) :precision binary64 (* 0.5 (sqrt (* 2.0 (+ re (hypot re im))))))
double code(double re, double im) {
return 0.5 * sqrt((2.0 * (re + hypot(re, im))));
}
public static double code(double re, double im) {
return 0.5 * Math.sqrt((2.0 * (re + Math.hypot(re, im))));
}
def code(re, im): return 0.5 * math.sqrt((2.0 * (re + math.hypot(re, im))))
function code(re, im) return Float64(0.5 * sqrt(Float64(2.0 * Float64(re + hypot(re, im))))) end
function tmp = code(re, im) tmp = 0.5 * sqrt((2.0 * (re + hypot(re, im)))); end
code[re_, im_] := N[(0.5 * N[Sqrt[N[(2.0 * N[(re + N[Sqrt[re ^ 2 + im ^ 2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
0.5 \cdot \sqrt{2 \cdot \left(re + \mathsf{hypot}\left(re, im\right)\right)}
\end{array}
Initial program 52.0%
sqr-neg52.0%
+-commutative52.0%
sqr-neg52.0%
+-commutative52.0%
distribute-rgt-in52.0%
cancel-sign-sub52.0%
distribute-rgt-out--52.0%
sub-neg52.0%
remove-double-neg52.0%
+-commutative52.0%
hypot-define98.6%
Simplified98.6%
Final simplification98.6%
(FPCore (re im) :precision binary64 (if (or (<= re 1.5e-25) (and (not (<= re 860000.0)) (<= re 1.75e+58))) (* 0.5 (sqrt (* 2.0 im))) (* 0.5 (sqrt (* 2.0 (* 2.0 re))))))
double code(double re, double im) {
double tmp;
if ((re <= 1.5e-25) || (!(re <= 860000.0) && (re <= 1.75e+58))) {
tmp = 0.5 * sqrt((2.0 * im));
} else {
tmp = 0.5 * sqrt((2.0 * (2.0 * re)));
}
return tmp;
}
real(8) function code(re, im)
real(8), intent (in) :: re
real(8), intent (in) :: im
real(8) :: tmp
if ((re <= 1.5d-25) .or. (.not. (re <= 860000.0d0)) .and. (re <= 1.75d+58)) then
tmp = 0.5d0 * sqrt((2.0d0 * im))
else
tmp = 0.5d0 * sqrt((2.0d0 * (2.0d0 * re)))
end if
code = tmp
end function
public static double code(double re, double im) {
double tmp;
if ((re <= 1.5e-25) || (!(re <= 860000.0) && (re <= 1.75e+58))) {
tmp = 0.5 * Math.sqrt((2.0 * im));
} else {
tmp = 0.5 * Math.sqrt((2.0 * (2.0 * re)));
}
return tmp;
}
def code(re, im): tmp = 0 if (re <= 1.5e-25) or (not (re <= 860000.0) and (re <= 1.75e+58)): tmp = 0.5 * math.sqrt((2.0 * im)) else: tmp = 0.5 * math.sqrt((2.0 * (2.0 * re))) return tmp
function code(re, im) tmp = 0.0 if ((re <= 1.5e-25) || (!(re <= 860000.0) && (re <= 1.75e+58))) tmp = Float64(0.5 * sqrt(Float64(2.0 * im))); else tmp = Float64(0.5 * sqrt(Float64(2.0 * Float64(2.0 * re)))); end return tmp end
function tmp_2 = code(re, im) tmp = 0.0; if ((re <= 1.5e-25) || (~((re <= 860000.0)) && (re <= 1.75e+58))) tmp = 0.5 * sqrt((2.0 * im)); else tmp = 0.5 * sqrt((2.0 * (2.0 * re))); end tmp_2 = tmp; end
code[re_, im_] := If[Or[LessEqual[re, 1.5e-25], And[N[Not[LessEqual[re, 860000.0]], $MachinePrecision], LessEqual[re, 1.75e+58]]], N[(0.5 * N[Sqrt[N[(2.0 * im), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(0.5 * N[Sqrt[N[(2.0 * N[(2.0 * re), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;re \leq 1.5 \cdot 10^{-25} \lor \neg \left(re \leq 860000\right) \land re \leq 1.75 \cdot 10^{+58}:\\
\;\;\;\;0.5 \cdot \sqrt{2 \cdot im}\\
\mathbf{else}:\\
\;\;\;\;0.5 \cdot \sqrt{2 \cdot \left(2 \cdot re\right)}\\
\end{array}
\end{array}
if re < 1.4999999999999999e-25 or 8.6e5 < re < 1.7499999999999999e58Initial program 55.7%
sqr-neg55.7%
+-commutative55.7%
sqr-neg55.7%
+-commutative55.7%
distribute-rgt-in55.7%
cancel-sign-sub55.7%
distribute-rgt-out--55.7%
sub-neg55.7%
remove-double-neg55.7%
+-commutative55.7%
hypot-define98.1%
Simplified98.1%
Taylor expanded in re around 0 29.1%
if 1.4999999999999999e-25 < re < 8.6e5 or 1.7499999999999999e58 < re Initial program 40.7%
sqr-neg40.7%
+-commutative40.7%
sqr-neg40.7%
+-commutative40.7%
distribute-rgt-in40.7%
cancel-sign-sub40.7%
distribute-rgt-out--40.7%
sub-neg40.7%
remove-double-neg40.7%
+-commutative40.7%
hypot-define100.0%
Simplified100.0%
Taylor expanded in re around inf 82.8%
*-commutative82.8%
Simplified82.8%
Final simplification42.3%
(FPCore (re im) :precision binary64 (if (<= re 6.5e+113) (* 0.5 (sqrt (* 2.0 (+ re im)))) (* 0.5 (sqrt (* 2.0 (* 2.0 re))))))
double code(double re, double im) {
double tmp;
if (re <= 6.5e+113) {
tmp = 0.5 * sqrt((2.0 * (re + im)));
} else {
tmp = 0.5 * sqrt((2.0 * (2.0 * re)));
}
return tmp;
}
real(8) function code(re, im)
real(8), intent (in) :: re
real(8), intent (in) :: im
real(8) :: tmp
if (re <= 6.5d+113) then
tmp = 0.5d0 * sqrt((2.0d0 * (re + im)))
else
tmp = 0.5d0 * sqrt((2.0d0 * (2.0d0 * re)))
end if
code = tmp
end function
public static double code(double re, double im) {
double tmp;
if (re <= 6.5e+113) {
tmp = 0.5 * Math.sqrt((2.0 * (re + im)));
} else {
tmp = 0.5 * Math.sqrt((2.0 * (2.0 * re)));
}
return tmp;
}
def code(re, im): tmp = 0 if re <= 6.5e+113: tmp = 0.5 * math.sqrt((2.0 * (re + im))) else: tmp = 0.5 * math.sqrt((2.0 * (2.0 * re))) return tmp
function code(re, im) tmp = 0.0 if (re <= 6.5e+113) tmp = Float64(0.5 * sqrt(Float64(2.0 * Float64(re + im)))); else tmp = Float64(0.5 * sqrt(Float64(2.0 * Float64(2.0 * re)))); end return tmp end
function tmp_2 = code(re, im) tmp = 0.0; if (re <= 6.5e+113) tmp = 0.5 * sqrt((2.0 * (re + im))); else tmp = 0.5 * sqrt((2.0 * (2.0 * re))); end tmp_2 = tmp; end
code[re_, im_] := If[LessEqual[re, 6.5e+113], N[(0.5 * N[Sqrt[N[(2.0 * N[(re + im), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(0.5 * N[Sqrt[N[(2.0 * N[(2.0 * re), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;re \leq 6.5 \cdot 10^{+113}:\\
\;\;\;\;0.5 \cdot \sqrt{2 \cdot \left(re + im\right)}\\
\mathbf{else}:\\
\;\;\;\;0.5 \cdot \sqrt{2 \cdot \left(2 \cdot re\right)}\\
\end{array}
\end{array}
if re < 6.5000000000000001e113Initial program 58.1%
sqr-neg58.1%
+-commutative58.1%
sqr-neg58.1%
+-commutative58.1%
distribute-rgt-in58.1%
cancel-sign-sub58.1%
distribute-rgt-out--58.1%
sub-neg58.1%
remove-double-neg58.1%
+-commutative58.1%
hypot-define98.3%
Simplified98.3%
Taylor expanded in re around 0 30.7%
if 6.5000000000000001e113 < re Initial program 22.1%
sqr-neg22.1%
+-commutative22.1%
sqr-neg22.1%
+-commutative22.1%
distribute-rgt-in22.1%
cancel-sign-sub22.1%
distribute-rgt-out--22.1%
sub-neg22.1%
remove-double-neg22.1%
+-commutative22.1%
hypot-define100.0%
Simplified100.0%
Taylor expanded in re around inf 91.7%
*-commutative91.7%
Simplified91.7%
Final simplification41.0%
(FPCore (re im) :precision binary64 (* 0.5 (sqrt (* 2.0 im))))
double code(double re, double im) {
return 0.5 * sqrt((2.0 * im));
}
real(8) function code(re, im)
real(8), intent (in) :: re
real(8), intent (in) :: im
code = 0.5d0 * sqrt((2.0d0 * im))
end function
public static double code(double re, double im) {
return 0.5 * Math.sqrt((2.0 * im));
}
def code(re, im): return 0.5 * math.sqrt((2.0 * im))
function code(re, im) return Float64(0.5 * sqrt(Float64(2.0 * im))) end
function tmp = code(re, im) tmp = 0.5 * sqrt((2.0 * im)); end
code[re_, im_] := N[(0.5 * N[Sqrt[N[(2.0 * im), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
0.5 \cdot \sqrt{2 \cdot im}
\end{array}
Initial program 52.0%
sqr-neg52.0%
+-commutative52.0%
sqr-neg52.0%
+-commutative52.0%
distribute-rgt-in52.0%
cancel-sign-sub52.0%
distribute-rgt-out--52.0%
sub-neg52.0%
remove-double-neg52.0%
+-commutative52.0%
hypot-define98.6%
Simplified98.6%
Taylor expanded in re around 0 25.2%
Final simplification25.2%
(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 2024052
(FPCore (re im)
:name "math.sqrt on complex, real part"
:precision binary64
:herbie-target
(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)))))