9.自由な座標

運動方程式を座標変換して、拘束に沿った自由な座標 θ で表すと、式(7)となる。例えば振り子の場合は式(9)。
この章のシミュレーション(クリックで計算開始)

自由な座標 θ での運動方程式(1)を求めたい

適当に曲げた針金を用意し、そこに物体を通した時、その物体の運動が知りたい。これは、これまで扱ってきた単純な拘束系であるが、1つ問題がある。これまでの章では、拘束条件が G(x)=0 という形の方程式で表されるとしていたわけだが、今の場合、適当に曲げた針金の形状を表す G を見つけることが困難なのである。
このような場合、右図のように、あるパラメータ θ によるパラメータ表示 x(θ) によって、曲線の形状を表すほうが簡単である。例えば2次元の振り子であれば、極座標の角度 θ が使えるし、一般の形状であれば、ベジェ曲線(以下の9.2節の【9.2-注1】)が使える。針金の場合、パラメータ θ は1つだけだが、曲面や多粒子の場合、このパラメータは複数個あるので、以後ベクトル表記 θ を用いることにする。 θ の成分の数を自由度という。 θ には拘束がかかっておらず自由な値をとれるので、 θ を、自由な座標と呼ぶことにする。(この用語は一般的ではない。一般化座標と呼ぶことが多いが、こちらは文字通りデカルト座標に限定しない一般の座標を考える際の用語であり、拘束条件とは無関係の文脈でも用いられる。)
自由な座標 θ から物体の位置 X (多粒子の場合も含めて大文字で表す)が求まるのだから、物体の運動を計算するには、 θ の時間変化が計算できれば十分である。即ち、前章までで扱ってきた運動方程式 M¨X=[] を、 ¨θ に対する方程式: ¨θ=[](1) の形に変形すればよい。これを、初期値 θ0,˙θ0 (任意の値をとれる)のもとで解けば、 θ(t) が求まり、 X(t) も分かることになる。
運動方程式を式(1)の形に変形するのは、純粋に数学的な問題である。この章では、「式(1)の導出」と「運動の計算」の2つの節に分けて議論を行う: 9.1:θ9.2:θ(t) この章からは、これまでのような、「力を求めれば運動方程式が確定して運動が計算できる」という素朴なニュートン力学の考え方から少し離れることになる。

9.1 自由な座標 θ での運動方程式

拘束条件 G(t,X)=0(2) のもとで、運動方程式は、以下の形で書ける: M¨X=F+GXTλFc(3) (拘束力 Fc については【8.2-注2】参照。 Fc が拘束面と「垂直」であることを示している。)この節では、上式を変形することで ¨θ=[] の形の運動方程式(7)を導く。

θ に対する運動方程式:式(7)

一般的に考えたいので、拘束が時間とともに変化する場合を含めることにする。その場合、座標 X は、時刻 t と自由な座標 θ の関数になる: X(t,θ)例えば、中心点が a(t) で動く振り子(第7章)の場合だと x(t,θ)=[lsinθlcosθ]a(t) となる。 θ が自由な座標であることにより、任意の t,θ に対して、 X(t,θ) は拘束条件(2)を満たす: G(t,X(t,θ))=0(4)
まずこの式を微分してみる。微小量 δt,δθ に対し、 G の変化 δG は、1次近似の範囲で 0 でなければならない: 0=δG(t,X(t,θ))Gtδt+GXδX(t,θ)Gtδt+GX(Xtδt+Xθδθ)=(Gt+GXXt)δt+GXXθδθ これが任意の δt,δθ で成り立つので、赤字部分と緑字部分はそれぞれ 0 となる: Gt+GXXt=0GXXθ=0(5) ただし、導出方法からも分かるように、第1式の Gt は、 Gt,X の関数と見た時の t での偏微分であるt,θ の関数と見た時ではない)
式(5)の第2式により、運動方程式(3)の両辺に Xθ の転置行列 XθT を左乗すれば、拘束力 Fc が消える: XθTM¨X=XθTF(6) 複雑な式となる拘束力 Fc (第8章の【8.2-注2】)が消えることは望ましい。しかも、この式の本数は、 θ の成分数と同じである。よって、上式は運動方程式 ¨θ=[] の形に変形することができるはずである。実際にそのような変形を行うと、以下の【9.1-注1】の式(7)のようになる。

【9.1-注1】自由な座標 θ に対する運動方程式:式(7)

座標 X が、時間 t と自由な座標 θ の関数 X(t,θ) で表せる時、 θ に関する運動方程式は以下のようになる: ¨θ=(XθTMXθ)1XθT[FM(ddtXθ)˙θMddtXt](7) ただし、 F は重力などの外力である。なんだか恐ろしい式だが、拘束条件が時間に依存しない場合、赤字部分は 0 である。

導出

運動方程式(3)と式(5)から、 ¨θ に対する式(6):(再掲) XθTM¨X=XθTF(8) が得られている。後は、 ¨X¨θ を用いたものに書き直せばよい: ¨X=ddtddtX(t,θ)=ddt(Xt+Xθ˙θ)=ddtXt+(ddtXθ)˙θ+Xθ¨θ これを、元の式(8)に代入して、 ¨θ=[] の形にすれば、式(7)になる。

9.2 物体の運動 θ(t) の計算

自由な座標 θ のもとで物体の運動 θ(t) を計算するには、任意の初期値 θ0,˙θ0 のもとで、運動方程式(7)を解けばよい。
この節では、まず、振り子と二重振り子の場合について、自由な座標 θ を使った形で書き下す。これらについては、以前の章で、デカルト座標での運動方程式を示している(振り子:第6章、二重振り子:第8章)。その後、冒頭で触れた任意形状の針金上の運動の計算を行う。

極座標での振り子の運動方程式:式(9)

2次元平面上の2次元振り子の自由な座標 θ を、右図のようにとる(おもりが垂れ下がった位置が θ=0。座標 x(θ) およびその微分は x(θ)=[xy]=[lsinθlcosθ]dxdθ=l[cosθsinθ]ddtdxdθ=l˙θ[sinθcosθ] となるので、これらを運動方程式(7):g=|g| ¨θ=(mdxdθ2)1dxdθT([0mg]m˙θddtdxdθ) に代入すれば、極座標での運動方程式が得られる: ¨θ=glsinθ(9) 数値計算は、既に示した第6章のデカルト座標でのものと、見かけ上変わらないので割愛する。
なお、3次元の場合の運動方程式を考えると、方位角を ϕ として、以下のようになる。 [¨θ¨ϕ]=⎢ ⎢glsinθ+˙ϕ2sinθcosθ2˙θ˙ϕcosθsinθ⎥ ⎥ 式(9)に比べて複雑な式になり、 θ=0 のところに特異性がある。一方、第6章のようにデカルト座標で表せば、2次元・3次元どちらでも同じ形になり、特異性もない。このように、自由な座標を使ったからといって必ずしも簡単になるわけではない。

極座標での二重振り子の運動方程式:式(10)

2次元平面上の二重振り子を考える。4次元2拘束なので、自由度は 42=2 である。その自由な座標 θ=[θ1θ2]T を、右図のようにとる。座標 X(θ) およびその微分は X(θ)=[x1x2]=⎢ ⎢ ⎢l1sinθ1l1cosθ1l1sinθ1+l2sinθ2l1cosθ1l2cosθ2⎥ ⎥ ⎥dXdθ=[Xθ1Xθ2]=⎢ ⎢ ⎢l1cosθ10l1sinθ10l1cosθ1l2cosθ2l1sinθ1l2sinθ2⎥ ⎥ ⎥ddtdXdθ=⎢ ⎢ ⎢ ⎢ ⎢l1˙θ1sinθ10l1˙θ1cosθ10l1˙θ1sinθ1l2˙θ2sinθ2l1˙θ1cosθ1l2˙θ2cosθ2⎥ ⎥ ⎥ ⎥ ⎥ となるので、これらを運動方程式(7):g=|g|
¨θ=(dXdθTMdXdθ)1dXdθT⎜ ⎜ ⎜⎢ ⎢ ⎢ ⎢0m1g0m2g⎥ ⎥ ⎥ ⎥M(ddtdXdθ)˙θ⎟ ⎟ ⎟ に代入すれば、運動方程式が得られる: ¨θ=[1cosθ12l1rcosθ12l1rmr]mrcos2θ12[mrsinθ1lr˙θ22sinθ2˙θ21][l11gsinθ12]lr=l1l2,mr=m1+m2m2,θ12=θ1θ2(10) 数値計算は、既に示した第8章のデカルト座標でのものと、見かけ上変わらないので割愛する。

ベジェ曲線上のおもりの運動方程式

任意の形状の曲線を拘束条件とする場合には、ベジェ曲線を用いることができる(以下の【9.2-注1】)。外力として重力がかかっているとする。運動方程式を求めるには、 x(θ) の微分:(同注の式(11)の微分) dxdθ=3(1θ)2p1+3(13θ)(1θ)h1+3θ(23θ)h2+3θ2p2ddtdxdθ=6˙θ[(1θ)p1+(2+3θ)h1+(13θ)h2+θp2] を、運動方程式(7):g=|g| ¨θ=dxdθT([0mg]m˙θddtdxdθ)mdxdθ2 に代入すればよい。特にきれいになるわけではないので、結果は割愛する。
数値計算を行うと、右図のようになる。グリッドの間隔は 1m である。これ以外の曲線の例については以下の数値計算を参照。

【9.2-注1】3次のベジェ曲線

ベジェ曲線は、右図のように、2つの任意の点 p1,p2 の間を結ぶ曲線を表すために用いられる。具体的には、曲線の形状をコントロールするためのハンドルと呼ばれる2つの点 h1,h2 を用いて、以下のように定義される:x(0)=p1,x(1)=p2 x(θ)=(1θ)3p1+3θ(1θ)2h1+3θ2(1θ)h2+θ3p2(11)

補足

複雑な曲線を補間する際には、曲線を分割し、各区間をベジェ曲線で補間するという方法がとられる。右図のように、2つのベジェ曲線が点 p2=p1 においてなめらかに接続する条件は、 h2h1 が反平行になることである。