calculation with qlib

1. `qlib` の概説

CPU に備わった最も高精度の浮動小数点数 (i386 で言えば long double) を越える精度で計算をしたい場合，より高精度な浮動小数点数を定義し，その演算をソフトウェアにてエミュレートするしか方策はありません． (浮動小数点演算のソフトウェア・エミュレートは，その機能のない CPU においては普通に行なわれていることです．) ただし，CPU に備わった浮動小数点数と違って FPU は使えないので，計算時間は激増します．

それを実現する library はいくつか配布されていまし，FreeBSD の ports になっているものもあります．私は拡張倍精度 (long double) に対応した算術関数の library (ports/math/ldouble) を通じてそれを配布している Cephes を知り，qlib の存在を知りました．要するに「芋蔓的な検索でたまたま見つけた」に過ぎず，いくつかの alternative を比較検討したわけではないので，「なぜ他でもなくこれを推奨するのか」の根拠はありません．ただ，私が使った限りでは問題なく動作していて，重宝しているので，ここに紹介致します．

なお，qlib にて新たに定義される浮動小数点数の型が qfloat であるのは c++ にて qfloat.h を利用する場合だけですが，ここではその呼称を qfloat で統一することにします．

qlib は現在でも bug fix は行なわれておりますが，少々不満や不便を感じている点もあります．すなわち，

演算に 64 bit の整数が使われていない．(16 or 32 bit です．)
64 bit 精度の long long (amd64 では long) と型変換する機能がない．

前者は assembler なので，私には歯が立ちません．後者は時間があれば作りたいです．

2. library の作成

2.1 `gcc` の確認

利用する gcc の path を確認します．例えば FreeBSD 標準の gcc は

% which gcc
/usr/bin/gcc

また FreeBSD 4 系で ports/lang/gcc30 (version 3 の gcc) を導入し，それを使う場合

% which gcc30
/usr/local/bin/gcc30

念のために version も確認します．

% /usr/bin/gcc -v

% /usr/local/bin/gcc30 -v

標準の gcc の version は FreeBSD の 4 系以前は 2，同 5, 6 系ならば 3，同 7 系以降ならば 4 のはずです．

2.2 配布ファイルの入手と展開

先ず Cephes から qlib.zip を download して下さい．それを次のように展開します．

% mkdir qlib          ← directory 名は自由に付けて下さい．
% cd qlib
% unzip ../qlib.zip   ← ports/archivers/unzip が必要です．

全部で 202 個のファイルが一つの directory に展開されます．それらを目的別に分類すれば次のようになります．

総括的な文書

qlib.doc : 配布ファイルの一覧 (必見です)

notes : 作者の連絡先
library に関係する source 等と文書

(1) header files, makefile, etc.

(2) 四則演算，型変換，平方根 (makefile では ARITH にて列挙される)

(3) 算術関数と定数 (makefile では LIBOBJS から ARITH を除いたもの)

(4) qlibdoc.html : 算術関数の使用法の一覧

なお，(1), (2) の中には，ここで対象としている環境 (i386 + gcc) では使わないものもある．
動作の点検や使用例の source 等と文書
機会があれば，これらの解説を書きます．
unix では使わないもの

*.rsp : Microsoft Linker script

*.zip : Borland C++ のためのファイル

2.3 ソースの修正

この節に記した作業は必須ではありません．修正，整理を必要としない方は，2.4 節に進んで下さい．

qlib 活用編に記したように，qlib には修正した方が良いと思えるところがいくつかあります．また 2.2 節で述べたように，ファイルが一つの directory に展開されてしまい，管理が面倒です．

そこで，修正のための patch とファイルを整理するための shell script を作りました．これをご利用下さい． → qlib.tgz
同封の文書 README に，2.4 節の代わりに行なう作業の手順が指示してあります．

2.4 library の make

i386 + gcc の環境にて日付が 2005-07-23 の qlib.zip で配布されたファイルを使う場合には，下記のような下準備なしにいきなり make を実行することができます． (この makefile は元々 Linux 用に書かれているので，Linux でも同様です．それ以外の環境においては下準備が必要となるかもしれません．それに関する情報は makefile の冒頭に書かれたコメントから得られます．)

library libmq.a は，上記のファイルを展開した directory にて次のような手順で作成します．

1. header file の選択

qhead32.h を qhead.h に cp する．

% diff qhead32.h qhead.h   ← 現在の配布ファイルではこれらは同一のはずだが，それを確認をする．
% cp qhead32.h qhead.h     ← 違いがなければ cp する必要はない．

2. makefile の修正 (gcc の指定)

CC= gcc

を変える必要があるなら修正する．例えば，

CC= /usr/bin/gcc

CC= /usr/local/bin/gcc30

3. makefile の修正 (CPU の指定等)

CPU の指定や最適化オプションの変更をしたい場合は，

CFLAGS= -g -O2 -Wall -fno-builtin

の行を変更する．例えば，

CFLAGS= -O3 -Wall -fno-builtin -march=cpu_type  ← "cpu_type" には具体的に CPU を規定する．

4. make

% make

添付された makefile により作成されるのは static library のみです．また make install のような機能はなく，「置き場」は自分で決めます． (他の directory に置きたいのならば，自分でそこに cp します．)

3. library の使用法

例えば自作のプログラム s1.c, s2.C は以下のようにして make します．

% gcc s1.c -o s1 -lmq -lm
% g++ s2.C -o s2 -lmq

qlib の関数の中には libm を使うものもあります．それを使う場合，option の -lm が必要になります． (g++ は標準で libm を link するので，-lm は不要です．)

自作のプログラムを library がある directory 以外のところで作るのであれば，libmq.a や使用する header file をそこに cp するか，または gcc に -I, -L option でそれらの所在を指示する必要があります． (header file は ../include に，libmq.a は ../lib にあるとして，)

% gcc -I../include s1.c -o s1 -L../lib -lmq -lm
% g++ -I../include s2.C -o s2 -L../lib -lmq

c 言語から使う場合，四則演算ですら関数の呼び出しで行なうので，演算を行なう部分がアセンブラみたいです．

#include <stdio.h>
#include "qhead.h"        /*  必ず include する．  */
#include "qprotos.h"      /*  必ず include する．  */

extern QELT qpi[];   /*  M_PI のような定数は変数として用意されている．  */
extern QELT qtwo[];  /*  定数 2.0 も変数として用意されている．          */

main(void)
{
   char c[128];
   QELT qfour[NQ],x[NQ],y[NQ];   /*  各々の配列が一つの qfloat の変数である．  */

   qadd(qtwo,qtwo,qfour);        /*  4 = 2 + 2 をつくる．  */
   qdiv(qfour,qpi,x);            /*  x = qpi/qfour         */
   qsin(x,y);                    /*  y = sin(x)            */

   qtoasc(y,c,30);               /*  qfloat を E 形式の文字列に変換する．  */
   printf("sin(pi/4) = %s\n",c);

   return 0;
}

c++ 言語から使う場合，演算子オーバーローディングの機能が使えるので，普通の c 言語のように書けます．

#include <stdio.h>
#include "qfloat.h"        // 必ず include する．

extern qfloat qpi;         // M_PI のような定数は変数として用意されている．

main(void)
{
   qfloat x,y;             // 型は qfloat である．

   x = 0.25*qpi;                        // 定数 0.25 の精度は double なので要注意．
   y = sin(x);                          // sin() は libmq.a の qsin() が使われる．

   printf("sin(pi/4) = %E\n",qtod(y));  // double に型変換して出力する．
   qprint(30,"sin(pi/4) =",y,"\n");     // qfloat として出力する．(E 形式)

   return 0;
}

c, c++ に備わった浮動小数点数との共存も可能です．

#include <stdio.h>
#include <math.h>
#include "qfloat.h"

extern qfloat qpi;

main(void)
{
   double x,y;
   qfloat qx,qy;

   x  = 0.25*M_PI;  y  = sin(x);   // libm  の sin()  が使われる．
   qx = 0.25*qpi;   qy = sin(qx);  // libmq の qsin() が使われる．

   return 0;
}

qlib の導入

1. `qlib` の概説

2. library の作成

2.1 `gcc` の確認

2.2 配布ファイルの入手と展開

2.3 ソースの修正

2.4 library の make

3. library の使用法

4. 精度の点検

qlib の導入

1. qlib の概説

2. library の作成

2.1 gcc の確認

2.2 配布ファイルの入手と展開

2.3 ソースの修正

2.4 library の make

3. library の使用法

4. 精度の点検

1. `qlib` の概説

2.1 `gcc` の確認