The algorithm

The central idea to numerically evaluate GPLs is to first map their parameters to the domain where the corresponding series representation is convergent (7) and to then use the series expansion up to some finite order. Thus, we will first look at how to remove trailing zeros in Section Removal of trailing zeros, and then how to make a GPL without trailing zeros convergent in Section Making GPLs convergent as presented in [7]. In Section Increase rate of convergence, we comment on accelerating the convergence of already convergent GPLs. Finally, in Section An example reduction we apply the algorithm to an explicit example.

Removal of trailing zeros

Consider a GPL of weight \(m\) with \(m\!-\!j\) trailing zeros

\[G(z_1,...,z_j,0_{m-j}\ ;y)\,.\]

We now shuffle \(\vec a = (z_1,...,z_j,0_{m-j-1})\) with \(\vec b = (0)\). This results in \(m\!-\!j\) times the original GPL as well as terms with less trailing zeros

(9)\[\begin{split} G(0\ ;y)\cdot G(z_1,...,z_j,0_{m-j-1}\ ;y) &= (m-j) G(z_1,...,z_j,0_{m-j}\ ;y) \\&\qquad + \sum_{\vec s} G(s_1,...,s_j,z_j,0_{m-j-1}\ ;y)\,,\end{split}\]

where the sum runs over all shuffle \(\vec s=(z_1,...,z_{j-1})\, \shuffle\, (0)\). We now solve (9) for \(G(z_1,...,z_j,0_{m-j};y)\) and obtain an expression with fewer trailing zeros. By applying this strategy recursively, we can remove all trailing zeros.

Making GPLs convergent

Reduction to pending integrals

Consider a GPL of the form

(10)\[ G(a_1,...,a_{i-1},s_r,a_{i+1},...,a_m\ ;y)\]

where \(s_r(=a_i)\) has the smallest absolute value among all the non-zero parameters in \(G\). If \(|s_r| < |y|\), (10) has no convergent series expansion. In order to remove the smallest weight \(s_r\), we apply the fundamental theorem of calculus to generate terms where \(s_r\) is either integrated over or not present anymore

\[\begin{split} G(a_1,...,a_{i-1},s_r,a_{i+1},...,a_m\ ;y) = G(a_1,...,a_{i-1}, 0 ,a_{i+1},...,a_m\ ;y) \\ + \int_0^{s_r} \D s_{r+1} \frac{\partial}{\partial s_{r+1}} G(a_1,...,a_{i-1},s_{r+1},a_{i+1},...,a_m\ ;y)\,.\end{split}\]

For the second term we use partial fraction decomposition and integration by parts. Then we obtain different results depending on where \(s_r\) is in the parameter list:

  • If \(s_r\) appears first in the list (i.e. \(i = 1\) and \(s_r=a_1\)) we find

    (11)\[\begin{split} G(s_r,a_{i+1},...,a_m\ ;y) = G(0 ,a_{i+1},...,a_m\ ;y) + \underbrace{\int_0^{s_r}\frac{\D s_{r+1}}{s_{r+1}-y }}_{G(y\ ;s_r)} G(a_{i+1},...,a_m\ ;y) \\ + \underbrace{\int_0^{s_r}\frac{\D s_{r+1}}{s_{r+1}-a_{i+1}} G(s_{r+1},a_{i+2},..,a_m\ ;y)}_{\text{pending integral}} - \underbrace{\int_0^{s_r}\frac{\D s_{r+1}}{s_{r+1}-a_{i+1}}}_{G(a_2\ ;s_r)} G(a_{i+1},...,a_m\ ;y) \,.\end{split}\]

    In the first term on the right-hand side, \(s_r\) is absent. Therefore the resulting GPL is simpler. It might still be non-convergent, but we can use this method recursively on the resulting GPLs until we end up with convergent GPLs.

    In the second and fourth terms the integration variable \(s_{r+1}\) does not appear in the parameters of the GPL, so that the integral can be solved (we write the solution as a GPL instead of a logarithm to be able to continue recursively).

    The third term does have the integration variable \(s_{r+1}\) among the weights and therefore yields what we refer to as a pending integral. This object can be written as a linear combination of simpler GPLs as we will see in Section Evaluation of pending integrals.

    Note that all GPLs on the right-hand side have depth reduced by one.

  • If \(s_r\) appears in the middle of the list, i.e. \(1 < i < m\), we find

    (12)\[\begin{split} G(a_1,...,a_{i-1},s_r, a_{i+1},...,a_m\ ;y) = & \\ + G(a_1,...,a_{i-1}, 0 ,&a_{i+1},...,a_m\ ;y) \\ - \int_0^{s_r} \frac{\D s_{r+1}}{s_{r+1}-a_{i-1}} & G(a_1,...,a_{i-2},s_{r+1},a_{i+1},...,a_m\ ;y) \\ + \underbrace{\int_0^{s_r} \frac{\D s_{r+1}}{s_{r+1}-a_{i-1}}}_{G(a_{i-1}\ ;s_r)} & G(a_1,...,a_{i-1}, a_{i+1},...,a_m\ ;y) \\ + \int_0^{s_r} \frac{\D s_{r+1}}{s_{r+1}-a_{i+1}} & G(a_1,...,a_{i-1},s_{r+1},a_{i+2},...,a_m\ ;y) \\ - \underbrace{\int_0^{s_r} \frac{\D s_{r+1}}{s_{r+1}-a_{i+1}}}_{G(a_{i+1}\ ;s_r)} & G(a_1,...,a_{i-1}, a_{i+1},...,a_m\ ;y) \,.\end{split}\]

    Again we obtain simpler GPLs (without \(s_r\) or lower depth) as well as pending integrals.

  • If \(s_r\) appears last in the list, i.e. \(i = m\), we use the shuffle algebra to remove \(s_r\) from the last place, just as we have done to remove trailing zeros.

We repeat these steps also for GPLs that are already under a pending integral.

Evaluation of pending integrals

The most general term created by the procedure of the last section is of the form

(13)\[\begin{split} \PI\Big(\vec p=(y',\vec b),i,\vec g = (\vec a,y)\Big) &\equiv \int_0^{y' } \frac{\D s_1}{s_1-b_1} \int_0^{s_1 } \frac{\D s_2}{s_2-b_2} \cdots \int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r} \\&\qquad G(a_1,...,a_{i-1},s_r,a_{i+1},...,a_m\ ;y) \, .\end{split}\]

Here we have adopted the convention that \(i=0\) implies that the integration variable does not appear inside the GPL. For example

\[\begin{split}\PI\Big(\vec p=(1,2,3),0,(4,5)\Big)&= \int_0^1 \frac{\D s_1}{s_1-2}\int_0^{s_1}\frac{\D s_2}{s_2-3} G(4;5) \, \\ \PI\Big(\vec p=(1,2,3),2,(4,5)\Big)&= \int_0^1 \frac{\D s_1}{s_1-2}\int_0^{s_1}\frac{\D s_2}{s_2-3} G(4,s_2;5)\,.\end{split}\]

As we use the algorithm, we need a way to collapse the pending integrals back down again. As an example, consider the case \(i=1\)

\[\begin{split} & \PI\Big(\vec p=(y',\vec b),1,\vec g=(\vec a,y)\Big) = \int_0^{y'} \frac{\D s_1}{s_1-b_1}\cdots\int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r} G(s_r,a_{i+1},...,a_m\ ;y) = \\ & \underbrace{\int_0^{y'} \frac{\D s_1}{s_1-b_1}\cdots\int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r}}_{\PI(\vec p,0,())} G( 0, a_{i+1},...,a_m\ ;y) \\ & + \underbrace{\int_0^{y'} \frac{\D s_1}{s_1-b_1}\cdots\int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r} \int_0^{s_r}\frac{\D s_{r+1}}{s_{r+1}-y}}_{\PI\Big((\vec p,y),0,()\Big)} G(a_{i+1},...,a_m\ ;y) \\ & + \underbrace{\int_0^{y'} \frac{\D s_1}{s_1-b_1}\cdots\int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r} \int_0^{s_r}\frac{\D s_{r+1}}{s_{r+1}-a_{i+1}} G(s_{r+1},a_{i+2},...,a_m\ ;y) }_{\PI\Big((\vec p,a_{i+1}),1,(a_{i+2},...,a_m;y)\Big)}\\ & - \underbrace{\int_0^{y'} \frac{\D s_1}{s_1-b_1}\cdots\int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r} \int_0^{s_r}\frac{\D s_{r+1}}{s_{r+1}-a_{i+1}}}_{\PI\Big((\vec p,a_{i+1}),0,()\Big)} G(a_{i+1},...,a_m\ ;y) \\&= \PI\Big(\vec p, 0, ()\Big) G(0,a_{i+1},...,a_m\ ; y) + \PI\Big((\vec p, y), 0, () \Big) G(a_{i+1}, ...., a_m\ ; y)\\& + \PI\Big((\vec p, a_{i+1}), 1, (a_{i+2}, ..., a_m\ ; y)\Big) - \PI\Big((\vec p, a_{i+1}), 0, ()\Big) G(a_{i+1}, ..., a_m\ ; y)\,.\end{split}\]

The other combinations follow similarly

\[\begin{split} \PI\Big( \vec p , i ,(\vec{a};y) \Big) &= + \PI\Big( \vec p , 0 ,() \Big) \, G(a_1,...,a_{i-1},0,a_{i+1},...,a_m\ ;y) \\& - \PI\Big((\vec p,a_{i-1}),i-1,(a_{i+1},...,a_m; y) \Big) \\& + \PI\Big((\vec p,a_{i-1}), 0 ,() \Big) \, G(a_1,...,a_{i-1}, a_{i+1},...,a_m\ ;y) \\& + \PI\Big((\vec p,a_{i+1}), i ,(a_1,...,a_{i-1},a_{i+2},...,a_m;y)\Big) \\& - \PI\Big((\vec p,a_{i+1}), 1 ,() \Big) \, G(a_1,...,a_{i-1}, a_{i+1},...,a_m\ ;y)\,.\end{split}\]

As we recursively apply the algorithm, we increase the number of pending integrals in front but decrease the depth of the \(G\)-functions by one unit in every recursion step. We do this until

  1. the only GPLs remaining under pending integrals are of depth one, i.e. \(G_m(s_r\; y)\),

  2. \(s_r\) is the argument, i.e. \(G(...\ ;s_r)\), or

  3. there are no GPLs under pending integrals.

We now discuss all these cases in turn:

  1. For GPLs of depth one, i.e. \(G_m(s_{r \pm}; y)\), we will be working with explicit logarithms. Hence, we need to indicate the infinitesimal imaginary part. We have to distinguish two cases: \(m=1\) and \(m>1\). For \(m = 1\) we have

    \[G_1(s_{r \pm};y) = G_1(y_{2 \mp};s_r) - G(0;\ s_r) + \log(-y)\,.\]

    Note that we will most likely have pending integrals in front, thus each term gives again a simpler pending integral

    \[ \PI\Big(\vec p=(y_\pm', \vec b), 1, (y) \Big) = G(\vec b, y_\mp\ ; y') - G(\vec b, 0\ ; y') + \log(-y_\mp) G\big(\vec b, y')\]

    The first and second terms have been reduced to case 2. and the third term to case 3.

    For \(m > 1\), we note

    (14)\[ G_m(s_{r \pm}\ ; y) = -\zeta(m) + \int_0^ {y} \frac{\D t}{t}G_{m-1}(t_{\pm}\ ; y) - \int_0^{s_r} \frac{\D t}{t}G_{m-1}(t_{\pm}\ ; y)\,.\]

    The second and third terms are now longer pending integrals, albeit with reduced weight

    \[\begin{split} \PI\Big( \vec{p} , m ,(0_{m-1},y)\Big) &= - \zeta(m)\PI\Big( \vec{p} , 0 ,( )\Big) \\&\qquad + \PI\Big(( y,0 ),m-1,(0_{m-2};y)\Big) \PI\Big( \vec{p} , 0 ,( )\Big) \\&\qquad - \PI\Big((\vec{p},0),m-1,(0_{m-2};y)\Big)\,.\end{split}\]
  2. In this case we end up simply with one large GPL

    \[\int_0^{y' } \frac{\D s_1}{s_1-b_1} \cdots \int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r} \, G(\vec a\ ;s_r) = G((\vec b,\vec a)\ ; y' )\,.\]

    In terms of pending integrals this is written as

    \[\PI\Big(\vec p=(y',\vec b), m+1, \vec g\Big) = G(\vec b, \vec g\ ; y')\,.\]
  3. If there is no GPL under the pending integral, the integral evaluates to a GPL

    \[\int_0^{y' } \frac{\D s_1}{s_1-b_1} \cdots \int_0^{s_{r-1}} \frac{\D s_r}{s_r-b_r} = G(b_1,...,b_r\ ;y' )\,.\]

In each case we end up with GPLs that are simpler in the sense that \(s_r\) has been eliminated. These might still be non-convergent due to other (non-zero) \(z_i\) elements being smaller in absolute value than \(y\). But applying the removal of \(s_r\) recursively we can eliminate all \(z_i\) for which \(|z_i| < |y|\). Therefore in the end we always obtain convergent GPLs.

Increase rate of convergence

Even though we have now only convergent GPLs, that does not imply that the convergence is fast enough for numerical applications. From now on we will only consider \(y=1\), as we can normalise any convergent GPL using (6). Convergence of such a GPL is slow if some \(z_i\) is close to the unit circle, i.e.

\[1 \le |z_i| \le \lambda < 2\,,\]

where \(\lambda\) is a parameter to be chosen.

Only for such \(z_i\) we apply the following strategy: to increase the rate of convergence we can use the fact that GPLs satisfy the Hölder convolution equation [1]

\[G(z_1,...,z_k\ ;1) = \sum_{j=0}^k (-1)^j G\Big(1-z_j,...,1-z_1\ ;1-\tfrac1p\Big) G\Big(z_{j+1},...,z_k\ ; \tfrac1p\Big)\,,\]

where \(p\) is an arbitrary non-zero complex number. Separating the first and the last term of this sum we obtain for \(p=2\) and again normalising the GPLs on the right-hand side

\[\begin{split}G(z_1,...,z_k\ ;1)& = G\big(2 z_1 ,...,2 z_k \ ; 1\big) + (-1)^k G\big(2(1-z_k),...,2(1-z_1)\ ; 1\big) \\& +\sum_{j=1}^{k-1}(-1)^j G\Big(2(1-z_j),...,2(1-z_1)\ ; 1\Big) G\Big(2 z_{j+1},...,2z_k \ ; 1\Big) \,.\end{split}\]

The first term has now better convergence as all parameters are twice as big. The GPL appearing in the sum all have reduced weight and are therefore not relevant for the present discussion.

The second term may or may not be convergent. If not, we repeat the algorithm outlined in Section Making GPLs convergent, including if necessary, Hölder convolution. At this stage it is not obvious why this recipe does indeed lead to a final answer and not to an infinite recursion. This can be shown by noting that the algorithm does only replace parameters with zero or permutes them; it does not introduce new non-trivial parameters. By carefully considering all possible behaviours under transformation \(z\mapsto 2(1-z)\), [7] proved that this method indeed works.

The choice of \(\lambda\) is a trade-off between accuracy and speed. A typical choice would be \(\lambda=1.1\) which is the default in handyG. \(\lambda\) can be changed using the hCircle option in set_options.

An example reduction

To illustrate the various aspects discussed so far, we include here an example of how the algorithm works in practice. For this purpose we reduce \(G(1,0,3;2)\) according to this algorithm until we end up with logarithms, polylogarithms and convergent MPLs. In our notation of a non-convergent GPL we have

(15)\[ G(\underbrace{1}_{s_r}, \underbrace{0}_{a_2}, \underbrace{3}_{a_3}\ ; \underbrace{2}_{y}) = G(0,0,3;2) + \int_0^1 \D s_1 \frac{\partial}{\partial s_1}G(s_1,0,3;2) \,.\]

The first term corresponds to \(G_3(3;2)\) and therefore it is a convergent trilogarithm. The second term has \(s_r\) appearing at the first place. Using (11) we obtain for the second term

(16)\[\begin{split} \int_0^1 \D s_1 \frac{\partial}{\partial s_1}G(s_1,0,3\ ;2) &= \int_0^1 \frac{\D s_1}{s_1-2}G(0,3\ ;2) + \int_0^1 \frac{\D s_1}{s_1-0}G(s_1,3\ ;2)\\& - \int_0^1 \frac{\D s_1}{s_1-0} G(0,3\ ;2)\,.\end{split}\]

The first and last terms are both conventional functions. Hence, we only need to worry about the second term which involves a pending integral. In order to evaluate it, we apply again (11) to the GPL under the pending integral to find

(17)\[\begin{split} G(s_1,3\ ;2) &= G(0, 3\ ;2) +\int_0^{s_1} \frac{\D s_2}{s_2-2}G( 3\ ;2) +\int_0^{s_1} \frac{\D s_2}{s_2-3}G(s_2\ ;2) \\& -\int_0^{s_1} \frac{\D s_2}{s_2-3}G( 3\ ;2) \, .\end{split}\]

Substituting this back into (16) gives

(18)\[\begin{split} \int_0^1 \frac{\D s_1}{s_1-0}&G(s_1,3\ ;2) =\int_0^1 \frac{\D s_1}{s_1 }G( 0,3\ ;2) +\int_0^1 \frac{\D s_1}{s_1 }\int_0^{s_1} \frac{\D s_2}{s_2-2}G( 3 \ ;2)\\& +\int_0^1 \frac{\D s_1}{s_1 }\int_0^{s_1} \frac{\D s_2}{s_2-3}G(s_2\ ;2) -\int_0^1 \frac{\D s_1}{s_1 }\int_0^{s_1} \frac{\D s_2}{s_2-3}G( 3 \ ;2 ) \,.\end{split}\]

Here only the third term is interesting, as the others are (poly)logarithms. The third term is a pending integral over a GPL of depth one. Thus,

(19)\[\begin{split} \int_0^1 \frac{\D s_1}{s_1} \int_0^{s_1}&\frac{\D s_2}{s_2-3} G(s_2\ ;2) \\&= \int_0^1 \frac{\D s_1}{s_1} \int_0^{s_1} \frac{\D s_2}{s_2-3} \Big(G(2\ ;s_2) - G(0\ ;s_2) + \log(-2) \Big)\end{split}\]

The first two terms have \(s_r\) as the argument and hence they are GPLs. The last term is independent of \(s_r\), making the integration trivial. Unfortunately, the second term \(G(0,3,0;1)\) has a trailing zero. To remove it, we shuffle \(G(0,3;1)\) with \(G(0;1)\) to find

(20)\[ G(0,3\ ;1)G(0\ ;1) = \sum_{\vec c = (0,3) \shuffle (0)} G(\vec c\ ;1) = G(0,3,0\ ;1) + 2\times G(0,0,3\ ;1)\,,\]

which we solve for \(G(0,3,0\ ;1)\).

Gathering all terms we obtain with \(G(0;\ 1)=\log1=0\)

\[\begin{split}G(1,0,3\ ;2) &= \underbrace{G(0,0,3\ ;2)}_{-\Li_3(2/3)} + \underbrace{G( 2 \ ;1)}_{\log(1/2)} \underbrace{G(0,3\ ;2)}_{-\Li_2(2/3)} - \cancel{ G( 0 \ ;1) G(0,3\ ;2)}\\& + \cancel{ G( 0 \ ;1) G(0,3\ ;2)} + \underbrace{G( 0,2\ ;1)}_{-\Li_2(1/3)} \underbrace{G(3\ ;2)}_{\log(1/3)} - \underbrace{G( 0,3\ ;1)}_{-\Li_2(2/3)} \underbrace{G(3\ ;2)}_{\log(1/3)} \\& + \underbrace{G(0,3,2\ ;1)}_{ \Li_{2,1}(1/3,3/2)} + \underbrace{G( 0,3\ ;1)}_{-\Li_2(1/3)} \log(-2) - \cancel{ G( 0 \ ;1) G( 0,3\ ;1)} \\& + 2 \underbrace{G(0,0,3;1)}_{-\Li_3(1/3)} = -0.81809 - 1.15049\I \,.\end{split}\]