扁平楕円体の重力ポテンシャルと重力の導出|地球の重力のモデル化

スポンサーリンク
ホーム » 天体力学 » 扁平楕円体の重力ポテンシャルと重力の導出|地球の重力のモデル化

今回は、扁平楕円体(断面が楕円となる回転体のこと)の重力ポテンシャル(=ある点での単位質量当たりのポテンシャルエネルギー)を導き、また、この扁平楕円体周辺の重力を導出します。

扁平楕円体の重力ポテンシャル

$r’$ を原点からの質量素片 $\diff M$ までの距離、$G$ を万有引力定数、$P_n$ をルジャンドル多項式とする。このとき、$r$ の位置における重力ポテンシャル $\psi$ を次のように表せる。

\begin{eqnarray}
\psi=-\ff{G}{r}\sum_{n=0}^{\infty}\int\left( \ff{r’}{r} \right)^nP_n(\cos\phi)\diff M\\
\end{eqnarray}

ただし、$\phi$ を $\B{r}$ と $\B{r’}$ の成す角とする。

扁平楕円体による重力

扁平楕円体が $z$ 軸について対称とする。扁平楕円体の質量を $M$、万有引力定数を $G$ として、位置 $\B{r}$ に置かれた単位質量に作用する重力 $\B{f}$ は次のように与えられる。

\begin{split}
\B{f}=-\ff{G}{r^2}\left\{M+\ff{3(I_z-I_e)}{2r^2}\left(1-\ff{5}{r^2} \right) \right\}\hat{\B{r}}-\ff{3G(I_z-I_e)z}{r^5}\B{k}
\end{split}

ただし、$\hat{\B{r}}$ を $\B{r}$ 方向の単位ベクトル、$\B{k}$ を $z$ 軸方向の単位ベクトル、$I_z$ を $z$ 軸周りの慣性モーメント、$I_e$ を $x,y$ 軸周りの慣性モーメントとする。

これらの導出を行うため、最初に扁平楕円体の重力ポテンシャルの導出を行います。

スポンサーリンク

扁平楕円体の重力ポテンシャルの導出

地球は自転による遠心力などにより、赤道がわずかに膨らんだ形状をしています。このように、その断面が楕円となる回転体のことを扁平楕円体と呼びます。

厳密には、地球はもっと複雑な形状をしていますが、近似的に扁平楕円体と考えて、この周囲に形成される重力ポテンシャル(=ある点での単位質量当たりのポテンシャルエネルギー)の導出を行います。

扁平楕円体の模式図

今、全質量が $M$ で不均一な質量分布を持つ扁平楕円体について考えます。そして、原点(質量中心)から距離 $r$ 離れた位置 $\RM{P}$ に単位質量が置かれているとします(=人工衛星などに相当)。

このとき、質量中心から $r’$ の距離にある、扁平楕円体内部の質量素片 $\diff M$ によって $\RM{P}$ に形成される重力ポテンシャル $\diff \psi$ を次のように表せます。ただし、$G$ を万有引力定数、$\rho$ を$\RM{P}$ と質量素片との距離とします。

\begin{split}
\diff\psi=-G\ff{\diff M}{\rho}
\end{split}

ここで、$\B{r}$ と $\B{r}’$ の成す角度を $\phi$ とします。(※太字はベクトルを表します)すると、余弦定理から、$\rho=\sqrt{r^2+r’^2-2rr’\cos\phi}$ とできるので、上式を

\begin{split}
\diff\psi=-G\ff{\diff M}{\sqrt{r^2+r’^2-2rr’\cos\phi}}
\end{split}

とできます。次に、分母を変形して、

\begin{split}
\ff{1}{\sqrt{r^2+r’^2-2rr’\cos\phi}}&=\ff{1}{r}\left\{1+\left(\ff{r’}{r} \right)^{2}-2\ff{r’}{r}\cos\phi \right\}^{-\ff{1}{2}}=(1+x)^{-\ff{1}{2}}
\end{split}

について考えます。ただし、$x=\DL{\left(\ff{r’}{r} \right)^{2}-2\ff{r’}{r}\cos\phi}$ として整理しました。これを二項定理を用いて展開すると、以下のようにできます。

\begin{split}
(1+x)^{-\ff{1}{2}}=1-\ff{1}{2}x+\ff{3}{8}x^2-\ff{15}{48}x^3+\cdots
\end{split}

これを考慮して $\DL{\ff{r’}{r}}$ を昇冪の順に展開すると

\begin{split}
\left\{1+\left(\ff{r’}{r} \right)^{2}-2\ff{r’}{r}\cos\phi \right\}^{-\ff{1}{2}}&=1+\ff{r’}{r}\cos\phi+\left( \ff{r’}{r} \right)^{2}\left(-\ff{1}{2}+\ff{3}{2}\cos^2\phi \right)\EE
&\quad+\left( \ff{r’}{r} \right)^{3}\left(-\ff{3}{2}\cos\phi+\ff{5}{2}\cos^3\phi \right)+\left( \ff{r’}{r} \right)^{4}\left(\ff{3}{8}-\ff{15}{4}\cos^2\phi+\ff{35}{8}\cos^4\phi \right)+\cdots
\end{split}

式の中身を良く見ると、右辺にルジャンドル多項式が含まれていることに気が付きます。したがって、

\begin{split}
\left\{1+\left(\ff{r’}{r} \right)^{2}-2\ff{r’}{r}\cos\phi \right\}^{-\ff{1}{2}}&=1+\ff{r’}{r}P_1(\cos\phi)+\left( \ff{r’}{r} \right)^2P_2(\cos\phi)\EE
&\quad+\left( \ff{r’}{r} \right)^3P_3(\cos\phi)+\left( \ff{r’}{r} \right)^4P_4(\cos\phi)\EE
&=\sum_{n=0}^{\infty}\left( \ff{r’}{r} \right)^nP_n(\cos\phi)
\end{split}

という関係が導けます。これらの結果より、最初に示した質量素片による微小な重力ポテンシャルを

\begin{split}
\diff\psi&=-G\ff{\diff M}{\sqrt{r^2+r’^2-2rr’\cos\phi}}\EE
&=-\ff{G}{r}\sum_{n=0}^{\infty}\left( \ff{r’}{r} \right)^nP_n(\cos\phi)\diff M
\end{split}

と表すことができます。これを以下のように積分すると、扁平楕円体全体が形成する重力ポテンシャル $\phi$ が得られます。

\begin{eqnarray}
\psi=-\ff{G}{r}\sum_{n=0}^{\infty}\int\left( \ff{r’}{r} \right)^nP_n(\cos\phi)\diff M\tag{1}
\end{eqnarray}

扁平楕円体の重力ポテンシャルと力学的形状係数

上で得られた扁平楕円体の重力ポテンシャルに登場する角度 $\phi$ は、単位球面上で図のような関係があるとします。

このとき、ルジャンドル多項式の加法定理より

\begin{split}
P_n(\cos\phi)&=P_n(\cos\q)P_n(\cos\q’)+2\sum_{k=1}^n\ff{(n-k)!}{(n+k)!}P_n^k(\cos\q)P_n^k(\cos\q’)\cos\{k(\lambda-\lambda’) \}
\end{split}

と変形できます。ただし、$P_n^k(s)$ をルジャンドル陪関数とします。これを式 $(1)$ に適用すると

\begin{split}
\psi&=-\sum_{n=0}^{\infty}\left[ \ff{P_n(\cos\q)}{r^{n+1}}\cdot G\int_Mr’^{n}P_n(\cos\q’)\diff m \right.\EE
&\quad+\sum_{k=1}^n\ff{P_n^k(\cos\q)}{r^{n+1}}\left\{2\ff{(n-k)!}{(n+k)!}G\int_M r’^{n}P_n^k(\cos\q’)\cos k\lambda’\cos k\lambda \diff m \right.\EE
&\quad+\left.\left. 2\ff{(n-k)!}{(n+k)!}G\int_M r’^{n}P_n^k(\cos\q’)\sin k\lambda’\sin k\lambda \diff m \right\}\right]
\end{split}

ここで、以下のように係数を定義すると、

$$
\left\{
\begin{split}
&A_{n0}=G\int_Mr’^{n}P_n(\cos\q’)\diff m \EE
&A_{nn}=2\ff{(n-k)!}{(n+k)!}G\int_M r’^{n}P_n^k(\cos\q’)\cos k\lambda’ \diff m \EE
&B_{nn}=2\ff{(n-k)!}{(n+k)!}G\int_M r’^{n}P_n^k(\cos\q’)\sin k\lambda’ \diff m
\end{split}
\right.
$$

上式を、

\begin{split}
\psi&=-\sum_{n=0}^{\infty}\left\{\ff{A_{n0}}{r^{n+1}}+\sum_{k=1}^n\ff{P_n^k(\cos\q)}{r^{n+1}}(A_{nk}\cos k\lambda+B_{nk}\sin k\lambda ) \right\}
\end{split}

とできます。

上の定数を用いてさらに、力学的形状係数と呼ばれる扁平楕円体の形状に関する定数、$J_n,C_{nm},S_{nm}$ を導入します。($r_e$ を扁平楕円体の赤道半径とします)

$$
\left\{
\begin{split}
&J_{n}=-\ff{A_{n0}}{GMr_e^n} \EE
&C_{nm}=\ff{A_{nm}}{GMr_e^n} \EE
&S_{nm}=\ff{B_{nm}}{GMr_e^n}
\end{split}
\right.
$$

さて、$A_{00}=GM,A_{11}=B_{11}=A_{21}=B_{21}=0$ となります。ゆえに、$\phi$ を

\begin{split}
\psi&=-\ff{GM}{r}\left\{1-\sum_{n=2}^{\infty}J_n\cdot \left( \ff{r_e}{r} \right)^nP_n(\cos\q)+ \sum_{n=2}^{\infty}\left( \ff{r_e}{r} \right)^n\sum_{k=1}^{\infty}P_n^k(\cos\q)(C_{nk}\cos k\lambda+S_{nk}\sin k\lambda) \right\}
\end{split}

と変形できます。こうして得られた重力ポテンシャルの式は、地球などの重力ポテンシャルの一般的な場合の表示として用いられます。

慣性モーメントによる扁平楕円体の重力ポテンシャルの表示

重力ポテンシャルのもう一つの表示方法として、慣性モーメントによる表示を導出します。まず、式 $(1)$ を第三項まで展開して、

\begin{split}
\psi&=-\ff{GM}{r}-\ff{G}{r^2}\int r’\cos\phi\diff M-\ff{G}{2r^3}\int r’^2(3\cos^2\phi-1)\diff M \EE
&\qquad-G\sum_{n=3}^{\infty}\ff{1}{r^{n+1}}\int r’^nP_n(\cos\phi)\diff M
\end{split}

$\cos\phi$ を具体的な形で表すことを考えます。ここで、$\B{r}=(x,y,z),\B{r}’=(x’,y’,z’)$ と置いて、各単位ベクトルを $\hat{\B{r}}=(\A,\beta,\gamma),\hat{\B{r}’}=(\A’,\beta’,\gamma’)$ とします。すると、

$$
\left\{
\begin{split}
&\A=\ff{x}{r},\beta=\ff{y}{r},\gamma=\ff{z}{r} \EE
&\A’=\ff{x’}{r’},\beta’=\ff{y’}{r’},\gamma’=\ff{z’}{r’}
\end{split}
\right.
$$

の関係で表せ、さらに、$\cos\phi=\hat{\B{r}}\cdot\hat{\B{r}’}$ であるので、

\begin{split}
\cos\phi&=\hat{\B{r}}\cdot\hat{\B{r}’}=\ff{1}{r’}(\A x’+\beta y’+\gamma z’)
\end{split}

と求められます。この結果を最初の式の第二項に適用すると、

\begin{split}
\ff{G}{2r^3}\int r’^2(3\cos^2\phi-1)&=-\ff{G}{r^2}\int (\A x’+\beta y’+\gamma z’)\diff M\EE
&=-\ff{G}{r^2}\left( \A \int x’\diff M+\beta \int y’\diff M+\gamma \int z’ \diff M \right)
\end{split}

とできます。ここで、連続体の重心の公式を思い出ましょう。今、重心の座標を $\bar{x},\bar{y},\bar{z}$ とすると、扁平楕円体の重心が原点と一致していることから、$(\bar{x},\bar{y},\bar{z})=(0,0,0)$ となります。ゆえに、

\begin{split}
\bar{x}=\ff{\int x’\diff M}{M}=0,\,\bar{y}=\ff{\int y’\diff M}{M}=0,\,\bar{z}=\ff{\int z’\diff M}{M}=0
\end{split}

となるので、

\begin{split}
-\ff{G}{r^2}\int (\A x’+\beta y’+\gamma z’)\diff M=0
\end{split}

となって、ゆえに第二項は $0$ となります。次に、第三項の計算を行います。第三項を展開すると、

\begin{split}
&-\ff{G}{2r^3}\left\{(3\A^2-1)\int x’^2\diff M+(3\beta^2-1)\int y’^2\diff M+(3\gamma^2-1)\int z’^2\diff M \right. \EE
&\qquad +\left.6\A\beta\int x’y’\diff M+6\beta\gamma\int y’z’\diff M+6\gamma\A\int z’x’\diff M \right\}
\end{split}

このとき、後半の項については

\begin{split}
\int x’y’\diff M=\int y’z’\diff M=\int z’x’\diff M=0
\end{split}

となります。したがって、$3\A^2-1=\DL{\ff{3x^2}{r^2}-1 }$ などの関係があることを用いて、

\begin{split}
&-\ff{G}{2r^3}\left[\ff{3}{r^2}\left\{\int (r^2-y^2-z^2)x’^2\diff M+\int (r^2-z^2-x^2)y’^2\diff M \right. \right.\EE
&\qquad+\left.\left.\int (r^2-x^2-y^2)z’^2\diff M\right\}-\int(x’^2+y’^2+z’^2)\diff M\right]\EE
=&-\ff{G}{2r^3}\left[\ff{3}{r^2}\left\{\int r^2(x’^2+y’^2+z’^2)\diff M-x^2\int (y’^2+z’^2)\diff M \right. \right.\EE
&\qquad-\left.\left.y^2\int (z’^2+x’^2)\diff M-z^2\int (x’^2+y’^2)\diff M\right\}-\int(x’^2+y’^2+z’^2)\diff M\right]
\end{split}

と変形できます。ところで、$x$ 軸、$y$ 軸、$z$ 軸周りでの扁平楕円体の慣性モーメントを考えると、それぞれ次のように表せます。

$$
\left\{
\begin{split}
I_x&=\int(y’^2+z’^2)\diff M\EE
I_y&=\int(z’^2+x’^2)\diff M\EE
I_z&=\int(x’^2+y’^2)\diff M
\end{split}
\right.
$$

これを用いると、

\begin{split}
&-\ff{G}{2r^3}\left\{2\int(x’^2+y’^2+z’^2)\diff M-\ff{3}{r^2}(I_xx^2+I_yy^2+I_zz^2) \right\}\EE
=&-\ff{G}{2r^3}\left\{(I_x+I_y+I_z)-3(I_x\A^2+I_y\beta^2+I_z\gamma^2) \right\}\EE
=&-\ff{G}{2r^3}(I_x+I_y+I_z-3I)
\end{split}

が得られます。ただし、$I=I_x\A^2+I_y\beta^2+I_z\gamma^2$ としました。以上より、慣性モーメントによる重力ポテンシャル $\phi$ の表示式を、以下のように表せます。

\begin{split}
\psi&=-\ff{GM}{r}-\ff{G}{2r^3}(I_x+I_y+I_z-3I)-G\sum_{n=3}^{\infty}\ff{1}{r^{n+1}}\int r’^nP_n(\cos\phi)\diff M
\end{split}

扁平楕円体による重力の導出

上で導いた、慣性モーメントにより表した扁平楕円体の重力ポテンシャルから、扁平楕円体の重力を導出します。

まず、$3$ 次以上の高次のポテンシャルを無視して、

\begin{eqnarray}
\psi&=-\ff{GM}{r}-\ff{G}{2r^3}(I_x+I_y+I_z-3I)\tag{2}
\end{eqnarray}

さらに、地球などの天体では、その形状が $z$ 軸について対称と見なせるため、$I_x=I_y$ と近似できることに注目します。すると、$I_e=I_x=I_y$ として、$I$ を再計算でき、

\begin{split}
I&=I_x\A^2+I_y\beta^2+I_z\gamma^2\EE
&=\ff{I_e(x^2+y^2)+I_zz^2}{r^2}\EE
&=\ff{I_e(x^2+y^2+z^2)-I_ez^2+I_zz^2}{r^2}\EE
&=I_e+(I_z-I_e)\ff{z^2}{r^2}
\end{split}

これを式$(2)$に戻して、

\begin{eqnarray}
\psi&=-\ff{GM}{r}-\ff{G(I_z-I_e)}{2r^3}+\ff{3G(I_z-I_e)}{2r^5}z^2
\end{eqnarray}

とできます。

準備が整ったので、単位質量当たりの重力 $\B{f}$ がようやく計算できます。すなわち、グラディエントを用いると、重力ポテンシャルと重力の関係が以下のように表せ、さらにグラディエントの性質を考慮しつつ計算を進めると、

\begin{split}
\B{f}&=\nabla \psi\EE
&=\nabla\left\{-\ff{GM}{r}-\ff{G(I_z-I_e)}{2r^3}+\ff{3G(I_z-I_e)}{2r^5}z^2 \right\}\EE
&=-\ff{G}{r^2}\left\{M+\ff{3(I_z-I_e)}{2r^2}\left(1-\ff{5}{r^2} \right) \right\}\hat{\B{r}}-\ff{3G(I_z-I_e)z}{r^5}\B{k}
\end{split}

が得られます。ただし、$\hat{\B{r}}$ を $\B{r}$ 方向の単位ベクトル、$\B{k}$ を $z$ 軸方向の単位ベクトルとします。

タイトルとURLをコピーしました