nptclのブログ

Common Lisp処理系nptの開発メモです。https://github.com/nptcl/npt

warnの出力抑制はどうやって実現しているのか

warnの出力を抑制する方法

(handler-bind ((warning #'muffle-warning))
  (warn "Hello"))

やった!

warn関数とその周辺の実装を考える

Common Lispwarningの話題です。
warn関数では警告を出力できますが、 その出力を抑止したい場合はどうしたらいいでしょうか。

先行して回答を示したように、handler-bindmuffle-warningを組み合わせて 設定することでwarnの出力を抑止することができます。

つまり、単純に

(warn "Hello")

と実行すると、

WARNING: Hello

みたいな出力が出るのですが、

(handler-bind ((warning #'muffle-warning))
  (warn "Hello"))

では何も出力されません。
これは一体どのような仕組みで実行されているのでしょうか?

handler-bindmuffle-warningという組み合わせから、 conditionrestartが手を組んでいることがわかると思います。 たかがwarning抑止だと舐めてかかれるようなモノじゃなく、 なんか知らんが結構複雑ですので、ちゃんと説明していきます。

muffle-warningとは何か

muffle-warningrestart関数と呼ばれています。
最後に使い方を説明しますが、warningを中断させる機能があります。
restart関数は他にもあり、abort関数やらcontinue関数がそれにあたります。 restart関数の内容は、多少の差異はあるものの大体決まっており、 find-restartinvoke-restartの2つにより構成されています。

muffle-warningの実装例である、muffle-warning!を下記に示します。

(defun muffle-warning! (&optional condition)
  (let ((restart (find-restart 'muffle-warning! condition)))
    (if restart
      (invoke-restart restart)
      (error (make-condition 'control-error)))))

たぶんどの処理系でもこんな感じに実装されていると思います。 find-restartで探して、見つかったならそれをinvokeするというもの。

つまりは、warningというconditionを受け取ったら、 muffle-warningより、再起動関数が呼び出されることになります。 逆にいうなら、warning conditionsignalで呼び出す前には、 必ずrestart-bindrestart-caseにてmuffle-warning restartを用意してあげなければ control-errorになってしまいます。 この要件が、まずはwarn関数を作成する際に必要になります。

注意点としては、restart関数は自分の関数名と同じ名前のrestartを探します。 今の場合、muffle-warning関数は、muffle-warningという名前のrestartを 探してinvokeしているのです。 関数名とrestartの名前はわざと同じにしているのですが、 同じsymbolでも無関係であることを覚えておいた方がいいと思います。

自作のwarning conditionを作る

テスト用にwarning conditionを作りましょう。
要件は1つ。formatの2要素である、format文字列と引数を受け取れることです。 slotを作るのが面倒なのでsimple-conditionを継承します。 名前はwarning!とします。

(define-condition warning! (simple-condition) ())

(defmethod print-object ((instance warning!) stream)
  (let ((format (simple-condition-format-control instance))
        (args (simple-condition-format-arguments instance)))
    (apply #'format stream format args)))

print-object methodは、例えばprincなんかで引数として与えられたときに、 どういう文字列を出力するかを指定するものです。 深く考えずに、warning!が受け取った引数から、formatをそのまま出力してやります。

例えばこんな感じ。

(let ((inst (make-condition
              'warning! :format-control "HELLO" :format-arguments nil)))
  (format t "~S~%" inst))
→ HELLO

自作のwarn関数を作る

では、上記の要件から、自作のwarn関数である、warn!関数を作成することを考えます。

要件は2つ。

  • muffle-warning! restartを用意する。
  • warning! conditionを実行する。

warn!関数の実装は下記のようになります。

(defun warn! (format &rest args)
  (restart-case
    (signal
      (make-condition
        'warning! :format-control format :format-arguments args))
    (muffle-warning! ())))

それではwarn!関数を次のように実行した場合はどうなるでしょうか。

(handler-bind ((warning! #'muffle-warning!))
  (warn! "Hello"))

warn!関数は、まずはmuffle-warning! restartを用意したあと、 warn! instanceに出力内容を設定して、signalconditionを呼び出します。
するとhandler-bindにより、muffle-warning!関数が呼び出されます。
muffle-warning!関数はinvokeによりwarn!関数に戻り、restart-caseが実行されます。
restart-caseの本体には何も記載がないものの、 bindではなくcaseであるため、warn!関数に制御がそのまま残り、 そして何もしないでwarn!関数が終了します。

流れは次のような感じになります。

  • warn!起動
  • signalによりhandler-bindwarning!起動
  • muffle-warning!関数実行
  • restart-casemuffle-warning!起動
  • restart-case終了
  • warn!終了

これは本来のwarn関数の挙動と同じになります。

それではhandler-bind無しでwarn!を呼び出した場合はどうなるでしょうか。
残念ながらこのままでは、warning! conditionhandlerが存在しないため、 signal関数が何もしないで終わってしまいます。 warn!関数が警告文字を出力するためには、 出力用のhandler-bindが必要となります。

出力用handler-bindを用意する

これはプログラマーが明示的に用意するものではなく、 Common Lisp処理系があらかじめ用意しておくものです。 warning!を出力するための、システムのコードを示します。

(handler-bind
  ((warning!
     (lambda (c)
       (format *error-output* "WARNING-TEST: ~A~%" c))))
  ;; code
  )

上記の;; codeの部分で、eval-loopであったり、 load関数だったりが動作するわけです。

この要件から、warn!関数は

(warn! "Hello")

としていたものは

(handler-bind
  ((warning!
     (lambda (c)
       (format *error-output* "WARNING-TEST: ~A~%" c))))
  (warn! "Hello"))

みたいな感じになります。 warn!関数がsignalを起動するため、handler-bindに制御が渡って、 format関数によりWARNING-TESTが出力されます。

一方、出力抑止である

(handler-bind ((warning! #'muffle-warning!))
  (warn! "Hello"))

は、

(handler-bind
  ((warning!
     (lambda (c)
       (format *error-output* "WARNING-TEST: ~A~%" c))))
  (handler-bind ((warning! #'muffle-warning!))
    (warn! "Hello")))

となります。

warn!関数はsignalによりhandler-bindを起動しますが、 内側のhandler-bindに制御が渡り、mufffle-warning!関数が呼び出されるため、 一番外側のhandler-bindには制御が渡らず、WARNING-TESTが出力されないのです。

この仕組みは、warning! conditionを途中で盗み見することができることを意味しています。 例えば次のコードを考えます。

(handler-bind
  ((warning!
     (lambda (c)
       (format *error-output* "<<<~A>>>~%" c))))
  (warn! "Hello"))

出力結果

<<<Hello>>>
WARNING-TEST: Hello

用意したhandler-bindで一度conditionを受け取ってから、 format<<<Hello>>>を出力させています。 conditionはそのまま外側に伝搬するので、 WARNING-TEST: Helloも出力されるわけです。

それでは、これ以上伝搬させたくない場合はどうしたらいいでしょうか?
そんな時に使用するのがmuffle-warning関数です。

(handler-bind
  ((warning!
     (lambda (c)
       (format *error-output* "<<<~A>>>~%" c)
       (muffle-warning! c))))
  (warn! "Hello"))

出力結果

<<<Hello>>>

テストコード

説明で作成した自作のwarn!関数を下記に示します。

(defun muffle-warning! (&optional condition)
  (let ((restart (find-restart 'muffle-warning! condition)))
    (if restart
      (invoke-restart restart)
      (error (make-condition 'control-error)))))

(define-condition warning! (simple-condition) ())

(defmethod print-object ((instance warning!) stream)
  (let ((format (simple-condition-format-control instance))
        (args (simple-condition-format-arguments instance)))
    (apply #'format stream format args)))

(defun warn! (format &rest args)
  (restart-case
    (signal
      (make-condition
        'warning! :format-control format :format-arguments args))
    (muffle-warning! ())))

(defmacro progn-warning (&body body)
  `(handler-bind
     ((warning! (lambda (c)
                  (format *error-output* "~&WARNING-TEST: ~A~%" c))))
     ,@body))


;;
;;  main
;;
(progn-warning
  (warn! "AAA")
  (warn! "BBB: ~A, ~A, ~A" 100 200 300)
  (handler-bind ((warning! #'muffle-warning!))
    (warn! "CCC")
    (warn! "DDD"))
  (warn! "EEE"))

出力例

WARNING-TEST: AAA
WARNING-TEST: BBB: 100, 200, 300
WARNING-TEST: EEE

clispだと2回出力されるのですが何故だろう。
無視します。

非正規化浮動小数の有効桁数を求める

C言語にてCommon Lispfloat-precision関数を作成したときのメモです。

非正規化の浮動小数とは、値が小さすぎるけど何とかして 浮動小数にて表現したいというギリギリの場合の方法です。

感覚としては、

1.23456 × 10の-103乗

みたいに表現したいけど指数-103は小さすぎて無理だから

0.00012 × 10の-99乗

のように表すのが非正規化です。

上記の例の場合、正規化の場合は1.23456という値から 6桁まで表現できているのがわかりますが、 非正規化の場合は0.000に押されて、12という2桁にまで減ってしまっています。 この「2桁」というのをC言語の世界で求める方法を説明します。

一応注意書きしておくと、上記の例は基数が10になっていますが、 実際はだいたい2なので注意してください。

まず普通の時、つまり正規化された浮動小数の桁数は、 C言語float.hにdefineが用意されています。

FLT_MANT_DIG     24
DBL_MANT_DIG     53
LDBL_MANT_DIG    53, 64, 113とか?

もし浮動小数が正規化の場合は、何も考えずに上記の値を返却します。

では正規化か非正規化かを調べるにはどうしたらいいでしょうか。 思いつく方法は大小を比較することですが、 C99にはfpclassifyという関数が用意されています。 例えばdoubleを対象にした場合、

fpclassifyの戻り値が
  FP_NORMAL なら DBL_MANT_DIG

とします。 浮動小数点の値がゼロの場合はどうしましょうか? ここでは有効桁数は0ということにします。

fpclassifyの戻り値が
  FP_ZEOR なら 0

さて問題はfpclassifyの戻り値がFP_SUBNORMALの場合、つまり非正規化数の場合です。
返却値から先に示すと、

fpclassifyの戻り値が
  FP_SUBNORMAL なら DBL_MANT_DIG - (DBL_MIN_EXP - 1) + ilogb(value);

となります。
DBL_MANT_DIGはすでに説明しました。
ilogb(value)は、指数部を整数で返却する関数です。

問題はDBL_MIN_EXPなのです。
DBL_MIN_EXPは何なのでしょうか。
検索すると2通りの意味が出てきます。

  • 基数のe乗で表現できる正規化数の最小値e
  • 基数のe-1乗で表現できる正規化数の最小値e

どっちなんだ。
つまり最小値そのものなのか、そうじゃなく最小値に1を足した数なのか? 現物を調査すれば早いので結果を言いますと、e-1です。

IEEE754では、double型の正規化数の指数部の最小値は-1022になります。 つまり-1023に相当する値(binary値だと0になる)が非正規化数となります。

正規化数の最小値-1022より、

e-1 = -1022
e = -1021

となるので、

DBL_MIN_EXP = -1021

なのです。

では本当にC言語ではそうなのか?
コンパイラの問題だということもあり得ますので、 ちゃんと調べるためには規格書を見るしかありません。

ここであっているのかな。

Page 31
minimum negative integer such that FLT_RADIX raised to one less than that power is
 a normalized floating-point number, emin
FLT_MIN_EXP
DBL_MIN_EXP
LDBL_MIN_EXP

英語わからないです。 でもone less thanとは、「1を引いた数」という意味のようです。 なるほど、なぜ1を引いているのかは知りませんが、規格がそうなっているんですね。

ということなので、

FP_SUBNORMAL なら DBL_MANT_DIG - (DBL_MIN_EXP - 1) + ilogb(value);

になります。

C言語にて、doubleのコードを示します。

#include <float.h>
#include <math.h>

int float_precision_double(double v)
{
    switch (fpclassify(v)) {
        case FP_NORMAL:
            return DBL_MANT_DIG;

        case FP_ZERO:
            return 0;

        case FP_SUBNORMAL:
            return DBL_MANT_DIG - (DBL_MIN_EXP - 1) + ilogb(v);

        default:
            fprintf(stderr, "error.\n");
            exit(1);
            return 0;
    }
}

実行の例を示します。

float_precision_double(123.456);
  →53
float_precision_double(0x01.0000p-1021);
  →53
float_precision_double(0x01.0000p-1022);
  →53
float_precision_double(0x01.0000p-1023);
  →52
float_precision_double(0x01.0000p-1024);
  →51
float_precision_double(0);
  →0

nptでは、ここで使っています。

最後に。
「非正規化」という名称はいろいろ種類があるように思えます。 デノーマル数という人もいれば、C言語では「SUBNORMAL」と表現されています。 Wikipediaの「非正規化数」が正しいのかもしれません。

rationalizeのアルゴリズムが全く分からなかった

Common Lispの関数rationalizeの話題です。

関数rationalizeとは、浮動小数を分数に変換する関数です。
例えば浮動小数0.1を分数に直すと次のようになります。

* (rationalize 0.1)
1/10

実は浮動小数0.11/10に変換するって難しい事なのです。
0.251/4に変換するなら簡単です。
0.25は基数2の浮動小数では自然に表せる数ですから。

0.1は違います。
浮動小数0.1というのは、基数が2の浮動小数では 正確に表すことができずに仮数が循環してしまうものとして有名です。 例えば十進数で103で割ると0.3333...になるのと同じです。

0.1を単純に分数へ変換するにはrational関数を使います。

* (rational 0.1)
13421773/134217728

すごい値ですね。
でもちゃんと0.1なんです。

* (coerce 13421773/134217728 'float)
0.1

関数rationalの実装は簡単です。 仮数×2の指数乗をそのまま分数形式にして通分するだけですから。

しかし関数rationalizeは違います。
どうも浮動小数0.1と誤差の範囲から適切なものを探し出しているようなのです。
ではどうやって実装されているのでしょうか。

考えてもわかりませんでした。
仕方がないのでcmuclsbclのソースを見てみました。

cmuclのソースを見ると、1618行目から下記のように始まっています。

;;; RATIONALIZE  --  Public
;;;
;;; The algorithm here is the method described in CLISP.  Bruno Haible has
;;; graciously given permission to use this algorithm.  He says, "You can use
;;; it, if you present the following explanation of the algorithm."

ヒュー!! カッコいい!!

たぶんcmuclさんがrationalizeアルゴリズムを使用したいってお願いして、 clispさんがいいよって言ったのかなって思います。

Bruno Haibleという名前には見覚えがあります。
なにせclispを起動するときに毎回出てくるのですから。
いつもお世話になっております。

$ clisp
  i i i i i i i       ooooo    o        ooooooo   ooooo   ooooo
  I I I I I I I      8     8   8           8     8     o  8    8
  I  \ `+' /  I      8         8           8     8        8    8
   \  `-+-'  /       8         8           8      ooooo   8oooo
    `-__|__-'        8         8           8           8  8
        |            8     o   8           8     o     8  8
  ------+------       ooooo    8oooooo  ooo8ooo   ooooo   8

Welcome to GNU CLISP 2.49 (2010-07-07) <http://clisp.cons.org/>

Copyright (c) Bruno Haible, Michael Stoll 1992, 1993
Copyright (c) Bruno Haible, Marcus Daniels 1994-1997
Copyright (c) Bruno Haible, Pierpaolo Bernardi, Sam Steingold 1998
Copyright (c) Bruno Haible, Sam Steingold 1999-2000
Copyright (c) Sam Steingold, Bruno Haible 2001-2010

Type :h and hit Enter for context help.

[1]>

こういうところで、Lisp仲間がつながっているのですね。

clispのソースはここで確認できます。

cmucl, sbcl, clispのソースには、コメントにアルゴリズムが載っています。

そのアルゴリズムを見ていたのですが、全く分からない。
sbclのソースには「An introduction to number theory」を見ろと書かれており、 さらにharvard.eduのリンク切れのサイトが紹介されている。 何となく思うんですが、 分数の割り算でつまずいている自分が 理解するのは非常に困難なんじゃないか。

理解できないなら既存のソースをたどって作るしかありません。 cmuclのソースはPublic Domainとの事なので、ありがたく参考にさせていただきます。

nptではcmucllispコードをもとにC言語で実装しました。 - https://github.com/nptcl/npt/blob/v0.1.6/src/real_decode.c#L746

他の関数はC言語に極力逆らわないように実装しているのですが、 関数rationalizeに関しては、C言語の皮をかぶったlispになっております。 しかしちゃんと動作しているようなので安心しました。

closureには何が保存されるのか

closureとは関数にデータを保存するための仕組みです。
Lisp大好きな人はlambdaと一緒に多用します。
例えばこんな使い方をします。

(let ((x 0))
  (setq *call* (lambda () (incf x))))

(funcall *call*)
 -> 1
(funcall *call*)
 -> 2
(funcall *call*)
 -> 3
(funcall *call*)
 -> 4

上記の例では、lambdaにより生成された無名関数が、 変数xをclosureにて保存しています。 letで宣言された変数xのスコープはとっくに終わっているにもかかわらず、 funcallで関数を呼び出すと、変数xは立派に役目を果たしているというわけです。

では本題ですが、closureは何を対象にどんなものを保存するのでしょうか。
lexical変数だよ!と思った人は正解です。 問題はそれだけなのかということ。

典型的な回答としては下記の通り。

  • lexical変数
  • flet, labels関数

flet, labelsは局所関数を作成するためのものです。
変数ではなく関数です。
例えばこんな感じ。

(flet ((call () :hello))
  (setq *call* (lambda () (call))))

(funcall *call*)
 -> :HELLO

面倒な話題として、setf関数もclosureの対象となります。
説明すると長くなるので手短にしますが、 defsetfdefine-setf-expanderで定義されたマクロ形式の方ではなく、 flet, labelsで定義した(setf 名前)形式の関数がclosureの保存対象となります。 何が違うんだって思われるかもしれませんが、なんか結構違うんです。 そのうち詳しく説明します。

例をあげます。

(flet (((setf aaa) (value cons)
          (rplaca cons (* value value))
          value))
  (setq *call* (lambda (c v)
                 (setf (aaa c) v))))

(let ((c (cons 10 20)))
  (funcall *call* c 999)
  c)(998001 . 20)

ではdefun宣言による関数はclosureに保存されないのでしょうか?
実はされません。
下記に例をあげます。

(defun aaa ()
  :hello)

(setq *call* (lambda () (aaa)))

(funcall *call*)
  →:HELLO

(defun aaa ()
  :zzz)

(funcall *call*)
  →:ZZZ

defunによる大域関数はspecial変数みたいな扱いなんですね。

ここまではどのCommon Lispでも通用するでしょう。
ここからは処理系に依存する話になるかもしれません。

話は大きく変わり、tagbodyの実装をしたときの話です。
tagbodyのスコープを脱出したgoについて考えていました。

tagbodyは動的エクステントなので、スコープの外で呼び出されたgoは無効となります。 問題は「無効」という意味が、放っておいてもいいオブジェクト という意味では「ない」ということです。 無効になったtagは、呼び出された時点でエラーになります。 処理系は「エラーになる」という所まで面倒見てやらなければいけないのです。

例えば次の通り。

(tagbody
  (format t "Hello~%")
  again
  (format t "Tagbody~%")
  (setq *call* (lambda () (go again))))

(funcall *call*)
  →エラー、againはすでに閉じられている

どうすればいいのか。
goに結びつけられたtagを専用のオブジェクトにするのです。 tagには戻り先のtagbodyに加えて、有効無効を表すフラグを持たせます。 tagはtagbodygoで共有されます。 もしtagbodyがスコープを抜けたら、持っているすべてのtagに無効フラグを立たせます。 もしgoが無効のtagを呼び出したらエラーにします。

では、そのtagオブジェクトは誰がいくつ保有するのか。 最初の考えだと、tagオブジェクトはtagbodyのコードが唯一保有すること考えていました。 つまり、(tagbody ...)という表記が現れたら、 そのtagbodyがtagオブジェクトを保有するのだという考えです。 でもこれだと、再帰呼出か、スレッドによる同時実行により破綻します。 オブジェクトが一個しかないのですから。

解決方法は、tagオブジェクトはtagbodyを実行するたびに生成するしかなかったのです。 感覚としてはこんな感じ。

(let ((again (make-tagobject%)))
  (tagbody%
    (format t "Hello~%")
    again
    (format t "Tagbody~%")
    (setq *call* (lambda () (go% again)))
  (close-tagobject% again))

となると、当然tagオブジェクトもclosureに保存する必要が出てきます。
話が長くなりましたが、結局はtagbodyのtagもclosureに入れるんだよっていうことです。

全く同じ話がblock/return-fromにも当てはまります。 しかしcatch/throwは対象外。 catchspecial変数みたいな扱いなので必要ありません。

結論は、closureには下記のオブジェクトが格納される。

  • lexical変数
  • flet, labels関数
  • flet, labelssetf関数
  • tagbodyのtag
  • blocksymbol

分数の割り算で余りを求める

分数の割り算で、余りを求める方法がわからなかったという内容です。
どうも小学校の算数の内容らしいですが、わからなかったのだから仕方がない。
あと「帯分数」というモノを初めて知りました。 習った記憶が全くない。

Common Lispでは、truncate関数を使うことによって求められます。
では計算方法はどうなるのでしょうか。

例えば「100/3 ÷ 7/3」の商と余りはいくつになるでしょう。
素直に割り算すると100/7になるので、 「14 余り 2/7」が答えになりそうですが違います。 これだと割り算の結果を、ただ整数と分数の部分に分けているだけです。 つまりはこんな感じ。

* (truncate (/ 100/3 7/3))
14 ; 2/7

正しい回答を出力してみます。

* (truncate 100/3 7/3)
14 ; 2/3

商の14までは合っていますが、余りが違います。 この余りをちゃんと求めるには、商と余りの関係を式で表す必要があります。

割り算

a ÷ b = 商...余り

のとき、書き換えると

余り = a - b×商

なので、

余り = 100/3 - (7/3)×14

計算をすると余り = 2/3となります。

この分数の余りの求め方はfloor, ceiling, truncate, round, mod, remのすべてで使っています。

当たり前と言えばそうなんですが、実装面だと割り算で分数を扱うときに 初めて上記の計算が必要になります。 分数以外の型であるfixnum, bignum, floatは別の方法で求めるので、 この計算を自分でしなくてもいいというのがなんとなく不思議な感じがします。

Common LispとC言語で虚数部のゼロ符号に違いが出る

ずっと悩んでいたことがわかったのでメモします。
例えば下記の関数

log(1-z)

のzに10を入れると、Common LispC言語では結果が違います。

Common Lispの例

(defun log1-z (z)
  (log (- 1.0 z)))
(log1-z 10)
  -> #C(2.1972246 3.1415927)

C言語の例

#include <stdio.h>
#include <complex.h>
#include <math.h>

complex double log1_z(complex double z)
{
    return clog(1.0 - z);
}

int main()
{
    complex double z = log1_z(10.0);
    printf("%f, %f\n", creal(z), cimag(z));
    return 0;
}
$ cc main.c -lm
$ ./a.out
2.197225, -3.141593

結果

Common Lisp
2.1972246  3.1415927
C言語
2.197225  -3.141593

虚数部の符号が違いますね。
それだけなんですが、なぜ?

確認していったら、どうも関数logに渡される引数の 虚数部ゼロの符号が違っているからのようです。 下記の式の部分です。

1.0-z

Common Lispでは

1.0 - 10.0 = -9.0

と当たり前の数を出しますが、 C言語では引数がそもそも複素数で指定されているため、

1.0 - 10.0 = -9.0 -0.0*I

のようになります。 このゼロにマイナスがあるかないかで、logの値が変わっていました。 数学に全然詳しくないんですが、これはどっちも正しいでいいんですかね。 ダメだとしても、差異はどうしようもなさそうです。

言語レベルでの話をしますと、Common Lispでは簡単にはC言語の様な挙動に できないと思います。
関数呼び出し時に

(log1-z 10)

(log1-z #c(10 0))

としたところで、

#c(10 0)

10

に変換されてから関数に渡されてしまいます。

【追記】間違いです。 #c(10 0)ではなく浮動小数#c(10.0 0.0)とすればcomplex型のままになります。 以降の内容にはそれほど影響がなかったため、内容はそのままにします。

ではlog1-zの関数定義を下記のように変更したらどうか。

(defun log1-z (z)
  (declare (type complex z))
  (log (- 1.0 z)))

ダメですね。 declare typeは型チェックを行う役目しかありません。 sbcl, clisp, cclで確認しましたが期待通りの動きはしません。

ただし、負ゼロの値を渡すことはできる場合はあります。 下記をsbclとcclにて確認しました。

* (log -10.0)
#C(2.3025851 3.1415927)
* (log #c(-10.0 -0.0))
#C(2.3025851 -3.1415927)

ただしclispではダメ。

まあ、本当に必要なら符号見て場合わけすりゃいい話なんですが、 標準関数内に上記のコードが埋め込まれている場合があります。 上記の差異はarcsinを確認しているときに見つかりました。

Common Lispでのasin(30)
  #C(1.5707964 -4.0940666)

C言語でのcasin(30.0)
  1.570796, 4.094067

arcsinの定義を追っていくと、

arcsin(z) = -i*log(i*z + sqrt(1-z*z)) 

となり、上記sqrt部が

sqrt(1-z*z) = exp(log(1-z*z) / 2.0)

となるため、log(1-z*z)が含まれるので、虚数部ゼロの符号に差異が生じるというわけです。

何が正しいのか全くわからん。

C言語にて関数csinlが使えない場合がある

数学のsin関数に複素数って渡せるんですね。 全然知りませんでした。

C言語のc99標準ライブラリにも用意されているようです。
具体的には下記の関数。

csinf, csin, csinl, sin(tgmath)

問題は、clangのlong double complexを引数とする関数がいくつか使えない事です。 例えばclangではcsinlを使おうとするとエラーになります。 clangというとFreeBSDの標準コンパイラですが、 manコマンドのマニュアルにも載ってなようなので、 少なくともFreeBSD開発側は認識できているようです。
確認した環境は下記の通り。

% freebsd-version
12.0-RELEASE-p3
% clang --version
FreeBSD clang version 6.0.1 (tags/RELEASE_601/final 335540) (based on LLVM 6.0.1)
Target: x86_64-unknown-freebsd12.0
Thread model: posix
InstalledDir: /usr/bin

モジュール自体が存在しないので、csinlがあるとコンパイルが通りません。 tgmath.hを利用してsin関数にlong double complexを渡しても、 あれは結局マクロでcsinlにするだけなのでエラーになります。

そうなんですか、無理ですか。
完全に他人任せで実装されるまで大人しく待ってようと思います。

しかし妥協すれば解決策はいくつかあると思います。
私の大のお気に入りは、下記の極めて最低な方法。

#ifdef __clang__
#define csinl csin
#endif

他には、C言語をあきらめて、C++を利用する方法。

#include <stdio.h>
#include <cmath>
#include <complex>

int main()
{
    std::complex<float> f(1.2f, 3.4f);
    std::complex<double> d(1.2, 3.4);
    std::complex<long double> l(1.2L, 3.4L);

    f = std::sin(f);
    d = std::sin(d);
    l = std::sin(l);

    printf("%30.25e, %30.25e\n", real(f), imag(f));
    printf("%30.25e, %30.25e\n", real(d), imag(d));
    printf("%30.25Le, %30.25Le\n", real(l), imag(l));

    return 0;
}

実行結果

$ clang++ test.cpp
$ ./a.out
1.3979410171508789062500000e+01, 5.4228153228759765625000000e+00
1.3979408806017994848502894e+01, 5.4228154724634016758955113e+00
1.3979408806017997938912767e+01, 5.4228154724634012448167275e+00

うまく行っているように見えます。
なぜC++でうまく行ってC99じゃダメなのか。

他には、clangではなくgccを利用する方法。
2月くらいのGentoo Linuxで確認しましたが、こちらは問題ありません。

そしてあまりお勧めできないのは自分でcsinlを実装する方法です。
正直csinl程度なら何とかなりそうなもんですけど、 arcsinやらになると、数値計算に自信がないならやめた方がいいかもしれません。

もしかして自分が勘違いしているだけでcsinlはc99標準じゃない? 関数clogなんて丸ごと存在しないんですよね。 わざわざ/usr/includeを確認してしまった。