2isqrt (example 3.6)

Average Error: 19.5 → 0.4

Time: 9.0s

Precision: binary64

Cost: 20288

Math

\[\frac{1}{\sqrt{x}} - \frac{1}{\sqrt{x + 1}} \]

↓

\[{\left(x \cdot \left({\left(1 + x\right)}^{-0.5} + \frac{1}{\sqrt{x}}\right)\right)}^{-1} \cdot \frac{1}{1 + x} \]

FPCore

(FPCore (x) :precision binary64 (- (/ 1.0 (sqrt x)) (/ 1.0 (sqrt (+ x 1.0)))))

↓

(FPCore (x)
 :precision binary64
 (*
  (pow (* x (+ (pow (+ 1.0 x) -0.5) (/ 1.0 (sqrt x)))) -1.0)
  (/ 1.0 (+ 1.0 x))))

C

double code(double x) {
	return (1.0 / sqrt(x)) - (1.0 / sqrt((x + 1.0)));
}

↓

double code(double x) {
	return pow((x * (pow((1.0 + x), -0.5) + (1.0 / sqrt(x)))), -1.0) * (1.0 / (1.0 + x));
}

Fortran

real(8) function code(x)
    real(8), intent (in) :: x
    code = (1.0d0 / sqrt(x)) - (1.0d0 / sqrt((x + 1.0d0)))
end function

↓

real(8) function code(x)
    real(8), intent (in) :: x
    code = ((x * (((1.0d0 + x) ** (-0.5d0)) + (1.0d0 / sqrt(x)))) ** (-1.0d0)) * (1.0d0 / (1.0d0 + x))
end function

Java

public static double code(double x) {
	return (1.0 / Math.sqrt(x)) - (1.0 / Math.sqrt((x + 1.0)));
}

↓

public static double code(double x) {
	return Math.pow((x * (Math.pow((1.0 + x), -0.5) + (1.0 / Math.sqrt(x)))), -1.0) * (1.0 / (1.0 + x));
}

Python

def code(x):
	return (1.0 / math.sqrt(x)) - (1.0 / math.sqrt((x + 1.0)))

↓

def code(x):
	return math.pow((x * (math.pow((1.0 + x), -0.5) + (1.0 / math.sqrt(x)))), -1.0) * (1.0 / (1.0 + x))

Julia

function code(x)
	return Float64(Float64(1.0 / sqrt(x)) - Float64(1.0 / sqrt(Float64(x + 1.0))))
end

↓

function code(x)
	return Float64((Float64(x * Float64((Float64(1.0 + x) ^ -0.5) + Float64(1.0 / sqrt(x)))) ^ -1.0) * Float64(1.0 / Float64(1.0 + x)))
end

MATLAB

function tmp = code(x)
	tmp = (1.0 / sqrt(x)) - (1.0 / sqrt((x + 1.0)));
end

↓

function tmp = code(x)
	tmp = ((x * (((1.0 + x) ^ -0.5) + (1.0 / sqrt(x)))) ^ -1.0) * (1.0 / (1.0 + x));
end

Wolfram

code[x_] := N[(N[(1.0 / N[Sqrt[x], $MachinePrecision]), $MachinePrecision] - N[(1.0 / N[Sqrt[N[(x + 1.0), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]

↓

code[x_] := N[(N[Power[N[(x * N[(N[Power[N[(1.0 + x), $MachinePrecision], -0.5], $MachinePrecision] + N[(1.0 / N[Sqrt[x], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], -1.0], $MachinePrecision] * N[(1.0 / N[(1.0 + x), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]

Target

Original	19.5
Target	0.6
Herbie	0.4

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

Derivation

1. Initial program 19.5

   \[\frac{1}{\sqrt{x}} - \frac{1}{\sqrt{x + 1}} \]

2. Applied egg-rr 19.6

   \[\leadsto \color{blue}{\frac{1}{\frac{{x}^{-0.5} + {\left(1 + x\right)}^{-0.5}}{\frac{1}{x} - \frac{1}{1 + x}}}} \]

3. Applied egg-rr 5.0

   \[\leadsto \frac{1}{\frac{{x}^{-0.5} + {\left(1 + x\right)}^{-0.5}}{\color{blue}{\frac{\frac{1 + \left(x - x\right)}{x}}{1 + x}}}} \]

4. Applied egg-rr 0.4

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

5. Applied egg-rr 0.5

   \[\leadsto {\left(\color{blue}{\mathsf{fma}\left({\left(\sqrt{x}\right)}^{-0.5}, {\left(\sqrt{x}\right)}^{-0.5}, {\left(x + 1\right)}^{-0.5}\right)} \cdot x\right)}^{-1} \cdot \frac{1}{x + 1} \]

6. Simplified 0.4

   \[\leadsto {\left(\color{blue}{\left({\left(1 + x\right)}^{-0.5} + \frac{1}{\sqrt{x}}\right)} \cdot x\right)}^{-1} \cdot \frac{1}{x + 1} \]
   Proof
   (+.f64 (pow.f64 (+.f64 1 x) -1/2) (/.f64 1 (sqrt.f64 x))): 0 points increase in error, 0 points decrease in error
   (+.f64 (pow.f64 (Rewrite<= +-commutative_binary64 (+.f64 x 1)) -1/2) (/.f64 1 (sqrt.f64 x))): 0 points increase in error, 0 points decrease in error
   (+.f64 (pow.f64 (+.f64 x 1) -1/2) (Rewrite<= unpow-1_binary64 (pow.f64 (sqrt.f64 x) -1))): 0 points increase in error, 0 points decrease in error
   (+.f64 (pow.f64 (+.f64 x 1) -1/2) (pow.f64 (sqrt.f64 x) (Rewrite<= metadata-eval (*.f64 2 -1/2)))): 0 points increase in error, 0 points decrease in error
   (+.f64 (pow.f64 (+.f64 x 1) -1/2) (Rewrite<= pow-sqr_binary64 (*.f64 (pow.f64 (sqrt.f64 x) -1/2) (pow.f64 (sqrt.f64 x) -1/2)))): 66 points increase in error, 30 points decrease in error
   (Rewrite<= +-commutative_binary64 (+.f64 (*.f64 (pow.f64 (sqrt.f64 x) -1/2) (pow.f64 (sqrt.f64 x) -1/2)) (pow.f64 (+.f64 x 1) -1/2))): 0 points increase in error, 0 points decrease in error
   (Rewrite<= fma-udef_binary64 (fma.f64 (pow.f64 (sqrt.f64 x) -1/2) (pow.f64 (sqrt.f64 x) -1/2) (pow.f64 (+.f64 x 1) -1/2))): 19 points increase in error, 21 points decrease in error

7. Final simplification 0.4

   \[\leadsto {\left(x \cdot \left({\left(1 + x\right)}^{-0.5} + \frac{1}{\sqrt{x}}\right)\right)}^{-1} \cdot \frac{1}{1 + x} \]

Alternatives

Alternative 1
Error	0.4
Cost	26884

\[\begin{array}{l} \mathbf{if}\;\frac{1}{\sqrt{x}} + \frac{-1}{\sqrt{1 + x}} \leq 10^{-12}:\\ \;\;\;\;\frac{1}{1 + x} \cdot {\left(\sqrt{x} \cdot 2\right)}^{-1}\\ \mathbf{else}:\\ \;\;\;\;{x}^{-0.5} - {\left(1 + x\right)}^{-0.5}\\ \end{array} \]

Alternative 2
Error	0.4
Cost	13760

\[\frac{\frac{1}{x \cdot \left({\left(1 + x\right)}^{-0.5} + {x}^{-0.5}\right)}}{1 + x} \]

Alternative 3
Error	0.6
Cost	13696

\[\frac{-1}{\left(\sqrt{x} + x \cdot {\left(1 + x\right)}^{-0.5}\right) \cdot \left(-1 - x\right)} \]

Alternative 4
Error	0.4
Cost	13380

\[\begin{array}{l} \mathbf{if}\;x \leq 37000000:\\ \;\;\;\;{x}^{-0.5} - {\left(1 + x\right)}^{-0.5}\\ \mathbf{else}:\\ \;\;\;\;\frac{1}{1 + x} \cdot \left(\sqrt{\frac{1}{x}} \cdot 0.5\right)\\ \end{array} \]

Alternative 5
Error	1.1
Cost	7236

\[\begin{array}{l} \mathbf{if}\;x \leq 0.66:\\ \;\;\;\;-1 + \left({x}^{-0.5} + x \cdot 0.5\right)\\ \mathbf{else}:\\ \;\;\;\;\frac{1}{1 + x} \cdot \left(\sqrt{\frac{1}{x}} \cdot 0.5\right)\\ \end{array} \]

Alternative 6
Error	1.1
Cost	7044

\[\begin{array}{l} \mathbf{if}\;x \leq 1:\\ \;\;\;\;-1 + \left({x}^{-0.5} + x \cdot 0.5\right)\\ \mathbf{else}:\\ \;\;\;\;0.5 \cdot {\left(\frac{1}{x}\right)}^{1.5}\\ \end{array} \]

Alternative 7
Error	1.3
Cost	6916

\[\begin{array}{l} \mathbf{if}\;x \leq 0.66:\\ \;\;\;\;-1 + {x}^{-0.5}\\ \mathbf{else}:\\ \;\;\;\;0.5 \cdot {\left(\frac{1}{x}\right)}^{1.5}\\ \end{array} \]

Alternative 8
Error	1.6
Cost	6788

\[\begin{array}{l} \mathbf{if}\;x \leq 0.66:\\ \;\;\;\;-1 + {x}^{-0.5}\\ \mathbf{else}:\\ \;\;\;\;\frac{0.5}{{x}^{1.5}}\\ \end{array} \]

Alternative 9
Error	31.4
Cost	6528

\[{x}^{-0.5} \]

Alternative 10
Error	62.8
Cost	64

\[-1 \]

Reproduce

herbie shell --seed 2022329 
(FPCore (x)
  :name "2isqrt (example 3.6)"
  :precision binary64

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

  (- (/ 1.0 (sqrt x)) (/ 1.0 (sqrt (+ x 1.0)))))