diffeq2nd_planet_verlet.py 技術ドキュメント

プログラムの動作

diffeq2nd_planet_verlet.py は、複数の天体間の万有引力相互作用による運動をシミュレーションするためのPythonスクリプトです。本来はVerlet法を用いて連立2階微分方程式を解き、天体の軌道を計算することを目的としていますが、現在のコードでは主要な時間発展計算部分がコメントアウトされているため、シミュレーションは実行されません。

プログラムの主な機能は以下の通りです。

  • planet_db.csv というCSVファイルから、太陽系内の天体（惑星など）の初期データ（質量、初期位置、初期速度）を読み込みます。

  • 読み込んだ天体の初期情報を標準出力に表示します。

  • 各天体に働く万有引力に基づく初期の合力（加速度）を計算します。

  • シミュレーション結果を保存するためのCSVファイル diffeq2nd_Planet_Verlet.csv を作成し、ヘッダー行を書き込みます。

現状では、時間発展シミュレーションは行われず、初期設定と出力ファイルの準備までで処理が終了します。

原理

このプログラムは、本来以下の物理法則と数値計算アルゴリズムに基づいて天体運動をシミュレートする設計となっています。

万有引力の法則

プログラムは、複数の天体間に働く万有引力を計算します。質量 \(M_i\) を持つ天体 \(i\) と質量 \(M_j\) を持つ天体 \(j\) の間に働く引力の大きさ \(F_{ij}\) は、万有引力定数 \(G\) と天体間の距離 \(r_{ij}\) を用いて以下のように表されます。

\[ F_{ij} = G \frac{M_i M_j}{r_{ij}^2} \]

各天体には、他のすべての天体からの引力が合力として作用し、それによって天体は加速度を得ます。プログラム内では、天体 \(i\) に働く他の天体 \(j\) からの力成分を Fij 関数で計算し、Fi 関数で合力を計算しています。

Verlet法 (Velocity Verlet法)

天体の運動方程式（2階微分方程式）を数値的に解くために、Verlet法（特にVelocity Verlet法）が用いられることを意図しています。Verlet法はエネルギー保存性が高く、長時間のシミュレーションに適しています。基本的な更新式は以下の通りです。

  1. 位置の更新: $\( \vec{r}(t+\Delta t) = \vec{r}(t) + \vec{v}(t)\Delta t + \frac{1}{2}\vec{a}(t)(\Delta t)^2 \)$

  2. 半ステップ速度の更新: $\( \vec{v}(t+\frac{1}{2}\Delta t) = \vec{v}(t) + \frac{1}{2}\vec{a}(t)\Delta t \)$

  3. 新しい加速度の計算: 更新された位置 \(\vec{r}(t+\Delta t)\) を用いて、新しい時刻 \(t+\Delta t\) での加速度 \(\vec{a}(t+\Delta t)\) を計算します。

    \[ \vec{a}(t+\Delta t) = \frac{\vec{F}(\vec{r}(t+\Delta t), t+\Delta t)}{m} \]

  4. 最終速度の更新: $\( \vec{v}(t+\Delta t) = \vec{v}(t+\frac{1}{2}\Delta t) + \frac{1}{2}\vec{a}(t+\Delta t)\Delta t \)$

ただし、現在の diffeq2nd_planet_verlet.py のコードでは、上記Verlet法による時間発展のループ計算部分が3重引用符でコメントアウトされており、実際に天体運動のシミュレーションは行われません。

単位系

プログラム内で使用される定数 G1 は、以下の単位系に合わせてスケーリングされています。

  • 距離: 天文単位 (Astronomical Unit, AU)

  • 時間: 日 (Day)

これにより、天体の位置や速度の数値が扱いやすい範囲になります。

必要な非標準ライブラリとインストール方法

このプログラムは、以下の非標準ライブラリを使用しています。

  • numpy: 数値計算を効率的に行うためのライブラリです。プログラム内でインポートされていますが、現在のコード（コメントアウトされた部分を除く）では直接使用されていません。

これらのライブラリは、Pythonのパッケージ管理ツール pip を使用してインストールできます。

pip install numpy

必要な入力ファイル

プログラムを実行するには、planet_db.csv という名前の入力ファイルが必要です。このファイルには、シミュレーション対象の各天体の初期情報がCSV（カンマ区切り）形式で記述されます。

ファイル名

planet_db.csv

期待されるファイル形式とデータ構造

CSVファイルはヘッダー行を含み、以下の列を持つ必要があります。

列名

データ型

説明

Name

文字列

天体の名前（例: Sun, Earth）

Mass

浮動小数点数

天体の質量 (kg)

Revolution Radius

浮動小数点数

初期公転半径 (AU)。中心天体からの初期距離。

Revolution Velocity

浮動小数点数

初期公転速度 (AU/日)。

planet_db.csv の例

Name,Mass,Revolution Radius,Revolution Velocity
Sun,1.989e30,0.0,0.0
Earth,5.972e24,1.0,0.0172
Mars,6.39e23,1.52,0.0139

生成される出力ファイル

プログラムは、シミュレーション結果を以下のCSVファイルに出力することを意図しています。

ファイル名

diffeq2nd_Planet_Verlet.csv

内容

  • ヘッダー行: シミュレーション時間 (t) と、各天体のX座標、Y座標のラベルが並びます。 例: t,x(Sun),y(Sun),x(Earth),y(Earth),x(Mars),y(Mars) ただし、nprint_planets の設定により、すべての天体が出力されるわけではありません。

  • データ行: 本来は、各時刻における各天体のX、Y座標が記録される予定です。

重要: 現在のプログラムでは、Verlet法による時間発展の計算部分が3重引用符でコメントアウトされているため、diffeq2nd_Planet_Verlet.csv にはヘッダー行のみが書き込まれ、それ以降のシミュレーションデータは記録されません。

コマンドラインでの使用例 (Usage)

diffeq2nd_planet_verlet.py はコマンドライン引数を取りません。以下の形式で実行します。

python diffeq2nd_planet_verlet.py

コマンドラインでの具体的な使用例

以下の入力ファイル planet_db.csv を用意します。

入力ファイル planet_db.csv の準備

Name,Mass,Revolution Radius,Revolution Velocity
Sun,1.989e30,0.0,0.0
Earth,5.972e24,1.0,0.0172
Mars,6.39e23,1.52,0.0139

実行コマンド

python diffeq2nd_planet_verlet.py

期待される標準出力

Planet simulator: Solve simulataneous second order diffrential equations by Verlet method
Planets:
   Sun
      Mass: 1.989e+30
      Revolution Radius: 0.0
      Revolution Velocity: 0.0
   Earth
      Mass: 5.972e+24
      Revolution Radius: 1.0
      Revolution Velocity: 0.0172
   Mars
      Mass: 6.39e+23
      Revolution Radius: 1.52
      Revolution Velocity: 0.0139

Write to [diffeq2nd_Planet_Verlet.csv]
  t        x(Sun)        y(Sun)      x(Earth)      y(Earth)

生成される diffeq2nd_Planet_Verlet.csv の内容

t,x(Sun),y(Sun),x(Earth),y(Earth)

実行結果の説明

プログラムはまず、planet_db.csv から読み込んだ太陽、地球、火星の初期情報を標準出力に表示します。次に、出力ファイル名 (diffeq2nd_Planet_Verlet.csv) が表示され、続けて出力されるCSVファイルのヘッダー行（時間、太陽のX/Y座標、地球のX/Y座標）が標準出力にも表示されます。

重要な注意点として、Verlet法による時間発展の計算コードがコメントアウトされているため、シミュレーションは初期化ステップまでしか実行されません。したがって、diffeq2nd_Planet_Verlet.csv ファイルにはヘッダー行のみが書き込まれ、それ以降の時刻における天体位置データは含まれません。