VisualC++ と OpenGL を利用した仮想物理実験室
【2-10-1】ばね弾性力＋重力がある場合の運動方程式とアルゴリズム

文責：遠藤理平（2010年8月30日）カテゴリ：仮想物理実験室(325)

【2-9-1】ばね弾性力による運動の運動方程式と（２次元）アルゴリズムでは、ばね弾性力がある場合の運動方程式の導出と、コンピュータでシミュレーションを行うためのアルゴリズムの導出を行ないました。そして、【2-9-2】ばね弾性力による単振動運動（２次元）のシミュレーションでは、ばねが x-y 平面内で変位するばねの単振動のシミュレーションを行ないました。

本節では、xyz 空間内でばねが変位する場合をシミュレーションするためのアルゴリズムの導出をおこないます。座標(0,0,z_{poll})に置かれた動かない支柱に、自然長 L_0[m] のばねの片側を固定し、反対側をボールつけ、xyz 空間内に伸ばします。時刻 t[s] のときのボールの位置を r(t)[m] とします。

ばね弾性力＋重力による運動の変数

■ 時刻 t[s] ■ ボールの位置 r(t)=(x(t),y(t),z(t))[m]
■ ばねの自然長 L_0[m]
■ ばねの長さ L(t)[m]
■ ばねの変位 ΔL(t)[m]
■ ばね弾性力 f_k = (f_{kx},f_{ky},f_{kz})[N]
■ 重力 f_g = (0,0,f_f_{gz})[N]
■ 角度 Φ[rad]、θ[rad]

ボールに加わる力がばね弾性力だけでなく重力も加わるので、ばね弾性力には添え字の k を、重力には、添え字に g をつけて区別することにします。

コンピュータアルゴリズム導出までのステップ

【2-9-2】ばね弾性力による単振動運動（２次元）のシミュレーションと同様にすすめます。

１．ばねの長さ L(t)
２．ばねの変位 ΔL(t)
３．ばね弾性力 f_k
４．ばね弾性力の x,y,z 各成分 f_{kx},f_{ky},f_{kz}
５．ボールへ加わる力の x,y,z 各成分 F_x, F_y, F_z
６．ニュートン運動方程式から加速度の x,y,z 各成分 a_x(t),a_y(t),a_z(t)

ばね弾性力 f_k の導出

【2-8-1】ばね弾性力による運動の運動方程式（１次元）とアルゴリズムで導出したとおり、ばね弾性力 f[N] とばねの変位 ΔL[m] との関係は、フックの法則と呼ばれる次の関係式を満たします。

(2.10.1)

ばね弾性力 f[N] を導出するために、ばねの変位 ΔL[m] を得る必要があります。時刻 t[s] のときのばねの変位 ΔL(t) は、ばねの長さ L(t) からばねの自然長 L_0[m] を引いたものになります。ばねの長さ L(t) は、次の図の関係から求めることができます。

L(t) は、三平方の定理から、

(2.10.2)

ばねの変位 ΔL(t) は次のように決まります。

(2.10.3)

上式を(2.10.1)に代入することで、ばね弾性力 f[N] を求めることができます。

(2.10.4)

ついでに、角度Φも求まります。

(2.10.5)

ばね弾性力の x,y,z 各成分（f_{kx},f_{ky},f_{kz}）の導出

３次元の場合では、x,y,z 成分の一括で求めることができず、２段階で行ないます。本節でははじめに、z 成分（f_{kz}）と xy平面成分（f_{kxy}）の２つに分割後、さらに、xy平面成分（f_{kxy}）を x 成分（f_{kx}）と y 成分（f_{ky}）とに分解することで、全成分を導出します。

z 成分（f_{kz}）と xy平面成分（f_{kxy}）への分割

ばね弾性力 f_k と角度Φ が求まっているので、z 成分（f_{kz}）と xy平面成分（f_{kxy}）はそれぞれ、sinΦ, cosΦ を掛けることで得られます。

x 成分（f_{kx}）と y 成分（f_{ky}）への分割

x 成分（f_{kx}）と y 成分（f_{ky}）への分割のために必要なのは、角度θです。

上図から、sinθ, cosθが求まります。

(2.10.6)

xy平面成分（f_{kxy}）に、sinθ, cosθをそれぞれ掛けることで、x 成分（f_{kx}）と y 成分（f_{ky}）が得られます。

ここまで得られた式を代入することで、ばね弾性力の x,y,z 各成分（f_{kx},f_{ky},f_{kz}）の導出することができました。

(2.10.7)

ばね弾性力＋重力がある場合の運動方程式

ボールに、ばね弾性力と重力が働く場合、ボールに作用するすべての力 F[N] は、ばね弾性力の各成分（f_{kx},f_{ky},f_{kz}）と重力（0,0,f_{gz}）となります。x,y,z 成分それぞれに作用する力は、

(2.10.8)

となります。式(2.10.8)を【５日目】重力加速度、重力、ニュートンの運動方程式で導出したニュートンの運動方程式

(2.10.9)

に代入することで、空気抵抗がある場合のボールの運動方程式が得られます。

(2.10.10)

ばね弾性力が働く場合、ボールの加速度の各成分は、質量 m[kg]に反比例することがわかります。つまり上式は、「ばね弾性力による加速度は重いほど小さく、軽いほど大きくなる」ことを意味しています。これは、直感的に成り立っていることがわかります。

ばね弾性力による運動のアルゴリズムの導出

式(2.10.10)として得られたニュートンの運動方程式を【13日目】ばね弾性力がある場合の運動方程式と同様に、コンピュータでシミュレーションするためのアルゴリズム導出します。式(2.10.10)を差分で表すと、

(2.10.11)

となります。L_n は、支柱とボールとの距離です。

(2.10.12)

式(2.10.11)は、t = t_n[s] のときの加速度 a_n[m/s^2] は、そのときの速度 v_n[m/s]で得られることを表しています。得られた加速度と等加速度直線運動のアルゴリズムを用いることで、速度 v[m/s] と位置 x[m]とを決定することができます。

(2.10.12)

(2.10.13)

本節では、ばねの変位を x軸方向に限定しているので、初期位置 x_0, 初速度 v_{x0} を指定することで、様々な運動のシミュレーションを行うことができます。次節では、具体的なシミュレーションを行ないます。