Jmat.Real.lambertw, newton loop step

Average Error: 13.6 → 0.3

Time: 1.3m

Precision: 64

• could not determine a ground truth for program body

  1. wj = -7.594672547429651e+73
  2. x = -2.4734580902544025e+122

Math

\[wj - \frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}\]

↓

\[\begin{array}{l} \mathbf{if}\;wj \le 7.273522667070596 \cdot 10^{-05}:\\ \;\;\;\;\sqrt{(\left(wj \cdot wj\right) \cdot \left(wj \cdot wj - wj\right) + \left(wj \cdot wj\right))_*} \cdot \sqrt{(\left(wj \cdot wj\right) \cdot \left(wj \cdot wj - wj\right) + \left(wj \cdot wj\right))_*} + \frac{x}{e^{wj} \cdot wj + e^{wj}}\\ \mathbf{else}:\\ \;\;\;\;wj - \frac{wj - \frac{x}{e^{wj}}}{1 + wj}\\ \end{array}\]

C

double f(double wj, double x) {
        double r36325569 = wj;
        double r36325570 = exp(r36325569);
        double r36325571 = r36325569 * r36325570;
        double r36325572 = x;
        double r36325573 = r36325571 - r36325572;
        double r36325574 = r36325570 + r36325571;
        double r36325575 = r36325573 / r36325574;
        double r36325576 = r36325569 - r36325575;
        return r36325576;
}

↓

double f(double wj, double x) {
        double r36325577 = wj;
        double r36325578 = 7.273522667070596e-05;
        bool r36325579 = r36325577 <= r36325578;
        double r36325580 = r36325577 * r36325577;
        double r36325581 = r36325580 - r36325577;
        double r36325582 = fma(r36325580, r36325581, r36325580);
        double r36325583 = sqrt(r36325582);
        double r36325584 = r36325583 * r36325583;
        double r36325585 = x;
        double r36325586 = exp(r36325577);
        double r36325587 = r36325586 * r36325577;
        double r36325588 = r36325587 + r36325586;
        double r36325589 = r36325585 / r36325588;
        double r36325590 = r36325584 + r36325589;
        double r36325591 = r36325585 / r36325586;
        double r36325592 = r36325577 - r36325591;
        double r36325593 = 1.0;
        double r36325594 = r36325593 + r36325577;
        double r36325595 = r36325592 / r36325594;
        double r36325596 = r36325577 - r36325595;
        double r36325597 = r36325579 ? r36325590 : r36325596;
        return r36325597;
}

Error

Bits error versus wj

Bits error versus x

Target

Original	13.6
Target	13.1
Herbie	0.3

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

Derivation

1. Split input into 2 regimes

2. if wj < 7.273522667070596e-05

   1. Initial program 13.4

      \[wj - \frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}\]
   2. Using strategy rm

   3. Applied div-sub 13.4

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

   4. Applied associate--r- 6.8

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

   5. Taylor expanded around 0 0.3

      \[\leadsto \color{blue}{\left(\left({wj}^{2} + {wj}^{4}\right) - {wj}^{3}\right)} + \frac{x}{e^{wj} + wj \cdot e^{wj}}\]

   6. Simplified 0.3

      \[\leadsto \color{blue}{(\left(wj \cdot wj\right) \cdot \left(wj \cdot wj - wj\right) + \left(wj \cdot wj\right))_*} + \frac{x}{e^{wj} + wj \cdot e^{wj}}\]
   7. Using strategy rm

   8. Applied add-sqr-sqrt 0.3

      \[\leadsto \color{blue}{\sqrt{(\left(wj \cdot wj\right) \cdot \left(wj \cdot wj - wj\right) + \left(wj \cdot wj\right))_*} \cdot \sqrt{(\left(wj \cdot wj\right) \cdot \left(wj \cdot wj - wj\right) + \left(wj \cdot wj\right))_*}} + \frac{x}{e^{wj} + wj \cdot e^{wj}}\]

   if 7.273522667070596e-05 < wj

   1. Initial program 23.5

      \[wj - \frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}\]
   2. Using strategy rm

   3. Applied add-cube-cbrt 24.2

      \[\leadsto wj - \color{blue}{\left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right) \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}}\]

   4. Applied add-cube-cbrt 24.4

      \[\leadsto \color{blue}{\left(\sqrt[3]{wj} \cdot \sqrt[3]{wj}\right) \cdot \sqrt[3]{wj}} - \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right) \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\]

   5. Applied prod-diff 24.4

      \[\leadsto \color{blue}{(\left(\sqrt[3]{wj} \cdot \sqrt[3]{wj}\right) \cdot \left(\sqrt[3]{wj}\right) + \left(-\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right)\right))_* + (\left(-\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right) \cdot \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right) + \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right)\right))_*}\]

   6. Simplified 23.5

      \[\leadsto \color{blue}{\left(wj - \frac{wj - \frac{x}{e^{wj}}}{1 + wj}\right)} + (\left(-\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right) \cdot \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right) + \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \left(\sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}} \cdot \sqrt[3]{\frac{wj \cdot e^{wj} - x}{e^{wj} + wj \cdot e^{wj}}}\right)\right))_*\]

   7. Simplified 1.2

      \[\leadsto \left(wj - \frac{wj - \frac{x}{e^{wj}}}{1 + wj}\right) + \color{blue}{0}\]
3. Recombined 2 regimes into one program.

4. Final simplification 0.3

   \[\leadsto \begin{array}{l} \mathbf{if}\;wj \le 7.273522667070596 \cdot 10^{-05}:\\ \;\;\;\;\sqrt{(\left(wj \cdot wj\right) \cdot \left(wj \cdot wj - wj\right) + \left(wj \cdot wj\right))_*} \cdot \sqrt{(\left(wj \cdot wj\right) \cdot \left(wj \cdot wj - wj\right) + \left(wj \cdot wj\right))_*} + \frac{x}{e^{wj} \cdot wj + e^{wj}}\\ \mathbf{else}:\\ \;\;\;\;wj - \frac{wj - \frac{x}{e^{wj}}}{1 + wj}\\ \end{array}\]

Reproduce

herbie shell --seed 2019107 +o rules:numerics
(FPCore (wj x)
  :name "Jmat.Real.lambertw, newton loop step"

  :herbie-target
  (- wj (- (/ wj (+ wj 1)) (/ x (+ (exp wj) (* wj (exp wj))))))

  (- wj (/ (- (* wj (exp wj)) x) (+ (exp wj) (* wj (exp wj))))))