最速降下曲線 (1回目)

2023/9/12(火)
最速降下曲線 (1回目)
 
(Brachistochrone curve)
 
 
■ 問題
重力加速度gのもとで(下向きを正とする)
点A(0, 0)から点B(w, h)に降下するとき
質点P(x, y)の最速降下曲線を求めよ
 
ベルヌーイ(Bernoulli)が問い
ニュートンはpm4:00に仕事が終わり
問題を見てam4:00に解き終わった
 
 
■ 解法
▼ 総時間t
距離
ds = √(dx² + dy²) = dx√(1 + dy²/dx²)
より
ds/dx = √(1 + dy²/dx²)
dsを進む時間をdtとすると
速度
v = ds/dt
= (ds/dx)(dx/dt) … (x方向の速度をs方向に補正)
= (dx/dt)√(1 + dy²/dx²)
よって
dt = √(1 + dy²/dx²) / v dx
= √{(1 + y'²) / v²}dx
 
総時間
t = ∫dt = ∫√{(1 + y'²) / v²}dx
 
高さy=0地点の速度v=0より
(1/2)mv² + (-mgy) = 0
v = √(2gy)  … ⓪
を代入
t = ∫√{(1 + y'²) / v²}dx
= ∫√{(1 + y'²)/(2gy)} dx
 
総時間
t = ∫₀^w √{(1 + y'²)/(2gy)} dx
 
▼ 変分法
最速降下曲線からのズレを
δy(x)とすると
δy(0) = δy(w) = 0
δT(y, y') = T(y+δy, y'+δy') - T(y, y')
≒ (∂T(y,y')/∂y)δy + (∂T(y,y')/∂y')δy'
 
tの式の積分の中を次の様に置くと
t = ∫₀^w √{(1 + y'²)/(2gy)} dx
T(y, y') = √{(1 + y'²)/(2gy)}  … ①
t = ∫₀^w T(y, y') dx
なので
δt = ∫₀^w δT(y, y') dx
= ∫₀^w{(∂T/∂y)δy + (∂T/∂y')δy'}dx
ここで
h = (∂T/∂y)δy , f = (∂T/∂y') , g' = δy'
と置くと
g = δy , f' = (∂T/∂y')'
[fg]₀^w = [(∂T/∂y')δy]₀^w = 0
(δy(0) = δy(w) = 0より)
と
(fg)' = f'g + fg' , ∫(fg)' = ∫f'g + ∫fg'
∫fg' = fg - ∫f'g  … 部分積分
をつかって
 
δt = ∫(h + fg')dx = ∫hdx +  fg - ∫f'gdx
= fg + ∫(h - f'g)dx = ∫(h - f'g)dx
つまり
 
δt = ∫₀^w {(∂T/∂y)δy - (∂T/∂y')'δy} dx
= ∫₀^w {(∂T/∂y) - (d/dx)(∂T/∂y')}δy dx
ここで
δt = 0にするには積分の中が0であればよい
(∂T/∂y) - (d/dx)(∂T/∂y') = 0  … ②
(オイラー･ラグランジュ方程式)
 
▼ 式の変形1
T(y, y') = √{(1 + y'²)/(2gy)}  … ①
(∂T/∂y) - (d/dx)(∂T/∂y') = 0  … ②
を解く
 
ここでT(y, y')の全微分を求め
dT = (∂T/∂y)dy + (∂T/∂y')dy'
dxで割る
dT/dx = (∂T/∂y)(dy/dx) + (∂T/∂y')(dy'/dx)
(∂T/∂y)y' = dT/dx - (∂T/∂y')y"
を
②式にy'を掛けた式
y'(∂T/∂y) - y'(d/dx)(∂T/∂y') = 0
に代入して(∂T/∂y)y'を消す
dT/dx - (∂T/∂y')(dy'/dx) - y'(d/dx)(∂T/∂y') = 0
dxで積分する
 
部分積分
∫f'g = fg - ∫fg'を使って
∫(d/dx)(∂T/∂y')y'dx
= (∂T/∂y')y' - ∫(∂T/∂y')y"dx
より
T - ∫(∂T/∂y')y"dx - (∂T/∂y')y'
+ ∫(∂T/∂y')y"dx = C  … C:積分定数
T - y'(∂T/∂y') = C  … ③
 
▼ 式の変形2
T(y, y') = √{(1 + y'²)/(2gy)}  … ①
T - y'(∂T/∂y') = C  … ③
を解く
 
①式を③式に代入
∂T/∂y' = (1/2)2y'(1 + y'²)^-1/2/(2gy)
= y'/√{2gy(1 + y'²)}
より
T - y'(∂T/∂y')
= √{(1 + y'²)/(2gy)} - y'y'/√{2gy(1 + y'²)}
= {(1 + y'²) - y'²} / √{2gy(1 + y'²)}
= 1 / √{2gy(1 + y'²)} = C
これを両辺を2乗し逆数にする
2gy(1 + y'²) = 1/C² 
 
変形して定数1/(gC²)をAに置き直す
y(1 + y'²) = 1/(2gC²) = 2A
y'² = 2A/y - 1 = (2A - y)/y
y' = √{(2A - y)/y}  (0 ＜ y ≦ 2A)  … ④
 
▼ 微分方程式を解く
y' = √{(2A - y)/y}  (0 ＜ y ≦ 2A)  … ④
を解く
 
cos2θ = cos²θ-sin²θ = 2cos²θ-1 = 1-2sin²θ
2cos²θ = 1 + cos2θ
2sin²θ = 1 - cos2θ
sin2θ = 2sinθcosθ
 
ここで
y = A(1 - cosθ)
と置く
(y = 0 のとき x = 0, θ = 0)
 
④式に代入して
y' = √{(2A - y)/y}
= √{(2 - 1 + cosθ)/(1 - 1cosθ)}
= √{(1 + cosθ)/(1 - cosθ)}
= cos(θ/2) / sin(θ/2) = dy/dx
dy/dx = cos(θ/2) / sin(θ/2)  … ⑤
dy = {cos(θ/2) / sin(θ/2)}dx
と
y = A - Acosθ
dy/dθ = Asinθ = 2Asin(θ/2)cos(θ/2)
dy = 2Asin(θ/2)cos(θ/2)dθ
より
 
{cos(θ/2) / sin(θ/2)}dx = 2Asin(θ/2)cos(θ/2)dθ
dx = 2Asin²(θ/2)dθ = A(1 - cosθ)dθ
∫dx = ∫A(1 - cosθ)dθ
x = A(θ - sinθ) + D
x = 0のときθ=0なので D = 0より
 
x = A(θ - sinθ)
y = A(1  - cosθ)
(0 ≦ θ ≦ 2π)
 
 
▼ tを求める
dx = 2Asin²(θ/2)dθ = A(1 - cosθ)dθ = ydθ
より
dx = ydθ
sin²(θ/2) = y/(2A)
 
点P(x, y)について
 
速度は
⓪式より
v = √(2gy)  … ⓪
傾きは
⑤式より
dy/dx = cos(θ/2) / sin(θ/2)  … ⑤
よって速度ベクトルは
v = ( vsin(θ/2) , vcos(θ/2) )
距離dxの移動にかかる時間dtは
dt = dx / vsin(θ/2)
= dx/√{2gy²/(2A)}
= (dx/y)√(A/g)
 
t = √(A/g)∫(1/y)dx = √(A/g)∫₀^θdθ = θ√(A/g)
 
t = θ√(A/g)
θ = t√(g/A)
 
▼ θ,Aを求める1
h = 0の時
x = A(θ - sinθ)
y = A(1  - cosθ)
 
y = A(1  - cosθ) = h = 0
1  - cosθ = 0
cosθ = 1
θ = 0, 2π
x = A(2π – sin2π) = 2πA = w
 
θ_max = 2π
A = w/2π
 
▼ θ,Aを求める2(ニュートン法)
0 ≦ θ ≦ 2π
x = A(θ - sinθ)
y = A(1  - cosθ)
 
x = A(θ - sinθ) = w
y = A(1  - cosθ) = h
A = w/(θ - sinθ)
A = h/(1  - cosθ)
w/(θ - sinθ) = h/(1 - cosθ)
h/w = (1 - cosθ) / (θ - sinθ)
(1 - cosθ) / (θ - sinθ) - h/w = 0
 
ニュートン法
f(x) = (1 - cosx) / (x - sinx) - h/w
f'(x) = sinx/(x-sinx) - (1-cosx)(1-cosx)/(x-sinx)² 
= {sinx(x-sinx) - (1-cosx)²}/(x-sinx)² 
= (xsinx - sin²x - 1 - cos²x + 2cosx)/(x-sinx)² 
= (xsinx + 2cosx - 2)/(x-sinx)² 
 
x_n ≠ 0, 2π
x_n+1 = x_n - f(x_n)/f'(x_n)
 
θ_max = xを次の式に代入
A = h/(1  - cosθ)
 
 
■ 結果
▼ 前提
重力加速度gのもとで(下向きを正とする)
点A(0, 0)から点B(w, h)に降下するとき
質点P(x, y)の最速降下曲線を求める
 
▼ 軌道
0 ≦ θ ≦ 2π
x = A(θ - sinθ)
y = A(1  - cosθ)
θ = t√(g/A)
t = θ√(A/g)
 
θ_max = 2π, A = w/2π  (if h = 0)
 
f(θ) = (θ - sinθ) / (1  - cosθ) - h/w
f(θ_max) = 0の解, A = h/(1  - cosθ)  (if h > 0)
 
▼ ニュートン法(θ = xと置いて書いている)
f(x) = 0の解を求める
 
f(x) = (1 - cosx) / (x - sinx) - h/w
f'(x) = (xsinx + 2cosx - 2)/(x-sinx)² 
 
x_n ≠ 0, 2π
Δx = f(x_n)/f'(x_n)
x_n+1 = x_n – Δx
x = x_n (if Δx < ε)