アップデート 2005.7.25
単振り子
振り子の長さをl、振り子の振れ角をxとすると、振り子は円周方向に沿って−m･g･sinxの力を受ける。したがって運動方程式は、連立微分方程式に直すと
m･l･x"＝−m･g･sinx
あるいは
x'＝y
y'＝−ω2sinx （ω2＝g/l）
となる。これを単振り子(Simple Pendulum)と言う。
微分方程式系
x'＝y
y'＝−ω2sinx
は保存系である。全エネルギーは
H(x,y)＝y2/2−ω2cos(x)＝C
である。y2について解くと
y2＝2{C＋ω2cos(x)}
である。
全エネルギーCを一定にしたとき、相平面に描かれる軌道をここでは等高線と言うことにする。
さてC＞ω2のときはC＋ω2cos(x)＝(Ｃ−ω2)＋ω2{1＋cos(x)}＞0であるので、yは0にはならない。このことより、yは常に正、または常に負である。相平面での等高線は無限に伸びた、x軸に関して対称な、2本の曲線である。y＞0のとき、x'＝y＞0であるから、回転角速度x'が常に正となる。この場合振り子は振動するのではなく、回転する。
−ω2＜C＜ω2のとき、y2＝2{C＋ω2cos(x)}=2ω2{C/ω2＋cos(x)}である。このことより、cos(x)＞−C/ω2の範囲でyが存在する。変数xは2nπを中心としたある範囲内にあり、相平面上では楕円に良く似た等高線(閉曲線)になる。振り子は回転ではなく、往復振動となる。
C＝ω2のとき、y2＝2ω2{1＋cos(x)}であるから、x＝(2n＋1)πでy＝0となる。この点P((2n＋1)π,0）を通る曲線をセパラトリックスという。この点Pを初期値とする解は、x'＝0、y'＝0になるため、変化しない。
C＝−ω2のとき、y2＝−2ω2{1−cos(x)}であるから、x＝2nπ、y＝0以外の解を持たない。すなわち等高線は、1点Q(2nπ,0）に退化する。この点Qを初期値とする解は、x'＝0、y'＝0になるため、変化しない。
C＜−ω2のとき解曲線は存在しない。
(3.12,0)、(1.5,0)、(5,0)、(−5,0)、(−3.5,0.2)、(−10,1.5)、(3.5,−0.2)、(10,−1.5)から出発する解を同じ相平面上に描いてみると、次の図の様になる。
微分方程式
x'＝f(t,x,y,p)
y'＝g(t,x,y,p)
の計算に、関数マクロを用いる方法もある。黄緑色で表したセル、C3からC7に名前を付ける。C3には、積分回数という意味で、stepsという名前を付け、計算数628を入力する。C3の左欄B3には付けた名前、C3の右の欄D3には、変数の意味を書いておく。セルB3とDとはメモであり、EXCELでの計算には無関係であるが、この様にしておくと、後からEXCELのシートを見るときに便利である。更に、後から見直すときに何を書いたか思い出すためには必要な情報である。C4には、計算終了時刻という意味で、t_endという名前を付け、計算の最終時刻37.699を入力する。この値は12πである。C5には、xの初期値という意味で、x0という名前を付け、変位xの初期値3.145を入力する。C6には、yの初期値という意味で、y0という名前を付け、速度ｙの初期値0を入力する。初期値として、x＝π、y＝0を代入すると、解曲線が得られないので、それに近い値を初期値とする。C7には、定数ωを表す、omegaという名前を付け、1を入力する。水色のセルB3からB7には付けた名前を書く。また水色のセルD3からD7には入力した数値の説明を書く。関数を用いたマクロにより計算すると、10行目以下の黄色のセルに結果が表示される。xとyについての散布図を描くと下の様な図形が得られる。セパラトリックスにとても近い曲線である。
名前の付け方は、「セルの選択」→「挿入」→「名前」→「定義」の順に操作する。
|
|
A
|
B
|
C
|
D
|
|
1
|
単振り子の方程式:x"=-ω*ω*sin(x)
|
|
2
|
pendulum
|
|
3
|
|
steps
|
628
|
計算数
|
|
4
|
|
t_end
|
37.699
|
最終時間
|
|
5
|
|
x0
|
3.1415
|
変位の初期値
|
|
6
|
|
y0
|
0
|
速度の初期値
|
|
7
|
|
omega
|
1
|
ω
|
|
8
|
|
|
|
|
|
9
|
|
時間t
|
変位x
|
速度y
|
|
10
|
|
0
|
3.1415
|
0
|
|
11
|
|
0.06003
|
3.1414998
|
-5.56534E-06
|
今までの方法と比べると、関数マクロを使う方法は、入力するセルの位置に気を使わなくてすむ。マクロも簡単で見やすくなる。一方で微分方程式を変える度にマクロの中にある関数を変更しなくてはならない。
問題
速度に比例する抵抗を受ける場合の単振り子
x"＋ax'＋ω2sinx＝0
の解曲線を描け。
参考文献：丹羽敏男、微分方程式と力学系の理論入門、遊星社