Imprecision of complex li2 with Float64 and BigFloat for small real parts

[barton-willis]

Here is a case where a binary64 and BigFloat (with default precision) evaluation of li2 have a relative difference of about 8 x 10^-8. This relative difference is greater than what I expected.

julia> using PolyLog

julia> x = 4.831285545908206e-6 + 0.004396919500211628im
4.831285545908206e-6 + 0.004396919500211628im

julia> y0 = li2(x)
-1.941665612914504e-9 + 0.004396920676572405im

# convert x to a complex BigFloat, evaluate li2(x), and convert back to a complex binary64
julia> ye = convert(Complex{Float64}, li2(convert(Complex{BigFloat},x)))
-1.941665782029937e-9 + 0.004396920676572405im

julia> abs(real(y0)-real(ye))/abs(real(ye))
8.709811678577144e-8

I found this case by searching about 10 million points in the unit circle.

[Expander]

Thanks a lot for finding this deviation!

I've just cross-checked with Mathematica 13.3.0 to have another reference:

z = 4831285545908206/10^21 + 4396919500211628/10^18 I
N[PolyLog[2,z], 20]

As output I get:

-1.94166578202937687444628936853`13.795540844218891*^-9 + 
 0.00439692067657240512726530759719387623`20.150514997831955*I

The imaginary parts of all three numbers y0, ye and from Mathematica are the same (at 64 bit floating point precision).

The real parts differ:

reldiff(a,b) = abs(real(a) - real(b))/abs(real(b));

y0 = -1.941665612914504e-9 + 0.004396920676572405im; # Float64
ye = -1.941665782029937e-9 + 0.004396920676572405im; # BigFloat
mm = -1.94166578202937687444628936853e-9 + 0.00439692067657240512726530759719387623im; # Mma

reldiff(y0, ye) # 8.709811678577144e-8
reldiff(y0, mm) # 8.709782837299149e-8
reldiff(ye, mm) # 2.8841280506735895e-13

So, given Mathematica is correct to the last digit, it seams that at least li2(::ComplexF64) is quite imprecise for the point you've shown. In addition, I'd expect li2(::BigFloat) to be more precise as well (given the default BigFloat precision).

I'll have a look at where exactly this deviation comes from.

[Expander]

Just for our reference: The imprecision that was observed originated from expressions of the form clog(1 + z). For small z some digits are lost in such expressions. PR #24 fixes this by using log1p(z) in such cases, which is more precise for small z.

Anyway, thanks a lot for reporting this issue!