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

Average Error: 1.9 → 2.2

Time: 28.4s

Precision: binary64

Cost: 26692

Math

\[\frac{x \cdot e^{\left(y \cdot \log z + \left(t - 1\right) \cdot \log a\right) - b}}{y} \]

↓

\[\begin{array}{l} \mathbf{if}\;\left(t - 1\right) \cdot \log a \leq -684.5:\\ \;\;\;\;x \cdot \frac{\frac{{a}^{t}}{a}}{y}\\ \mathbf{else}:\\ \;\;\;\;\frac{x \cdot e^{\left(y \cdot \log z - \log a\right) - b}}{y}\\ \end{array} \]

FPCore

(FPCore (x y z t a b)
 :precision binary64
 (/ (* x (exp (- (+ (* y (log z)) (* (- t 1.0) (log a))) b))) y))

↓

(FPCore (x y z t a b)
 :precision binary64
 (if (<= (* (- t 1.0) (log a)) -684.5)
   (* x (/ (/ (pow a t) a) y))
   (/ (* x (exp (- (- (* y (log z)) (log a)) b))) y)))

C

double code(double x, double y, double z, double t, double a, double b) {
	return (x * exp((((y * log(z)) + ((t - 1.0) * log(a))) - b))) / y;
}

↓

double code(double x, double y, double z, double t, double a, double b) {
	double tmp;
	if (((t - 1.0) * log(a)) <= -684.5) {
		tmp = x * ((pow(a, t) / a) / y);
	} else {
		tmp = (x * exp((((y * log(z)) - log(a)) - b))) / y;
	}
	return tmp;
}

Fortran

real(8) function code(x, y, z, t, a, b)
    real(8), intent (in) :: x
    real(8), intent (in) :: y
    real(8), intent (in) :: z
    real(8), intent (in) :: t
    real(8), intent (in) :: a
    real(8), intent (in) :: b
    code = (x * exp((((y * log(z)) + ((t - 1.0d0) * log(a))) - b))) / y
end function

↓

real(8) function code(x, y, z, t, a, b)
    real(8), intent (in) :: x
    real(8), intent (in) :: y
    real(8), intent (in) :: z
    real(8), intent (in) :: t
    real(8), intent (in) :: a
    real(8), intent (in) :: b
    real(8) :: tmp
    if (((t - 1.0d0) * log(a)) <= (-684.5d0)) then
        tmp = x * (((a ** t) / a) / y)
    else
        tmp = (x * exp((((y * log(z)) - log(a)) - b))) / y
    end if
    code = tmp
end function

Java

public static double code(double x, double y, double z, double t, double a, double b) {
	return (x * Math.exp((((y * Math.log(z)) + ((t - 1.0) * Math.log(a))) - b))) / y;
}

↓

public static double code(double x, double y, double z, double t, double a, double b) {
	double tmp;
	if (((t - 1.0) * Math.log(a)) <= -684.5) {
		tmp = x * ((Math.pow(a, t) / a) / y);
	} else {
		tmp = (x * Math.exp((((y * Math.log(z)) - Math.log(a)) - b))) / y;
	}
	return tmp;
}

Python

def code(x, y, z, t, a, b):
	return (x * math.exp((((y * math.log(z)) + ((t - 1.0) * math.log(a))) - b))) / y

↓

def code(x, y, z, t, a, b):
	tmp = 0
	if ((t - 1.0) * math.log(a)) <= -684.5:
		tmp = x * ((math.pow(a, t) / a) / y)
	else:
		tmp = (x * math.exp((((y * math.log(z)) - math.log(a)) - b))) / y
	return tmp

Julia

function code(x, y, z, t, a, b)
	return Float64(Float64(x * exp(Float64(Float64(Float64(y * log(z)) + Float64(Float64(t - 1.0) * log(a))) - b))) / y)
end

↓

function code(x, y, z, t, a, b)
	tmp = 0.0
	if (Float64(Float64(t - 1.0) * log(a)) <= -684.5)
		tmp = Float64(x * Float64(Float64((a ^ t) / a) / y));
	else
		tmp = Float64(Float64(x * exp(Float64(Float64(Float64(y * log(z)) - log(a)) - b))) / y);
	end
	return tmp
end

MATLAB

function tmp = code(x, y, z, t, a, b)
	tmp = (x * exp((((y * log(z)) + ((t - 1.0) * log(a))) - b))) / y;
end

↓

function tmp_2 = code(x, y, z, t, a, b)
	tmp = 0.0;
	if (((t - 1.0) * log(a)) <= -684.5)
		tmp = x * (((a ^ t) / a) / y);
	else
		tmp = (x * exp((((y * log(z)) - log(a)) - b))) / y;
	end
	tmp_2 = tmp;
end

Wolfram

code[x_, y_, z_, t_, a_, b_] := N[(N[(x * N[Exp[N[(N[(N[(y * N[Log[z], $MachinePrecision]), $MachinePrecision] + N[(N[(t - 1.0), $MachinePrecision] * N[Log[a], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - b), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] / y), $MachinePrecision]

↓

code[x_, y_, z_, t_, a_, b_] := If[LessEqual[N[(N[(t - 1.0), $MachinePrecision] * N[Log[a], $MachinePrecision]), $MachinePrecision], -684.5], N[(x * N[(N[(N[Power[a, t], $MachinePrecision] / a), $MachinePrecision] / y), $MachinePrecision]), $MachinePrecision], N[(N[(x * N[Exp[N[(N[(N[(y * N[Log[z], $MachinePrecision]), $MachinePrecision] - N[Log[a], $MachinePrecision]), $MachinePrecision] - b), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] / y), $MachinePrecision]]

Target

Original	1.9
Target	11.4
Herbie	2.2

\[\begin{array}{l} \mathbf{if}\;t < -0.8845848504127471:\\ \;\;\;\;\frac{x \cdot \frac{{a}^{\left(t - 1\right)}}{y}}{\left(b + 1\right) - y \cdot \log z}\\ \mathbf{elif}\;t < 852031.2288374073:\\ \;\;\;\;\frac{\frac{x}{y} \cdot {a}^{\left(t - 1\right)}}{e^{b - \log z \cdot y}}\\ \mathbf{else}:\\ \;\;\;\;\frac{x \cdot \frac{{a}^{\left(t - 1\right)}}{y}}{\left(b + 1\right) - y \cdot \log z}\\ \end{array} \]

Derivation

1. Split input into 2 regimes

2. if (*.f64 (-.f64 t 1) (log.f64 a)) < -684.5

   1. Initial program 0.1

      \[\frac{x \cdot e^{\left(y \cdot \log z + \left(t - 1\right) \cdot \log a\right) - b}}{y} \]

   2. Simplified 18.4

      \[\leadsto \color{blue}{x \cdot \frac{\frac{{a}^{t}}{a}}{\frac{y \cdot e^{b}}{{z}^{y}}}} \]
      Proof
      (*.f64 x (/.f64 (/.f64 (pow.f64 a t) a) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (/.f64 (Rewrite<= exp-to-pow_binary64 (exp.f64 (*.f64 (log.f64 a) t))) a) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (/.f64 (exp.f64 (*.f64 (log.f64 a) t)) (Rewrite<= rem-exp-log_binary64 (exp.f64 (log.f64 a)))) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 36 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (Rewrite<= exp-diff_binary64 (exp.f64 (-.f64 (*.f64 (log.f64 a) t) (log.f64 a)))) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 5 points increase in error, 6 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (-.f64 (Rewrite=> *-commutative_binary64 (*.f64 t (log.f64 a))) (log.f64 a))) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (-.f64 (*.f64 t (log.f64 a)) (Rewrite<= *-lft-identity_binary64 (*.f64 1 (log.f64 a))))) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (Rewrite=> distribute-rgt-out--_binary64 (*.f64 (log.f64 a) (-.f64 t 1)))) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (Rewrite<= *-commutative_binary64 (*.f64 (-.f64 t 1) (log.f64 a)))) (/.f64 (*.f64 y (exp.f64 b)) (pow.f64 z y)))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (*.f64 (-.f64 t 1) (log.f64 a))) (/.f64 (*.f64 y (exp.f64 b)) (Rewrite<= exp-to-pow_binary64 (exp.f64 (*.f64 (log.f64 z) y)))))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (*.f64 (-.f64 t 1) (log.f64 a))) (/.f64 (*.f64 y (exp.f64 b)) (exp.f64 (Rewrite<= *-commutative_binary64 (*.f64 y (log.f64 z))))))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (*.f64 (-.f64 t 1) (log.f64 a))) (Rewrite<= associate-*l/_binary64 (*.f64 (/.f64 y (exp.f64 (*.f64 y (log.f64 z)))) (exp.f64 b))))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (Rewrite<= associate-/l/_binary64 (/.f64 (/.f64 (exp.f64 (*.f64 (-.f64 t 1) (log.f64 a))) (exp.f64 b)) (/.f64 y (exp.f64 (*.f64 y (log.f64 z))))))): 1 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (Rewrite<= exp-diff_binary64 (exp.f64 (-.f64 (*.f64 (-.f64 t 1) (log.f64 a)) b))) (/.f64 y (exp.f64 (*.f64 y (log.f64 z)))))): 1 points increase in error, 31 points decrease in error
      (*.f64 x (Rewrite<= associate-/l*_binary64 (/.f64 (*.f64 (exp.f64 (-.f64 (*.f64 (-.f64 t 1) (log.f64 a)) b)) (exp.f64 (*.f64 y (log.f64 z)))) y))): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (Rewrite<= exp-sum_binary64 (exp.f64 (+.f64 (-.f64 (*.f64 (-.f64 t 1) (log.f64 a)) b) (*.f64 y (log.f64 z))))) y)): 3 points increase in error, 48 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (Rewrite<= +-commutative_binary64 (+.f64 (*.f64 y (log.f64 z)) (-.f64 (*.f64 (-.f64 t 1) (log.f64 a)) b)))) y)): 0 points increase in error, 0 points decrease in error
      (*.f64 x (/.f64 (exp.f64 (Rewrite<= associate--l+_binary64 (-.f64 (+.f64 (*.f64 y (log.f64 z)) (*.f64 (-.f64 t 1) (log.f64 a))) b))) y)): 0 points increase in error, 0 points decrease in error
      (Rewrite=> associate-*r/_binary64 (/.f64 (*.f64 x (exp.f64 (-.f64 (+.f64 (*.f64 y (log.f64 z)) (*.f64 (-.f64 t 1) (log.f64 a))) b))) y)): 11 points increase in error, 8 points decrease in error

   3. Taylor expanded in y around 0 9.6

      \[\leadsto x \cdot \frac{\frac{{a}^{t}}{a}}{\color{blue}{y \cdot e^{b}}} \]

   4. Taylor expanded in b around 0 0.2

      \[\leadsto x \cdot \frac{\frac{{a}^{t}}{a}}{\color{blue}{y}} \]

   if -684.5 < (*.f64 (-.f64 t 1) (log.f64 a))

   1. Initial program 2.9

      \[\frac{x \cdot e^{\left(y \cdot \log z + \left(t - 1\right) \cdot \log a\right) - b}}{y} \]

   2. Taylor expanded in t around 0 3.3

      \[\leadsto \frac{x \cdot e^{\left(y \cdot \log z + \color{blue}{-1 \cdot \log a}\right) - b}}{y} \]

   3. Simplified 3.3

      \[\leadsto \frac{x \cdot e^{\left(y \cdot \log z + \color{blue}{\left(-\log a\right)}\right) - b}}{y} \]
      Proof
      (neg.f64 (log.f64 a)): 0 points increase in error, 0 points decrease in error
      (Rewrite<= mul-1-neg_binary64 (*.f64 -1 (log.f64 a))): 0 points increase in error, 0 points decrease in error
3. Recombined 2 regimes into one program.

4. Final simplification 2.2

   \[\leadsto \begin{array}{l} \mathbf{if}\;\left(t - 1\right) \cdot \log a \leq -684.5:\\ \;\;\;\;x \cdot \frac{\frac{{a}^{t}}{a}}{y}\\ \mathbf{else}:\\ \;\;\;\;\frac{x \cdot e^{\left(y \cdot \log z - \log a\right) - b}}{y}\\ \end{array} \]

Alternatives

Alternative 1
Error	1.9
Cost	20160

\[\frac{x \cdot e^{\left(y \cdot \log z + \left(t - 1\right) \cdot \log a\right) - b}}{y} \]

Alternative 2
Error	8.2
Cost	13768

\[\begin{array}{l} t_1 := x \cdot \frac{{z}^{y}}{y \cdot a}\\ \mathbf{if}\;y \leq -4.4 \cdot 10^{+177}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;y \leq 1.2 \cdot 10^{+169}:\\ \;\;\;\;\frac{x \cdot e^{\left(t - 1\right) \cdot \log a - b}}{y}\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

Alternative 3
Error	11.4
Cost	7440

\[\begin{array}{l} t_1 := x \cdot \frac{{z}^{y}}{y \cdot a}\\ t_2 := x \cdot \frac{\frac{{a}^{t}}{a}}{y}\\ \mathbf{if}\;b \leq -1 \cdot 10^{-229}:\\ \;\;\;\;t_2\\ \mathbf{elif}\;b \leq -2.1 \cdot 10^{-306}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;b \leq 7.5 \cdot 10^{-230}:\\ \;\;\;\;t_2\\ \mathbf{elif}\;b \leq 1.35 \cdot 10^{-8}:\\ \;\;\;\;t_1\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{a \cdot \left(y \cdot e^{b}\right)}\\ \end{array} \]

Alternative 4
Error	11.4
Cost	7440

\[\begin{array}{l} t_1 := x \cdot \frac{\frac{{a}^{t}}{a}}{y}\\ \mathbf{if}\;b \leq -6.6 \cdot 10^{-229}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;b \leq -1.6 \cdot 10^{-306}:\\ \;\;\;\;\frac{x}{\frac{y \cdot a}{{z}^{y}}}\\ \mathbf{elif}\;b \leq 1.45 \cdot 10^{-230}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;b \leq 1.3 \cdot 10^{-8}:\\ \;\;\;\;x \cdot \frac{{z}^{y}}{y \cdot a}\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{a \cdot \left(y \cdot e^{b}\right)}\\ \end{array} \]

Alternative 5
Error	11.9
Cost	7440

\[\begin{array}{l} t_1 := \frac{{a}^{t}}{a}\\ \mathbf{if}\;b \leq -1.8 \cdot 10^{-230}:\\ \;\;\;\;\frac{x \cdot t_1}{y}\\ \mathbf{elif}\;b \leq -9.6 \cdot 10^{-307}:\\ \;\;\;\;\frac{x}{\frac{y \cdot a}{{z}^{y}}}\\ \mathbf{elif}\;b \leq 4.8 \cdot 10^{-229}:\\ \;\;\;\;x \cdot \frac{t_1}{y}\\ \mathbf{elif}\;b \leq 1.3 \cdot 10^{-8}:\\ \;\;\;\;x \cdot \frac{{z}^{y}}{y \cdot a}\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{a \cdot \left(y \cdot e^{b}\right)}\\ \end{array} \]

Alternative 6
Error	11.8
Cost	7440

\[\begin{array}{l} t_1 := \frac{{a}^{t}}{a}\\ \mathbf{if}\;b \leq -1.1 \cdot 10^{-227}:\\ \;\;\;\;\frac{x \cdot t_1}{y}\\ \mathbf{elif}\;b \leq -2.1 \cdot 10^{-306}:\\ \;\;\;\;\frac{\frac{x \cdot {z}^{y}}{y}}{a}\\ \mathbf{elif}\;b \leq 1.12 \cdot 10^{-228}:\\ \;\;\;\;x \cdot \frac{t_1}{y}\\ \mathbf{elif}\;b \leq 1.4 \cdot 10^{-8}:\\ \;\;\;\;x \cdot \frac{{z}^{y}}{y \cdot a}\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{a \cdot \left(y \cdot e^{b}\right)}\\ \end{array} \]

Alternative 7
Error	10.8
Cost	7044

\[\begin{array}{l} \mathbf{if}\;b \leq 8 \cdot 10^{-40}:\\ \;\;\;\;x \cdot \frac{\frac{{a}^{t}}{a}}{y}\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{a \cdot \left(y \cdot e^{b}\right)}\\ \end{array} \]

Alternative 8
Error	16.8
Cost	6980

\[\begin{array}{l} \mathbf{if}\;b \leq 8 \cdot 10^{-40}:\\ \;\;\;\;x \cdot \left(\left(1 + \frac{1}{y \cdot a}\right) + -1\right)\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{a \cdot \left(y \cdot e^{b}\right)}\\ \end{array} \]

Alternative 9
Error	29.9
Cost	1104

\[\begin{array}{l} t_1 := \left(1 + \frac{x}{y \cdot a}\right) + -1\\ \mathbf{if}\;y \leq -1.65 \cdot 10^{-88}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;y \leq 5 \cdot 10^{-253}:\\ \;\;\;\;\frac{1}{y} \cdot \frac{x}{a}\\ \mathbf{elif}\;y \leq 6.5 \cdot 10^{-205}:\\ \;\;\;\;\frac{b}{a} \cdot \frac{-x}{y}\\ \mathbf{elif}\;y \leq 3.3 \cdot 10^{-41}:\\ \;\;\;\;\frac{x}{a \cdot \left(y + y \cdot b\right)}\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

Alternative 10
Error	25.4
Cost	968

\[\begin{array}{l} t_1 := x \cdot \left(\left(1 + \frac{1}{y \cdot a}\right) + -1\right)\\ \mathbf{if}\;x \leq -1.5 \cdot 10^{+89}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;x \leq 5 \cdot 10^{-56}:\\ \;\;\;\;\left(1 + \frac{x}{y \cdot a}\right) + -1\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

Alternative 11
Error	29.2
Cost	840

\[\begin{array}{l} t_1 := \left(1 + \frac{x}{y \cdot a}\right) + -1\\ \mathbf{if}\;y \leq -3.4 \cdot 10^{-88}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;y \leq 3.1 \cdot 10^{-80}:\\ \;\;\;\;\frac{1}{y} \cdot \frac{x}{a}\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

Alternative 12
Error	26.5
Cost	840

\[\begin{array}{l} t_1 := \left(1 + \frac{x}{y \cdot a}\right) + -1\\ \mathbf{if}\;y \leq -3.25 \cdot 10^{-93}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;y \leq 1.35 \cdot 10^{-40}:\\ \;\;\;\;\frac{\frac{x}{a + a \cdot b}}{y}\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

Alternative 13
Error	39.7
Cost	776

\[\begin{array}{l} t_1 := \frac{b}{a} \cdot \frac{-x}{y}\\ \mathbf{if}\;y \leq -1.45 \cdot 10^{+19}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;y \leq 1.2:\\ \;\;\;\;\frac{1}{y} \cdot \frac{x}{a}\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

Alternative 14
Error	38.4
Cost	712

\[\begin{array}{l} \mathbf{if}\;a \leq 10^{+24}:\\ \;\;\;\;\frac{\frac{x}{y}}{a}\\ \mathbf{elif}\;a \leq 5.3 \cdot 10^{+97}:\\ \;\;\;\;\frac{1}{y} \cdot \frac{x}{a}\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{y \cdot a}\\ \end{array} \]

Alternative 15
Error	39.9
Cost	584

\[\begin{array}{l} t_1 := \frac{x}{y \cdot a}\\ \mathbf{if}\;x \leq -5.6 \cdot 10^{-151}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;x \leq 2 \cdot 10^{-33}:\\ \;\;\;\;\frac{\frac{x}{a}}{y}\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

Alternative 16
Error	38.4
Cost	584

\[\begin{array}{l} \mathbf{if}\;a \leq 2 \cdot 10^{+27}:\\ \;\;\;\;\frac{\frac{x}{y}}{a}\\ \mathbf{elif}\;a \leq 5 \cdot 10^{+97}:\\ \;\;\;\;\frac{\frac{x}{a}}{y}\\ \mathbf{else}:\\ \;\;\;\;\frac{x}{y \cdot a}\\ \end{array} \]

Alternative 17
Error	41.4
Cost	320

\[\frac{x}{y \cdot a} \]

Reproduce

herbie shell --seed 2022325 
(FPCore (x y z t a b)
  :name "Numeric.SpecFunctions:incompleteBetaWorker from math-functions-0.1.5.2, A"
  :precision binary64

  :herbie-target
  (if (< t -0.8845848504127471) (/ (* x (/ (pow a (- t 1.0)) y)) (- (+ b 1.0) (* y (log z)))) (if (< t 852031.2288374073) (/ (* (/ x y) (pow a (- t 1.0))) (exp (- b (* (log z) y)))) (/ (* x (/ (pow a (- t 1.0)) y)) (- (+ b 1.0) (* y (log z))))))

  (/ (* x (exp (- (+ (* y (log z)) (* (- t 1.0) (log a))) b))) y))