# Asymptotic inversion of ExpIntegralEi function

I'm looking at the small-x and large-x asymptotic expansions of the inverse of exponential integral $$E_1$$ (https://dlmf.nist.gov/6.2#E1)

$$\begin{array}{lll} E_1 & = & \int_z^\infty \frac{e^{-t}}{t}\mathrm{d}t\\ E_{1\text{small}}^{ - 1} (x) &\sim& - \log x - \log ( - \log x) - \frac{{\log ( - \log x) - 1}}{{\log x}} + \mathcal{O}\!\left( {\frac{{(\log ( - \log x))^2 }}{{\log ^2 x}}} \right)\\ E_{1\text{large}}^{ - 1} (x) &\sim& e^{ - x - \gamma } + e^{ - 2x - 2\gamma } + \frac{5}{4}e^{ - 3x - 3\gamma } + \cdots \end{array}$$

Gergo Nemes derived them by hand here:

I expressed −𝛾−𝐸1(𝑧) via dlmf.nist.gov/6.6.E2, took the exponential of both sides, expanded the exponential of the power series and finally employed series reversion to solve for 𝑧

I'm interested in improving these expansions by adding terms and generally learning how to automate the process. What is the most elegant way to do this in Mathematica?

E1[z_] := -ExpIntegralEi[-z]; (* https://dlmf.nist.gov/6.2#E1 *)

E1inv[y_] := x /. First@Flatten@Solve[{E1[x] == y, x > 0}, x];
E1invSmall[x_] := -Log[x] - Log[-Log[x]] - (Log[-Log[x]] - 1)/
Log[x];(* + O(Log[-Log[x]]^2/Log[x]^2)*)
E1invLarge[x_] :=
Exp[-x - EulerGamma] + Exp[-2 x - 2 EulerGamma] +
5/4 Exp[-3 x - 3 EulerGamma];

Plot[{E1inv[y], E1invSmall[y], E1invLarge[y]}, {y, 0, 2},
PlotLabel ->
"Approximating \!$$\*SubscriptBox[SuperscriptBox[\(E$$, $$-1$$], \
$$1$$]\)",
PlotLegends -> {"\!$$\*SubscriptBox[SuperscriptBox[\(E$$, $$-1$$], \
$$1$$]\)",
"\!$$\*SubscriptBox[SubscriptBox[SuperscriptBox[\(E$$, $$-1$$], \
$$1$$], $$small$$]\)",
"\!$$\*SubscriptBox[SubscriptBox[SuperscriptBox[\(E$$, $$-1$$], \
$$1$$], $$large$$]\)"}] edit July 21

Large x expansion:

E1[z_] := EulerGamma - ExpIntegralEi[-z];
series1 = Series[E1[z], {z, 0, 4}, Assumptions -> z > 0] // Normal;
result = InverseSeries@Series[Exp[series1], {z, 0, 4}] // Normal;
result /. (z -> Exp[z + EulerGamma]) // TraditionalForm


$$\frac{907}{240} e^{-6 z-6 \gamma }+\frac{361}{144} e^{-5 z-5 \gamma }+\frac{31}{18} e^{-4 z-4 \gamma }+\frac{5}{4} e^{-3 z-3 \gamma }+e^{-2 z-2 \gamma }+e^{-z-\gamma }$$

• On the obvious track E1[z_] = -ExpIntegralEi[-z]; E1inv = InverseFunction[E1]; Series[E1inv[x], {x, 0, 3}] I get a result that is definitely wrong (hitting a wrong branch?). Jul 16, 2021 at 20:06
• Also plotting InverseFunction[E1] gives a weird discontinuous plot Jul 16, 2021 at 21:01
• May be variable substitution of the original integral definition of E1 helps a litte:........................................ t ->Log[x]  yields Integrate[1/(x Log[x]) 1/x, {x, E^z, Infinity}, Assumptions -> z > 0]  and then x- -> Log[r] yields Integrate[1/(Log[r]^2 Log[Log[r]]) 1/r, {r, E^E^z, Infinity}, Assumptions -> z > 0]  both forms are equal to -ExpIntegralEi[-z]  Jul 19, 2021 at 6:50
• Your function fulfills f'[x] == f[x] Exp[-f[x]]. Maybe it would be helpful somehow Jul 20, 2021 at 8:31
• The recipe is known (I added in the edit), I was more interested in how to do the algebraic steps in Mathematica, how to 1) expand exponential of the power series and 2) use series reversion to solve for z Jul 20, 2021 at 14:59

You got correct hints on Math.SE. Define:

E1[z_] := EulerGamma - ExpIntegralEi[-z];
E2[z_] := -ExpIntegralEi[-z] Exp[z];


They can be inverted like this

Series[E1[z], {z, 0, 4}, Assumptions -> z > 0]
(*-Log[z]+z-z^2/4+z^3/18-z^4/96+O[z]^5*)
InverseSeries[Series[Exp[z - z^2/4 + z^3/18 - z^4/96], {z, 0, 4}]/z]
(*1/z+(1/z)^2+5/(4 z^3)+31/(18 z^4)+361/(144 z^5)+O[1/z]^7*)


Some plot for the confirmation

ParametricPlot[{{E1[z], z}, {z, E1[z]},
{z, Log[Exp[z - z^2/4 + z^3/18 - z^4/96]/z]},
{Log[z],1/z + (1/z)^2 + 5/(4 z^3) + 31/(18 z^4) + 361/(144 z^5)},
{z, (361 E^(-5 z))/144 + (31 E^(-4 z))/18 + (5 E^(-3 z))/4 + E^(-2 z) +
E^-z}},
{z, 0.01, 5}, AspectRatio -> 1,
PlotRange -> {{0, 5}, {0, 5}}]


Now let us focus on the other limit. E2 has a nice series expansion at infinity

s2 = Series[E2[x], {x, Infinity, 1}]


$$\frac{1}{x}+O\left(\frac{1}{x^2}\right)$$

The dominant decay of E1 at large values of z comes therefore from the exponent---seek the solution in the form

$$z= -\log(y)+u(y).$$

Substitute this in the defining series expansion of E1. We get an equation of the form

$$e^u=\text{rational function of u}.$$

For instance, to the leading order we have

$$e^{u_0}=\frac{1}{\log(z)}.$$

Determine $$u_0$$, substitute back, determine the next term $$u_1$$, iterate.

• Hm, I wonder why that expansion looks so different from the expression in math.SE post Jul 21, 2021 at 17:45
• @YaroslavBulatov As far as I can see it looks the same (361 E^(-5 z))/144 + (31 E^(-4 z))/18 + (5 E^(-3 z))/4 + E^(-2 z) + E^-z and even more terms... I skipped some trivial substitutions :) Jul 21, 2021 at 18:02
• thanks, it works for Einvlarge....what do you mean by treated similarly for E2? InverseSeries command returns unevaluated if I just follow through the same steps Jul 21, 2021 at 20:31