fredrikj.net / blog /

Arb 0.7: getting serious with power series

August 7, 2013

I’ve tagged version 0.7 of Arb, as the changelog was starting to grow bit long. Most importantly, Arb now supports power series much better (providing operations such as reversion and evaluation of various special functions, for power series over both the real and complex numbers).

Here’s a demo: the following program computes the terms in the Taylor expansion of $\tan(\log(1+\operatorname{atan}(\cos(\exp(\sqrt{1+x})))))$ at $x = 0$ to order $O(x^{500})$ using 1000-digit working precision.

#include "fmprb_poly.h"

int main()
{
    fmprb_poly_t a, b, c, d, e, f, g;
    long n, prec;

    n = 500;
    prec = 1000 * 3.323;

    fmprb_poly_init(a); fmprb_poly_init(b); fmprb_poly_init(c);
    fmprb_poly_init(d); fmprb_poly_init(e); fmprb_poly_init(f);
    fmprb_poly_init(g);

    fmprb_poly_zero(a);
    fmprb_poly_set_coeff_si(a, 0, 1);
    fmprb_poly_set_coeff_si(a, 1, 1);

    fmprb_poly_sqrt_series(b, a, n, prec);
    fmprb_poly_exp_series(c, b, n, prec);
    fmprb_poly_cos_series(d, c, n, prec);
    fmprb_poly_atan_series(e, d, n, prec);
    fmprb_add_ui(e->coeffs, e->coeffs, 1, prec);
    fmprb_poly_log_series(f, e, n, prec);
    fmprb_poly_tan_series(g, f, n, prec);

    fmprb_poly_printd(g, 10); printf("\n\n");

    fmprb_poly_clear(a); fmprb_poly_clear(b); fmprb_poly_clear(c);
    fmprb_poly_clear(d); fmprb_poly_clear(e); fmprb_poly_clear(f);
    fmprb_poly_clear(g);

    flint_cleanup();
    return 0;
}

This program takes 0.5 seconds to run (64-bit Xeon X5675 3.07 GHz), and the last coefficients are accurate to about 900 digits (setting the precision to about 1100 digits would give 1000 accurate digits while only slightly increasing the running time).

In Mathematica 9.0.1, evaluating

Timing[Series[Tan[Log[1+ArcTan[
    Cos[Exp[Sqrt[1`1000 + x]]]]]], {x, 0, 500 - 1}];]

on the same computer takes 752 seconds. Mathematica’s slower algorithm is more numerically stable, giving about 990 accurate digits (it does not guarantee that they are provably correct).

Also the gamma, log gamma and zeta functions of power series are supported, and are much faster than in any other software I’m aware of (with the added bonus of using rigorous error bounds!). For example, one can compute the Taylor series of $\Gamma(\exp(1+x))$ at $x = 0$ with:

    fmprb_poly_set_coeff_si(a, 0, 1);
    fmprb_poly_set_coeff_si(a, 1, 1);
    fmprb_poly_exp_series(b, a, n, prec);
    fmprb_poly_gamma_series(c, b, n, prec);

Computing 100 terms to 100 accurate digits takes 0.03 seconds, and computing 1000 terms to 1000 accurate digits takes 17 seconds. Meanwhile, Mathematica takes one minute to compute just 13 terms to 100 digit accuracy with the command Series[Gamma[Exp[1`100+x]], {x,0,12}] (the running time appears to increase exponentially with the number of terms).

Regarding the zeta function, the code for computing many terms in the Taylor series has been available for some time (I wrote about it in the earlier post about Arb 0.4). The news in 0.7 is that user-friendly functions are provided for the fmprb_poly and fmpcb_poly types, and they use the reflection formula in the negative half-plane. Computing no more than a few terms has also gotten much faster, as mentioned in the recent post about computing zeros of the zeta function.

Here is the overview of changes from 0.6 to 0.7, copied from the documentation:

I did not post about the release of version 0.6 to this blog, so here is the overview of changes from 0.5 to 0.6 as well:

On a final note, I’ve made available the slides from my ISSAC ’13 presentation about Arb (for which I received the best software demo award, winning an awesome gömböc courtesy of Maplesoft!).