ENA, Section 1.4, Exercise 4d

Average Error: 24.5 → 0.3

Time: 5.6s

Precision: binary64

Cost: 13572

\[\left(0 \leq x \land x \leq 1000000000\right) \land \left(-1 \leq \varepsilon \land \varepsilon \leq 1\right)\]

Math

\[x - \sqrt{x \cdot x - \varepsilon} \]

↓

\[\begin{array}{l} t_0 := \sqrt{x \cdot x - \varepsilon}\\ \mathbf{if}\;\sqrt{-\varepsilon} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + t_0}\\ \mathbf{else}:\\ \;\;\;\;x - t_0\\ \end{array} \]

FPCore

(FPCore (x eps) :precision binary64 (- x (sqrt (- (* x x) eps))))

↓

(FPCore (x eps)
 :precision binary64
 (let* ((t_0 (sqrt (- (* x x) eps))))
   (if (!= (sqrt (- eps)) 0.0) (/ eps (+ x t_0)) (- x t_0))))

C

double code(double x, double eps) {
	return x - sqrt(((x * x) - eps));
}

↓

double code(double x, double eps) {
	double t_0 = sqrt(((x * x) - eps));
	double tmp;
	if (sqrt(-eps) != 0.0) {
		tmp = eps / (x + t_0);
	} else {
		tmp = x - t_0;
	}
	return tmp;
}

Fortran

real(8) function code(x, eps)
    real(8), intent (in) :: x
    real(8), intent (in) :: eps
    code = x - sqrt(((x * x) - eps))
end function

↓

real(8) function code(x, eps)
    real(8), intent (in) :: x
    real(8), intent (in) :: eps
    real(8) :: t_0
    real(8) :: tmp
    t_0 = sqrt(((x * x) - eps))
    if (sqrt(-eps) /= 0.0d0) then
        tmp = eps / (x + t_0)
    else
        tmp = x - t_0
    end if
    code = tmp
end function

Java

public static double code(double x, double eps) {
	return x - Math.sqrt(((x * x) - eps));
}

↓

public static double code(double x, double eps) {
	double t_0 = Math.sqrt(((x * x) - eps));
	double tmp;
	if (Math.sqrt(-eps) != 0.0) {
		tmp = eps / (x + t_0);
	} else {
		tmp = x - t_0;
	}
	return tmp;
}

Python

def code(x, eps):
	return x - math.sqrt(((x * x) - eps))

↓

def code(x, eps):
	t_0 = math.sqrt(((x * x) - eps))
	tmp = 0
	if math.sqrt(-eps) != 0.0:
		tmp = eps / (x + t_0)
	else:
		tmp = x - t_0
	return tmp

Julia

function code(x, eps)
	return Float64(x - sqrt(Float64(Float64(x * x) - eps)))
end

↓

function code(x, eps)
	t_0 = sqrt(Float64(Float64(x * x) - eps))
	tmp = 0.0
	if (sqrt(Float64(-eps)) != 0.0)
		tmp = Float64(eps / Float64(x + t_0));
	else
		tmp = Float64(x - t_0);
	end
	return tmp
end

MATLAB

function tmp = code(x, eps)
	tmp = x - sqrt(((x * x) - eps));
end

↓

function tmp_2 = code(x, eps)
	t_0 = sqrt(((x * x) - eps));
	tmp = 0.0;
	if (sqrt(-eps) ~= 0.0)
		tmp = eps / (x + t_0);
	else
		tmp = x - t_0;
	end
	tmp_2 = tmp;
end

Wolfram

code[x_, eps_] := N[(x - N[Sqrt[N[(N[(x * x), $MachinePrecision] - eps), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]

↓

code[x_, eps_] := Block[{t$95$0 = N[Sqrt[N[(N[(x * x), $MachinePrecision] - eps), $MachinePrecision]], $MachinePrecision]}, If[Unequal[N[Sqrt[(-eps)], $MachinePrecision], 0.0], N[(eps / N[(x + t$95$0), $MachinePrecision]), $MachinePrecision], N[(x - t$95$0), $MachinePrecision]]]

Target

Original	24.5
Target	0.3
Herbie	0.3

\[\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}} \]

Derivation

1. Initial program 24.5

   \[x - \sqrt{x \cdot x - \varepsilon} \]

2. Applied egg-rr 39.5

   \[\leadsto \color{blue}{\begin{array}{l} \color{blue}{\mathbf{if}\;x + \sqrt{x \cdot x - \varepsilon} \ne 0:\\ \;\;\;\;\frac{x \cdot x - \sqrt{\left(x \cdot x - \varepsilon\right) \cdot \left(x \cdot x - \varepsilon\right)}}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ } \end{array}} \]

3. Taylor expanded in x around 0 0.3

   \[\leadsto \begin{array}{l} \mathbf{if}\;x + \sqrt{x \cdot x - \varepsilon} \ne 0:\\ \;\;\;\;\frac{\color{blue}{\varepsilon}}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]

4. Taylor expanded in x around 0 0.3

   \[\leadsto \begin{array}{l} \mathbf{if}\;\color{blue}{\sqrt{-1} \cdot \sqrt{\varepsilon}} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]

5. Simplified 0.3

   \[\leadsto \begin{array}{l} \mathbf{if}\;\color{blue}{\sqrt{-\varepsilon}} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]
   Proof

   [Start] 0.3	\[ \begin{array}{l} \mathbf{if}\;\sqrt{-1} \cdot \sqrt{\varepsilon} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]
   exponential-simplify-20 [=>] 0.3	\[ \begin{array}{l} \mathbf{if}\;\color{blue}{\sqrt{-1 \cdot \varepsilon}} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]
   rational_best-simplify-3 [=>] 0.3	\[ \begin{array}{l} \mathbf{if}\;\sqrt{\color{blue}{\varepsilon \cdot -1}} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]
   rational_best-simplify-18 [<=] 0.3	\[ \begin{array}{l} \mathbf{if}\;\sqrt{\color{blue}{-\varepsilon}} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]

6. Final simplification 0.3

   \[\leadsto \begin{array}{l} \mathbf{if}\;\sqrt{-\varepsilon} \ne 0:\\ \;\;\;\;\frac{\varepsilon}{x + \sqrt{x \cdot x - \varepsilon}}\\ \mathbf{else}:\\ \;\;\;\;x - \sqrt{x \cdot x - \varepsilon}\\ \end{array} \]

Alternatives

Alternative 1
Error	1.0
Cost	20360

\[\begin{array}{l} t_0 := x - \sqrt{x \cdot x - \varepsilon}\\ \mathbf{if}\;t_0 \leq -1 \cdot 10^{-154}:\\ \;\;\;\;t_0\\ \mathbf{elif}\;\sqrt{-\varepsilon} \ne 0:\\ \;\;\;\;\frac{\varepsilon \cdot 0.5}{x - 0.5 \cdot \frac{\varepsilon}{x}}\\ \mathbf{else}:\\ \;\;\;\;t_0\\ \end{array} \]

Alternative 2
Error	1.0
Cost	13764

\[\begin{array}{l} t_0 := x - \sqrt{x \cdot x - \varepsilon}\\ \mathbf{if}\;t_0 \leq -1 \cdot 10^{-154}:\\ \;\;\;\;t_0\\ \mathbf{else}:\\ \;\;\;\;0.5 \cdot \frac{\varepsilon}{x}\\ \end{array} \]

Alternative 3
Error	9.0
Cost	6788

\[\begin{array}{l} \mathbf{if}\;x \leq 1.1 \cdot 10^{-104}:\\ \;\;\;\;x - \sqrt{-\varepsilon}\\ \mathbf{else}:\\ \;\;\;\;0.5 \cdot \frac{\varepsilon}{x}\\ \end{array} \]

Alternative 4
Error	9.1
Cost	6724

\[\begin{array}{l} \mathbf{if}\;x \leq 1.06 \cdot 10^{-104}:\\ \;\;\;\;-\sqrt{-\varepsilon}\\ \mathbf{else}:\\ \;\;\;\;0.5 \cdot \frac{\varepsilon}{x}\\ \end{array} \]

Alternative 5
Error	35.4
Cost	320

\[0.5 \cdot \frac{\varepsilon}{x} \]

Alternative 6
Error	61.2
Cost	192

\[x - x \]

Alternative 7
Error	61.7
Cost	64

\[x \]

Reproduce

herbie shell --seed 2023104 
(FPCore (x eps)
  :name "ENA, Section 1.4, Exercise 4d"
  :precision binary64
  :pre (and (and (<= 0.0 x) (<= x 1000000000.0)) (and (<= -1.0 eps) (<= eps 1.0)))

  :herbie-target
  (/ eps (+ x (sqrt (- (* x x) eps))))

  (- x (sqrt (- (* x x) eps))))