MPFRについて：任意精度数値計算ライブラリ

以前、１億桁のπの計算のためGMPを使いました。しかし、これには四則演算しか存在していなくて、通常の科学計算にはとても使えませんでした。
その後、調べていたらGMPを元に三角関数等をできるようにした、MPFR、MPCがありました。これら３つのライブラリは何とgccをビルドするのに必要なパッケージだとかで（何に使ってるんだ？）、最新のgccを使いたければインストールしておかないといけません。依存関係があり、GMP←MPFR←MPCという感じになるので、その順番にインストールしていかないといけません。大雑把に概要をいうと、以下のようになります。

・GMP
　任意精度算術演算ライブラリ

・MPFR
　こちらも任意精度の浮動小数点数演算ライブラリ（各種関数あり）

・MPC
　複素数演算ライブラリ

MPFRの簡単なサンプルを以下に示します。（これらのライブラリはCなのに注意）

/*                                                                                                     
mpfrの使い方を調べてみる                                                                        
 */

#include <stdio .h>
#include <stdlib .h>
#include <gmp .h>
#include <mpfr .h>

int main()
{
  mpfr_t work;

  mpfr_init2(work, 200);  // 有効桁数200bitで初期化                                                  
  mpfr_set_d(work, 1.0, MPFR_RNDD);  // double変数を代入                                              

  printf("work is ");
  mpfr_out_str(stdout, 10, 10, work, MPFR_RNDD);  // 出力　2個目はn進数での出力指定　                  
  // 3個目は何桁出力するか(size_t) '0'を指定すると当該変数の有効桁数すべて出力する                     
  putchar('\n');

  mpfr_const_pi(work, MPFR_RNDD);
  printf("pi is ");
  mpfr_out_str(stdout, 10, 0, work, MPFR_RNDD);
  putchar('\n');

  mpfr_clear(work);  // メモリの解放                                                                   

  return 0;
}

２回目の処理でπをセットしています。当然、任意の桁数を指定できますが、おそらく以前やったようなアルゴリズムで毎度必要な桁数を計算しているんでしょう。
ということは、桁数が大きいとそれに応じて処理時間がかかることになる？ちょっと計測してみました。（マシンはMacBookPro Retinaです。現在の我が家ではこいつが最速マシン）

①πの場合
以下の関数の実行前後でclockを取得し、関数の処理時間を調べてみました。

  gettimeofday( &tb, NULL );
  mpfr_const_pi(work, MPFR_RNDD);  // こいつの時間計測を行う
  gettimeofday( &ta, NULL );

変数workの精度を変えて計測してみました。（それぞれ数回計測し、5msec単位で丸めてます）

100,000bit　15msec
500,000bit　140msec
1,000,000bit　340msec
2,000,000bit　800msec

案の定、桁数と処理時間の関係はリニアではありませんね。

②三角関数の場合（atan）
三角関数（時間のかかりそうなatan）で計測してみました。

  gettimeofday( &tb, NULL );
  mpfr_atan(rslt, work, MPFR_RNDD);  // こいつの時間計測を行う
  gettimeofday( &ta, NULL );

100,000bit　15msec
500,000bit　140msec
1,000,000bit　340msec
2,000,000bit　795msec

う〜ん、ほとんどπと変わりませんね。たまたま演算量が同じ程度だったのかもしれません。
ただ何れにしても通常の倍精度浮動小数点の計算なんかとは比べ物にならないくらいの処理時間がかかっています。こういうライブラリがあるから安易に「100桁の精度で計算してよ」（１桁、約4bit弱必要なので400bitくらいですか）と依頼されても事前に注意しとかないと大変なことになります。
あとプログラム中では"mpfr_t"と変数を定義していますが、この実体はポインタで、必要な桁数のエリアをmallocしてそのポインタを保持しているだけです。従って、この変数自体はthread safeではありません。マルチスレッドのプログラムでこれを使おうって時には注意が必要です。