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.

Answer 2

Here's another way to characterize the polynomials $P_n(x)$ in addition to the reps earlier presented here and even earlier in A055356.

Let

$$y=\mathrm{Ei}(x)=f(x)$$

and its inverse in some neighborhood of the $xy-$ plane be

$$x=\mathrm{Ei}^{<-1>}(y)= h(y) = f^{<-1>}(y).$$

The differential of

$$x =f^{<-1>}(f(x))= h(f(x))$$

is

$$dx = h'(f(x))f'(x)dx,$$

implying for $(x,y) = (f^{<-1>}(y),f(x))$, that

$$h'(y) = \frac{1}{f'(x)}=g(x) = g(h(y)),$$

the inverse function theorem, with the first string of equalities in Enta's answer as a specific example, and the autonomous differential equation for the inverse function $h=f^{<-1>}$ of $f$ (see OEIS A145271).

But also, as observed by Charles Graves in the 1850s, for sufficiently analytic functions $W(x)$,

$$e^{tg(x)\partial_x}W(x)=e^{tg(x)\frac{\partial}{\partial x}} W(x) = e^{t\frac{\partial}{\partial f(x)}}W(x)$$

$$ = e^{t\frac{\partial}{\partial y}} W(f^{<-1>}(y)) = W(f^{<-1>}(y+t))= W(f^{<-1>}(f(x)+t)).$$

Then for $f(x_0) = 0$ and $W(x) = x$,

$$e^{tg(x)\partial_x} x |_{x=x_0} = f^{<-1>}(t), $$

so

$$f^{<-1>}(t) = \sum_{n \geq 0} (g(x)\partial_x)^n x|_{x=x_0} \frac{t^n}{n!} = \sum_{n \geq 0}c_n \frac{t^n}{n!}.$$

For $f(x) = \mathrm{Ei}(x)$, we have $g(x) = x e^{-x}$, and

$$c_n = (x e^{-x}\partial_x)^nx|_{x=x_0} = (x e^{-x}\partial_x)^{n-1} xe^{-x}|_{x=x_0}.$$

Now

$$\partial_x^n g(x)= \partial_x^n x e^{-x} = \sum_{k=0}^n \binom{n}{k}(\partial_x^{k}x)(\partial_x^{n-k}e^{-x}) = x\partial_x^{n}e^{-x}+n\partial_x^{n-1}e^{-x}$$

$$ = ((-1)^nx+(-1)^{n-1} n)e^{-x}= (-1)^n (x-n)e^{-x},$$

and the general expressions of $R^n g(x)= (g(x)\partial_x)^n g(x)$ in OEIS A145271 for the refined Eulerian partition polynomials can be applied to determine the polynomials $P_m(x)$, which are also characterized in A055356. For example, with $(k')^j = (\partial_x^kg(x))^j$,

$P_1(x)e^{-x} = R^0 g(x) = 1 (0')^1 = x e^{-x},$

$P_2(x)e^{-2x}= R^1 g(x) = 1 (0')^1 (1')^1$

$ = x e^{-x} (1-x) e^{-x} =[x(1-x)] e^{-2x},$

$P_3(x)e^{-3x} = R^2 g(x) = 1 (0')^1 (1')^2 + 1 (0')^2 (2')^1$

$= [x(1-x)^2 + x^2 (x-2)] e^{-3x} = (2x^3 -4 x^2 + x) e^{-3x},$

$P_4(x)e^{-4x} =R^3 g(x) = 1 (0')^1 (1')^3 + 4 (0')^2 (1')^1 (2')^1 + 1 (0')^3 (3')^1$

$= [x (1-x)^3 + 4 x^2 (1-x) (x-2) + x^3 (3-x)] e^{-4x}$

$= (-6 x^4 + 18 x^3 - 11 x^2 + x) e^{-4x}.$