vendredi 31 octobre 2008

OSPE 2.0

一応古在が見えるけれども、I が高い時に、e_min がどんどん上がっていってしまう。
OSPE1.0 と比べて、quadrupole terms は全く同じ。quadrupole の dei/dt, dgi/dt を OSPE1.0 から 2.0 に持ってきて貼り付けても、やはりこのドリフトが出てくる。と言うことは、GSL の Runge-Kutta の誤差が溜まってきているんではないだろうか。そういえば Fabrcyky や Wu のシミュレーションでも、この e_min が上がってくる現象が見られる。もしかして潮汐のせいではなくて積分の誤差のせい?

実験1:
ー タイムステップを、GSL_ODEIV_STEP に任せないで、自分で設定してやってみる。
解答1:
ー 随分遅くなるが、精度は上がる様。1億年の計算で dt = 1000 yr でも誤差が目立つ。100年だとドリフトは出ないが、計算にやたらと時間がかかる。

そして。

stepping algorithm によって結果が全然違う。2nd-order embedded Runge-Kutta だと、quadrupole だけの場合は rkf45 と比べても全くドリフトが出ない!ところが Octupole を入れるとどんどん上がっていく。
rk4.......quadropole,僅かな誤差累積。octupole, 全然ダメ。
(これらはどれも厳密にタイム・ステップを1000年に設定した場合)


........解決。eps を忘れていた。1e-6 のままだった。これではあかん。1e-8 にすると、quadrupole の場合はきちんと計算される様になった。Octupole はどこかで書き間違えているかも知れない。これからチェックする。

Fabrycky は 1e-8 でやっている。これでは充分ではないかも知れない。後で検査してみる。

------

間違い発見。最後の最後で e1 と e2 を間違っていた。何回もチェックしたのに見逃していた。
とにもかくにも、一応新しい OSPE は完成。ファイル1つに簡潔に収まったのは良かったけど、 GSL は格別速くはないみたい。次は潮汐など他の摂動。

本日ハロウィーン。IV は大変な事になっている模様。今週ももうお終い。取り敢えず OSPE が書き直せたのは成果と言って良い。学振の結果はどうなったか。

Aucun commentaire: