Colloquium

NO.263
Weekend Mathematicsコロキウム室/NO.263

NO.1894    ダイクストラ法と領域認識    2010.12.26.  DDT

皆さん、ダイクストラ法ってご存知ですか?。知りませんよね。知ってたらマニア過ぎます^^。ダイクストラ法とは、道路ネットワークが与えられた時に、ある交差点から別の交差点にいたる最短経路を探索する標準的(古典的)方法です。主に交通計画分野にいる学生か講師陣しか知らないと思います。ダイクストラさんが1959年に提出しました。
道路ネットワークとは、図-1にようなものを指します。図-1は、札幌や京都のような碁盤の目の市街地に対応した、道路ネットワークです。

道路ネットワークでは、考える範囲の全ての交差点と、一つの交差点に隣接する全ての交差点のデータ、および隣り合う交差点間の道路距離が与えられます。
図-1のように規則性のある場合は、例えば点(交差点)Aから点Bにいたる最短経路(赤)は、(北に3ブロック,東に2ブロック)進む任意の経路とすぐわかりますが、ネットワークに規則性がないと、そうは行きません。いきおいトライandエラーの探索という事になります。
すぐ思いつくのは、道路ネットワークの全ての可能な通り方を多分木で全探索して経路長順に並べ、AからBにいたる経路だけを抽出し最短を選ぶというやり方です。一つの交差点に隣接する全ての交差点のデータと、隣り合う交差点間の道路距離が与えられるので、これは原理的には可能です。ところが少し概算すると、このやり方は現実的に実行不可能なのがわかります。
多分木全探索を図-1に試すとします。交差点は概ね3分木になっています。すなわち、ある道路から交差点に進入した場合、3通りの分かれ道というわけです。3通りの内の一つを選ぶと、また3通りの分かれ道が現れます。もし交差点が100個あるとすれば(縦横に道路10本)、超概算で3×3×・・・×3=3100の探索数になります。 3100は、ほぼ1050のオーダーです。前もやりましたが、一回の探索を1GHzのCPU一回の動作で完了できるとしても、1030年ほどかかる計算になります。もちろん1050個の中には、同じブロックをぐるぐる巡回したり、Aから出発しなかったり、Bに行かなかったり、通過経路に戻ったりするのもあるので、制約条件をつけて大幅に数を減らす事は可能です。例えば、Aから出発するものだけに限れば、1/100になるのは明らかです。
それでもその結果、例え探索数が100億分の1になったとしても焼け石に水です。1020年です。図-1に示したように、わざわざ遠回りする経路(青)も全て含まれるからです。交差点100個(縦横に道路10本程度)なんて、ちょっと大きな町なら持ってるはずです。
という訳でダイクストラさんは、最短経路探索を効率化するダイクストラ法を考えました。ダイクストラ法は帰納法で証明できますが、探索効率を上げるための手順が追加されていて、わかりにくいものになっています。そこで、ダイクストラ法ほど効率は良くない、次の基本ダイクストラ(仮称です)を考えます。


1.基本ダイクストラの証明
最初に出発点Aと目的点Bを決め、Aから行ける隣接点C1(1),C2(1),・・・,Ca(1)への経路長を計算します。それらをL1(1),L2(1),・・・,La(1)とすれば、これはAで分岐する道路の道路長そのものです。また任意のCi(1)が、Aからの最短経路上にあるのは明らかです。
 経路長Li(1)を、経路長リストL(1)={Li(1)}として記録します。
 経路Pi(1)=A→Ci(1)を、 経路リストP(1)={Pi(1)}として記録します。
次に、各Ci(1)から行けるCj(2)を全て集めて C1(2),C2(2),・・・, Cb(2)とし、それらへの経路長をLj(2)とします。 Lj(2)は、Li(1)にそれぞれの道路長を足すだけです。 経路Pj(2)=A→Ci(1)→Cj(2)も作り、 それぞれリストに加えます。

  L(2)={Li(1),Li(2)}
  P(2)={Pi(1),Pi(2)}

上記のL,Pは省略して書いてますが、 本当はL(2)={Li(1)}∪{Li(2)}などの意味です。
{Ci(2)}には、Ci(1)へ逆行するものも含まれますし、 逆行しなくても他のCj(2)に一致するものも出る可能性もあります。 ありますが、とりあえず気にしません。 重要なのは{Ci(2)}が、2回のステップでAから行ける全ての点を網羅している事です。
(2)に含まれる一つの経路Pi(k)に注目します(k=1,2)。 Pi(k)はCi(k)への、最短経路とは限りません。何故なら、k回以下のステップmでCi(k)達するPj(m)があって、 Li(k)よりLj(m)の方が短いかもしれないからです。
また(k+1)回以上のステップnでCi(k)達する、より短いPj(n)が、この先出るかもしれないからです。しかしL(k)の最小値Lp(k)に対応する、経路Pp(k)とその端点Cp(k)については違います。
Lp(k)はk回以下のステップで生成される経路長の最小なので、k回以下のステップに限れば、経路Pp(k)は、端点Cp(k)への最短経路です。 それより小さな経路長が、k回以下のステップにはないからです。 さらに、(k+1)回以上のステップでCp(k)に達する経路長は、 全てLp(k)より大きい事が言えます。 (k+1)回以上のステップで到達できる点は全て、k回以下で行ける点を通過しなければならないからです。 ここで、{Ci(k)}が、k回のステップでAから行ける全ての点を網羅している事と、Lp(k)がL(k)の最小値である事が効きます。 以上のステップを、Bに到達するまで繰り返します。 いつか必ず達するのは明らかです。
ステップの繰り返しの終了条件は、 経路Pp(k)の端点Cp(k)が、Cp(k)=Bとなる事です。

もう帰納法に持ち込む準備は整いましたよね?(整いました!・・・何か某謎掛け芸人みたいですね^^)。なので証明抜きで、次の基本定理(だと思います)をあげます。

[ダイクストラ法の基本定理]
 L(k)={Li(1),Li(2),・・・,Li(k)}の最小値Lp(k)に対応する経路Pp(k)は、 その端点Cp(k)への最短経路。ただし(k-1)以下のステップで最小値となったものは、最小値評価から除外する。

[証明]
 帰納法による。k=1の時は、P(1)={Pi(1)}が全て最短経路なので、帰納法は成立する。
[証明終]

この基本定理さえやってしまえば、後は探索をもっと効率化するだけです。基本ダイクストラの不経済な点は以下です。
  (B-1)逆行する。
  (B-2)ループする(以前のステップで探索した点に戻る)。
  (B-3)既に最短経路がわかった点も探索・評価する。
  (B-4)より短い経路があるのに、そちらに探索を絞らない。
  (B-5)一回のステップで、不要に遠い点まで探索する。

(B-1)〜(B-5)を解消すると、次のわかりにくいダイクストラ法になります。

[0]初期状態
  (1)出発点Aと目的点Bを決める。
  (2)Aから到達できる隣接点C1,C2,・・・,Csへの経路長と経路をリストする。
  (3)リストの最小値となる点Cmを見出す。経路A→Cmが、Cmへの最短経路なのは明らか。
  (4)Cmを探索点から除外する。
   (最短経路のわかった点は探索しない)
  (5)Cmを経路長と経路リストから除外し、別途用意した最短経路リストに、A→Cmを記録する。
   (最短経路のわかった点は評価しない)

[1]以後の手順
  (1)出発点をCmとする。
   (一回のステップで、不要に遠い点まで探索しない)
  (2)Cmから到達できる隣接点D1,D2,・・・,Dtへの、Aからの経路長と、Aからの経路をリストする。
   (差分を足すだけ)
    ただしD1,D2,・・・,Dtの中に、経路A→・・・→Cmの点は含まない。
   (逆行とループ禁止)
  (3)D1,D2,・・・,Dtの中に、Ciと一致するものがあれば、既存の経路長と現在の経路長を比較し、短かいほうで、Aからの経路長リストとAからの経路リストを置き換える。そうでないものは、たんにリストに追加する。
   (より短い経路を優先)
  (4)C1,C2,・・・,Cs,D1,D2,・・・,Duの経路長の中から、最小のDmを選ぶ。
   (基本定理)
  (5)Dmを探索点から除外する。
   (最短経路のわかった点は探索しない)
  (6)Dmを経路長と経路リストから除外し、最短経路リストに、A→・・・→Cm→Dmを記録する。
   (最短経路のわかった点は評価しない)

[2]終了条件
  (1)DmがBなら終了する(目的点に達した)。
  (2)探索点がないなら終了する(Bが最後に探索された)。
  (3)(1)でも(2)でもないなら、Cm=Dmとして[1]-(1)に戻る。
 基本ダイクストラは、一回の探索ステップごとにBへの最短経路か否かを判定できる点を除いて、じつは多分木全探索と同じです。しかし判定アルゴリズムがあるので、多分木探査の最小ステップでBへの最短経路が得られ、そこで計算を打ち切れます。これが、ダイクストラ法が効率的である基本的な理由です。  (B-1),(B-2)の解消は、最短経路探索ならば必ず満たすべき制約なので、これらで本質的に速くなる訳ではありませんが、多少速くなるのは確実です。  (B-3)の解消は、ダイクストラ法を本質的に速くします。(B-3)を解消すれば、一回の探索ステップごとに探索点が一個づつ減っていくので、必要探索量の上限が1/nづつ複利で減る事態を招きます。ここでnは、各点の分岐数の平均程度の整数です。多分木探索で探索数が、nm(m:探索点数)の複利で増えるのと逆の状況です。ダイクストラ法は、後半に行くほどスピードアップします。  (B-4),(B-5)の解消も。セットで本質的な加速です。アルゴリズムとしては、[0]-(3)または[2]-(3)から[1]-(1)への移行と、[1]-(3)のリスト更新として実装されます。  [1]で選択されたDm点への経路は常に、現ステップでの「最最短経路」です。そこに次の探索ステップを絞ります。そうするとダイクストラ法では、k回目のステップで、k回のステップで到達できる点を「全て網羅はしない」事になります。ただし、k回以下のステップでは、必ず網羅されているステップが存在します。例えばk=1は、明らかに網羅です。  現ステップでの「最最短経路」を選ぶので、「最最短経路」は、網羅されたステップの点の側から決まっていきます。結果として、A点からの近さの順に最短経路が決まります。これによって不要に遠い点まで探索する事が防がれます。BがAから10番目に近い点なら、10ステップで到達できる訳です。  このようなダイクストラ法の計算量はm2に比例します(証明できません)。ですが、個人的に比較できるアルゴリズムとしては、m元連立一次方程式を解く、ガウスの消去法があります。これの計算量は、m3に比例です。現行の標準市販PCのCPU、1GHz×quad Coreなら、m=10000であっても1000本の最短経路決定に、長くて数秒と予想できます。  ダイクストラ法は、「言われてみれば当たり前」の事をやっています。しかしそれを、ここまで整理して定式化したダイクストラさんは、やっぱりすごいと思います。

2.有限要素法
自分の本当の専門は、構造屋で橋梁屋です。それがどうしてダイクストラ法なんかを持ち出したかと言うと、久しぶりに有限要素法解析を頼まれたからです。指定された解析ソフトは馴染み深いもので、最初はそんなに気にしてませんでした。ところがその解析ソフトを久しぶりに使用したところ、恐ろしく不便なのに気づきました(切れかけました)。昔は不便さに、気づいていなかったのです。という訳で、有限要素法を概説させて頂きます(←全く理由になってない^^;)。

有限要素法が適用される典型的な状況は、例えば正方形の鉄の塊があって、それが自重でどのように変形するのかを計算したい場合です。図-2の碁盤の目で区切られた正方形が、今回の鉄の塊です。碁盤の目の事をメッシュと言い、目の一つ一つを要素と言います。
有限要素法の前提は、じつは一般的な物理的仮定と数学的原理にあります。いま正方形の変形状況を知りたいので、図-2にx-y座標を入れ、正方形の各点のx方向,y方向の変位dx,dyが、(x,y)の関数としてわかればOKです。

   

つまりd=(dx,dy)は、変形前と変形後のある点の座標の差d=(x+dx,y+dy)−(x,y)です。
有限要素法の最初の前提は、次の経験事実と数学的原理です。

  (1)物理量は、みな微分可能である。
  (2)微分可能なら、局所的に線形近似(比例配分)可能。これは微分の意味そのもの。

(1)より(dx,dy)は、微分可能です。そして「局所的に」なので、図-2の要素が十分小さいとして、(dx,dy)を要素上で考えると結論として、要素の変形に対して次の歪み分解が得られます。



局所的には要素の変形は、x,y方向の伸びに対応する変位、せん断変形に対応する変位、回転に対応する変位の、たんなる足し算になります。足し算なるのは、(dx,dy)が微分可能で局所的に線形近似可能だからです。
材料は変形に対して抵抗します。例えばx方向に伸びたら、ある力でx方向に縮もうとします(バネを想像して下さい)。通常の有限要素法では、ここで2つ目の物理的仮定の登場です。

  (3)変形は、微小である。

これは大抵成り立ちます。例えば何10万tというタンカーでも、自重で微小にしか変形せず自立できるので、建造中に乾ドックに置いておけます。
変形に応答する抵抗も物理量です。(1)より微分できます。しかも変形は微小です。(2)を言い換えると「微分可能なら、微小入力に対する応答は、入力に比例する」です。つまり歪み分解の各変形モードに対して、本当にバネを想定できるという事です。歪み分解の変形モードのうち、要素を本当に変形させているのは、最初の3つです。それで要素には、次のようなバネが内蔵されていると考えます。



歪み分解の図からわかるように、これらのバネの伸び縮みは、要素の四つ角の点の変位が支配します。バネ定数は、鉄なら鉄の変形にふさわしいように材質に合わせて調整し、実際の未知数は、要素の四つ角の変位という離散化量です。
以上の手順で、図-2の碁盤の目の目を全てバネモデルでおきかえ、物体への外力を与えると、バネを縮ませる(伸ばす)四つ角の変位に関する連立一次方程式が得られます。それを解いて全ての要素の四つ角変位を得たら、再び要素が小さい事を利用して((1),(2)を利用して)、四つ角変位の値から要素内部の点の変位を、線形近似で求めます。これで終了です。
有限要素法は(1),(2),(3)という一般的原理に基づいているので、ほとんどの場合で現象に相当するバネモデルが存在し、ほとんどあらゆる分野に適用可能な近似計算法です(宇宙論にも)。とは言え「メッシュを切らなきゃ始まらない」・・・事は、わかって頂けると思います。
物体などの形状(解析領域)が単純なら、図-2のような碁盤の目メッシュで済みますが、形状が複雑だとそうはいきません。じっさい有限要素法で一番手がかかるのは、実はこの「メッシュ切り作業」なのです。

3.領域認識

図-3は、とある山の断面です。山の形はラジコン機による航空測量の結果と思われます。図中のA,B,Cの領域では、土や岩の固さが違います。A,B,Cの境のラインは、ボーリング調査をもとに、地質屋さんが予想したものです(固さも)。真ん中辺りにある半円状の領域は、この山に掘削予定のトンネル範囲です。
トンネルを掘削したら、山の応力状態(内部の力)がどうなるか?を確認するために、有限要素法を行います。つまり落盤などが心配なのです。応力状態は変形量からバネモデルのバネの伸び縮みを計算し、伸び縮み量にバネ定数をかければ得られるので、結局やる事はさっきといっしょです。
ただし今度は、解析領域の形がやや複雑な事、材料定数(バネ定数)が領域ごとに異なる事から、最低限A,B,Cの領域を、解析ソフトに入力する必要があります。トンネル範囲がたくさんの領域に分かれているのは、地盤改良などの工法を併用し、地山の性質(バネ定数)が変わってくるからで、これらも入力する必要があります。その他のラインは、綺麗なメッシュを得るための補助ラインですが、それらによってさらに領域が増えます。

解析ソフトに領域を認識させた後は、領域ごとに、通常はAdvancing Front法を用いて半自動のメッシュ生成を行います。得られたメッシュは、図-4ですが、ここまで来るのが大変でした。
図-3の丸は幾何学的な点を表し、図-3の形状は点をつなぐラインの集積です。領域を定義するとは、そのラインの集積に位相情報(点のつながり方)を与えるという事です。基本は、領域境界の全ての点を、境界に沿う順番で手入力です。そうすると例えば、右側のBの領域の認識などは、けっこうな面倒地道な作業になります。しかも使用した馴染みのソフトは、点の認識が甘いという欠点を持ち、領域に目に見えない隙間ができる事も、しばしばでした。その上、馴染みのこの古いソフトは、既にWindowsの古典標準と化している「Undo機能」を持たない始末です。つまり隙間をなくする作業は、全て「最初からやり直し」なのです。
昔の領域認識作業は、紙の図面上でデシタイザーという機械を用い、虫眼鏡のターゲット内に一点一点を捉えて位置を記録し、位置の数値データをPCに食わせるという、今思えば、とても人間がやるような作業ではありませんでした。それを画面上のクリック一発で処理する、このソフトが現れた時には画期的で、不便さに気づいていなかったのです。

で、領域の自動認識はできないか?、と考えました。ダイクストラ法の出番です。図-5に示すように、点Aから出発して、点Aに戻るような最短経路を探せばいいじゃないか?、という発想です。それはダイクストラ法で、出発点と目的点を一致させれば可能です。
ただこの方法には、2つの欠点があるのは、最初からわかっていました。一つは島が残る事です。
各領域の分岐点(三叉路以上になっている交差点)から出発点に戻る探索を順次行うとすると、外周が短い領域から順番に認識されます。またある分岐点の分岐辺が属する領域が確定したなら、その点から出発する探索においては、確定辺は二度と探索するべきでありません。同じ領域が二重に生じるからです。
そうすると領域の大きさによっては、その周囲の領域が全て先に認識され、図-6の「島」と図示した領域については、探索が行われない事になります。
これへの解決策は、たんに無視すればOKだと気づきました。初期の状態でやれるところまでダイクストラをやった後で、2つの領域に属するラインと点を全て削除すれば、孤立した島だけが残ります。孤立した島では探索不要なので、それらの認識は一瞬です。
もう一つの欠点は、複合領域が出来うる事です。出発点に戻る最短経路が、必ずしも単一領域に対応しない事です。図-7はその例です。赤のルートは青のルートより短い、という結果になります。それで初回は赤のルートが認識され、複合領域になります。

可能性としては、赤と青が正確に等長という場合もあり得ます。その場合は、ダイクストラ法の最小値判定で、偶然に赤が選ばれた訳です。
複合領域になったら、後述の繰り返しに持ち込んで処理できるので、まずは無視してOKとわかりました。ただし複合領域判定は必要です。それには生成された位相情報を利用できます。
単一でも複合でも、とにかく領域として認識されたら、領域境界上の点の順序集合が得られます(右周りか左周り)。その中から三叉路以上になっている点を選び、分岐が境界上にない事を確認して、分岐に沿って1ライン進みます。図-7のC1とC2です。全ての三叉路以上の交差点で、C1のような内部点がなければ、単一領域です。
点の領域内部/外部/境界上の判定にも、位相情報を利用できます。点と境界点を結ぶ直線の角度を、順序集合に従って右回りか左回りかで順次計算し、角度増分を積算すれば結果は、外部:0,内部:±2π,境界:±πです。これは、順序集合に従って自分が境界点を目視していった場合、自分の首がどのように回るか想像すれば、すぐわかります。0,π,2πの分離は非常に明確なので、コンピューターの数値誤差を考慮しても、まず間違う事はありません。
以上の処理を完了すれば、残るのは孤立した単一島か複合島です。複合島があれば単一島は全て無視し、複合島の内部探索を行います。

複合島の外周は、その外側の(消してしまった)領域の確定に伴って、確定辺になっています。一方内部の辺は、まわりの領域が確定していないので、全て未確定辺です。よって今までのアルゴリズムをたんに繰り返し適用するだけで、図-8の赤で示したような探索が自動で始まり、複合領域は単一領域にいつかは分解されて最後に残るのは、もともとの孤立島と、複合島から発生した孤立島のみです(図-9)。
図-9の状態まで来たら、後は孤立島を全て認識させて終了・・・と思ったのですが、最後の問題がありました。

島の定義はこうです。
  「境界上の全ての点が、二分木である事」。
これは全ての点が交差点でない事を言っていますが、この条件に、解析領域の最外周が該当してしまいます。
解析領域の最外周を領域として認識させると、内部の全ての単一領域と重複が生じて、領域ごとにメッシュを切るという作業ができなくなります。
現在は削除したラインと点を復旧した後、全ての点を含む領域は不可、という方法で対処していますが、もっと完全な方法が必要だと、その後思いました。
というのは、ここまでの作業を実行するプログラムは、既にApplicationの態を成しています。このまでで終わらせるのは勿体ないので、領域ごとの完全自動メッシュまでやろうと思っています。その価値は、図-3にようにやたらと領域を増やさずに、必要最小限の内部領域だけで解析領域を定義し、何の手もかけずに自動メッシュできる事です。
図-3のような「図面」は、いくら製図ソフトが発達したとは言え、結局は手で描きます。また領域入力後に図-4のメッシュを発生させるため、領域境界の分割数を手動で調整するという、もう一山作業もありました。
完全自動の領域認識−メッシュ生成が遭遇する典型的状況は、図-10のようなものだと思います。つまり孤立領域が基本です。という事は、今まで外周しか持たない単一領域で考えて来ましたが、内周を持つ領域への対処が是非必要です。孤立が基本だからです。それには、それぞれを外周しか持たない単一領域として認識させ、重複を起こさせないために領域の包含関係を整理して、重複を持つ領域に内周を追加する必要があります。つまり領域の包含順序の階層化です。位相情報を利用して、領域間の含む/含まない関係はすぐ判定できますが、単純にやると多分木を使った経路探索のような目に合うような気がして、現在は方針を検討中です。

とは言え、ここまでの結果にはいちおう満足しています。最後に領域認識過程の図をあげます。




NO.1893     e と π について       2010.12.25.  水の流れ

第251回数学的な応募問題


先日同僚との会合で次のような問題を頂きました。 その場ですぐに発想が出てきませんでしたが、・・・

e(超越数)とπ(円周率)について

  eπとπ の大小を調べよ。

<水の流れ:ある関数のグラフを利用して解きましたが・・・他の解法があれば・・・>

注:この記事に関する投稿の掲載は、2011年1月10日以降とします。

NO.1892    変曲点を通る直線(2)    2010.12.25.  夜ふかしのつらいおじさん

この問題の要点は、関数の偶奇性だと思います。
偶関数のグラフはy軸対称、奇関数のグラフは原点対称になります。
問題1は2つの変曲点なので求める直線はx軸に平行、問題2は3つの変曲点なので求める直線は原点を通ります。
また、関数 のグラフをx軸方向にp、y軸方向にq平行移動すると、グラフの式は

  y−q=f(x−p)

となります。

問題1
y=x−4x−2x+12x+3 をx軸方向に−1平行移動すると、



となります。この関数は偶関数なのでグラフはy軸対称になります。





問題2



となります。
この関数は奇関数なのでグラフは原点対称になります。





NO.1891    電磁気学Minimum-9(エネルギー保存則-2)    2010.12.5.  DDT

 以後、NO.1744 ベクトル解析Minimum(1)NO.1749 ベクトル解析Minimum(2)NO.1752 ベクトル解析Minimum(3)NO.1757 ベクトル解析Minimum(4)の式は、無条件に使います。 今回も、荷電粒子と電磁場を含む系を扱い、エネルギー保存則を導く事が目標です。ここでエネルギーには、荷電粒子の(質点の)運動エネルギーと電磁場のエネルギーの両方が含まれます。前回は流体力学に寄り道して、エネルギー保存則の導き方に「当たり」を付けました。そこでの結論は以下です。

 (1) Maxwell方程式を、荷電粒子へ仕事をするように変形する。つまり相互作用項を特定する。
 (2) 荷電粒子を全て含む体積Vで、上記の変形を体積積分する。
 (3) divを含む式を考え出す(ガウスの発散定理の部分積分版と思われる)。
 (4) (3)を用いて、表面積分に直せるところは、全て直す。
 (5) 相互作用項を、荷電粒子の運動方程式を利用して、荷電粒子の物理量で書き直す。
 (6) エネルギー保存則の成立条件を考える。
 (7) 得られた関係式が、仕事の受給関係になっている事を確認する。
 (8) それが時間的に不変であれば、エネルギー保存則。

 今回は(1)〜(8)を実行しようと思うのですが、その前にやはり前回のように、荷電粒子−電磁場系の決定方程式を観察した方が良いような気がします。「そんな物いらない!」という声が聞こえてきそうですが、これは自習ノートなので、許して下さい。

1.荷電粒子−電磁場系の決定方程式
ここは前回の引用です。



ここにδは、荷電粒子の位置に特異点を持つデルタ関数。
   ε0:真空の誘電率
   μ0:真空の透磁率

2.決定方程式の観察
 冒頭の方針(1)より、相互作用項を特定するのが重要です。それは(N-1)の右辺の形をしています。 (N-1)と(M-1)〜(M-6)を結びつけるために、(N-1)の右辺をデルタ関数表示に戻し、全てのiについて和をとります。 Vは、全てのmを含む体積です。

     (1)

(1)の右辺を見ると、eEという項と、ev×Bという項が、Maxwell方程式側に必要です。
(M-4)より、

  

が得られますが、(1)の右辺は体積積分になっているので、これもそうします。

     (2)

同様に(M-3)から、

     (3)

(2),(3)より、

     (4)

となり、(4)の右辺=(1)の右辺です。これで、力レベルでの方針(1),(2)はやった事になります。 ここで前回の手順を思い出すと、(4)の左辺で、時間変動を示す量, Vの外に流れ出す量を見分ければ、運動量の関係式になるはずです。
 時間変動を示す量には、∂/∂tが付くはずです。そこで(4)の左辺3項目に注目すると、

     (5)

です。何故、

  

が成り立つかと言うと、ベクトルのスカラー微分の定義、

  

に戻って、外積に線形性がある事を思い出せば、すぐわかります。
 以後しばらくは、(4)の左辺だけ書きます。
 (5)を(4)に投入し、∂/∂tの部分を分けると、

     (6)

が得られます。(6)の形を整えるために、(M-1)を使います。(M-1)からは、

  

が導けるので(6)は、

     (7)

です。
 次は、Vの外に流れ出す量です。
divに関する関係式(ガウスの発散定理の部分積分版と思われる)をひねり出せ、という方針(3)の部分です。

     (8)

を、どう処理したものかと、けっこう悩みました。重要な関係式を忘れていたからです。  (8)の中身は、tを転置として、

     (9)

と書けます。No.1825を参照して下さい。そして(9)の行列部分を成分で書き下せば、 (9)の成立はあっと言う間です。  (∇H)Bに関する体積積分は、(M-5)よりB=μ0Hなので、 前回やったu(div・u)の体積積分に帰着されます。

     (10)

 (∇H)tBの体積積分も前回やりました。(∇u)tuに関する結果です。

     (11)

(9),(10),(11)より(8)は、

     (12)

となります。(rot×E)×Dの部分も、(M-6)よりD=ε0Eなので同じです。

     (13)

(12),(13)を(7)へ投入し、(M-2)からdiv・B=0を使うと、信じられない事に、

     (14)

が得られます。あれだけ不揃いだった、相互作用項(4)の左辺に関して、Maxwell方程式(M-1)〜(M-6)を使い切った瞬間に、電場と磁場に関して完全に対称(反対称)な(14)が現れます。はっきり言って、びっくりしました。  以上で力レベルでの目標は達成できたのですが、後の事を考えて、やり直します。ここまでは、運動方程式側から運動量関係式を導いた事になります。逆も可能なはずです。 その際、方針(1)を実行するのにポイントになる式は、(12),(13)だろうと、もう当りは付いています。  実際(12),(13)の辺々を足せば、

  

ですが、(M-2)よりdiv・B=0であり、右辺の対称性を保つために最後の項を移項すれば、

     (15)

です。
 方針で言えば、(15)の左辺が相互作用項だと仮定して、方針(4)まで一気にやった形です。
 次に方針(5)に進みます。(M-1)と(M-3)と(M-4)より、

  

なので、

  

が得られますが、外積の時間に関する微分公式(5)と、運動方程式(N-1)を用いて、

     (16)

     (17)

となります。
 最後に方針(1)に戻り、(15)の左辺が、相互作用項であった事を確認します。それは明らかです。(16)でLorents力が出てきたので、(15)左辺は、少なくとも相互作用項を含むものだったと言えます。関係式を導くのが目的なので、「少なくとも含めば」良いわけです。
 以上の過程では、数学的恒等式(5),(12),(13)と、Maxwellの方程式(M-1)〜(M-4),運動方程式(N-1)より導かれる、代入関係しか使っていません。従って、運動方程式(N-1)から出発しても同じ結果が得られるのは、保証されます。
 ところで、厳密な議論は次回にまわしますが、(17)って運動量保存則みたいなものですよね?。

3.荷電粒子−電磁場系のエネルギー保存則
 2.に習って、まず相互作用項を特定します。じつは前回、これはすでにやってあります。前回の2.(4)です。

     (1)

 (1)を眺めていると、今回の2.(4)と似ている事に気づきます。

     2.(4)

 2.(4)の左辺の中身とviとの内積を合計したものが、(1)の左辺です。2.(4)の右辺2項目と3項目が、 (1)の右辺に対応します。
(1)に2.(4)の右辺1項目に対応するものがないのは、 (vi×B)とviが直交するため、0になるからです。ところで2.(17)は、運動量保存則みたいなものでした。
 一般に、運動量とエネルギーの関係式は、項構成がパラレルに対応します。質点系の場合で言えば、

  

ですし、前回の流体の例で言えば、

  

です。何故かと言うと、エネルギーも運動量も、同じ決定方程式から(上述の例は、いずれも運動方程式から)出発して導くものだからです。ここに注意すると、(1)と2.(4)の対応より、「少なくとも」エネルギーの相互作用項を含むものは、2.(12),(13)からの類推で、

     (2)

みたいなものであろう事は、あっけないくらい簡単にわかります。(2)を方針(1)を満たすものと仮定すると、方針(2)は既に満たされているので、方針(3)に進みます。divの関係式をひねり出すところです。
 それがあるんですよ。奇跡のように・・・。

     (3)

が成り立ちます。(3)もまた、あっけないくらい簡単に示せます。(3)の右辺とは、

     (4)

の事です。Aを普通のベクトルとして、内積と外積の循環公式、

     (5)

もあります。(3),(4)より(4)は、∇・に関する積の微分です。 積の微分公式は、微分される側をそれぞれ定数として微分し、足せば得られます。(4)を成分に書き下し、成分ごとに積の微分公式を適用すると、 これはすぐわかります。内積「・」にも線形性があるからです。
 (4)でHを定数とすれば、もう積の微分の影響を考慮する必要はありません。 ∇は普通のベクトルのように扱えます。(5)より(4)は、

     (6)

です。逆にEを定数とすれば、

     (7)

が得られます。HもEも定数でない場合は、(6)と(7)を足せば良いので、結局(3)です。 従って(2)は、ガウスの発散定理を使って、

  

すなわち、

     (8)

となります。
 ここまでで思うところはないでしょうか?。ベクトル解析の公式群は、余りにも電磁気学と相性が良すぎると。

   ・まるで、電磁気学のために造られたみたいだ.

 そうなんです。現在のベクトル解析の公式群は、連続体力学(弾性学と流体力学)を扱う目的で、ガウスやコーシーやベルヌーイらが導入し、電磁場もその発想で扱かおうとして、グラスマン,ストークス,グリーンといった人達が発展させたものです。奇跡のように、電磁気学と相性の良いベクトル解析の公式が存在するのは、じつは当たり前なんです。
 ともあれ(8)で表面積分に直ってしまったので、方針(4)まで来た状態です。方針(5)に進みます。(M-3),(M-1)より、

  

が得られるので、これらを(8)に投入すると、

     (9)

が得られます。(M-5),(M-6)を使うと、

     (10)

     (11)

なので、(10),(11)から(9)は、

     (12)

となります。さらに(12)の左辺2項目は、運動方程式(N-1)のvi・を取れば、

  

を導けるので(vi・(vi×B)=0)、

     (13)

であり、(13)を(12)へ代入すれば、

     (14)

が得られます。
 (14)を導くのに本質的だったのは、(M-1),(M-3),(N-1)でした。今回も、電場,磁場,荷電粒子に対して、それぞれの決定方程式を使い切った瞬間に、電場,磁場,荷電粒子に関して完全に対称(反対称)な(14)が現れます。Maxwellの方程式って、本当に綺麗な関係だと思います。
 (14)右辺のベクトルE×Hを、Poynting(ポインティング)ベクトルと言います。ポインティングさんが、(14)を導いたからです。(14)の「意味」からPoyntingベクトルは、荷電粒子−電磁場系のエネルギーの流出速度を表すものだと考えられます(直後参照)。
 後は、方針(6),(7),(8)をやり、方針(1)を確認するだけです。エネルギー保存則は、孤立系に関して導けるものだという事に注意すると、荷電粒子miと電磁場(E,H)のみを対象とすべきだとなります(方針(6))。そうすると電磁場を発生し得るものは荷電粒子miしかないので、積分体積Vを無限に広く取ると、その境界では、荷電粒子の影響は微小になり、E=H=0と考えるのが自然です。従って(14)は、

     (15)

になります。
 (15)は仕事の受給関係を表します(方針(7))。何故なら(15)左辺第1項は、(12)の

     (16)

から生じたもので(従って、方針(1)は満たされる)、生じたものは、電磁場が荷電粒子に対して成した仕事だからです。(15)より、その関係は時間的に不変です(方針(8))。
 (15)は、荷電粒子−電磁場系のエネルギー保存則です。今度も使用したのは、数学的恒等式と、決定方程式から導かれる代入関係だけです。よって、運動方程式から出発しても、同じ結果が得られます。
 とにもかくにもここに到って、ようやっと次のように、カッコつけて言えるようになりました。

   ・古来より、質点系のエネルギーと供に保存される量は、みなエネルギーと呼ばれて来た.

 (15)より、

     (17)

が電磁場のエネルギーだ!、と言う事です。エネルギーと言って良い理由は、

   ・仕事の受給関係が時間的に不変である事を述べたものが、エネルギー保存則.

に尽きると思えます。

4.静電場,静磁場における確認
 3.(14)の形だと、静電場や静磁場でもエネルギー流出が起きているように見えるので、それは無い事を確認しておきます。静電場,静磁場の典型例は、電気スタンドの横に磁石を置いたような状況です。こんな系からエネルギーの漏れ出しが起きるようなら、電気スタンドも磁石もすぐに使えなくなりますよね?。
 ところで今まではどちらかと言うと、電荷密度ρと電流密度jを与えて、EとHはどうなるか?、という方向でやってきたと思いますが、今は逆です。静電場,静磁場という条件をEとHに与えて、どうなるか見たい訳です。これは電荷密度と電流密度の状態を、電場,磁場側からの条件によって、逆に決める事を意味します。
 電荷密度と電流密度は、ρ=覇δi,j=覇viδiと荷電粒子に起因するので、 それは荷電粒子の走行状態を決定する事になります。荷電粒子の走行状態を決めないと答えは出ないよ、と気づくまでしばらくかかった次第です。

 荷電粒子−電磁場系の方程式を再記します。Maxwellの方程式系の計算は、常に完全な組を手元に置いておかないと、失敗する危険大です。

  

 最初に、電荷保存則に関する注意です。これは数学的関係式なので、 Maxwellの方程式系に含めていませんが、電磁気学内部での位置付けを考えると、電荷保存則は荷電粒子(電子)の物性を決めます。 そういう意味では、やはり物理の式です。

     (1)

が電荷保存則ですが、

  

です。何故なら電子の素電荷eもδ関数も、時間に陽に依存しない関数だからです。2個の電子が合体して、2倍の電荷密度なんかになるわけないだろう!、電子が圧縮される訳ないだろう!、という「物性」です。よって荷電粒子だけを扱っている時は常に、

     (2)

が成り立ちます。
 エネルギー流出が0である事を言うには、3.(14)より、Poytingベクトルのdivの体積積分が0、を言えば良いですが。4.(3)より、

     (3)

なので、計算には右辺を利用できます。最初に、少なくとも静磁場(磁石のまわり)の場合。
 Maxwell方程式の第1の組(M-1)と(M-2)は、

     (4)

     (5)

とおくと、自動的に満たされます(Minimum-4)。いま静磁場なので、磁束密度Bの時間変動はありません。よって(4)のベクトルポテンシャルAも、時間に陽に依存しないと考えるのが妥当です。そうすると(5)から、

     (6)

となります。φはスカラーポテンシャルです。(6)を(M-4)に代入し、時間で偏微分すると、(1),(2)より、

  

が得られます。div・∇=Δなので、

     (7)

です。ラプラス方程式(7)が、電場を決定します。偏微分方程式(7)の解を完全に決定するには、境界条件が必要です。目的は(3)右辺の計算なので、Vの表面Sでの条件をとりたくなりますが、電場EはVを越えてどこまでも広がっています。しかもVは、全ての荷電粒子miを含むのでした。これは宇宙じゅうの荷電粒子がVの中にあると考えるのと同じです。よって、境界条件を考える表面は、Vを超える体積であれば、完全に自由です。無限遠をとります。−∇φとはEの事なので、無限遠ではE=0で、∂φ/∂tも0と考えるのが自然です。
 有界領域に閉じ込められた、渦なし流体のところで言いましたが、湧き出しも吸い込みもなくて、境界上に流れがなければ、そのラプラス方程式の解は0です。いま、φはスカラーポテンシャルなので、∂φ/∂tも渦なしです。△(∂φ/∂t)=div・∇(∂φ/∂t)=0なので、湧き出しも吸い込みもありません。無限遠の境界上ではとにかく何でも0なので、境界上の流れに相当するものも0のはずです。よってφは場所によらず、∂φ/∂t=0です。

     (8)

 (8)を(M-3)に持ち込み、HをベクトルポテンシャルAにおきかえて、時間で偏微分すると、

     (9)

が得られます。(9)が、磁場の状態を決定します。
 ところが静磁場なので、Aの時間変動はないのでした。従って、(9)の左辺は0です。

     (10)

 (10)の左辺の意味は、dj/dt(x,y,z,t)で、(x,y,z,t)は全ての場所と時刻です。従って、デルタ関数の特異点でも、dvi/dtの任意の時刻でも、dj/dt (x,y,z,t)=0となり、 dvi/dtは個々に0とおけます。荷電粒子は全て、等速直線運動します。電線ですよ、電線!。定常電流です。ビオ・サバールの実験ですよ!。
 Maxwell方程式を完全に使い切ったところで、荷電粒子の走行状態が決まりました。だからMaxwell方程式は、全部必要なんです。走行状態が決まったので、(3)の計算を行います。

  

 Vからのエネルギー流出なしです。副産物は、静磁場⇒静電場。

 少なくとも静電場(電気スタンドのまわり)の場合。
 さっきと同じく、スカラーポテンシャルφとベクトルポテンシャルAを使います。(5)より、

     (12)

ですが、静電場なので、(12)のEに時間変動はありません。よって∂D/∂t=0ですが、この条件を(M-4)に使っても、恒等式0=0が得られるだけです((2)のためです)。(12)をもう一回時間で微分すれば、Eの時間変動がないので、

     (13)

ですが、(13)の右辺はスカラーポテンシャル∂φ/∂tの勾配なので、渦なし場を定義します(rot×右辺=0)。 よって∂2A/∂t2は渦なしで、∂A/∂tも渦なしです。  この結果を、残った(M-3)に使います。∂D/∂t=0なので、

  

ですが、HをAにおきかえて時間で偏微分すれば、

     (13)

となりますが、∂A/∂tは渦なしです。rot×で0になります。  再び(10)が得られ、荷電粒子は等速直線運動します。電場,磁場の発生原因は、荷電粒子とその運動です。従って、今度も静磁場だと結論できます。結局最後に得られるものは、さっきと同じです。

  

 (14)〜(18)を使えば、(11)の計算が、そっくりそのまま再現されます。Vからのエネルギー流出なしです。副産物は、静電場⇒静磁場。

 結局、静電場⇔静磁場なので、磁石と電気スタンドは共存できるわけですが、こうなるのには理由があります。Minimum-6で述べたように、電磁場の時間発展は、Eの勾配=Hの時間変動→次の瞬間のHの勾配=Eの時間変動→次の次の瞬間のEの勾配、という連鎖になっているので、どこかで磁束か電束の時間変動が0になると、連鎖が途切れて、以後はずっと止まったまま、という訳です。

5.まとめ
 今回の感想は次の一言です。「答えを知ってたから、出来たようなもので・・・」。それで「ポインティングさん、よくこんな事考えつきましたね・・・」です。
 このまま運動量保存則へ行こうと思っていましたが、色々新しいベクトル解析の式も出しましたし、流体力学をやり直したおかげで、ランダウの言う対偶テンソルの価値もわかった気がします。また運動量保存則⇔作用反作用の法則⇔Maxwellの応力テンソルなので、テンソル記法がないと、運動量保存則は難しい気もします。
 次はベクトル解析Minimum-5として、今回の公式に関連し、テンソル記法と対偶テンソルをまとめさせて下さい。スローガンを言うと、「テンソル記法とは、添字と闘う技術だよ!」です。

[参考文献]
[1]理論電磁気学,砂川重信,紀伊国屋書店,1999年.
[2]ランダウ理論物理学教程2,場の古典論,ランダウ・リフシッツ,東京図書,1978年.


戻る