Short calculation of the dilogarithm?

Question

Is there a nice way to implement the dilogarithm function for real values, without actually performing the integration?

A series solution would have been nice, but the series around $0$ has a convergence radius of $1$, so it doesn't work for larger $x$. Ideally, I'm looking for an "elegant" method, rather than the "fastest" method.

For reference, the gsl library uses over 650 lines of code for the implementation, but I'm looking for somthing a bit more compact.

Answer 1

The cited Wikipedia page gives the expansion $$ \mathrm{Li}_2(x)=\frac{\pi^2}{3}-\frac12\log(x)^2-\sum_{k=1}^\infty\frac1{k^2x^k}-i\pi\log(x)\tag{1} $$ for $x\ge1$. Combined with $$ \mathrm{Li}_2(x)=\sum_{k=1}^\infty\frac{x^k}{k^2}\tag{2} $$ for $|x|\le1$, you should get what you need.

---

Equation $(1)$ also works for $x\le-1$ if we use $\log(x)=\log(-x)-\pi i$: $$ \mathrm{Li}_2(x)=-\frac{\pi^2}{6}-\frac12\log(-x)^2-\sum_{k=1}^\infty\frac1{k^2x^k}\tag{3} $$ for $x\le-1$.

---

Inversion Formula $$ \begin{align} \mathrm{Li}_2(x) &=-\int_0^x\log(1-t)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{6}-\int_1^x\log(1-t)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{6}-\pi i\log(x)-\int_1^x\log(t-1)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{6}-\pi i\log(x)-\int_{1/x}^1\log(1/t-1)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{6}-\pi i\log(x)-\int_{1/x}^1\Big(\log(1-t)-\log(t)\Big)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{6}-\pi i\log(x)+\frac{\pi^2}{6}+\int_0^{1/x}\log(1-t)\frac{\mathrm{d}t}{t}+\int_{1/x}^1\log(t)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{3}-\pi i\log(x)-\mathrm{Li}_2(1/x)-\frac12\log(x)^2\tag{4}\\ &=-\frac{\pi^2}{6}-\mathrm{Li}_2(1/x)-\frac12\log(-x)^2\tag{5} \end{align} $$ $(2)$ and $(4)$ prove expansion $(1)$. $(2)$ and $(5)$ prove expansion $(3)$.

---

Duplication Formula $$ \begin{align} \mathrm{Li}_2(x) &=-\int_0^x\log(1-t)\frac{\mathrm{d}t}{t}\\ &=-\int_0^1\log(1-t)\frac{\mathrm{d}t}{t}+\int_x^1\log(1-t)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{6}+\int_0^{1-x}\log(t)\frac{\mathrm{d}t}{1-t}\\ &=\frac{\pi^2}{6}-\log(x)\log(1-x)+\int_0^{1-x}\log(1-t)\frac{\mathrm{d}t}{t}\\ &=\frac{\pi^2}{6}-\log(x)\log(1-x)-\mathrm{Li}_2(1-x)\tag{6}\\ \mathrm{Li}_2(x)+\mathrm{Li}_2(1-x) &=\frac{\pi^2}{6}-\log(x)\log(1-x)\tag{7} \end{align} $$

Answer 2

I do not know how well this compares to other methods, but it is one I have devised for the purpose of computing the polylogarithm for values of z near the unit circle. In the case of the dilogarithm, the group of anharmonic ratios allows one to reduce the computation in the general case to a fundamental region for that group. However that still appears to involve computation for values near (or on) the unit circle, in particular the cube roots of -1 in the right half plane are fixed points for that group, so one needs to compute those. And nearby values.

A rational approximation to the polylogarithm can be had with poles in the branch cut (where they ought to be for any decent approximation). Runge's theorem guarantees there will be many.

We arrived at this one:

$$\text{Li}_{s}(z)\approx \frac{\pi z}{Γ(s)\sqrt{N}} \sum_{k=-N}^{N}\frac{e^{k\pi/\sqrt{N}}}{(e^{k\pi/\sqrt{N}}+1-z)(1+e^{k\pi/\sqrt{N}})} \log^{s-1}(1+e^{k\pi/\sqrt{N}})$$ $$= z\sum_{k=-N}^{N}\frac{w_{k}}{λ_{k}-z}$$

which is a Riemann sum for the integral representation

$$\text{Li}_{s}(z)=\frac{z}{Γ(s)}\int_{-\infty}^{\infty} \frac{e^{\tau}}{(e^{\tau}+1-z)(1+e^{\tau})} \log^{s-1}(1+e^{\tau})d\tau$$

which follows from a simple change of variables in the well known integral representation

$$\text{Li}_{s}(z)=\frac{1}{Γ(s)}\int_0^\infty \frac{ze^{-t}}{1-ze^{-t}}t^{s-1}dt$$

which can be obtained easily from the Laplace transform.

Answer 3

For further reference, based on @robjohn's answer, I came up with the following code for the calculation of the Spence function: $$\Phi(x)=-\int_0^x \frac{\ln|1-y|}{y}dy$$

double spenceFunction(double x, double precision = 1.0e-13) {
    if (x >= 0.5 && x < 1.0) {
        return  M_PI*M_PI/6.0 - spenceFunction(1.0 - x) - std::log( x)*std::log(1-x);
    }
    if (x >= 1.0) {
        return  M_PI*M_PI/3.0 - spenceFunction(1.0 / x) - std::log( x)*std::log(x) / 2.0;
    }
    if (x < -1.0) {
        return -M_PI*M_PI/6.0 - spenceFunction(1.0 / x) - std::log(-x)*std::log(-x) / 2.0;
    }

    double sum = 0.0, i = 1.0;
    for (; i < INT_MAX; i += 1.0) {
        double delta = std::pow(x,i) / (i*i);
        sum += delta;
        if (std::abs(delta/sum) < precision) break;
    }
    return sum;
}