gmpで除算でmpfで

C GMP/MPFR

こういう検索語群もあったのだが、あれ？
GMPを使ったコードについて言及した記憶はあるけれど、
それは整数型mpz_tに関するものであって、
浮動小数点数型mpf_tを使ったコードについて書いたことあったかな？
と思って記事のローカルバックアップ、というかローカルオリジナルをgrepした。
見つかった。
見つかったけれど、それはGMPのmpf_t型でなくMPFRのmpfr_t型のことだった。
コイン投げの確率計算で最終的に確率を小数で表記するために使っている。
「GMPのmpf系関数ではなくMPFRを使用している」という一節があるので、
この唯一単独で使ったmpfという語に引っかかったのかもしれない。

GMPにおけるmpf_t型の除算はGMPのブランチであるMPFRのmpfr_t型の除算と変わらない。

#include <stddef.h>
#include <gmp.h>

int main(void)
{
    mpf_t a, b;
    mpf_init_set_ui(a, 1);
    mpf_init_set_ui(b, 29);
    mpf_div(a, a, b);
    gmp_printf("%.60Ff\n", a);
    gmp_printf("%lu\n", mpf_get_default_prec());
    mpf_clears(a, b, NULL);
    return 0;
}

これは、

...snip
    double a = 1., b = 29.;
    a /= b;
    printf("%d\n", a);
...snip

と同じような処理で、29の逆数を小数で表示しようとしている。
実行すると、

0.034482758620689655172400000000000000000000000000000000000000
64

出力の2行目のようにmpf_t型のデフォルトの精度が手元の環境では64ビットであるため、
29の逆数の表示は20桁程度（\log_{10}{2^{64}}\approx 19.3）で途切れてしまっている。
60桁程度が必要なら、\log_{2}{10^{60}}\approx 199.3なので、

...snip
    mpf_t a, b;
    mpf_set_default_prec(200);
    mpf_init_set_ui(a, 1);
...snip

のように、デフォルトで200ビットの精度を要求しておくと、

0.034482758620689655172413793103448275862068965517241379310345
224

28桁の循環節 0344827586206896551724137931 を持つ小数点以下60桁の小数表示を得る。
最後の桁の5は次の桁が8であるので4が繰り上がっている。
また、mpf_t型の新たなデフォルトの精度は要求した200ビットを満たす224ビットとなっている。