ページビューの合計

2017年7月23日日曜日

最速降下曲線の導き方

阪大に毎週水曜日に非常勤で通っている。先週は学生さんから変分法の質問を受けた。オイラー方程式から、汎関数 f が位置 x に無関係な場合にベルトラミの公式を求めよ、という問題でした。熱心な小栗旬似の学生さんでサイト物理のかぎしっぽ変分法の項を勉強していた。講義終了後だったので、時間のないことも手伝って暫し考えてみたが、彼の疑問点には答えられなかった。それで、自宅に帰って上のサイトを見て、もう一度考えてみると簡単で、単に積の微分公式にすぎなかった。テキストの説明不足でしかない。私の頭の悪さを再確認した次第であった。学生さんごめんね。


それを解説して、ついでに最速降下曲線を変分法を用いて解くというのが今回の目的である。

関数 f(x,y,y') によって決まる汎関数

I[y]\equiv \int _{a}^{b}f(x,y,y')dx

の変分問題(停留値問題とも云う)を解くためには、つぎのオイラー方程式を解けば良い.

\displaystyle {\partial f\over \partial y}-{d\over dx}\Big({\partial f\over \partial y'}\Big)=0 \tag{0}


上の方程式(0)だね。今 f =f(y,y')  つまり f  が x に関係しない場合を考える。

y=y(x),   y'=y'(x)  と考えられるので 合成関数の微分公式より、

df/dx = (∂f/∂y)(dy/dx) + (∂f/∂y')(dy'/dx)   (1)

がいえる。これから

 (∂f/∂y)(dy/dx) = df/dx - (∂f/∂y')(dy'/dx)   (2)

オイラー方程式(0)に y'=(dy/dx) を掛けると


y' \Big({\partial f\over \partial y}-\frac{d}{dx}{\partial f\over \partial y'} \Big)= y'{\partial f\over \partial y}-y'\frac{d}{dx}{\partial f\over \partial y'}=0 \tag{3}


(2)  から (3)  を引くと (ていねいに掛け算記号×を書いて)


df/dx =  (∂f/∂y')×(dy'/dx) + y'×(d/dx)[∂f/∂y']   (a)

ところで積の微分公式より、

(d/dx)[y'×(∂f/∂y')] = (dy'/dx)×(∂f/∂y') + y'×(d/dx)[∂f/∂y'] (b)

なので、(a)と (b)から

(d/dx)[f - y'×(∂f/∂y')]  =  0 

となり、

f-y'\displaystyle {\partial f\over \partial y'}=\mathrm{const.} \tag{4}

がしたがう。


さて上の(4)を用いて、最速降下曲線がサイクロイドになることを証明しよう。

Wikipedia から最速降下曲線の定義を思い出そう。


Brachistochrone.png
最速降下曲線(さいそくこうかきょくせん 英: Brachistochrone curve

任意の2点間を結ぶ全ての曲線のうちで、曲線上に軌道を束縛された物体に対して重力 (に代表される保存力) のみが作用する仮定の下、物体が速度0でポテンシャルが高い方の点を出発してからもう一方の点に達するまでの所要時間がもっとも短いような曲線である。

最速降下曲線はサイクロイドである。AとBが与えられAがBよりも高いとき、Aを無限斜面で通り、またBも通りAとBの間で最大値をとらない上下逆のサイクロイドがひとつだけある。これが最速降下曲線である。

したがって最速降下曲線は物体の重さと重力定数の強さにはよらない。この問題は変分法の道具を使って解くことが出来る。

歴史的には次のような経緯がある。

ガリレオは1638年に著書"Two New Sciences"で、最速降下曲線は円弧であるとしたが誤りであった。
ヨハン・ベルヌーイは(以前に解析した等時曲線を参照して)この問題を解いた後、1696年6月に著書"Acta Eruditorum"で読者に対して問題を提示した。4人の数学者がこれに応じて解答した。アイザック・ニュートンヤコブ・ベルヌーイ(ヨハンの兄)、ゴットフリート・ライプニッツギヨーム・ド・ロピタルである。ロピタルを除く3人の解答は1697年の同じ版で出版された。
弟に対抗してヤコブ・ベルヌーイはより難しい最速降下曲線問題を作った。それを解いている間に新しい手法を開発し、それがレオンハルト・オイラーによって改良され後に変分法と呼ばれるものになった。

点 P1 から別の点 P2 に移動する時間 12 は、s をその移動する弧の長さ、

ν をその速度とすれば、(時間×速さ=移動距離)なので、

微小時間で考えると dt=ds/ν  と考えられるので、その時間和を考えることにより

積分


 t_(12)=int_(P_1)^(P_2)(ds)/v,
(1)

で与えられる。

ところで移動物体(上図の黒玉●)には重力以外の外力はかかっていないので,

落ちた距離  y  と速度  v  の間にはエネルギー保存則

 1/2mv^2=mgy,

が成り立っている。これから


 v=sqrt(2gy).


(3)

となる。さらにピタゴラスの定理より、

 ds=sqrt(dx^2+dy^2)=sqrt(1+y^('2))dx
(4)

がいえるので、代入して

t_(12)=int_(P_1)^(P_2)(sqrt(1+y^('2)))/(sqrt(2gy))dx
(5)
=int_(P_1)^(P_2)sqrt((1+y^('2))/(2gy))dx.
(6)


これは汎関数 f  を

 f=(1+y^('2))^(1/2)(2gy)^(-1/2).
(7)

とした時の変分問題に他ならない。良く知られているように、オイラー方程式は

 (partialf)/(partialy)-d/(dx)((partialf)/(partialy^'))=0.
(8)

の形にかける。物理のかぎしっぽ 変分法1 (Joh著)参照。

今の場合 f(y,y^') は、独立変数xを含んでいないので、

上で述べたことよりベルトラミの公式


 f-y^'(partialf)/(partialy^')=C.
(9)

が導かれる。 (7)の f を y'  で偏微分してやると、

 (partialf)/(partialy^')=y^'(1+y^('2))^(-1/2)(2gy)^(-1/2),
(10)

となる。(9)に(10)を代入して計算すると、(注:コピペ上次式も番号(9)です)

\sqrt {1+y'^{2}\over 2gy} - y' \Big({y' \over \sqrt{2gy(1+y'^{2})}}  \Big) =\frac{1}{\sqrt{2gy(1+y'^{2})}}= C  \tag{9}



となる。つまり


 1/(sqrt(2gy)sqrt(1+y^('2)))=C.
(11)

となり両辺二乗して整理すると


[1+((dy)/(dx))^2]y=1/(2gC^2)
(12)



なる微分方程式が得られる。(12)式の 定数を 2A と置き直して、

式 (12)を次のように変形する.(式番号は以下(11)からに変更する)

y' =\sqrt{ { 2A - y \over y } }  \tag{11}


この微分方程式を解けばサイクロイドが現われるのだが、

それを示すのは少しテクニックがいる。

ここで単純化を行う。いま,曲線の定義域としては y \geq 0 を考えよう.

実際は -y を考えることになっているのだが、それは後で注意する。

また y' は実数でないと困るので,これより 2A \geq y \geq 0 という条件が必要になる.

もう一つ,初期条件として \theta =0 のとき, y=0 という条件も加える.

出発点 P1 を原点とするわけだ。

このとき y を次のようにパラメーター表示することができる.

y=A-A\cos \theta \tag{12}


点 P1 と P2 が同じ高さならば  \theta =0 から θ=2π まで動くと考えるわけである。

式(12)を両 辺 θ で微分すれば次のようになる.

dy=A\sin \theta d\theta = 2A\cos \frac{\theta}{2} \sin \frac{\theta}{2} d\theta \tag{13}


この y のパラメーター表示を使うと,式(11)は次のように書き直すことが出来る.


y' &=\sqrt{ { 2A - y \over y } }  \\&=\sqrt{ { A + A\cos \theta \over A-A\cos \theta} }  \\&=\sqrt{ { \cos^{2} \frac{\theta}{2} \over \sin^{2} \frac{\theta}{2}} }  \\&=\frac{\cos \frac{\theta}{2}}{\sin \frac{\theta}{2}} \tag{14}


両辺に dx を掛けて,次のように書くことができる.
dy= \frac{\cos \frac{\theta}{2}}{\sin \frac{\theta}{2}} dx \tag{15}


(13)と(15)を連立して, dy を消去すれば x\theta の関係式を求めることが出来る.

つまり、

dx=2A\sin^{2} \frac{\theta}{2} d\theta =A(1-\cos \theta) d\theta

が成り立つ。

この式の両辺積分を積分して, x のパラメーター表示

x=A(\theta -\sin \theta) + D

がでてくる。ここで D は積分定数だが,いま初期条件として \theta=0 のとき x=0 を考えているので, D=0 が決まる.

結局,最速降下曲線の方程式はパラメーター表示で次のように求まった.

x= A(\theta - \sin \theta)  \tag{12}
y= A(1 - \cos \theta)  \tag{16}


co-cycloid-fig02.gif


重力は負の方向に働くので、実は初期条件の設定から -y を考えることになる。

したがって、実際の x , y によって表されるのは,次のような下向きのサイクロイドとなる.


cycloid-joh.gif


すなわち,最速降下曲線サイクロイドだということが分かりました.

面白いシュミレーション動画があります。赤線が最速降下曲線です。







いずれもサイクロイドの一部で、常に最短時間なのである。

以上、物理のかぎしっぽ最速降下曲線 と Weisstein, Eric W. "Brachistochrone Problem". MathWorld 、および www.mathcurve.com  からのパクリ編集でした。

今回はこれでおしまい。