任意精度実数計算用にPARI/GPを導入してみる。
最近、35段14次の高次ルンゲクッタ法をExcel VBAで計算するシリーズをやってみたけれど、VBAが倍精度実数しか使えないので精度ではいまいちだった。そこで多倍長実数(任意精度実数)か、そこまで言わなくても4倍精度実数が計算できるプログラミング言語がないかな、と少し探す。CのライブラリならすぐGMPとか
出るんですが、言語のほうがおもしろそう(うちのブログでやっているのは単に私の暇つぶしなので、覚えることが多いほうがいいということで。。。)
昔から知っているのはUBASIC
http://www.rkmath.rikkyo.ac.jp/~kida/ubasic.htm
なんですが、Windows版がないしちょっと古くなった感はどうしても否めない。同じくBASIC系で、(仮称)十進BASIC
http://hp.vector.co.jp/authors/VA008683/
も249桁できるということ。でも日本語の文献があり過ぎる!(いや、暇つぶしだから、文献はないほうが。。。)
で、見つけたのがPARI/GP。
http://pari.math.u-bordeaux.fr/
ダウンロードしてインストールはすぐにできた。
\p 100
とかすると、100桁の精度になります。
(表示が長すぎるんで、表示のみ短くするのは、default(format, "g0.28");とか。)
でPiと打つと、
3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068
と表示。簡単。
繰り返しがCと似てる分混乱するんだが、
for (i=1, 5, print(i))
見たいな感じ。最初for(i=1, i<=5,...)とか書いて全然動かなかった。ネストも()で行うんで、
for (i=1,5,
for (j=1,5,
...
)
)
みたいな。関数の定義も特にreturnで返すんじゃなくて
f(x) = {
x;
}
なんて感じ。この関数内だけで使う変数はlocal(a,b);とか。
ベクトルと行列は
v=vector(5,i,0.); m=matrix(5,5,i.j,0.);
で定義して、v[1]やm[1,2]とか。これも微妙にCと違うので混乱。また添え字は1から始まるのも微妙に違う。
でもこれだけあればもうルンゲクッタくらい書けそうだ。やってみよう(続く)。
« 「Yの悲劇」を読んだ。 | トップページ | 35段14次のルンゲクッタ法をPARI/GPで計算(Lorenz方程式) »
「学問・資格」カテゴリの記事
- 高周波(RF・マイクロ波・ミリ波・5G)関連ニュース2021年2月16日 IEEE Microwave Magazineの特集はオールデジタルのRFID、Microwave JournalはEバンド ミリ波通信に衛星や気球を使う話、アメリカの半導体企業がバイデンに投資を迫る、(2021.02.17)
- カオスを生じる電気回路、Chua’s circuitをLTspiceで回路シミュレーションしてみる。(2021.02.19)
- Labyrinth Chaos(迷宮カオス)を生むThomas-Rössler方程式のパラメータbを色々変えて、Python+Scipyでルンゲクッタ8次のDOP853(Dormand&Prince)を使って計算してGIFアニメ(2021.02.16)
- フィッツヒュー・南雲 (FitzHugh-Nagumo) 方程式をPython+Scipyでルンゲクッタ8次のDOP853(Dormand Prince)で計算。(2021.02.23)
- 「水晶振動子の等価回路計算」をカシオの高精度計算サイトkeisan.casio.jpの自作式としてUP! インピーダンスの大きさと位相がグラフ化できる。(2021.02.12)
« 「Yの悲劇」を読んだ。 | トップページ | 35段14次のルンゲクッタ法をPARI/GPで計算(Lorenz方程式) »
コメント