Converting between radix B and binary

radix_set_mpn converts a binary limb array to radix $B = b^e$ and radix_get_mpn converts back. For $B = 10^{19}$ the first is the expensive half of printing a number in decimal, and FLINT's flint_mpn_get_str (hence fmpz_get_str) uses it for base 10 when fft_small is available; once the number is in radix $10^{19}$, writing the digits takes linear time (radix_get_str_decimal).

Both directions use divide and conquer with multiplications only. To convert from radix $B$ to binary, the input is split as $a = a_{\mathrm{lo}} + a_{\mathrm{hi}} B^k$ with $k$ about half the length, both halves are converted recursively, and the result is $a_{\mathrm{lo}} + a_{\mathrm{hi}} \cdot B^k$, evaluated with binary arithmetic and a precomputed binary power $B^k$. To convert from binary to radix $B$, the binary input is read as a number in radix $c = 2^s$ with $c < B$ (cutting it into $s$-bit chunks, which is free), and the same recursion is run in reverse: the halves are converted recursively and combined with radix-$B$ arithmetic (radix_mul, which uses fft_small) and a precomputed radix-$B$ power $c^k$. No division by a large power is needed, unlike in the classical subquadratic conversion.

Sizes and thresholds#

• radix_get_mpn uses Horner's rule (radix_get_mpn_basecase) below 50 limbs, also at the leaves of the recursion.
• radix_set_mpn uses repeated division by $B$ (radix_set_mpn_basecase) below 40 input limbs, and the recursion switches to it on pieces of fewer than 25 chunks.
• The tables of powers ($c^{k_0}, c^{k_1}, \ldots$ with $k_{i+1} = \lceil k_i / 2 \rceil$) are built by repeated squaring, with one exact division by $c$ whenever $k_i$ is odd. Trailing zero limbs of the powers (the powers of $B = 10^{19} = 2^{19} 5^{19}$ have many zero bits) are skipped in the multiplications. radix_set_mpn caches the most recent table per thread.
• Above 500 limbs, the two recursive halves are converted in parallel when threads are available (see Threads).

Example#

Convert $3^{200}$ from binary to radix $10^{19}$, print it, and convert it back:

#include <flint/flint.h>
#include <flint/fmpz.h>
#include <flint/radix.h>

int main(void)
{
    radix_t radix;
    fmpz_t x;
    nn_ptr d, t;
    slong xn, dn, tn;
    char * s;

    radix_init(radix, 10, 0);   /* B = 10^19 on 64-bit machines */

    fmpz_init(x);
    fmpz_ui_pow_ui(x, 3, 200);
    xn = fmpz_size(x);

    d = flint_malloc(radix_set_mpn_need_alloc(xn, radix) * sizeof(ulong));
    dn = radix_set_mpn(d, COEFF_TO_PTR(*x)->_mp_d, xn, radix);

    s = radix_get_str_decimal(NULL, d, dn, 0, radix);
    flint_printf("%s\n", s);
    flint_printf("%wd digits, %wd decimal limbs\n",
        (slong) _radix_size_digits(d, dn, radix), dn);

    t = flint_malloc((xn + 1) * sizeof(ulong));
    tn = radix_get_mpn(t, d, dn, radix);
    flint_printf("round trip: %d\n",
        tn == xn && mpn_cmp(t, COEFF_TO_PTR(*x)->_mp_d, xn) == 0);

    flint_free(s);
    flint_free(d);
    flint_free(t);
    fmpz_clear(x);
    radix_clear(radix);
    return 0;
}

265613988875874769338781322035779626829233452653394495974574961739092490901302182994384699044001
96 digits, 6 decimal limbs
round trip: 1