Numeric.SpecFunctions:incompleteBetaApprox from math-functions-0.1.5.2, A

Average Error: 19.8 → 7.5

Time: 5.4s

Precision: binary64

\[ \begin{array}{c}[x, y] = \mathsf{sort}([x, y])\\ \end{array} \]

Math

\[\frac{x \cdot y}{\left(\left(x + y\right) \cdot \left(x + y\right)\right) \cdot \left(\left(x + y\right) + 1\right)} \]

↓

\[\begin{array}{l} t_0 := \frac{x}{\frac{\mathsf{fma}\left(y + x, y + x, {\left(y + x\right)}^{3}\right)}{y}}\\ t_1 := \frac{\mathsf{hypot}\left(y + x, {\left(y + x\right)}^{1.5}\right)}{\sqrt{y}}\\ \mathbf{if}\;y \leq -2.002618035851781 \cdot 10^{-307}:\\ \;\;\;\;t_0\\ \mathbf{elif}\;y \leq 7.100604378657969 \cdot 10^{-174}:\\ \;\;\;\;\frac{x}{\mathsf{fma}\left(3, x \cdot x, \mathsf{fma}\left(x, 2, \mathsf{fma}\left(y, \mathsf{fma}\left(x, 3, 1\right), \mathsf{fma}\left(y, y, \frac{x}{\frac{y}{x}}\right)\right)\right)\right)}\\ \mathbf{elif}\;y \leq 5.82500637883476 \cdot 10^{+90}:\\ \;\;\;\;t_0\\ \mathbf{elif}\;y \leq 1.1768391999675483 \cdot 10^{+153}:\\ \;\;\;\;\frac{x}{y \cdot \left(y + \mathsf{fma}\left(3, x, 1\right)\right)}\\ \mathbf{else}:\\ \;\;\;\;\frac{1}{t_1} \cdot \frac{x}{t_1}\\ \end{array} \]

FPCore

(FPCore (x y)
 :precision binary64
 (/ (* x y) (* (* (+ x y) (+ x y)) (+ (+ x y) 1.0))))

↓

(FPCore (x y)
 :precision binary64
 (let* ((t_0 (/ x (/ (fma (+ y x) (+ y x) (pow (+ y x) 3.0)) y)))
        (t_1 (/ (hypot (+ y x) (pow (+ y x) 1.5)) (sqrt y))))
   (if (<= y -2.002618035851781e-307)
     t_0
     (if (<= y 7.100604378657969e-174)
       (/
        x
        (fma
         3.0
         (* x x)
         (fma x 2.0 (fma y (fma x 3.0 1.0) (fma y y (/ x (/ y x)))))))
       (if (<= y 5.82500637883476e+90)
         t_0
         (if (<= y 1.1768391999675483e+153)
           (/ x (* y (+ y (fma 3.0 x 1.0))))
           (* (/ 1.0 t_1) (/ x t_1))))))))

C

double code(double x, double y) {
	return (x * y) / (((x + y) * (x + y)) * ((x + y) + 1.0));
}

↓

double code(double x, double y) {
	double t_0 = x / (fma((y + x), (y + x), pow((y + x), 3.0)) / y);
	double t_1 = hypot((y + x), pow((y + x), 1.5)) / sqrt(y);
	double tmp;
	if (y <= -2.002618035851781e-307) {
		tmp = t_0;
	} else if (y <= 7.100604378657969e-174) {
		tmp = x / fma(3.0, (x * x), fma(x, 2.0, fma(y, fma(x, 3.0, 1.0), fma(y, y, (x / (y / x))))));
	} else if (y <= 5.82500637883476e+90) {
		tmp = t_0;
	} else if (y <= 1.1768391999675483e+153) {
		tmp = x / (y * (y + fma(3.0, x, 1.0)));
	} else {
		tmp = (1.0 / t_1) * (x / t_1);
	}
	return tmp;
}

Julia

function code(x, y)
	return Float64(Float64(x * y) / Float64(Float64(Float64(x + y) * Float64(x + y)) * Float64(Float64(x + y) + 1.0)))
end

↓

function code(x, y)
	t_0 = Float64(x / Float64(fma(Float64(y + x), Float64(y + x), (Float64(y + x) ^ 3.0)) / y))
	t_1 = Float64(hypot(Float64(y + x), (Float64(y + x) ^ 1.5)) / sqrt(y))
	tmp = 0.0
	if (y <= -2.002618035851781e-307)
		tmp = t_0;
	elseif (y <= 7.100604378657969e-174)
		tmp = Float64(x / fma(3.0, Float64(x * x), fma(x, 2.0, fma(y, fma(x, 3.0, 1.0), fma(y, y, Float64(x / Float64(y / x)))))));
	elseif (y <= 5.82500637883476e+90)
		tmp = t_0;
	elseif (y <= 1.1768391999675483e+153)
		tmp = Float64(x / Float64(y * Float64(y + fma(3.0, x, 1.0))));
	else
		tmp = Float64(Float64(1.0 / t_1) * Float64(x / t_1));
	end
	return tmp
end

Wolfram

code[x_, y_] := N[(N[(x * y), $MachinePrecision] / N[(N[(N[(x + y), $MachinePrecision] * N[(x + y), $MachinePrecision]), $MachinePrecision] * N[(N[(x + y), $MachinePrecision] + 1.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]

↓

code[x_, y_] := Block[{t$95$0 = N[(x / N[(N[(N[(y + x), $MachinePrecision] * N[(y + x), $MachinePrecision] + N[Power[N[(y + x), $MachinePrecision], 3.0], $MachinePrecision]), $MachinePrecision] / y), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$1 = N[(N[Sqrt[N[(y + x), $MachinePrecision] ^ 2 + N[Power[N[(y + x), $MachinePrecision], 1.5], $MachinePrecision] ^ 2], $MachinePrecision] / N[Sqrt[y], $MachinePrecision]), $MachinePrecision]}, If[LessEqual[y, -2.002618035851781e-307], t$95$0, If[LessEqual[y, 7.100604378657969e-174], N[(x / N[(3.0 * N[(x * x), $MachinePrecision] + N[(x * 2.0 + N[(y * N[(x * 3.0 + 1.0), $MachinePrecision] + N[(y * y + N[(x / N[(y / x), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], If[LessEqual[y, 5.82500637883476e+90], t$95$0, If[LessEqual[y, 1.1768391999675483e+153], N[(x / N[(y * N[(y + N[(3.0 * x + 1.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(N[(1.0 / t$95$1), $MachinePrecision] * N[(x / t$95$1), $MachinePrecision]), $MachinePrecision]]]]]]]

Target

Original	19.8
Target	0.1
Herbie	7.5

\[\frac{\frac{\frac{x}{\left(y + 1\right) + x}}{y + x}}{\frac{1}{\frac{y}{y + x}}} \]

Derivation

1. Split input into 4 regimes

2. if y < -2.00261803585178102e-307 or 7.10060437865796928e-174 < y < 5.82500637883476e90

   1. Initial program 14.9

      \[\frac{x \cdot y}{\left(\left(x + y\right) \cdot \left(x + y\right)\right) \cdot \left(\left(x + y\right) + 1\right)} \]

   2. Simplified 8.4

      \[\leadsto \color{blue}{\frac{x}{\frac{\mathsf{fma}\left(x + y, x + y, {\left(x + y\right)}^{3}\right)}{y}}} \]

   if -2.00261803585178102e-307 < y < 7.10060437865796928e-174

   1. Initial program 27.5

      \[\frac{x \cdot y}{\left(\left(x + y\right) \cdot \left(x + y\right)\right) \cdot \left(\left(x + y\right) + 1\right)} \]

   2. Simplified 21.3

      \[\leadsto \color{blue}{\frac{x}{\frac{\mathsf{fma}\left(x + y, x + y, {\left(x + y\right)}^{3}\right)}{y}}} \]

   3. Taylor expanded in x around 0 11.7

      \[\leadsto \frac{x}{\color{blue}{\left(\frac{1}{y} + 3\right) \cdot {x}^{2} + \left(y + \left({y}^{2} + \left(2 + \left(2 \cdot y + y\right)\right) \cdot x\right)\right)}} \]

   4. Simplified 11.7

      \[\leadsto \frac{x}{\color{blue}{\mathsf{fma}\left(3 + \frac{1}{y}, x \cdot x, y\right) + \mathsf{fma}\left(y, y, \mathsf{fma}\left(3, y, 2\right) \cdot x\right)}} \]

   5. Taylor expanded in y around 0 11.6

      \[\leadsto \frac{x}{\color{blue}{3 \cdot {x}^{2} + \left(2 \cdot x + \left(\left(3 \cdot x + 1\right) \cdot y + \left({y}^{2} + \frac{{x}^{2}}{y}\right)\right)\right)}} \]

   6. Simplified 4.1

      \[\leadsto \frac{x}{\color{blue}{\mathsf{fma}\left(3, x \cdot x, \mathsf{fma}\left(x, 2, \mathsf{fma}\left(y, \mathsf{fma}\left(x, 3, 1\right), \mathsf{fma}\left(y, y, \frac{x}{\frac{y}{x}}\right)\right)\right)\right)}} \]

   if 5.82500637883476e90 < y < 1.176839199967548e153

   1. Initial program 28.0

      \[\frac{x \cdot y}{\left(\left(x + y\right) \cdot \left(x + y\right)\right) \cdot \left(\left(x + y\right) + 1\right)} \]

   2. Simplified 22.5

      \[\leadsto \color{blue}{\frac{x}{\frac{\mathsf{fma}\left(x + y, x + y, {\left(x + y\right)}^{3}\right)}{y}}} \]

   3. Taylor expanded in y around inf 5.3

      \[\leadsto \frac{x}{\color{blue}{y \cdot \left(1 + \left(2 \cdot x + x\right)\right) + {y}^{2}}} \]

   4. Simplified 5.3

      \[\leadsto \frac{x}{\color{blue}{y \cdot \left(\mathsf{fma}\left(3, x, 1\right) + y\right)}} \]

   if 1.176839199967548e153 < y

   1. Initial program 24.2

      \[\frac{x \cdot y}{\left(\left(x + y\right) \cdot \left(x + y\right)\right) \cdot \left(\left(x + y\right) + 1\right)} \]

   2. Simplified 11.3

      \[\leadsto \color{blue}{\frac{x}{\frac{\mathsf{fma}\left(x + y, x + y, {\left(x + y\right)}^{3}\right)}{y}}} \]

   3. Applied egg-rr 8.0

      \[\leadsto \color{blue}{\frac{1}{\frac{\mathsf{hypot}\left(x + y, {\left(x + y\right)}^{1.5}\right)}{\sqrt{y}}} \cdot \frac{x}{\frac{\mathsf{hypot}\left(x + y, {\left(x + y\right)}^{1.5}\right)}{\sqrt{y}}}} \]
3. Recombined 4 regimes into one program.

4. Final simplification 7.5

   \[\leadsto \begin{array}{l} \mathbf{if}\;y \leq -2.002618035851781 \cdot 10^{-307}:\\ \;\;\;\;\frac{x}{\frac{\mathsf{fma}\left(y + x, y + x, {\left(y + x\right)}^{3}\right)}{y}}\\ \mathbf{elif}\;y \leq 7.100604378657969 \cdot 10^{-174}:\\ \;\;\;\;\frac{x}{\mathsf{fma}\left(3, x \cdot x, \mathsf{fma}\left(x, 2, \mathsf{fma}\left(y, \mathsf{fma}\left(x, 3, 1\right), \mathsf{fma}\left(y, y, \frac{x}{\frac{y}{x}}\right)\right)\right)\right)}\\ \mathbf{elif}\;y \leq 5.82500637883476 \cdot 10^{+90}:\\ \;\;\;\;\frac{x}{\frac{\mathsf{fma}\left(y + x, y + x, {\left(y + x\right)}^{3}\right)}{y}}\\ \mathbf{elif}\;y \leq 1.1768391999675483 \cdot 10^{+153}:\\ \;\;\;\;\frac{x}{y \cdot \left(y + \mathsf{fma}\left(3, x, 1\right)\right)}\\ \mathbf{else}:\\ \;\;\;\;\frac{1}{\frac{\mathsf{hypot}\left(y + x, {\left(y + x\right)}^{1.5}\right)}{\sqrt{y}}} \cdot \frac{x}{\frac{\mathsf{hypot}\left(y + x, {\left(y + x\right)}^{1.5}\right)}{\sqrt{y}}}\\ \end{array} \]

Reproduce

herbie shell --seed 2022210 
(FPCore (x y)
  :name "Numeric.SpecFunctions:incompleteBetaApprox from math-functions-0.1.5.2, A"
  :precision binary64

  :herbie-target
  (/ (/ (/ x (+ (+ y 1.0) x)) (+ y x)) (/ 1.0 (/ y (+ y x))))

  (/ (* x y) (* (* (+ x y) (+ x y)) (+ (+ x y) 1.0))))