The $n$th prime in terms of Lambert-W

Question

Cipolla showed in his 1902 paper that there exists a sequence $(P_m), \ m\geq 1,$ of polynomials with integer coefficients such that, for any integer $m$,

$$p_n=n\bigg[\log n + \log\log n -1+\sum_{j=1}^{m}\frac{(-1)^{j-1}P_j (\log \log n)}{\log ^j n}+o\left(\frac{1}{\log^m n}\right)\bigg]$$

where $p_n$ is the $n$th prime, where the first few polynomials $P_j$ are:

$ \begin{align} P_2(X)=&\ (X^2-6X+11)\color{red}{/2!}\\ P_3(X)=&\ (2X^3 − 21X^2 + 84X − 131)\color{red}{/3!}\\ P_4(X)=&\ (6X^4 − 92X^3 + 588X^2 − 1908X + 2666)\color{red}{/4!}\\ \vdots& \\ \end{align} $

as reiterated in section $5$ in this paper $($though I believe the first equation should read $P^\prime_k=k(k-1)P_{k-1}\color{red}{\bf{-}} P^\prime_{k-1}k,\quad $ as opposed to $\color{red}{\bf{+}}.$ Also note that in the above paper, the polynomials don't include the $\color{red}{/j!}$ term $).$

---

Now, if $y=n/\log n, \ n=-y\ W_{-1}(-1/y)$, and if $y=n/(\log n-1), \ n=-y\ W_{-1}(-e/y).$ A heuristic argument might follow that since $n/(\log n-1)$ is a better approximation of $\pi (n)$ than $n/\log n$ ,$-n\ W_{-1}(-e/n)$ should be a better approximation of $p_n$ than $n \log n$.

---

Following on from my previous question, it appears that there exists a sequence $Q_m$ of rational coefficients such that, for any integer $m$,

$$p_n=n \left[-W_{-1}\left(-e/n\right)-\sum _{j=1}^m Q_j(W_{-1}\left(-e/n\right){}^{-j})+o\left(\frac{1}{\log^m n}\right)\right]$$

It turns out that the sum of the coefficients of the Cipolla polynomials divided by the corresponding factorial are precisely the coefficients $Q_j.$

eg

for $P_2(X)=\ X^2-6X+11,$ the corresponsing coefficient $Q_2=(1-6+11)/2!=3;$ and for $P_3(X)=\ 2X^3 − 21X^2 + 84X − 131,$ the corresponsing coefficient $Q_3=(2-21+84-131)/3!=-11.$

The first few terms then are $-1, 3, -11, 105/2, -613/2 \dots$

---

Calculating the coefficient $Q_j:$

Clear[a, b, c, q, r, k, i, j, n, x, cipolla, pnw]
fL[n_, wp_] := Quiet[x /. FindRoot[LogIntegral@x == n, {x, N[n Log[n], wp]}, WorkingPrecision -> wp]]; fL[n_] := fL[n, 200];
a[1, x_] := a[1,x]=x - 2;
a[k_, x_] :=a[k,x]= With[{t = Integrate[Expand[k (k - 1) a[k - 1, x] - k D[a[k - 1, x], x]], x]}, t + First[List @@ (Quiet[-(k - 1) (Sum[Binomial[(k - 1) - 1, j] a[j, x] (a[k - 1 - j, x] + D[a[k - 1 - j, x], x]), {j, 1, k - 2}] + (a[k - 1, x] + D[a[k - 1, x], x])) - k a[k - 1, x]]  //  Expand)] - (CoefficientList[t, x])[[2]]];
b[k_, x_] := a[k, x]/k!;
c[k_] := Tr[CoefficientList[b[k, x], x]];
cipolla[r_, n_] := n Log[n] + n Log[Log[n]] - n + n Sum[(-1)^(i - 1) b[i, q]/Log[n]^i, {i, 1, r}] /. q -> Log[Log[n]];cipolla[ n_] :=cipolla[Floor[Log[n]-2], n] -3E/2;
pnw[t_, n_] := n (Sum[-(c /@ Range@t)[[k]] ProductLog[-1, -E/n]^(-k), {k, 1, t}] - ProductLog[-1, -E/n]);pnw[ n_] :=pnw[Floor[-ProductLog[-1,-E/n]-2], n]-E ;
compare[n_]:= Quiet@Module[{a,b,c,d,e,f},n=N[n,200];a=cipolla[n];b=pnw[n];c=Floor[Log10@a];d=Floor[Log10@b];e=fL@n;f=Floor[Log10[e]];Column[{N[a,c],N[b,d],N[e,f]}]];
asymCoeff[x_]:=x !(x+2)/E;

where fL is a numerical FindRoot for $\operatorname{li}^{-1}(n),$ a calculates the Cipolla polynomials and c calculates the coefficients for the Lambert-W expansion. cipolla calculates $\operatorname{li}^{-1}(n)$ as given in the first display, and pnw calculates the full Lambert-W expansion.

eg c /@ Range@7 calculates $Q_1$ to $Q_7.$

As a side note, it seems that $\lfloor W_{-1}(-e/n) -2\rfloor$ terms $-e$ appears to be close to optimal for the Lambert-W expansion, and $\lfloor \log n -2\rfloor$ terms $-3e/2$ is close to optimal for the Cipolla expansion (see below). Truncated in this way, both approximations should $=\operatorname{li}^{-1}(n)+O (\log n).$

example implementations:

a[#,x]&/@Range@6//Column
Plot[{pnw@n,cipolla@n,fL@n},{n,5,100}]
Plot[{pnw@n-fL@n,cipolla@n-fL@n},{n,5,10^12}]

(*
-2+x
11-6 x+x^2
-131+84 x-21 x^2+2 x^3
2666-1908 x+588 x^2-92 x^3+6 x^4
-81534+62860 x-22020 x^2+4380 x^3-490 x^4+24 x^5
3478014-2823180 x+1075020 x^2-246480 x^3+35790 x^4-3084 x^5+120 x^6
*)

---

I can't explain why this expansion translates to a Lambert-W expansion in this way - why does this work? Is this peculiar to the expansion of $\operatorname{li}^{-1}(n)?$ Is there a simple proof for this, or am I being too optimistic?

---

Observations so far

The Lambert-W expansion is equivalent to the Cipolla expansion if $\log\log n\rightarrow 1$ and $\log n\rightarrow -W_{-1}-e/n.$

ie it can be rewritten as

$$p_n=n\bigg[\log n + \color{red}{1} -1+\sum_{j=1}^{m}\frac{(-1)^{j-1}P_j (\color{red}{1})}{\color{red}{-W_{-1}}^j\color{red}{-e/n}}+o\left(\frac{1}{\log^m n}\right)\bigg]$$

Also, $\log\log n=1$ at $n=e^e$ and $\log n= -W_{-1}-e/n$ at $n=e^e.$ I don't know whether this helps at all, but it might be relevant.