Inverse function of the Exponential Integral $\mathrm{Ei^{-1}}(x)$

Question

The Exponential integral is defined by

$$ \mathrm{Ei}(x) = \int_{-\infty}^x \frac{e^t}{t} \mathrm dt, $$ and has the following expansion

$$ \mathrm{Ei}(x) = \log x + \gamma + \sum_{k=1}^\infty \frac{x^k}{k \cdot k!}. $$ Someone recently made a contribution on the french Wikipedia page for it and he claimed the following expansion for its compositional inverse function: $$ \forall |x| < \frac{\mu}{\log \mu},\quad \mathrm{Ei}^{-1}(x) = \sum_{n=0}^\infty \frac{x^n}{n!} \frac{P_n(\log\mu)}{\mu^n}, $$ with $\mu$ being the Soldner constant ($\log\mu$ being the positive zero of $\mathrm{Ei}$) and $(P_n)$ is defined by

$$ P_0(x) = x, \quad P_{n+1}(x) = x(P_n'(x) - nP_n(x)). $$

So $P_0(x) = P_1(x) = x, P_2(x) = x - x^2, P_3(x) = 2x^3 - 4x^2 + x, \ldots$

This is quite an incredible result, but his source did not mention any of this and I couldn't find the result anywhere on the web. I told him about it and he responded that he had probably found it somewhere on a web page dedicated to maths, but couldn't anymore, so he undid his contribution.

I think this is quite a shame considering how incredible the formula looks, so I wanted to prove it to add it back. My first instinct was to use Lagrange's inversion formula at the point $\log\mu$:

$$ \mathrm{Ei}^{-1}(x) = \log\mu + \sum_{n=1}^\infty \frac{x^n}{n!} \left.\left(\frac{\mathrm d}{\mathrm dt}\right)^{n-1} \left(\frac{t-\log\mu}{\mathrm{Ei}(t)}\right)^n\right|_{t=\log\mu}. $$

With this expression, we immediatly see that the formula holds for the constant term, and we can compute the terms one by one with L'Hôpital's rule. However proving the general formula and making the polynomial sequence appear seems to be much more difficult.

The radius of convergence can easily be computed as for $n\geq1$, $\deg P_n = n$.

Evaluating the polynomial sequence at $1$ and $-1$ yields the sequences A089963 and A029768, which further supports the claim.

---

EDIT:

Numerical evidence shows that the formula given by Tyma Gaidash for $P_n$ is true:

$$ P_n(x) = \left.\left(\frac{\mathrm d}{\mathrm dt}\right)^{n-1} \left(\frac{te^x}{\mathrm{Ei}(t+x) - \mathrm{Ei}(x)}\right)^n\right|_{t=0}. $$ The current goal is to show that the above satisifies the recurrence relation, after which it is easy to show the rest.

Answer 1

This result can be obtained directly from a Maclaurin expansion of the function. By denoting \begin{equation} y=\mathrm{Ei}^{-1}(x) \end{equation} the integral representation of $\mathrm{Ei}(y)$ gives \begin{align} \frac{dy}{dx}&=\left(\frac{dx}{dy}\right)^{-1}\\ &=\left( \frac{d\mathrm{Ei}(y)}{dy} \right)^{-1}\\ &=ye^{-y} \end{align} Now, \begin{align} \frac{d^2y}{dx^2}&=\frac{d}{dx}\left( \frac{dy}{dx} \right)\\ &=\frac{dy}{dx}.\frac{d}{dy}\left( \frac{dy}{dx} \right)\\ &=ye^{-y}.\frac{d}{dy}\left(ye^{-y}\right)\\ &=(y-y^2)e^{-2y} \end{align} Assuming \begin{equation} \frac{d^ny}{dx^n}=P_n(y)e^{-ny} \end{equation} we have \begin{align} \frac{d^{n+1}y}{dx^{n+1}}&=\frac{dy}{dx}.\frac{d}{dy}\left( \frac{d^ny}{dx^n} \right)\\ &=ye^{-y}\left[P'_n(y)e^{-ny}-nP_n(y)e^{-ny}\right]\\ &=ye^{-(n+1)y}\left( P'_n(y)-nP_n(y)\right) \end{align} and thus this term has the same form as $\frac{d^ny}{dx^n}$ with \begin{equation} P_{n+1}(y)=y\left[ P'_n(y)-nP_n(y)\right] \end{equation}

For $x=0$, we have $y=\log\mu$, and then \begin{align} \left.\frac{d^ny}{dx^n}\right|_{x=0}&=\left.\frac{d^ny}{dx^n}\right|_{y=\log\mu}\\ &=\mu^{-n}P_n(\log\mu) \end{align} Maclaurin expansion writes \begin{equation} \mathrm{Ei}^{-1}(x)=\sum_{n=0}^\infty\frac{x^n}{n!}\mu^{-n}P_n(\log\mu) \end{equation} with the given recurrence relation and $P_0(y)=y$, as expected.