float and double

切れるOS
検算
インタープリタ対コンパイラ対 …
mine.so 観光
many scheme
float to double
- マシン固有のFPU
- ちょっとびっくり

切れるOS

キレるOSじゃなくて、コネクションが切れてしまうOSの事ね。

古いThinkpadを Debian + FreeBSD のディアルブートにして、１月半程経過した。朝にdebianを使えば、午後はFreeBSDをブートって具合に、ほぼ平等に使ってる。

で、Debianの方なんだけど、時々セッションが重くなり、そのうちにWiFiの接続が切れてしまうんだ。そんな時は、再びsshで接続しなおして、tmux a とかして難を逃れている。

なんでや？重くなるって事を考えると、通信ロスが頻発して、再送を繰り返しているんではなかろうか？何とか証拠を掴みたいものだ。

その点、FreeBSDの方は極めて安定。まだ一度もそのような現象にみまわれた事は無い。ただ一つ難があるとすれば、起動時にWiFiのセッションが確立せず、30秒ぐらいの時間を於いての再接続で繋る事。一度電波を掴まえれば後は問題無いんで、安心していられる。

やっぱり、安定、安心のFreeBSD だな。

あれ？ OpenBSDはどうした？

obusb$ doas syspatch
doas (sakae@obusb.hoge.jp) password:
Get/Verify syspatch69-016_sshd.tgz 100% |***************|  1938 KB    00:01
Installing patch 016_sshd
Get/Verify syspatch69-017_x509.tgz 100% |***************| 15582 KB    00:05
Installing patch 017_x509
Errata can be reviewed under /var/syspatch

おお、patchが降ってきた。そろそろ次の版が出る季節なんで、先行のバックポート・パッチの提供なんだな。

x509の方は兎も角、sshdが気になる。

OpenBSD 6.9 errata 016, Sep 27, 2021:

sshd(8) failed to clear supplemental groups when executing an
AuthorizedUsersCommand or AuthorizedPrincipalsCommand helper program.

Apply by doing:
    signify -Vep /etc/signify/openbsd-69-base.pub -x 016_sshd.patch.sig \
        -m - | (cd /usr/src && patch -p0)

And then rebuild and install sshd(8)
    cd /usr/src/usr.bin/ssh
    make obj
    make
    make install

WiFiも安定しているし、これぞ本当に安心して使えるOSだ、な。

ああ、忘れていた事がある。このdual-bootなマシン、言わば敷地一杯までののよくばり２世帯住宅なんだけど、世帯間を結ぶ通路が無いんだ(完全な設計ミスって認めるよ)。お届け物は、Windowsを介してしか出来無い。すこぶる不便。

なんとかしたいんで、リナのfs(5)を見たんだけど、そんなの知らないって態度。大人のFreeBSD側からなら、何とかなるかな。

検算

前回作った、LC回路の同調周波数を計算するやつ。それぞれのOSで実施

sakae@pen:~$ gosh        ;; Debian at amd64
gosh> (use mine)
gosh> (resf 1e-6 500e-12)
7117625.5

[sakae@fb ~]$ gosh        ;; FreeBSD at i386
gosh> (use mine)
gosh> (resf 1e-6 500e-12)
7117625.300323632

結果が違うじゃん、って事で、Gauche本のTESTの章を参考に、何も考えずに、Hz以下は無視していいでねぇって事にしちゃったんだ。

(define (nearly-equal? a b)
  (< (abs (- a b)) 1.0))

(test* "test-mine-resf" 7117625 (resf 1e-6 500e-12) nearly-equal?)

でも、真値はどこにあるの？素直にC言語で確認。

[sakae@fb ~]$ ./a.out
   7117625.5000

これ、FreeBSDでC言語の結果だ。そうなると、FreeBSDのgauche-packageの挙動がおかしい？一応、第３者委員会の意見も聞いてみましょうかね。なるべくC言語のコンパイラーから遠い所に居る、独立機関の方がいいな。

[sakae@fb ~]$ maxima
;;; Loading #P"/usr/local/lib/ecl-20.4.24/sb-bsd-sockets.fas"
;;; Loading #P"/usr/local/lib/ecl-20.4.24/sockets.fas"
Maxima 5.44.0 http://maxima.sourceforge.net
using Lisp ECL 20.4.24
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
(%i1) 1.0 / (6.283185307179586 * sqrt(1e-6 * 500e-12)) ;
(%o1)                          7117625.434171771

FreeBSDではeclて言う組み込みCLを土台にしてた。

sakae@pen:~$ maxima

Maxima 5.44.0 http://maxima.sourceforge.net
using Lisp GNU Common Lisp (GCL) GCL 2.6.12
Distributed under the GNU Public License. See the file COPYING.
Dedicated to the memory of William Schelter.
The function bug_report() provides bug reporting information.
(%i1) 1.0 / (6.283185307179586 * sqrt(1e-6 * 500e-12)) ;
(%o1)                          7117625.434171771
(%i2) fpprec:50;
(%o2)                                 50
(%i3)  bfloat(1.0 / (6.283185307179586 * sqrt(1e-6 * 500e-12)));
(%o3)              7.117625434171770699322223663330078125b6
(%i4)  bfloat(1.0 / (6.283185307179586 * sqrt(1b-6 * 500b-12)));
(%o4)        7.1176254341717708622248965831569037806983783971197b6

入力式 %i3 は、少数点の変換で精度落ちしてるから、間違いだ。%i4式が、余す所なく、50桁まで、信用出来る。

Debianでは、GNU一族になったGCLを使ってるな。なんか系列企業な小会社みたいだな。これぞ、社会の縮図だ。

それは、さておき、両maximaでは同じ結果になった。値をFreeBSD,Debianと比べると、まあ、真ん中っぽい。果たして何が真値？

参考まで、Debian(i386)に固定し、次の比較をした

S式での定義         モジュール版(so)      C言語(a.out)
7117625.434171771   7117625.5           7117625.300324

インタープリタ対コンパイラ対 …

gauche-packageを使って用意した関数って、gccでコンパイルされた結果だよな。結果的に、みんなのあこがれchez風になってる訳だ。じゃ、普通にインタープリタ上で定義したやつと比べてスピードどうよってのが、頭をよぎる。

まずは、コンパイラ版

(use gauche.time)
(use mine)
(time-this 10000000 (lambda () (resf 1e-6 500e-12)))
#<time-result 10000000 times/  1.191 real/  1.190 user/  0.000 sys>
(time-this 10000000 (lambda () (resf 1e-6 500e-12)))
#<time-result 10000000 times/  1.218 real/  1.210 user/  0.000 sys>

時間計測用のモジュールが用意されてたので、それを使って１千万回実行。

(use gauche.time)
(use math.const)

(define (resf L C)
  (/ (* 2pi (sqrt (* L C)))))
resf
(time-this 10000000 (lambda () (resf 1e-6 500e-12)))
#<time-result 10000000 times/  3.172 real/  3.170 user/  0.000 sys>
(time-this 10000000 (lambda () (resf 1e-6 500e-12)))
#<time-result 10000000 times/  3.198 real/  3.200 user/  0.000 sys>

インタープリタ上で定義して実行。約 2.6倍の差がついた。

それじゃ、もう一つ getpid をmypidって言うscheme名で呼び出すやつ。前回は偏屈な理由でって煙に卷くような理由を付けたけど、これが出来るか確認したかったのさ。

(time-this 10000000 (lambda () (mypid)))
#<time-result 10000000 times/  9.826 real/  8.210 user/  1.600 sys>
(time-this 10000000 (lambda () (mypid)))
#<time-result 10000000 times/  9.825 real/  8.050 user/  1.760 sys>

えらい時間がかかってますなあ。と思ったら、根底はシステムコールなんだね。とっても重い呼出って事が、火を見るより明らかだ。

そんじゃ、お手軽に生getpidを使ってみる。gauche-c-wrapperね。これが本当にruby-ffiの対抗馬です。gauche-packageの方は、スタティックな挙動。それに対してc-wrapperはダイナミックな挙動になる。両者にどれぐらいの差が出るか、興味津々ですよ。

sakae@deb:/tmp/mine$ gosh
(use gauche.time)
(use c-wrapper)
(c-load-library "libc")
#<undef>
(getpid)
 *** ERROR: unbound variable: getpid
Stack Trace:
____________________________
        at "(standard input)":4
  1  (eval expr env)
        at "/usr/share/gauche-0.97/0.9.10/lib/gauche/interactive.scm":303
(c-include "sys/types.h")
(c-include "unistd.h")
(getpid)
1635

使うライブラリィの宣言だけでは、当然怒られた。ちゃんとmanして、インクルードファイルも教えてあげましょう（幾つもあると、うんざりする事間違いなし）。

(time-this 10000000 (lambda () (getpid)))
#<time-result 10000000 times/ 15.797 real/ 11.830 user/  3.910 sys>
(time-this 10000000 (lambda () (getpid)))
#<time-result 10000000 times/ 14.627 real/ 10.560 user/  4.020 sys>

５割増しで遅くなった(カーネル側の消費時間が大幅に増加してるけど、何に油売ってるのだろう)。常用するなら、やっぱりドラエモンのポケットだね(意味不なら、前回の記事を参照)。

mine.so 観光

お決まりの観光をしてみる。その前に、

sakae@deb:/tmp/mine$ nm mine.so
00001330 t minelib_mypid
00004160 d minelib_mypid__STUB
00001350 t minelib_resf
000041a0 d minelib_resf__STUB
00001300 t minelib_test_mine
00004120 d minelib_test_mine__STUB
00001230 T res
00001280 T test_mine
         U getpid@GLIBC_2.0

自分で書き下した関係者がこんな風に、反映されてた。

sakae@deb:/tmp$ gdb -q gosh
Reading symbols from gosh...
(No debugging symbols found in gosh)
(gdb) b resf
Function "resf" not defined.
Make breakpoint pending on future shared library load? (y or [n]) y
Breakpoint 1 (resf) pending.
(gdb) r
Starting program: /usr/bin/gosh
[Thread debugging using libthread_db enabled]
Using host libthread_db library "/lib/i386-linux-gnu/libthread_db.so.1".
gosh> (use mine)
gosh> (resf 1e-6 500e-12)

Breakpoint 1, resf (l=9.99999997e-07, c=4.99999986e-10) at ./mine.c:9
9         float tmp =  6.283185307179586 * sqrt(l * c);

関数に渡される引数の精度ってこんなものなの？ goshが何かやってないかな？後で確認。

(gdb) bt
#0  resf (l=9.99999997e-07, c=4.99999986e-10) at ./mine.c:9
#1  0xb7fc83c6 in minelib_resf (SCM_FP=0x9519c, SCM_ARGCNT=2, data_=0x0)
    at ./minelib.stub:9
#2  0xb7bd0b76 in ?? () from /lib/libgauche-0.97.so.0
#3  0xb7be47a5 in ?? () from /lib/libgauche-0.97.so.0
#4  0xb7be559a in Scm_ApplyRec () from /lib/libgauche-0.97.so.0
#5  0xb7be5ae6 in ?? () from /lib/libgauche-0.97.so.0
#6  0x0040445d in enter_repl ()
#7  0x00402dba in main ()

そして呼出履歴。aptから入れた本体はdebug情報が落されているんで、物足りないぞ。呼出のためのスタブは、何か呪文っぽいな。disassemble して読んでみろって助言は、却下します。

(gdb) finish
Run till exit from #0  resf (l=9.99999997e-07, c=4.99999986e-10) at ./mine.c:10
0xb7fc83c6 in minelib_resf (SCM_FP=0x9519c, SCM_ARGCNT=2, data_=0x0)
    at ./minelib.stub:9
9         (result (resf l c)))
Value returned is $1 = 7117625.5
(gdb) c
Continuing.
7117625.5

sakae@deb:/tmp$ gdb -q a.out
Reading symbols from a.out...
(gdb) b resf
Breakpoint 1 at 0x11bb: file z.c, line 5.
(gdb) r
Starting program: /tmp/a.out

Breakpoint 1, resf (l=9.99999997e-07, c=4.99999986e-10) at z.c:5
5         float tmp =  6.283185307179586 * sqrt(l * c);
(gdb) finish
Run till exit from #0  resf (l=9.99999997e-07, c=4.99999986e-10) at z.c:5
0x00401224 in main () at z.c:10
10        printf("%15.6f\n", resf(1e-6, 500e-12));
Value returned is $1 = 7117625.5
(gdb) c
Continuing.
 7117625.300324
[Inferior 1 (process 1916) exited normally]

出力の変換の所で、うまく丸めこまれている？？？

many scheme

色々なschemeが、どんな回答を返すか、下記のコードで確認してみる。

sakae@deb:/tmp$ cat resf.scm
(define (resf L C)
  (/ (* 6.283185307179586 (sqrt (* L C)))))

(display (resf 1e-6 500e-12))
(newline)
(exit)

chez scheme

sakae@deb:/tmp$ scheme
Chez Scheme Version 9.5.4
Copyright 1984-2020 Cisco Systems, Inc.

> (compile-file "resf.scm")
compiling resf.scm with output to resf.so

コンパイルすると、思わせぶりなsoファイルが出来るけど、unixとは全く関係ない。

sakae@deb:/tmp$ scheme resf.so
Chez Scheme Version 9.5.4
Copyright 1984-2020 Cisco Systems, Inc.

7117625.434171771

こんな風に使うしかないのさ。

petite

chez schemeのインタープリタ版

sakae@deb:/tmp$ petite resf.scm
Petite Chez Scheme Version 9.5.4
Copyright 1984-2020 Cisco Systems, Inc.

7117625.434171771

ちゃんと一致してる。

gsi

gambitのインタープリタ版

sakae@deb:/tmp$ gsi resf.scm
7117625.434171771

Gambit v4.9.3と言う版です。これのコンパイラ版が用意されてて、簡単にシングルバイナリーが作成出来る。

gsc

sakae@deb:/tmp$ gsc -exe resf.scm
sakae@deb:/tmp$ ./resf
7117625.434171771
sakae@deb:/tmp$ ls -l resf
-rwxr-xr-x 1 sakae sakae 57648 Sep 30 14:15 resf*

なかなか良い感じ！！

scm

sakae@deb:/tmp$ scm resf.scm
7.117625434171771e6

何故これらを選ぶのか？

schemeなら他にも色々あるだろう。guileとかracket(drscheme)とかchickenとか、星の数程ある。とても全部を試すって訳にはいかない。

オイラーの選択基準は、なるべく多数のOSでも動く事。gaucheは残念ながらOpenBSDでコンパイル出来無い。基準から外れる訳だけど、それに目をつぶれば、余りにもぬくぬく出来るから。今回みたいな事が出来るんで、益々楽しいな。

chez系は、世界的に有名だから。昔はBSD系でも動いたけど、シスコのお達しによりリナ系中心。かろうじてpetiteはFreeBSDでも動く。まあ合格の範疇かな。オイラーのつたない記事も各国からアクセスがあるしね。

gambitは、OSを選ばずって所が良い。ipadでも動くのは、もしもの時に遊び道具になってもらえそう。もしもって何だ？

scmは、オイラーのschemeへの入口だったから。初恋の人は忘れられないっていうじゃないですか。極めて情緒的な理由でした。

それより、論理だ論理だ。

float to double

と言う事で、重い腰をあげて、resf周りの宣言をfloatからdoubleに変更した。まあ、最初に気付くべき事だけど、みんなに指摘されて、脳味噌が回転すると言うおそまつさ。

ちょっと弁明しとくと、とっかかりのGauche-glが、f32vectorを多用してた。そして、参考にした記事でもfloatを使ってた、と言う脳内CPU停止状態だったから。

Openglは昔の事だから、floatがデフォでメモリー節約設計だったんだね。それに追従すれば、f32vectorを使わざるを得ない。

浮動小数点数型と誤差

type     min                       max
float    1.175494351e-38           3.402823466e+38
double   2.2250738585072014e-308   1.7976931348623158e+308

黙って double を使っておけば間違い無し。多数桁が欲かったら、maximaの多倍長浮動少数点数型を使え。

マシン固有のFPU

浮動少数点演算の同値判定って事で調べていたら、イプシロンなんていう定数があるそうな。マシンに固有の設定になってる。OpenBSDだと

sakae@deb:/usr/src/OpenBSD_6.9/sys$ grep  __FLT_EPSILON -rIl .
./arch/powerpc64/include/_float.h
./arch/arm64/include/_float.h
./arch/alpha/include/_float.h
./arch/mips64/include/_float.h
./arch/m88k/include/_float.h
./arch/sh/include/_float.h
./arch/i386/include/_float.h
./arch/amd64/include/_float.h
./arch/powerpc/include/_float.h
./arch/hppa/include/_float.h
./arch/arm/include/_float.h
./arch/sparc64/include/_float.h

なるほど、マシン毎にある。インテル入ってるなら、i387のコプロなハードな。

#define __FLT_RADIX             2               /* b */
#define __FLT_ROUNDS            __flt_rounds()
#define __FLT_EVAL_METHOD       2               /* long double */

#define __FLT_MANT_DIG          24              /* p */
#define __FLT_EPSILON           1.19209290E-07F /* b**(1-p) */
#define __FLT_DIG               6               /* floor((p-1)*log10(b))+(b == 10) */

1Bitの重みって事なのか。そして、仮数部24Bit,10進数だと6桁精度って読めばいいのかな。 80Bit仕様の定数もこのファイルに定義されてた。

sparc64だと、仮数部 113Bit,指数部33Bitなんていう機構もサポートしてるみたい。

F社はこれを捨てて、16BitバージョンのAI向けFPUを載せた、スパコン富嶽なんてのを作っているな。多数のFPUを同時に動かして、AI頑張りましょって作戦だな。

ちょっとびっくり

C言語でのresfをOpenBSDのgdbにかけた。

ob$ gdb -q a.out
Reading symbols from a.out...done.
(gdb) b resf
Breakpoint 1 at 0x19f0: file resf.c, line 4.
(gdb) r
Starting program: /tmp/a.out

Breakpoint 1, resf (l=8.7304488599225978e-61, c=0) at resf.c:4
4       double resf(double l, double c){
(gdb) bt
#0  resf (l=8.7304488599225978e-61, c=0) at resf.c:4
#1  0x000003d6d7809a87 in main () at resf.c:10

えっ、こんな引数なの？びっくり仰天。

(gdb) n
5         double tmp =  6.283185307179586 * sqrt(l * c);
(gdb) p l
$1 = 9.9999999999999995e-07
(gdb) p c
$2 = 5.0000000000000003e-10

ああ、引数が生えてきた。これ、gdbの個性？？

This year's Index

Home