Excelで解くケプラーの方程式 | 池袋駅北口の「ぐんまのやぼう」

池袋駅北口の「ぐんまのやぼう」

「ぐんまのやぼう」攻略法の研究と実践

黒点の移動速度の求め方」で太陽の自転軸の傾いている方向を表す式として

  θi = θ + π(ti - t1)/1[year]

と書きました。

θは太陽の「自転軸の傾いている方向」なのですが、これは太陽から見た地球の方位角(の逆方向)(に一定値を足したもの)でもあります。

上の記事は「地球の軌道は円軌道」という前提ですからこれでいいのですが、じっさいの地球の軌道は楕円なのでこれは正確ではありません。

正しくθに相当する平均近点角lに対する太陽から見た地球の方位角f(真近点角)を求める必要があります。

それにはケプラーの方程式

  u - e sin u = l

を解いて離心近点角 u を求める必要があります(e:離心率です。真近点角fは離心近点角uから求めることができます)

この方程式は(非線形で)解析的には解けないので数値的に解を求める必要があります。

離心近点角uの初期値を適当に決めて上式の残差

  R= l - u + e sin u

を求め、残差Rが0に(近く)なるように計算を進めていきます。

残差 R を u で微分すると

  dR/du = -1 + e cos u

となりますから、まず

  R1 = l - u1 + e sin u1

を求め

  Δu1 = - R1 / ( dR/du ) = - R1 / ( -1 + e cos u1 )
  .
  .
  Δui = - Ri / ( dR/du ) = - Ri / ( -1 + e cos ui )

から

  u2 = u1 + Δu1
  .
  .
  ui = ui-1 + Δui-1

として新たな R を求めるということを繰り返しRが十分小さくなるu2,..,uiを求めればいいはずです。

u1をどう決めるか問題ですが

  平均近点角≒離心近点角

と考え

  u1 = l + e sin l

とします。

Excelでこれを計算してみます。

各セルに次のように計算式をセットします。
離心率eは0.01672とします。これは地球のじっさいの離心率です。
平均近点角l=1(57.29...度)の場合を計算します。
池袋駅北口のぐんまのやぼう

じっさいに計算するとあっけなく収束してしまい離心近点角 u=1.01420 が求まります。
池袋駅北口のぐんまのやぼう

-----

ところでExcelにはソルバーという機能があります。これを使うともっと安直です。

「ソルバー」はアドインなので最初にインストールする必要があるのですが、一度インストールしてしまえばExcelの一つの機能としてふつうに使えます。

次のようなシートを作ります。
池袋駅北口のぐんまのやぼう

データ=>分析=>ソルバーからソルバーを起動します。

uつまりA5を変化させてRつまりB5セルを0にすることが目標ですから
  目的セル: $B$5
  目標値: 値: 0
  変化させるセル: $A$5
と設定します。
池袋駅北口のぐんまのやぼう

そして実行ボタンをクリックすると次の結果が一瞬で得られます。
池袋駅北口のぐんまのやぼう

  u=1.01420

となりもちろん上の繰り返し計算を行ったときと同じ結果が得られます。

------

ここからペタしないでくださいね (^^) (自動ペタ回避・収集用です)
ペタ
ペタ
---------
PC・スマホの方はここからペタお願いします
携帯の方はここからペタお願いします