Why do these identical limits give different results?

Question

Bug introduced in 12.0.0.0. and persisting through 13.3.0.0.

---

I wanted to calculate this limit:

 Limit[E^(E^ProductLog[n] (-1 + E^(1/2 Sqrt[ProductLog[n]/n]))) (1 + 1/2 Sqrt[1/(n*ProductLog[n])])^-n, n -> Infinity]

The correct result is E^(1/8)

Furthermore, it is true that

 Limit[(E^(-(1/2) E^(ProductLog[n]/2)))/ (1 + 1/2 Sqrt[1/(n*ProductLog[n])])^-n, n -> Infinity]

is equal to 1

Therefore, it must be

 Limit[E^(E^ProductLog[n] (-1 + E^(Sqrt[ProductLog[n]]/(2 Sqrt[n]))))*(E^(-(1/2) E^(ProductLog[n]/2))) , n -> Infinity]

also equal to E^(1/8) But the result shocked me. Mathematica gives Infinity !!!

The graph confirms my original assumption

 Plot[E^(E^ProductLog[n] (-1 + E^(Sqrt[ProductLog[n]]/(2 Sqrt[n]))))*(E^(-(1/2) E^(ProductLog[n]/2))) , {n, 1, 1000000}]

On deeper analysis, I was able to better isolate the problem and have the following simpler example.

 Limit[-(1/2) Sqrt[n/ProductLog[n]] + ((-1 + E^(1/(2 Sqrt[n/ProductLog[n]]))) n)/ ProductLog[n], n -> Infinity]

again gives the wrong result of Infinity

The limit is equal to 1/8, which is supported by the following graph

 Plot[-(1/2) Sqrt[n/ProductLog[n]] + ((-1 + E^(1/(2 Sqrt[n/ProductLog[n]]))) n)/ ProductLog[n], {n, 1, 1000000}]

Answer 1

The easiest way, by a limited approach corresponding to approaching positive infinity:

expr = -(1/2) Sqrt[
     n/ProductLog[n]] + ((-1 + E^(1/(2 Sqrt[n/ProductLog[n]]))) n)/
    ProductLog[n];

Limit[expr /. n -> 1/n, n -> 0, Direction -> "FromAbove"]

(*  1/8  *)

Another workaround:

Limit[SeriesCoefficient[expr, {n, Infinity, 0}], n -> Infinity]

(*  1/8  *)

Possible cause of trouble:

The series at infinity seems divergent:

sd = Series[expr, {n, Infinity, 2}];
seriescoeffs = sd[[3]];
firstnumerator = sd[[4]];
denominator = sd[[6]];

kk = firstnumerator - 1;
(++kk; a[kk] (n^(kk/denominator)) -> 
    Limit[Simplify[# (n^(kk/denominator))], 
     n -> Infinity]) & /@ seriescoeffs
(*
{a[-1]/Sqrt[n] -> 0,
 a[0]          -> 1/8,
 Sqrt[n] a[1]  -> \[Infinity], 
 n a[2]        -> \[Infinity],
 n^(3/2) a[3]  -> \[Infinity], 
 n^2 a[4]      -> \[Infinity]}
*)

Answer 2

A partial workaround.

$Version

(* "13.1.0 for Mac OS X x86 (64-bit) (June 16, 2022)" *)

Clear["Global`*"]

f[n_] = -(1/2) Sqrt[
     n/ProductLog[n]] + ((-1 + E^(1/(2 Sqrt[n/ProductLog[n]]))) n)/
    ProductLog[n] // FullSimplify

(* -(1/2) Sqrt[E^ProductLog[n]] + 
 E^ProductLog[n] (-1 + E^(1/(2 Sqrt[E^ProductLog[n]]))) *)

Using Asymptotic

asymp = Asymptotic[f[n], n -> Infinity] // Simplify[#, n > 1] &;

Limit[asymp, n -> Infinity]

(* 1/8 *)

Or using Series

ser = Series[f[n], {n, Infinity, 1}] //
    Normal // Simplify[#, n > 1] &;

Limit[ser, n -> Infinity]

(* 1/8 *)