2

次のアルゴリズムで集中的に使用される係数の計算に依存するモデルチェッカーを作成しています。

![代替テキスト][1]

qは double であり、tdouble もkint です。e指数関数の略です。この係数は、(そのステップの) 前のすべての係数の合計が 1 に達するまで、常に 0 から始まり、変更さqれtないステップで使用されます。k

私の最初の実装は文字通りのものでした：

let rec fact k =
  match k with
    0 | 1 -> 1
    | n -> n * (fact (k - 1))

let coeff q t k = exp(-. q *. t) *. ((q *. t) ** (float k)) /. float (fact k)

kもちろん、小さなしきい値 (15-20) を超えた場合、階乗全体を計算することは不可能だったので、これはそれほど長くは続きませんでした: 明らかに結果がおかしくなり始めました。そこで、インクリメンタル除算を行って全体を再配置しました。

let rec div_by_fact v d =
  match d with
    1. | 0. -> v
    | d -> div_by_fact (v /. d) (d -. 1.)

let coeff q t k = div_by_fact (exp(-. q *. t) *. ((q *. t) ** (float k))) (float k)

qこのバージョンは、とtで十分な「通常」の場合は非常にうまく機能しますが、たとえばq = 50.0、一連の 0t = 100.0から計算を開始しk = 0 to 100、特定の数から最後まで NaN が続きます。

もちろん、これは、数値が 0 に近づき始める操作または同様の問題によって引き起こされます。

式を最適化して、幅広い入力に対して十分に正確な結果を得る方法について何か考えはありますか?

すべてがすでに 64 ビットである必要があります (デフォルトで double を使用する OCaml を使用しているため)。128 ビットの double を使用する方法もあるかもしれませんが、その方法はわかりません。

私は OCaml を使用していますが、C、C++、Java など、どの言語でもアイデアを提供できます。私はそれらすべてをかなり使用しました。

4

2 に答える 2

3 
qt^k/k! = e^[log[qt^k/k!]]
log[qt^k/k!] = log[qt^k] - log[k!] // log[k!] ~ klnk - k  by stirling
             ~ k ln(qt) - (k lnk - k)
             ~ k ln(qt/k) - k

kの値が小さい場合、スターリング近似は正確ではありません。ただし、既知の有限範囲を実行しているように見えるためlog[k!]、エラーを回避して、計算して配列に入れることができます。

もちろん、さらに実行できるバリエーションは複数あります。

于 2010-09-03T19:45:01.503 に答える
1

これは答えではありませんが（私は信じています）、おそらく単なる説明です。何か誤解した場合は、コメントの後に削除します。

私が理解しているように、次の合計が1に等しいなど、nを計算しようとしています。

代替テキスト

漸近的に 1 に近づくことがわかるかもしれませんが、決して 1 に等しいことはありません。質問を誤解していた場合は、訂正してください。

于 2010-09-04T00:50:08.763 に答える