3

Mathematica のボックス オプション式 (https://github.com/barrycarter/bcapps/blob/master/box-option-value.m) を Maxima-fy しようとしていますが、かなり単純な統合で Maxima がクラッシュします。

load(distrib); 
pdflp(x, p0, v, p1, p2, t1, t2) := pdf_normal(x,log(p0),sqrt(t1)*v); 
cdfmaxlp(x, p0, v, p1, p2, t1, t2) := 1-erf(x/(v*sqrt(t2-t1)/sqrt(2))); 

upandin(p0, v, p1, p2, t1, t2) :=  
 integrate( 
 float( 
 pdflp(x, p0, v, p1, p2, t1, t2)* 
 cdfmaxlp(log(p1)-x, p0, v, p1, p2, t1, t2) 
 ), 
 x, minf, log(p1));

特定の値で upandin を評価するとクラッシュします:

upandin(1, .15, 1.01, 1.02, 1/365.2425, 2/365.2425); 

rat: replaced -.00995033085316809 by -603/60601 = -.00995033085262619 

rat: replaced 2.718281828459045 by 23225/8544 = 2.718281835205993 

rat: replaced 8116.5 by 16233/2 = 8116.5 

rat: replaced 2.718281828459045 by 23225/8544 = 2.718281835205993 

rat: replaced -8116.5 by -16233/2 = -8116.5 

rat: replaced 1.0 by 1/1 = 1.0 

rat: replaced 1.792882852833688 by 4484/2501 = 1.792882846861255 

rat: replaced 180.1832400641081 by 126849/704 = 180.1832386363636 

rat: replaced 2.718281828459045 by 23225/8544 = 2.718281835205993 

rat: replaced -8116.5 by -16233/2 = -8116.5 

rat: replaced -1.0 by -1/1 = -1.0 

rat: replaced 1.792882852833688 by 4484/2501 = 1.792882846861255 

rat: replaced 180.1832400641081 by 126849/704 = 180.1832386363636 

rat: replaced 2.718281828459045 by 23225/8544 = 2.718281835205993 

rat: replaced -8116.5 by -16233/2 = -8116.5 

rat: replaced 1.0 by 1/1 = 1.0 

rat: replaced -1.0 by -1/1 = -1.0 
Maxima encountered a Lisp error: 

 The value 16090668801 is not of type FIXNUM.

upandin に float() がない場合、Maxima は積分を元の形式のままにします。

誰か助けてくれませんか？Mathematica を Maxima に変換するのは簡単だと思っていましたが、今は確信が持てません。

Mathematica のバージョンは問題なく動作します。

pdflp[x_, p0_, v_, p1_, p2_, t1_, t2_] :=  
 PDF[NormalDistribution[Log[p0],Sqrt[t1]*v]][x] 

cdfmaxlp[x_, p0_, v_, p1_, p2_, t1_, t2_] := 1-Erf[x/(v*Sqrt[t2-t1]/Sqrt[2])]; 

(* NIntegrate below "equivalent" to Maximas float(); no closed form *) 

upandin[p0_, v_, p1_, p2_, t1_, t2_] :=  
 NIntegrate[pdflp[x, p0, v, p1, p2, t1, t2]* 
           cdfmaxlp[Log[p1]-x, p0, v, p1, p2, t1, t2], 
{x, -Infinity, Log[p1]}] 

upandin[1, .15, 1.01, 1.02, 1/365.2425, 2/365.2425] 

0.0998337

編集:この関数を数値的に近似するオープンソースの Mathematica のようなプログラムはありますか? オープン ソース コードをオープン ソース プラットフォームにリリースしたいと思っています。

4

4 に答える 4

6

(たぶん、これに答えるビジネスはありませんが...)

単なる推測ですが、integrate は入力を再び正確にしたいと考えているようで、有理数演算を含む困難な bignum 計算を行っている可能性があります。おおよその e (オイラー数) を合理化するため、正確な入力を使用した integrate(0) とは異なる動作をする可能性があります。

確認したいことがあります

http://eagle.cs.kent.edu/MAXIMA/maxima_21.html

また

http://www.delorie.com/gnu/docs/maxima/maxima_62.html

Quadpack などの専用数値コード用。

(なぜ私がこれに答えようとしているのか、いまだに疑問に思っています。Stack Overflow のどこかに Maxima の専門家がいるに違いありません。)

Daniel Lichtblau Wolfram Research

于 2011-02-09T23:09:56.503 に答える
4

quad_qagi を使用して、無限区間にわたる積分を数値的に近似します。?? quad_ は、Quadpack 関数に関する情報を示します。

load (distrib);
pdflp (x, p0, v, p1, p2, t1, t2) := pdf_normal (x, log(p0), sqrt(t1)*v); 
cdfmaxlp (x, p0, v, p1, p2, t1, t2) := 1 - erf(x/(v * sqrt(t2 - t1)/sqrt(2))); 

upandin (p0, v, p1, p2, t1, t2) := block ([integrand],
   integrand : pdflp (x, p0, v, p1, p2, t1, t2) * cdfmaxlp (log(p1) - x, p0, v, p1, p2, t1, t2),
   quad_qagi (integrand, x, minf, log(p1))); 

upandin (1, .15, 1.01, 1.02, 1/365.2425, 2/365.2425);
 => [.09983372557898755, 2.839204848435967E-10, 225, 0]

返事が遅れて申し訳ありません。誰かが検索して見つけた場合に備えて、これをここに残します。

于 2012-05-31T19:25:02.017 に答える
3

私は Maxima がフロートを回避しようと懸命に努力していることを知っており、それがここでやろうとしていることだと思いますが、それを防ぐ方法を説明するには Maxima の第一人者ではありません。間隔を壊したり、被積分関数を手動で変換したりする必要があるかもしれませんが、ほとんどすべての数値がこれを処理できます。あなたはそれがかなり単純だと言っていることに注意してください。しかし、それは非常に急勾配です: これらのパラメータの場合、被積分関数は 0.1 で ~6*10^(-34)、-0.1 で ~3*10^(-206) です。これは、多くの単純な統合アルゴリズムを適合させるのに十分な範囲です。

とにかく、舞台裏で scipy と gsl のツールを使用して、Sage で十分に簡単に実行できます。

import scipy.stats

def pdflp(x,p0,v,t1):
    return scipy.stats.norm(log(p0), sqrt(t1)*v).pdf(x)

def cdfmaxlp(x,v,t1,t2):
    return (1-erf(x/(v*sqrt(t2-t1)/sqrt(2.))))

def upandin(p0, v, p1, p2, t1, t2):
    integrand = lambda x: (pdflp(x,p0,v,t1) * cdfmaxlp(log(p1)-x,v,t1,t2))
    return numerical_integral(integrand, -Infinity, log(p1))

sage: upandin(1, .15, 1.01, 1.02, 1/365.2425, 2/365.2425)
(0.099833725578983457, 7.5174412058308382e-07)

任意の精度が必要な場合は、mpmath の quad を使用します。[ここで「正しい」値を推測しましたが、最初はそれほど精度がないため、ばかげています。]

于 2011-03-27T07:15:57.423 に答える
3

Maxima の 'integrate' 関数は、数値ではなくシンボリックな積分を行います。積分から名詞形を返す場合、それは (シンボリック) 積分を実行できないことを意味します。エクスプレッションのパラメーターを正確な値から浮動小数点値に (「float」を使用して) 変更しても、それは変わりません。

あなたが探しているのは数値積分ルーチンだと思います.Maximaは、非常に基本的なrombergからさまざまなQuadpackメソッドまで、さまざまなものを提供しています(ドキュメントについては?? quadを試してください)。

    -s

PS「このくだらない「オープンソース」のもの」については、何がそれをもたらしたのですか? ウィキペディアの記事で Macsyma/Maxima の歴史を調べてみるとよいでしょう。

于 2011-04-01T18:18:28.420 に答える