Q. gccで "long double" を使用して数値計算を行ないましたが、同じ CPU での
Linux の結果と比べて計算精度が悪く、計算結果が異なります。
A. FreeBSD の標準カーネルでは浮動小数点の仮数部の精度が (IEEE754 と異なり)
53bit に設定されているために起きる現象です。
long double の仮数部の精度は64bitです。(double の仮数部の精度は 53bit です)
Intel x86 系の FPU の仮数部の精度は 24/53/64bit に設定可能です。
FreeBSD では仮数部の精度が 53bit に設定されているために、long double
を使用しても仮数部の精度が 53bit (つまり、double と同じ) しかありません。
精度を上げて数値計算を行なうには:
1. 関数 fpsetprec()/fpgetprec() を使用して、FPU の仮数部の精度を 64bit
に変更する事が出来ます。fpsetprec(3) を御覧下さい。
2. GNU Multiple Precision Arithmetic Library (GMP)、PARI 等 FreeBSD 標
準以外の数値演算ライブラリを利用してもよいでしょう。
1. fpsetprec()/fpgetprec() について
FreeBSD users MLの以下のメールから始まるスレッドが参考になるでしょう。
Subject: [FreeBSD-users-jp 24352] gcc bug
Subject: [FreeBSD-users-jp 46648] long double
その他、オンラインマニュアル fpsetround(3)、math(3)、ieee(3)などが参考
になります。
ただし、fpsetprec は double の内部演算も 64bit 精度にしてしまうという
副作用があります。double での演算結果が変わってしまう可能性がありますの
で注意して下さい。
long double 使用上の注意点:
* FreeBSD の printf(3) は 5.1-RELEASE から long double に対応しました。
それ以前 (少なくとも 4.8-RELEASE までと、5.0-RELEASE) は、書式で "%Lf"
等を使用してもエラーにはなりませんが、double の精度でしか表示されません。
* long double 用の算術関数 sqrtl(), sinl() 等が標準で用意されていません。
([管理番号 588] 参照)
プログラム例:
----------------------------
#ifdef __FreeBSD__
#include <ieeefp.h>
#endif
...
int main()
{
long double tmp1, tmp2;
...
#ifdef __FreeBSD__
fpsetprec(FP_PE); /* 仮数部の精度を 64bit にする。 */
#endif
...
}
----------------------------
2. 非標準の数値演算ライブラリについて
GMP (ports/math/libgmp4) は、任意多倍長の整数、整数を法とした計算
(modular arithmetic)、有理数、浮動小数点数の四則などが出来ます。但し、
三角関数、対数関数などの初等的な関数もないので、用途が限られるかも知れ
ません (自分で頑張って書けば、勿論問題ないですが)。
PARI (ports/math/pari) は、GMP 同様の数の他、p 進数や、それぞれの数
を係数とした多項式、羃級数、行列やベクトルの計算ができ、初等超越関数
(三角関数や指数対数関数) の他、これでもかというくらい様々な関数が実装
されています。本来は、整数論用のライブラリ集 (と、対話的に計算するため
の gp というコマンド (bc の様な感じ)) で、ちょっとした用途には大げさ
かも知れません。
いずれも C から使うことになるでしょう。(PARIは、Fortran や Pascal から
も使えるようです)。
他に参考になるものとしては
SAL: <URL:http://SAL.linet.gr.jp/index.shtml>
があります。
グループ名: FreeBSD-math