
(FPCore (a b) :precision binary64 (- (+ (pow (+ (* a a) (* b b)) 2.0) (* 4.0 (+ (* (* a a) (+ 1.0 a)) (* (* b b) (- 1.0 (* 3.0 a)))))) 1.0))
double code(double a, double b) {
return (pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0;
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
code = ((((a * a) + (b * b)) ** 2.0d0) + (4.0d0 * (((a * a) * (1.0d0 + a)) + ((b * b) * (1.0d0 - (3.0d0 * a)))))) - 1.0d0
end function
public static double code(double a, double b) {
return (Math.pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0;
}
def code(a, b): return (math.pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0
function code(a, b) return Float64(Float64((Float64(Float64(a * a) + Float64(b * b)) ^ 2.0) + Float64(4.0 * Float64(Float64(Float64(a * a) * Float64(1.0 + a)) + Float64(Float64(b * b) * Float64(1.0 - Float64(3.0 * a)))))) - 1.0) end
function tmp = code(a, b) tmp = ((((a * a) + (b * b)) ^ 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0; end
code[a_, b_] := N[(N[(N[Power[N[(N[(a * a), $MachinePrecision] + N[(b * b), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] + N[(4.0 * N[(N[(N[(a * a), $MachinePrecision] * N[(1.0 + a), $MachinePrecision]), $MachinePrecision] + N[(N[(b * b), $MachinePrecision] * N[(1.0 - N[(3.0 * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - 1.0), $MachinePrecision]
\begin{array}{l}
\\
\left({\left(a \cdot a + b \cdot b\right)}^{2} + 4 \cdot \left(\left(a \cdot a\right) \cdot \left(1 + a\right) + \left(b \cdot b\right) \cdot \left(1 - 3 \cdot a\right)\right)\right) - 1
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 4 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (a b) :precision binary64 (- (+ (pow (+ (* a a) (* b b)) 2.0) (* 4.0 (+ (* (* a a) (+ 1.0 a)) (* (* b b) (- 1.0 (* 3.0 a)))))) 1.0))
double code(double a, double b) {
return (pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0;
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
code = ((((a * a) + (b * b)) ** 2.0d0) + (4.0d0 * (((a * a) * (1.0d0 + a)) + ((b * b) * (1.0d0 - (3.0d0 * a)))))) - 1.0d0
end function
public static double code(double a, double b) {
return (Math.pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0;
}
def code(a, b): return (math.pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0
function code(a, b) return Float64(Float64((Float64(Float64(a * a) + Float64(b * b)) ^ 2.0) + Float64(4.0 * Float64(Float64(Float64(a * a) * Float64(1.0 + a)) + Float64(Float64(b * b) * Float64(1.0 - Float64(3.0 * a)))))) - 1.0) end
function tmp = code(a, b) tmp = ((((a * a) + (b * b)) ^ 2.0) + (4.0 * (((a * a) * (1.0 + a)) + ((b * b) * (1.0 - (3.0 * a)))))) - 1.0; end
code[a_, b_] := N[(N[(N[Power[N[(N[(a * a), $MachinePrecision] + N[(b * b), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] + N[(4.0 * N[(N[(N[(a * a), $MachinePrecision] * N[(1.0 + a), $MachinePrecision]), $MachinePrecision] + N[(N[(b * b), $MachinePrecision] * N[(1.0 - N[(3.0 * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - 1.0), $MachinePrecision]
\begin{array}{l}
\\
\left({\left(a \cdot a + b \cdot b\right)}^{2} + 4 \cdot \left(\left(a \cdot a\right) \cdot \left(1 + a\right) + \left(b \cdot b\right) \cdot \left(1 - 3 \cdot a\right)\right)\right) - 1
\end{array}
(FPCore (a b)
:precision binary64
(if (<=
(+
(pow (+ (* a a) (* b b)) 2.0)
(* 4.0 (+ (* (* a a) (+ a 1.0)) (* (* b b) (- 1.0 (* a 3.0))))))
INFINITY)
(+
(fma
4.0
(fma (* b b) (fma a -3.0 1.0) (fma a a (pow a 3.0)))
(pow (fma a a (* b b)) 2.0))
-1.0)
(+ -1.0 (pow a 4.0))))
double code(double a, double b) {
double tmp;
if ((pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (a + 1.0)) + ((b * b) * (1.0 - (a * 3.0)))))) <= ((double) INFINITY)) {
tmp = fma(4.0, fma((b * b), fma(a, -3.0, 1.0), fma(a, a, pow(a, 3.0))), pow(fma(a, a, (b * b)), 2.0)) + -1.0;
} else {
tmp = -1.0 + pow(a, 4.0);
}
return tmp;
}
function code(a, b) tmp = 0.0 if (Float64((Float64(Float64(a * a) + Float64(b * b)) ^ 2.0) + Float64(4.0 * Float64(Float64(Float64(a * a) * Float64(a + 1.0)) + Float64(Float64(b * b) * Float64(1.0 - Float64(a * 3.0)))))) <= Inf) tmp = Float64(fma(4.0, fma(Float64(b * b), fma(a, -3.0, 1.0), fma(a, a, (a ^ 3.0))), (fma(a, a, Float64(b * b)) ^ 2.0)) + -1.0); else tmp = Float64(-1.0 + (a ^ 4.0)); end return tmp end
code[a_, b_] := If[LessEqual[N[(N[Power[N[(N[(a * a), $MachinePrecision] + N[(b * b), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] + N[(4.0 * N[(N[(N[(a * a), $MachinePrecision] * N[(a + 1.0), $MachinePrecision]), $MachinePrecision] + N[(N[(b * b), $MachinePrecision] * N[(1.0 - N[(a * 3.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], Infinity], N[(N[(4.0 * N[(N[(b * b), $MachinePrecision] * N[(a * -3.0 + 1.0), $MachinePrecision] + N[(a * a + N[Power[a, 3.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] + N[Power[N[(a * a + N[(b * b), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] + -1.0), $MachinePrecision], N[(-1.0 + N[Power[a, 4.0], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;{\left(a \cdot a + b \cdot b\right)}^{2} + 4 \cdot \left(\left(a \cdot a\right) \cdot \left(a + 1\right) + \left(b \cdot b\right) \cdot \left(1 - a \cdot 3\right)\right) \leq \infty:\\
\;\;\;\;\mathsf{fma}\left(4, \mathsf{fma}\left(b \cdot b, \mathsf{fma}\left(a, -3, 1\right), \mathsf{fma}\left(a, a, {a}^{3}\right)\right), {\left(\mathsf{fma}\left(a, a, b \cdot b\right)\right)}^{2}\right) + -1\\
\mathbf{else}:\\
\;\;\;\;-1 + {a}^{4}\\
\end{array}
\end{array}
if (+.f64 (pow.f64 (+.f64 (*.f64 a a) (*.f64 b b)) 2) (*.f64 4 (+.f64 (*.f64 (*.f64 a a) (+.f64 1 a)) (*.f64 (*.f64 b b) (-.f64 1 (*.f64 3 a)))))) < +inf.0Initial program 99.8%
associate--l+99.8%
+-commutative99.8%
+-commutative99.8%
sub-neg99.8%
associate-+l+99.8%
+-commutative99.8%
associate-+l+99.8%
Simplified99.8%
if +inf.0 < (+.f64 (pow.f64 (+.f64 (*.f64 a a) (*.f64 b b)) 2) (*.f64 4 (+.f64 (*.f64 (*.f64 a a) (+.f64 1 a)) (*.f64 (*.f64 b b) (-.f64 1 (*.f64 3 a)))))) Initial program 0.0%
sub-neg0.0%
Simplified5.6%
Taylor expanded in a around inf 96.1%
Final simplification98.8%
(FPCore (a b)
:precision binary64
(let* ((t_0
(+
(pow (+ (* a a) (* b b)) 2.0)
(* 4.0 (+ (* (* a a) (+ a 1.0)) (* (* b b) (- 1.0 (* a 3.0))))))))
(if (<= t_0 INFINITY) (+ t_0 -1.0) (+ -1.0 (pow a 4.0)))))
double code(double a, double b) {
double t_0 = pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (a + 1.0)) + ((b * b) * (1.0 - (a * 3.0)))));
double tmp;
if (t_0 <= ((double) INFINITY)) {
tmp = t_0 + -1.0;
} else {
tmp = -1.0 + pow(a, 4.0);
}
return tmp;
}
public static double code(double a, double b) {
double t_0 = Math.pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (a + 1.0)) + ((b * b) * (1.0 - (a * 3.0)))));
double tmp;
if (t_0 <= Double.POSITIVE_INFINITY) {
tmp = t_0 + -1.0;
} else {
tmp = -1.0 + Math.pow(a, 4.0);
}
return tmp;
}
def code(a, b): t_0 = math.pow(((a * a) + (b * b)), 2.0) + (4.0 * (((a * a) * (a + 1.0)) + ((b * b) * (1.0 - (a * 3.0))))) tmp = 0 if t_0 <= math.inf: tmp = t_0 + -1.0 else: tmp = -1.0 + math.pow(a, 4.0) return tmp
function code(a, b) t_0 = Float64((Float64(Float64(a * a) + Float64(b * b)) ^ 2.0) + Float64(4.0 * Float64(Float64(Float64(a * a) * Float64(a + 1.0)) + Float64(Float64(b * b) * Float64(1.0 - Float64(a * 3.0)))))) tmp = 0.0 if (t_0 <= Inf) tmp = Float64(t_0 + -1.0); else tmp = Float64(-1.0 + (a ^ 4.0)); end return tmp end
function tmp_2 = code(a, b) t_0 = (((a * a) + (b * b)) ^ 2.0) + (4.0 * (((a * a) * (a + 1.0)) + ((b * b) * (1.0 - (a * 3.0))))); tmp = 0.0; if (t_0 <= Inf) tmp = t_0 + -1.0; else tmp = -1.0 + (a ^ 4.0); end tmp_2 = tmp; end
code[a_, b_] := Block[{t$95$0 = N[(N[Power[N[(N[(a * a), $MachinePrecision] + N[(b * b), $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision] + N[(4.0 * N[(N[(N[(a * a), $MachinePrecision] * N[(a + 1.0), $MachinePrecision]), $MachinePrecision] + N[(N[(b * b), $MachinePrecision] * N[(1.0 - N[(a * 3.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]}, If[LessEqual[t$95$0, Infinity], N[(t$95$0 + -1.0), $MachinePrecision], N[(-1.0 + N[Power[a, 4.0], $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := {\left(a \cdot a + b \cdot b\right)}^{2} + 4 \cdot \left(\left(a \cdot a\right) \cdot \left(a + 1\right) + \left(b \cdot b\right) \cdot \left(1 - a \cdot 3\right)\right)\\
\mathbf{if}\;t\_0 \leq \infty:\\
\;\;\;\;t\_0 + -1\\
\mathbf{else}:\\
\;\;\;\;-1 + {a}^{4}\\
\end{array}
\end{array}
if (+.f64 (pow.f64 (+.f64 (*.f64 a a) (*.f64 b b)) 2) (*.f64 4 (+.f64 (*.f64 (*.f64 a a) (+.f64 1 a)) (*.f64 (*.f64 b b) (-.f64 1 (*.f64 3 a)))))) < +inf.0Initial program 99.8%
if +inf.0 < (+.f64 (pow.f64 (+.f64 (*.f64 a a) (*.f64 b b)) 2) (*.f64 4 (+.f64 (*.f64 (*.f64 a a) (+.f64 1 a)) (*.f64 (*.f64 b b) (-.f64 1 (*.f64 3 a)))))) Initial program 0.0%
sub-neg0.0%
Simplified5.6%
Taylor expanded in a around inf 96.1%
Final simplification98.8%
(FPCore (a b) :precision binary64 (if (or (<= a -5.4e+75) (not (<= a 1.26e+47))) (+ -1.0 (pow a 4.0)) (+ -1.0 (pow b 4.0))))
double code(double a, double b) {
double tmp;
if ((a <= -5.4e+75) || !(a <= 1.26e+47)) {
tmp = -1.0 + pow(a, 4.0);
} else {
tmp = -1.0 + pow(b, 4.0);
}
return tmp;
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8) :: tmp
if ((a <= (-5.4d+75)) .or. (.not. (a <= 1.26d+47))) then
tmp = (-1.0d0) + (a ** 4.0d0)
else
tmp = (-1.0d0) + (b ** 4.0d0)
end if
code = tmp
end function
public static double code(double a, double b) {
double tmp;
if ((a <= -5.4e+75) || !(a <= 1.26e+47)) {
tmp = -1.0 + Math.pow(a, 4.0);
} else {
tmp = -1.0 + Math.pow(b, 4.0);
}
return tmp;
}
def code(a, b): tmp = 0 if (a <= -5.4e+75) or not (a <= 1.26e+47): tmp = -1.0 + math.pow(a, 4.0) else: tmp = -1.0 + math.pow(b, 4.0) return tmp
function code(a, b) tmp = 0.0 if ((a <= -5.4e+75) || !(a <= 1.26e+47)) tmp = Float64(-1.0 + (a ^ 4.0)); else tmp = Float64(-1.0 + (b ^ 4.0)); end return tmp end
function tmp_2 = code(a, b) tmp = 0.0; if ((a <= -5.4e+75) || ~((a <= 1.26e+47))) tmp = -1.0 + (a ^ 4.0); else tmp = -1.0 + (b ^ 4.0); end tmp_2 = tmp; end
code[a_, b_] := If[Or[LessEqual[a, -5.4e+75], N[Not[LessEqual[a, 1.26e+47]], $MachinePrecision]], N[(-1.0 + N[Power[a, 4.0], $MachinePrecision]), $MachinePrecision], N[(-1.0 + N[Power[b, 4.0], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;a \leq -5.4 \cdot 10^{+75} \lor \neg \left(a \leq 1.26 \cdot 10^{+47}\right):\\
\;\;\;\;-1 + {a}^{4}\\
\mathbf{else}:\\
\;\;\;\;-1 + {b}^{4}\\
\end{array}
\end{array}
if a < -5.39999999999999996e75 or 1.26e47 < a Initial program 34.6%
sub-neg34.6%
Simplified38.3%
Taylor expanded in a around inf 99.1%
if -5.39999999999999996e75 < a < 1.26e47Initial program 98.4%
sub-neg98.4%
Simplified98.4%
Taylor expanded in b around inf 93.9%
Final simplification96.1%
(FPCore (a b) :precision binary64 (+ -1.0 (pow a 4.0)))
double code(double a, double b) {
return -1.0 + pow(a, 4.0);
}
real(8) function code(a, b)
real(8), intent (in) :: a
real(8), intent (in) :: b
code = (-1.0d0) + (a ** 4.0d0)
end function
public static double code(double a, double b) {
return -1.0 + Math.pow(a, 4.0);
}
def code(a, b): return -1.0 + math.pow(a, 4.0)
function code(a, b) return Float64(-1.0 + (a ^ 4.0)) end
function tmp = code(a, b) tmp = -1.0 + (a ^ 4.0); end
code[a_, b_] := N[(-1.0 + N[Power[a, 4.0], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
-1 + {a}^{4}
\end{array}
Initial program 71.7%
sub-neg71.7%
Simplified73.3%
Taylor expanded in a around inf 65.8%
Final simplification65.8%
herbie shell --seed 2024043
(FPCore (a b)
:name "Bouland and Aaronson, Equation (25)"
:precision binary64
(- (+ (pow (+ (* a a) (* b b)) 2.0) (* 4.0 (+ (* (* a a) (+ 1.0 a)) (* (* b b) (- 1.0 (* 3.0 a)))))) 1.0))