地磁気双極子の地軸からのずれをもとめてみる

地磁気双極子(地球の中心に棒磁石があると仮定した場合の棒磁石にあたるもの)は地軸(地球の回転軸)から若干ずれています。今回はそれからどれだけずれているのかを求めてみようという記事。といっても完全に理論値から求められるわけではなくIGRF(国際標準地球磁場)の2010年度版のデータを使って求めてみます。(というのを授業でやったのでまとめもかねて...)
地球の磁場はダイナモ理論によれば地球の外核の運動によって電磁誘導により地球磁場が発生しています。運動によって磁場が引き起こされているので当然年々変化します。そこでIGRFは5年ごとに地球磁場の観測データからもっとも地球磁場と合うようなデータを提供しています。ではどういう形式でデータが与えられているのかということですが、磁気ポテンシャルWを与えています。磁場を与えてもいいんでしょうけど、ベクトルになるとデータ量が多くて不便ということなのかな?具体的には磁気ポテンシャルWは球面調和関数を用いて
W=a\sum^{\infty}_{n=1}\sum^{m}_{m=0}(\frac{a}{r})^{n+1}(g^{m}_{n}\cos m\phi+h^{m}_{n}\sin m\phi)P^{m}{n}(\cos\theta)
と表示されるので、この関数の係数(ガウス係数と呼ばれる)であるg^{m}_{n}およびh^{m}_{n}をデータとして提供しています。といってもこれは理論式なので実際nは13で打ちとめてあるみたいです。ここでaは地球の平均半径、φは緯度、rは地球の中心からの距離、θは余緯度(赤道では無く極からの角度),Pはルジャンドル陪関数となっています。
といってもわかりにくいですね...球面調和関数というのは下図(ウィキペディアより引用)のようにちょっと変形した球体のようなもののことでn,mにおおじていろいろな形をとります。
f:id:ikaro1192:20120705220137p:plain
これの大きさや割合を変えて粘土のようにくっつけてあげればどんな複雑な球に近いものでも表せそうですよね。その割合を決めるのがガウス係数というわけです。このように球面調和関数の和で近似してあげることで磁場をより少ないデータで表現しようというわけです。ちなみにjpgはこれと似たような考え方でsinとcosと係数の組み合わせで画像を表現して圧縮するという方法(フーリエ級数)を利用しています。他にもいろいろな関数で展開するということは頻繁に行われます。対象におおじてよりよく近似できる関数を用いるわけですね。地球磁場はだいたい中心からの距離と角度に依存して決まるので今回は球面調和関数が用いられているというわけです。ルジャンドル陪関数というのはルジャンドル多項式という微分方程式の解です。説明すると長くなるので今回は省略します。水素原子とかみたいな球対称の物理現象を扱うときによく出てくるものなので今回使われるものそういう意味では納得できるでしょう。
 さて、前置きが長くなりました。そろそろ本題に入りましょう。上のWから磁気双極子の成分だけを取り出せればよさそうです。そして、双極子磁場というのはWのうちn=1の部分であらわされます。(図を見れば明らかでしょう)。そこで上の式を具体的に書き下して見ましょう。
W=\frac{a^3}{r^3}(g^{0}_{1}P^{0}_{1}(\cos\theta)+(g^{0}_{1}(\cos\phi)+h^{0}_{1}(\sin\phi))P^{1}_{1}(\cos\theta))
またルジャンドル陪関数は計算することができ
P^{0}_{1}(cos\theta)=\cos\theta
P^{1}_{1}(cos\theta)=\sin\theta
となります。
また後々のことを考えていったん極座標から直交座標に直したいので関係式
x=r\sin\theta \cos\phi
y=r\sin\theta \sin\phi
z=r\cos\theta
を使いましょう。これとルジャンドル陪関数を代入すると
W=\frac{a^3}{r^3}((g^{0}_{1}\frac{z}{r}+g^{0}_{1}\frac{x}{r}+h^{0}_{1}\frac{y}{r})
を得ます。ポテンシャルがこれで求まったので積分すれば磁場がもとまるわけですが、それだとある任意の地点での磁場なので都合が悪いです。そこで磁気モーメントを求めることを考えます。磁気モーメントとは負極から正極に向かう向きのベクトルで大きさは極子の磁化の大きさと距離をかけたものです。
 そんなことを考えつつ双極子の磁気ポテンシャルを理論的に求めてみましょう。
まず、磁場の大きさを定義しましょう。これは単位磁荷から受ける力として定義されます。またに関するクーロンの法則によれば距離rだけ離れた磁荷m_1m_2に働く力は真空の透磁率を\mu_0}として
F=\frac{\mu_{0}}{4\pi}\cdot \frac{m_1 m_2}{r^2}
と表せられるのでm_2=Mm_2=1とすると磁荷Mが感じる磁場Bは
B=\frac{\mu_{0}}{4\pi}\cdot \frac{M}{r^2}
と表示できることになります。
さて、磁気ポテンシャルと磁場(ここでは大きさではなくベクトル量であることに注意!)の関係は
\bf{B}=\bigtriangledown W
と表示されます。\bigtriangledown極座標で表示すると結構めんどくさい計算(微分のチェインルールと基底の関係式を使って地道に計算すれば出ます。がここでは書くのが大変なのでやりません。)をおこなって
\bigtriangledown=(\frac{\partial}{\partial r},\frac{1}{r\sin\theta}\frac{\partial}{\partial \phi},\frac{1}{r}\frac{\partial}{\partial \theta})
みたいにあらわされます。これをWに作用させてあげるとWはrのみにしか依存しないので結局
\bf{B}=(\frac{\partial W}{\partial r},0,0)
となります。そうすると大きさもr成分のみなので
B=\frac{\partial W}{\partial r}
となります。Wを求めるためにBをrで∞からrまで積分してあげれば
W=\frac{\mu_{0} M}{4\pi r}
となります。
さて今もとめたのは単極子(モノポール)のポテンシャルです。これを双極子のポテンシャルに拡張しましょう。図のようにとってあげると
f:id:ikaro1192:20120706161725p:plain
W=\frac{\mu_{0} M}{4\pi r_a}+\frac{\mu_{0} -M}{4\pi r_b}
=\frac{\mu_{0} M}{4\pi}(\frac{r_b-r_a}{r_b r_a})
となります。ここでm-M間の距離をrとし、dを微小量とするとd<<r,r_a,r_bなのでr//r_a//r_bと考えることができ
\angle Mma=Mbm=\theta
とし、二次の微小量以上の項を無視してこれを先ほどのWの式に代入すると
W=\frac{\mu_{0} M}{4\pi}(\frac{d\cos\theta}{r^2})
をえます。ここで磁気モーメント\bf{I}の大きさIはI=dMを導入すると
W=\frac{\mu_{0} I\cos\theta}{4\pi r^2}
となります。磁気モーメントの向きは負極(S極)から正極(N極)への向きであることに注意して、とm-M方向の単位ベクトル\bf{n}を用いて内積の形でWを表現すれば
W=\frac{\mu_{0} \bf{I}\cdot{\bf{n}}}{4\pi r^2}
となります。
\bf{I}=(I_x,I_y,I_z)
とおき、
\bf{n}=(\frac{x}{r},\frac{y}{r},\frac{z}{r})
なので結局
W=\frac{\mu_{0}}{4\pi}(I_x\frac{x}{r}+I_y\frac{y}{r}+I_z\frac{z}{r})
となります。
さて球面調和関数で表示した双極子磁場と理論的に求めた双極子磁場のポテンシャルは一致するはずです。
\frac{a^3}{r^3}((g^{0}_{1}\frac{z}{r}+g^{1}_{1}\frac{x}{r}+h^{1}_{1}\frac{y}{r})=\frac{\mu_{0}}{4\pi}(I_x\frac{x}{r}+I_y\frac{y}{r}+I_z\frac{z}{r})
あとは係数比較することにより
(I_x,I_y,I_x)=\bf{I}==\frac{\mu_{0} a^3}{4\pi}(g^{1}_{1},h^{1}_{1},g^{0}_{1})
を得ることができます。磁気モーメントが求まったのであとは計算するだけ。ということで
(g^{1}_{1}=-1585.9,h^{1}_{1}=4945.1,g^{0}_{1}=-29496.5)
を代入して
\theta=cos^{-1}(\frac{I_z}{I})
からθをもとめて170度と求めることができます。つまり北から170度の位置にN極があるということですね。ということで何とかできました。