
(FPCore (a b) :precision binary64 (sqrt (- (* a a) (* b b))))
double code(double a, double b) {
return sqrt(((a * a) - (b * b)));
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
code = sqrt(((a * a) - (b * b)))
end function
public static double code(double a, double b) {
return Math.sqrt(((a * a) - (b * b)));
}
def code(a, b): return math.sqrt(((a * a) - (b * b)))
function code(a, b) return sqrt(Float64(Float64(a * a) - Float64(b * b))) end
function tmp = code(a, b) tmp = sqrt(((a * a) - (b * b))); end
code[a_, b_] := N[Sqrt[N[(N[(a * a), $MachinePrecision] - N[(b * b), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\sqrt{a \cdot a - b \cdot b}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 5 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (a b) :precision binary64 (sqrt (- (* a a) (* b b))))
double code(double a, double b) {
return sqrt(((a * a) - (b * b)));
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
code = sqrt(((a * a) - (b * b)))
end function
public static double code(double a, double b) {
return Math.sqrt(((a * a) - (b * b)));
}
def code(a, b): return math.sqrt(((a * a) - (b * b)))
function code(a, b) return sqrt(Float64(Float64(a * a) - Float64(b * b))) end
function tmp = code(a, b) tmp = sqrt(((a * a) - (b * b))); end
code[a_, b_] := N[Sqrt[N[(N[(a * a), $MachinePrecision] - N[(b * b), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\sqrt{a \cdot a - b \cdot b}
\end{array}
(FPCore (a b) :precision binary64 (if (<= a -2e-282) (- (/ 0.5 (/ (/ a b) b)) a) (* (sqrt (- a b)) (sqrt (+ a b)))))
double code(double a, double b) {
double tmp;
if (a <= -2e-282) {
tmp = (0.5 / ((a / b) / b)) - a;
} else {
tmp = sqrt((a - b)) * sqrt((a + b));
}
return tmp;
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8) :: tmp
if (a <= (-2d-282)) then
tmp = (0.5d0 / ((a / b) / b)) - a
else
tmp = sqrt((a - b)) * sqrt((a + b))
end if
code = tmp
end function
public static double code(double a, double b) {
double tmp;
if (a <= -2e-282) {
tmp = (0.5 / ((a / b) / b)) - a;
} else {
tmp = Math.sqrt((a - b)) * Math.sqrt((a + b));
}
return tmp;
}
def code(a, b): tmp = 0 if a <= -2e-282: tmp = (0.5 / ((a / b) / b)) - a else: tmp = math.sqrt((a - b)) * math.sqrt((a + b)) return tmp
function code(a, b) tmp = 0.0 if (a <= -2e-282) tmp = Float64(Float64(0.5 / Float64(Float64(a / b) / b)) - a); else tmp = Float64(sqrt(Float64(a - b)) * sqrt(Float64(a + b))); end return tmp end
function tmp_2 = code(a, b) tmp = 0.0; if (a <= -2e-282) tmp = (0.5 / ((a / b) / b)) - a; else tmp = sqrt((a - b)) * sqrt((a + b)); end tmp_2 = tmp; end
code[a_, b_] := If[LessEqual[a, -2e-282], N[(N[(0.5 / N[(N[(a / b), $MachinePrecision] / b), $MachinePrecision]), $MachinePrecision] - a), $MachinePrecision], N[(N[Sqrt[N[(a - b), $MachinePrecision]], $MachinePrecision] * N[Sqrt[N[(a + b), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;a \leq -2 \cdot 10^{-282}:\\
\;\;\;\;\frac{0.5}{\frac{\frac{a}{b}}{b}} - a\\
\mathbf{else}:\\
\;\;\;\;\sqrt{a - b} \cdot \sqrt{a + b}\\
\end{array}
\end{array}
if a < -2e-282Initial program 56.7%
Taylor expanded in a around -inf 96.3%
mul-1-neg96.3%
unsub-neg96.3%
unpow296.3%
associate-*r/96.3%
associate-/l*96.3%
associate-/r*99.7%
Simplified99.7%
if -2e-282 < a Initial program 49.9%
difference-of-squares50.2%
Simplified50.2%
*-commutative50.2%
sqrt-prod99.2%
Applied egg-rr99.2%
Final simplification99.4%
(FPCore (a b) :precision binary64 (if (<= a -2e-282) (- (/ 0.5 (/ (/ a b) b)) a) (fma -0.5 (/ b (/ a b)) a)))
double code(double a, double b) {
double tmp;
if (a <= -2e-282) {
tmp = (0.5 / ((a / b) / b)) - a;
} else {
tmp = fma(-0.5, (b / (a / b)), a);
}
return tmp;
}
function code(a, b) tmp = 0.0 if (a <= -2e-282) tmp = Float64(Float64(0.5 / Float64(Float64(a / b) / b)) - a); else tmp = fma(-0.5, Float64(b / Float64(a / b)), a); end return tmp end
code[a_, b_] := If[LessEqual[a, -2e-282], N[(N[(0.5 / N[(N[(a / b), $MachinePrecision] / b), $MachinePrecision]), $MachinePrecision] - a), $MachinePrecision], N[(-0.5 * N[(b / N[(a / b), $MachinePrecision]), $MachinePrecision] + a), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;a \leq -2 \cdot 10^{-282}:\\
\;\;\;\;\frac{0.5}{\frac{\frac{a}{b}}{b}} - a\\
\mathbf{else}:\\
\;\;\;\;\mathsf{fma}\left(-0.5, \frac{b}{\frac{a}{b}}, a\right)\\
\end{array}
\end{array}
if a < -2e-282Initial program 56.7%
Taylor expanded in a around -inf 96.3%
mul-1-neg96.3%
unsub-neg96.3%
unpow296.3%
associate-*r/96.3%
associate-/l*96.3%
associate-/r*99.7%
Simplified99.7%
if -2e-282 < a Initial program 49.9%
Taylor expanded in a around inf 95.3%
+-commutative95.3%
fma-def95.3%
unpow295.3%
associate-/l*99.0%
Simplified99.0%
Final simplification99.4%
(FPCore (a b) :precision binary64 (if (<= a -2e-282) (- (/ 0.5 (/ (/ a b) b)) a) a))
double code(double a, double b) {
double tmp;
if (a <= -2e-282) {
tmp = (0.5 / ((a / b) / b)) - a;
} else {
tmp = a;
}
return tmp;
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8) :: tmp
if (a <= (-2d-282)) then
tmp = (0.5d0 / ((a / b) / b)) - a
else
tmp = a
end if
code = tmp
end function
public static double code(double a, double b) {
double tmp;
if (a <= -2e-282) {
tmp = (0.5 / ((a / b) / b)) - a;
} else {
tmp = a;
}
return tmp;
}
def code(a, b): tmp = 0 if a <= -2e-282: tmp = (0.5 / ((a / b) / b)) - a else: tmp = a return tmp
function code(a, b) tmp = 0.0 if (a <= -2e-282) tmp = Float64(Float64(0.5 / Float64(Float64(a / b) / b)) - a); else tmp = a; end return tmp end
function tmp_2 = code(a, b) tmp = 0.0; if (a <= -2e-282) tmp = (0.5 / ((a / b) / b)) - a; else tmp = a; end tmp_2 = tmp; end
code[a_, b_] := If[LessEqual[a, -2e-282], N[(N[(0.5 / N[(N[(a / b), $MachinePrecision] / b), $MachinePrecision]), $MachinePrecision] - a), $MachinePrecision], a]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;a \leq -2 \cdot 10^{-282}:\\
\;\;\;\;\frac{0.5}{\frac{\frac{a}{b}}{b}} - a\\
\mathbf{else}:\\
\;\;\;\;a\\
\end{array}
\end{array}
if a < -2e-282Initial program 56.7%
Taylor expanded in a around -inf 96.3%
mul-1-neg96.3%
unsub-neg96.3%
unpow296.3%
associate-*r/96.3%
associate-/l*96.3%
associate-/r*99.7%
Simplified99.7%
if -2e-282 < a Initial program 49.9%
difference-of-squares50.2%
Simplified50.2%
Taylor expanded in a around inf 98.3%
Final simplification99.0%
(FPCore (a b) :precision binary64 (if (<= a -2e-282) (- a) a))
double code(double a, double b) {
double tmp;
if (a <= -2e-282) {
tmp = -a;
} else {
tmp = a;
}
return tmp;
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8) :: tmp
if (a <= (-2d-282)) then
tmp = -a
else
tmp = a
end if
code = tmp
end function
public static double code(double a, double b) {
double tmp;
if (a <= -2e-282) {
tmp = -a;
} else {
tmp = a;
}
return tmp;
}
def code(a, b): tmp = 0 if a <= -2e-282: tmp = -a else: tmp = a return tmp
function code(a, b) tmp = 0.0 if (a <= -2e-282) tmp = Float64(-a); else tmp = a; end return tmp end
function tmp_2 = code(a, b) tmp = 0.0; if (a <= -2e-282) tmp = -a; else tmp = a; end tmp_2 = tmp; end
code[a_, b_] := If[LessEqual[a, -2e-282], (-a), a]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;a \leq -2 \cdot 10^{-282}:\\
\;\;\;\;-a\\
\mathbf{else}:\\
\;\;\;\;a\\
\end{array}
\end{array}
if a < -2e-282Initial program 56.7%
Taylor expanded in a around -inf 99.4%
mul-1-neg99.4%
Simplified99.4%
if -2e-282 < a Initial program 49.9%
difference-of-squares50.2%
Simplified50.2%
Taylor expanded in a around inf 98.3%
Final simplification98.8%
(FPCore (a b) :precision binary64 a)
double code(double a, double b) {
return a;
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
code = a
end function
public static double code(double a, double b) {
return a;
}
def code(a, b): return a
function code(a, b) return a end
function tmp = code(a, b) tmp = a; end
code[a_, b_] := a
\begin{array}{l}
\\
a
\end{array}
Initial program 53.1%
difference-of-squares53.4%
Simplified53.4%
Taylor expanded in a around inf 52.9%
Final simplification52.9%
(FPCore (a b) :precision binary64 (* (sqrt (+ (fabs a) (fabs b))) (sqrt (- (fabs a) (fabs b)))))
double code(double a, double b) {
return sqrt((fabs(a) + fabs(b))) * sqrt((fabs(a) - fabs(b)));
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
code = sqrt((abs(a) + abs(b))) * sqrt((abs(a) - abs(b)))
end function
public static double code(double a, double b) {
return Math.sqrt((Math.abs(a) + Math.abs(b))) * Math.sqrt((Math.abs(a) - Math.abs(b)));
}
def code(a, b): return math.sqrt((math.fabs(a) + math.fabs(b))) * math.sqrt((math.fabs(a) - math.fabs(b)))
function code(a, b) return Float64(sqrt(Float64(abs(a) + abs(b))) * sqrt(Float64(abs(a) - abs(b)))) end
function tmp = code(a, b) tmp = sqrt((abs(a) + abs(b))) * sqrt((abs(a) - abs(b))); end
code[a_, b_] := N[(N[Sqrt[N[(N[Abs[a], $MachinePrecision] + N[Abs[b], $MachinePrecision]), $MachinePrecision]], $MachinePrecision] * N[Sqrt[N[(N[Abs[a], $MachinePrecision] - N[Abs[b], $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\sqrt{\left|a\right| + \left|b\right|} \cdot \sqrt{\left|a\right| - \left|b\right|}
\end{array}
herbie shell --seed 2023182
(FPCore (a b)
:name "bug366, discussion (missed optimization)"
:precision binary64
:herbie-target
(* (sqrt (+ (fabs a) (fabs b))) (sqrt (- (fabs a) (fabs b))))
(sqrt (- (* a a) (* b b))))