Unityで万有引力による運動の確認(2) ～月も回してみた

前回、地球が太陽から受ける万有引力によって楕円軌道で運動することをUnityで再現しましたが、今回は、さらに月を加えてみたいと思います。そんなの、月と地球の間に万有引力を発生させればいいじゃないかと思うかもしれませんが（自分も最初はそうしました）、実は太陽による万有引力も重要な役割を持っています。今回もパラメータやC#スクリプトを公開します。

月の公転の条件
月の公転を実装
実行結果
実行結果～太陽による万有引力を省略した場合

月の公転の条件

冒頭で述べた通り、太陽の周りを公転している地球の周りを月が公転するためには、月に対して地球による万有引力と太陽による万有引力の2つを加える必要があります（前回同様に単純化したいので、地球と太陽には万有引力は与えません）。

月も、地球と一緒に太陽の周りを公転していると考えると、太陽による万有引力を合成するのは自然な考えですね。

f:id:mzmlab:20191116000132p:plain

月の公転を実装

太陽、地球、月それぞれのパラメータと、C#スクリプトを示します。

パラメータ

■万有引力定数
・太陽　Gs = 6.67E-06 Nm^2/kg^2 (実際の1E+05倍)
・地球　Ge = 6.67 Nm^2/kg^2 (実際の1E+11倍)

　※スケールがでたらめなので、万有引力定数は太陽に引かれる場合と
　　地球に引かれる場合とで分ける。

■質量
・太陽　Ms = 2E+07 kg
・地球　Me = 50 kg
・月　　Mm = 1 kg

■距離
・太陽と地球の距離 Rse = 20 m
・地球と月の距離　Rem = 5 m

■初速
・地球 Ve = sqrt(Gs×Ms / Rse) = 2.582634 m/s (太陽に対して垂直方向に与える)
・月　 Vm = sqrt(Ge×Me / Rem) = 8.167007 m/s (地球に対して垂直に与える)

　※等速円運動となるvを求める。(式の導出は前回記事参照)

Unityコンポーネント設定

■太陽

コンポーネント			設定	備考
Transform	Position	X	0 m
		Y	0 m
		Z	0 m
Rigidbody	Mass		2E+07 kg
Rigidbody	Use Gravity		チェックなし

■地球

コンポーネント			設定	備考
Transform	Position	X	20 m
		Y	0 m
		Z	0 m
Rigidbody	Mass		50 kg
	Use Gravity		チェックなし
C# Script	初速（Z方向）		2.582634 m/s (=Ve)	C#スクリプトのStart（）で設定する。
	力		太陽から受ける万有引力を計算して設定	C#スクリプトのFixed Update（）で設定する。

■月

コンポーネント			設定	備考
Transform	Position	X	25 m
		Y	0 m
		Z	0 m
Rigidbody	Mass		1 kg
	Use Gravity		チェックなし
C# Script	初速（Z方向）		10.74964 m/s (=Ve+Vm)	C#スクリプトのStart（）で設定する。地球から見た相対速度がVmとなるようにする。
	力		太陽と地球それぞれから受ける万有引力を計算し、合成した結果を設定	C#スクリプトのFixed Update（）で設定する。

C#スクリプト

月のオブジェクトにアタッチするC#スクリプトを示します。
地球については前回記事を参照ください。

public GameObject のsunとearthはUnity画面のInspectorビューに現れるので、太陽と地球のオブジェクトをドラッグ＆ドロップしてください。万有引力や初速もpublicにするとInspectorビューから変更できるので便利かもしれません。

using System.Collections;
using System.Collections.Generic;
using UnityEngine;

public class MoonBehaviourScript : MonoBehaviour {
    //オブジェクト
    public GameObject sun;
    public GameObject earth;
    //パラメータ
    private float G_s = 6.67E-6f;    //万有引力定数(太陽)
    private float G_e = 6.67f;       //万有引力定数(地球)
    private float initVelocityZ = 2.582634f + 8.167007f;    //初速
    private float M_m;      //月の質量
    private float M_s;    //太陽の質量
    private float M_e;    //地球の質量
    private float r_s;    //太陽と月の距離
    private float r_e;    //地球と月の距離
    private float f_s;    //太陽から受ける万有引力
    private float f_e;    //地球から受ける万有引力

    void Start()
    {
        //質量の情報を取得
        M_m = GetComponent<Rigidbody>().mass;
        M_s = sun.GetComponent<Rigidbody>().mass;
        M_e = earth.GetComponent<Rigidbody>().mass;

        //初速を与える
        Vector3 initVelocity = new Vector3(0f, 0f, initVelocityZ); ;
        GetComponent<Rigidbody>().velocity = initVelocity;
    }

    void FixedUpdate()
    {
        //月から太陽に向かうベクトルを計算
        Vector3 direction_s = sun.transform.position - GetComponent<Transform>().position;
        //月から地球に向かうベクトルを計算
        Vector3 direction_e = earth.transform.position - GetComponent<Transform>().position;

        //距離rを計算
        r_s = direction_s.magnitude;
        r_e = direction_e.magnitude;

        //単位ベクトルに変換（後で万有引力の方向として使う）
        direction_s.Normalize();
        direction_e.Normalize();

        //万有引力を計算
        f_s = G_s * M_m * M_s / (r_s * r_s);
        f_e = G_e * M_m * M_e / (r_e * r_e);

        //加速度を与える
        GetComponent<Rigidbody>().AddForce((f_s * direction_s) + (f_e * direction_e), ForceMode.Force);
    }
}