Fast switch and spline scheme for accurate inversion of nonlinear functions: The new first choice solution to Kepler’s equation - Tommasini, Daniele; Olivieri, David N. (Applied mathematics and computation, 2020)
•A universal algorithm for function inversion that does not require any initial guess.•A novel theoretical analysis of the errors that agrees with the numerical results.•Smaller errors than expected from general spline analysis, by a ∼10−22 factor for W ∈ [0; 10].•For Kepler problem, its python code is  ∼ 2000 times faster than Newton methods for large N.•Method faster than Newton’s even for N  ≥  2 for Kepler problem with tolerance 10−6rad. Numerically obtaining the inverse of a function is a common task for many scientific problems, often solved using a Newton iteration method. Here we describe an alternative scheme, based on switching variables followed by spline interpolation, which can be applied to monotonic functions under very general conditions. To optimize the algorithm, we designed a specific ultra-fast spline routine. We also derive analytically the theoretical errors of the method and test it on examples that are of interest in physics. In particular, we compute the real branch of Lambert’s W(y) function, which is defined as the inverse of xexp (x), and we solve Kepler’s equation. In all cases, our predictions for the theoretical errors are in excellent agreement with our numerical results, and are smaller than what could be expected from the general error analysis of spline interpolation by many orders of magnitude, namely by an astonishing 3×10−22 factor for the computation of W in the range W(y) ∈ [0, 10], and by a factor 2×10−4 for Kepler’s problem. In our tests, this scheme is much faster than Newton-Raphson method, by a factor in the range 10−4–10−3 for the execution time in the examples, when the values of the inverse function over an entire interval or for a large number of points are requested. For Kepler’s equation and tolerance 10−6 rad, the algorithm outperforms Newton’s method for all values of the number of points N ≥ 2.