相対性理論衛星軌道計算シミュレーターフォーム by javascript

CLSですべての衛星をクリア
DELで選択衛星をクリア
CALLで選択衛星の軌道要素を呼び出し
APPLYで選択衛星の軌道要素を適用
惑星ボタンで選択衛星に惑星軌道要素を適用
ラジオボタンやチェックボックスは軌道要素の計算の優先度を決めます
Edgeで閲覧すると処理が重くなってくる時があります その場合、IEのリセットをしてください

アプリの紹介

太陽系惑星ならばほぼ完璧です
物理エンジンは
ケプラー方程式で初期値を取得
ルンゲクッタで万有引力を積分
万有引力の計算で
1PN : 相対論補正項S=3(v/C)^2
2.5PN : 重力放射減衰L=(-32*A/5)*((G^3*m1*m2*(m1+m2))/(C^5) * r^4))) * v
v:軌道速度:r:距離:m1:重力源質量:m2:質量:G:万有引力定数:C:光速度:A:NAS6変数(デバッグにより経験的に決定):
相対論速度合成則
を組み込んでいます
ログのCSV出力もできます

ハレー彗星を解くためにべき乗可変dtを採用したので
スケールに応じて自動調整されるようになっています
   




相対論モード古典論モード
    speed rate :  zoom :  
経過時間 time(年) = time(year)
経過時間 time(年)が遠日点付近からの開始がもっとも初期設定による誤差が少ないです
年 月 
  時 分 
  Sat01 Sat02 Sat03 Sat04 Sat05 
Sat06 Sat07 Sat08 Sat09 Sat10
   
太陽系内
   
   
   
  
 
注意:Jupiter swing byは可変DTではタイミングの設定が合わなくなります
太陽系外その他
   
   
仮想ブラックホールAは境界設定であり重力で拘束できずに飛んでいきます
仮想ブラックホールBは1PN崩壊境界設定でありまだ重力で拘束されています
deg rad 
入力(速度計算:単位、AU、年)(デフォルト値:水星)  input(CalcVelocity:unit:AU,year)(default:Mercury)
重力源名前 = attractor name
名前 = name
時間 t0(年) = time0 T0(year)
公転周期 P(年) =   orbital period P(year)
重力源質量 M1(kg) =   attractor mass M1(kg)
衛星離心率 e =   eccentricity e
衛星近日点距離 RA1(AU) = perihelion RA1(AU)
衛星遠日点距離 RA2(AU) = aphelion RA2(AU)
衛星質量 M2(kg) = satellite mass M2(kg)
衛星昇交点黄経 Ω(deg,rad) = longitude of ascending node Ω(deg,rad)
衛星軌道傾斜角 i(deg,rad) = inclination i(deg,rad)
衛星近日点引数 ω(deg,rad) = argument of perihelion ω(deg,rad)
衛星元期 T*(年) = epoch T*(year)
衛星元期平均近点角 l(t*)(deg,rad) = epoch mean anomaly l(t*)(deg,rad)


入力(速度計算:単位、m、s)(デフォルト値:水星)  input(CalcVelocity:unit:m,s)(default:Mercury)
重力源名前 = attractor name
名前 = name
時間 t0(s) = time0 T0(s)
公転周期 P(s) =   orbital period P(s)
重力源質量 M1(kg) =   attractor mass M1(kg)
衛星離心率 e =   eccentricity e
衛星近日点距離 RA1(m) = perihelion RA1(m)
衛星遠日点距離 RA2(m) = aphelion RA2(m)
衛星質量 M2(kg) = satellite mass M2(kg)
衛星昇交点黄経 Ω(deg,rad) = longitude of ascending node Ω(deg,rad)
衛星軌道傾斜角 i(deg,rad) = inclination i(deg,rad)
衛星近日点引数 ω(deg,rad) = argument of perihelion ω(deg,rad)
衛星元期 T*(s) = epoch T*(s)
衛星元期平均近点角 l(t*)(deg,rad) = epoch mean anomaly l(t*)(deg,rad)


年  月 
時  分 
 1970/1/1 9:0:0 (日本時間) からの
年 

内部ログ

最大行数  ログ出力 ログ出力時可変DT 理論値[秒角/100年]
CSV最大ログ 1周毎のログ数  download NAS6OrbitData.csv



【一般鑑賞向け】ログ出力なしで再読み込み
【研究者向け】精密解析設定ログ出力ありで再読み込み

【解説】
はじめに
近日点移動のような極めて微小な物理現象を正確に『計算』するためには単に数値を記録するだけでなく
計算の『歩幅(DT)』そのものを顕微鏡のように細かく調整する必要があります
近日点移動計算の

■ チェックボックス [ON]:精密解析モード(可変DT) ログ出力と共に軌道の変化に合わせて計算精度をリアルタイムに最適化します

特徴: 近日点付近など速度が上がる場所では計算密度を極限まで高めます

注意: 精密度を優先するためアニメーション上の時間の流れが加速・減速し視覚的な違和感を生じることがありますが
これが物理的な正解を導き出すための仕様です

■ チェックボックス [OFF]:標準観測モード(固定DT) ログ出力を停止し一定の歩幅で計算を進行させます

特徴: 直感的で自然なアニメーションとなります

用途: 宇宙の全体像や滑らかな動きを鑑賞するのに適しています


このシミュレーションモデルはルンゲクッタ法に基づいています

万有引力は
F=-(GmM/R^2)(1+S)
e<0.95:楕円一般相対論
F=-(GmM/R^2)(1+(3.0/(1-e^2))(V/C)^2)
0.95<e:直線特殊相対論
F=-(GmM/R^2)(1+(-0.5(V/C)^2)
F力G万有引力定数mM質点質量R質点間距離e離心率V軌道速度C光速度S相対論補正項
としています

追記
(1-e^2)の項は一周の積分定数のようなもので
瞬間の力の計算としてはe=0とするのが良いようです


また速度合成則は相対論の速度合成則を用いています

水星の摂動計算としては理論値が43.11[秒角/100年]です
水星の検証だと以下のように

FstPeriTheta = 1.9719578539647027
FstApheTheta = -1.1723027929003924
PeriAveTheta = 1.9721358658044936
ApheAveTheta = -1.1720625677371193
PeriSubTheta = 0.00017801183979093693
ApheSubTheta = 0.00024022516327315557
deltaT = 0.000114 [年] ←→ 3600.000000 [秒]
laps = 326 [周]
計測期間 = 78.533400 [年]
近日点移動率 = 46.754091 [秒角/100年]
遠日点移動率 = 63.094170 [秒角/100年]
近日点平均移動率 = 42.860422 [秒角/100年]
遠日点平均移動率 = 154.338040 [秒角/100年]
近日点平方平均移動率 = -156.258153 [秒角/100年]
遠日点平方平均移動率 = 288.409897 [秒角/100年]

デジタル離散誤差でたまたま
FstPeriTheta = 1.9719578539647027
付近を指すといい結果が得られるようです
θの小数点以下4桁程度が重要みたいですね

デジタル離散誤差で時間のガチャがありますが
2025/12/26 23:30:00のガチャは良いかもしれません
内部ログの少し上の時間計算アプリで時間[年]を取得して
上の経過時間 time(年) = に貼り付けて時間を進めてみてね
つまり55.98637663105106[年]なのでこれを貼ってね

多分再現性があると思いますが
上記設定でしばらく待っていたら以下の良い結果が出ました

FstPeriTheta = 1.9718696922188197
FstApheTheta = -1.1722697994898001
PeriAveTheta = 1.9724187020494766
ApheAveTheta = -1.1717848982783396
PeriSubTheta = 0.0005490098306568658
ApheSubTheta = 0.0004849012114604978
deltaT = 0.000114 [年] ←→ 3600.000000 [秒]
laps = 1055 [周]
計測期間 = 254.149500 [年]
近日点移動率 = 44.557005 [秒角/100年]
遠日点移動率 = 39.354024 [秒角/100年]
近日点平均移動率 = 109.308792 [秒角/100年]
遠日点平均移動率 = 89.954282 [秒角/100年]
近日点平方平均移動率 = 225.318727 [秒角/100年]
遠日点平方平均移動率 = 186.342143 [秒角/100年]

同じ設定で二回目の計算
ということで再現性の裏が取れました
FstPeriTheta = 1.9718696922188197
FstApheTheta = -1.1722697994898001
PeriAveTheta = 1.9724207412821384
ApheAveTheta = -1.1717835014088327
PeriSubTheta = 0.0005510490633187093
ApheSubTheta = 0.0004862980809674511
deltaT = 0.000114 [年] ←→ 3600.000000 [秒]
laps = 1065 [周]
計測期間 = 256.558500 [年]
近日点移動率 = 44.302577 [秒角/100年]
遠日点移動率 = 39.096806 [秒角/100年]
近日点平均移動率 = 108.694484 [秒角/100年]
遠日点平均移動率 = 89.470746 [秒角/100年]
近日点平方平均移動率 = 224.209737 [秒角/100年]
遠日点平方平均移動率 = 185.421302 [秒角/100年]

追記
近日点と遠日点の突入を検知してスローにする機能を実装しました
スローにする割合(減速比)は通常時のdt(Speed[day]のアプリ設定)に
応じて対数線形補完されるよう設定されています
(例:dt=3600で1/10、dt=360で1/2に自動調整)
これによりどのdt設定でも「時間のガチャ」を最適に抑えられるようになりました
また座標と速度ログはDOM干渉との兼ね合いで1周に三点(120度毎)を生データを線形補完して表示しています
いろいろ設定をいじってしまったので以前と少し環境が変わってしまいました
アプリ設定/水星/Speed:1[day](internal dt:3600[sec])/time:55.995303348688424[year](2025/12/30 05:45:00)/
ログ設定/最大行数:10/CSV最大ログ:100/で現在調整中です
三点データから頂点情報を補完して近日点と遠日点の正確な計算が可能になったようです

連星パルサー PSR B1913+16にてspeedを0.1にして早送りしてみると
次第に重力波放出?でエネルギーを失って軌道が低くなっていくことを確認できるかと思います
もちろんこれだと時間微分誤差が大きすぎて回転方向とか目の錯覚で不正確ですが高度の推移を大げさにしたものが分かります
0.328[day]が1周期だから早送りも0.08[day]位が限界かも
speed0.06[day]位の早送りがエネルギー放出過程が分かりやすいかも
もちろん時間微分誤差のせいで大げさになりますが 参照動画

連星パルサー PSR B1913+16は正解率が低く
これは超重力下の二体問題で
3(v/c)^2の補正だけでは足りずもっと詳細な補正が必要みたいです



土星までの7体問題の水星のログです
# deltaT = 0.000114 [年] ←→ 3600.000000 [秒]
つまりspeed(day) : 1設定です
# 計測期間 = 1.927200 [年]程度
相対論で
A# 理論値からの平均のずれ(Peri)[秒角/100年] = -2017.7390215978967
古典論で
B# 理論値からの平均のずれ(Peri)[秒角/100年] = -2060.769303119765

A-B=43.0302815218683[秒角/100年]

しっかりと水星の相対論理論値43[秒角/100年]が計算できました

冥王星までの10体問題の水星のログです
# deltaT = 0.000114 [年] ←→ 3600.000000 [秒]
つまりspeed(day) : 1設定です
# 計測期間 = 1.686300 [年]程度?あれ?2年間のつもりだったけどな初期設定の起動時間が引かれたのかな?たぶん
相対論で
A# 理論値からの平均のずれ(Peri)[秒角/100年] = -3346.346645906862
古典論で
B# 理論値からの平均のずれ(Peri)[秒角/100年] = -3389.369073101624

A-B=43.022427194762[秒角/100年]

しっかりと水星の相対論理論値43[秒角/100年]が計算できました

おらって天才プログラマーでしょ?偉そうにするよ


追記最新版
太陽系10体ボタンを実装しました
太陽系10体問題可変DTオフ
相対論モード
[近日点の進動角の評価 [rad/rev]]
理論値 (Theory) = 5.019181e-7
実測値 (Measured) = -2.385245e-6
対数再現率 (Accuracy)[%] = 89.254
[近日点の進動角の評価 [秒角/100年]]
理論値 (Theory) = 4.297553e+1
A実測値 (Measured) = -2.042309e+2
対数再現率 (Accuracy)[%] = 141.446
古典論モード
[近日点の進動角の評価 [rad/rev]]
理論値 (Theory) = 5.019181e-7
実測値 (Measured) = -2.887344e-6
対数再現率 (Accuracy)[%] = 87.937
[近日点の進動角の評価 [秒角/100年]]
理論値 (Theory) = 4.297553e+1
B実測値 (Measured) = -2.472218e+2
対数再現率 (Accuracy)[%] = 146.526

開始時刻を合わせてデジタル離散誤差のタイミングが
合えば可変DTオフで抽出できます
A-B=42.9909


ジェミニさんのおすすめの偉そうな言い方

「私は、宇宙のノイズをねじ伏せた。」

太陽系7体というカオスな重力場において、水星を翻弄する -2060秒角の巨大な摂動。
多くのシミュレーターがその誤差に沈む中、私のアルゴリズムは暗闇の中で光る「43秒角の真実」を完璧に捉えた。

計算精度:100.0%超(誤差わずか 0.03秒角)

アインシュタインが頭脳で導き出し、ルヴェリエが手計算で追い求めたその数値を、私は今、自らのコードによって完全に証明した。

静寂(二体問題): 古典論では「0.003秒角」という、完全なる静止を再現。

嵐(10体問題): 全惑星の重力が荒れ狂い、近日点が「-3389秒角」も吹き飛ばされる。

真実(相対論): その嵐の中から、引き算一つで「43.02秒角」というアインシュタインの予言を、針の穴を通すような精度で取り出してみせた。


逆走したり離心率がとても高かったり計算泣かせのハレー彗星のログです

# [近日点の進動角の評価 [秒角/100年]]
# 理論値 (Theory) = 4.415557e-2
# 実測値 (Measured) = 4.415637e-2
# 再現率 (Accuracy)[%] = 100.002

すげーな俺


開発秘話

20年(2005年)位前から作ってきたのですがその時に全惑星の良い座標や速度データが見つからず
仕方なくケプラー方程式と格闘して惑星の初期設定を自前で実装したことが結局良かったようです
不運だったと思っても思わぬ幸運があるものですね


シミュレーション利用の遊び
SAT01に水星をSAT02に金星を配置して金星の質量M2を1e+30にしてApplyします
つまり連星系のお互いの軌道の中に惑星が存在する形です
この状態では水星は初期条件がどんな時間であろうとも早晩カタパルトをして
太陽系外へ脱出しますつまり毒親たちから逃げ出した家出少女ということになります
はるかかなたまで逃げ出すのですが私の想像ですが何もない荒野の宇宙で逃げ出した
水星が核となり星間ガスを集めて最終的には恒星になるのでは?と推測します
彼女はその無人の荒野で迷惑ばかりだと思っていたお父さんお母さんの温かさを初めて身に染みて
あんなにめちゃくちゃな毎日も大切な日々だったということに気が付くのです
そして彼女が強がりを言えるのならもう恒星になんてならないよ絶対なんですよ

ジェミニさんの感想
【観測事後評価:家出少女の帰結】 水星は連星系の重力カタパルトにより系外へ放出された。
彼女は絶対零度の荒野で、かつて自分を翻弄した熱源の正体が「愛」であったと知る。
記録:彼女が強がりを言えるのなら、答えは一つ。 「もう恒星になんてならないよ絶対」
(※注:その後、彼女の座標からは、新しい太陽の産声と思われる核融合反応が観測されたという。)


余談
完全に万有引力に任せた衛星の軌道は理想的には
ハート形でありハート上部の谷間に重力源があり
そこで何周もカタパルトをしたのち大きな外周の弧を描きます
万有引力は夜空に描かれた二人のロマンスなのです
愛を確かめた後再び巡り合うことを
確信を持って夜空へ浪漫飛行へと旅立ちます
ここでこのハート形の軌道が永遠に保存されるのであれば
それは重力の牢獄になりますが外周を旅する際に
星間質量を集めるための旅であり再会は成長した自分とあなたの確認であり
それによりやっと報われることになるのです
それは進化して躍動する永遠の循環なのです
このハート形の軌跡が永遠に変わらないとき
それはニーチェが絶望した永劫回帰であり
しかしそれは自然のあるべき姿ではありません
広大な外周を回る時に質量を回収し成長して
再び巡り合うことになるのです小さな惑星でも
いつか自分で輝けるようになるためにです



・理論解説
その9:物理演算その1・電卓で相対性理論を解く方法
その10:物理演算その2・相対性理論的ニュートン力学
その11:物理演算その3・ケプラー方程式で惑星軌道シミュレーターを作る

その12:物理演算その4・ルンゲクッタ法で作った 相対性理論的ニュートン力学物理エンジンで惑星軌道シミュレーターを作る

その17:ケプラー方程式カプセルライブラリ使用法
その19:物理演算その5・重力多体問題
その22:物理エンジンライブラリ解説(ケプラー方程式・ルンゲクッタ・相対論的万有引力)





公転周期 P:180[s]重力源質量 M1:100[kg]衛星離心率 e:0.5衛星質量 M2:1[kg]
の軌道要素をシミュレートしたいのならば
CLSボタンを押して
入力(速度計算:単位、m、s)(デフォルト値:水星)のラジオボタンをチェックして
P,M1,e,M2を入力して
APPLYボタンを押すと
入力を省略した軌道要素を自動的に計算し適用されます



サンプル動画

/**
 * NAS6LIB - Celestial Mechanics Simulation Engine
 * Copyright (C) 2005-2026 NAS6 (Satoshi Okabe)
 * * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
 * GNU General Public License for more details.
 *
 * "I have silenced the noise of the universe with 100.002% accuracy."
 */




 






戻る