VisualC++ と OpenGL を利用した仮想物理実験室
【2-10-2】ばね弾性力+重力による単振動運動(3次元)のシミュレーション
(2-10-1)【2-10-1】ばね弾性力+重力による運動の運動方程式(3次元)とアルゴリズムでは、 ばねの変位が xyz平面内の場合で、ばね弾性力と重力がボールに作用する場合の運動方程式の導出と、コンピュータでシミュレーションを行うためのアルゴリズムの導出を行ないました。 t=t_n[s]のときの加速度 a_n[m/s^2]は、
となり、そのときのボールの位置 r_n = (x_n, y_n, z_n) によって決まる値となります。 L_n は、原点からボールまでの距離、
です。コンピュータシミュレーションでは、 初期条件として、l_0 = 20.0[m]とし、10個のボールの質量 m [kg]を変化させて運動の違いを確認します。
ばね弾性力+重力による運動のための計算アルゴリズム
支柱の位置とボールの位置・速度・加速度の初期値
ばね定数を k、自然長 L[m] のばねの片側を固定した動かない支柱を原点に置き、ばねの反対側をボールつけ、x-y 平面上に伸ばします。時刻 t=t_n[s] のときのボールの位置を r_n とします。
k = 10.0; //ばね定数 L = 30.0; //ばねの自然長 r = 20.0; //ボールの初期半径 g =9.8 ball_r = 4.0;//ボールの半径 //支柱の定義 poll_x = 0.0; poll_y = 0.0; poll_z = 60.0; poll_r = 5.0; //半径 poll_h = 20.0; //高さ
速度・位置を逐次計算するアルゴリズム
//位置の計算 x = x + vx * dt; y = y + vy * dt; z = z + vz * dt; l = sqrt(pow(ball[i].x,2)+pow(ball[i].y,2)+pow(poll_z-ball[i].z,2)); //<-----支柱からボールまでの距離 //加速度の計算 ax = - k/ball[i].m * (1.0 - L/l) * ball[i].x; //<-----ばね弾性力からの寄与 ay = - k/ball[i].m * (1.0 - L/l) * ball[i].y; //<-----ばね弾性力からの寄与 az = - k/ball[i].m * (1.0 - L/l) * (ball[i].z-poll_z) -g; //速度の計算 vx = vx + ax * dt; vy = vy + ay * dt; vz = vz + az * dt;
このアルゴリズムでは誤差が大きく、時間と共に、単振動の振幅が大きくなってしまいます。 アルゴリズムの改善は、次に行ないます。
ばね弾性力による単振動のシミュレーション
動かない支柱にばねで繋がったボールが運動するシミュレーションです。運動は、ばね弾性力と重力による運動となります。
VisualC++ + OpenGL プログラミング
「【0.2日目】仮想物理実験室の構築 (ver1.1)」をベースにプログラミングします。上記で解説した逐次計算を行うアルゴリズムをプログラミングします。
ばねの定義
ばねの定義については、ばね弾性力による運動の運動方程式(1次元)とアルゴリズムをご覧下さい。
仮想物理実験室変数の定義
////////////////////////////////////////////////////////////////////////// // 変数の定義 ////////////////////////////////////////////////////////////////////////// #define _BITMAP 0 //アニメーション作成用ビットマップの保存 0:しない 1:する ofstream ofs1( "ball1.data" ); //-------------------------------------------------------- // 仮想物理実験室変数の定義 //-------------------------------------------------------- double t = 0.0; //時刻 double dt= 0.05; //時間刻み int tn = 0; //ステップ数 double k = 10.0; //ばね定数 double L = 30.0; //ばねの自然長 double l = 20.0; //ばねの初期長 double g =9.8; //重力加速度 double radian, radian_phi, radian_theta; double theta; //支柱の定義 double poll_x = 0.0; double poll_y = 0.0; double poll_z = 60.0; double poll_r = 5.0; //半径 double poll_h = 20.0; //高さ // ボールの定義 double ball_r = 4.0;//ボールの半径 struct BALL { double m; double x, y, z; double vx, vy, vz; double ax, ay, az; }; const int N = 10; BALL ball[N]; void SetUp(void){ for(int i=0; i<N; i++ ){ radian = 2.0 * PI * double(i)/double(N); ball[i].m = 2.0*double(i+1)/double(N); //ボールの質量 ball[i].x = l * cos(radian); ball[i].y = l * sin(radian); ball[i].z = poll_z; //ボールがめり込まないように ball[i].vx = 0.0; ball[i].vy = 0.0; ball[i].vz = 0.0; ball[i].ax = 0.0; ball[i].ay = 0.0; ball[i].az = 0.0; } }
ボールの位置「x,y,z」、速度「vx,vy,vz」、加速度「ax,ay,az」を 構造体「BALL」のメンバとします。「BALL」型の変数「ball1」を宣言と同時に初期値を設定しています。
速度・位置を逐次計算するアルゴリズム
//-------------------------------------------------------- // 計算と物体の描画 //-------------------------------------------------------- void Calculate(){ t = dt * double(tn); //ofs1 << dt*tn << " " << ball1.x << " " << ball1.vx << " " << ball1.ax <<endl; for(int i=0; i<N; i++){ //位置の算出 ball[i].x = ball[i].x + ball[i].vx * dt; ball[i].y = ball[i].y + ball[i].vy * dt; ball[i].z = ball[i].z + ball[i].vz * dt; l = sqrt(pow(ball[i].x,2)+pow(ball[i].y,2)+pow(poll_z-ball[i].z,2)); // 支柱からボールまでの距離 //加速度の算出 ball[i].ax = - k/ball[i].m * (1.0 - L/l) * ball[i].x; ball[i].ay = - k/ball[i].m * (1.0 - L/l) * ball[i].y; ball[i].az = - k/ball[i].m * (1.0 - L/l) * (ball[i].z-poll_z) -g; //速度の算出 ball[i].vx = ball[i].vx + ball[i].ax * dt; ball[i].vy = ball[i].vy + ball[i].ay * dt; ball[i].vz = ball[i].vz + ball[i].az * dt; } tn++; } void DrawStructure(){ for(int i=0; i<N; i++){ glPushMatrix(); glMaterialfv(GL_FRONT, GL_AMBIENT, ms_silver.ambient); glMaterialfv(GL_FRONT, GL_DIFFUSE, ms_silver.diffuse); glMaterialfv(GL_FRONT, GL_SPECULAR, ms_silver.specular); glMaterialfv(GL_FRONT, GL_SHININESS, &ms_silver.shininess); glTranslated(ball[i].x, ball[i].y, ball[i].z); //平行移動値の設定 glutSolidSphere(ball_r, 20, 20); //引数:(半径, Z軸まわりの分割数, Z軸に沿った分割数) glPopMatrix(); //バネ glPushMatrix(); glMaterialfv(GL_FRONT, GL_AMBIENT, ms_yellow_rubber.ambient); glMaterialfv(GL_FRONT, GL_DIFFUSE, ms_yellow_rubber.diffuse); glMaterialfv(GL_FRONT, GL_SPECULAR, ms_yellow_rubber.specular); glMaterialfv(GL_FRONT, GL_SHININESS, &ms_yellow_rubber.shininess); drowSolidSpring(box_x, box_y, box_z, (ball[i].x, (ball[i].y,(ball[i].z); //引数:(始点x, 始点y, 始点z, 終点x, 終点y, 終点z ) glPopMatrix(); } //軸(円錐) glPushMatrix(); glMaterialfv(GL_FRONT, GL_AMBIENT, ms_copper.ambient); glMaterialfv(GL_FRONT, GL_DIFFUSE, ms_copper.diffuse); glMaterialfv(GL_FRONT, GL_SPECULAR, ms_copper.specular); glMaterialfv(GL_FRONT, GL_SHININESS, &ms_copper.shininess); glTranslated(poll_x, poll_y, poll_z); glutSolidCone(5.0,20.0,20,20);//引数:(半径, 高さ, Z軸まわりの分割数, Z軸に沿った分割数) glPopMatrix(); }
位置「ball1.x」「ball1.y」「ball1.z」
速度「ball1.vx」「ball1.vy」「ball1.vz」
加速度「ball1.ax」「ball1.ay」「ball1.az」
を逐次計算します。
先にも述べましたが、このアルゴリズムでは誤差が大きく、
時間と共に単振動の振幅が大きくなってしまいます。
アルゴリズムの改善は、次回行ないます。
VisualC++ と OpenGL を利用した仮想物理実験室
第0章 仮想物理実験室の構築
- 【0-1】OpenGL と Visual C++ 2008 Express Edition の準備
- 【0-2】仮想物理実験室の構築 ・(0-2-1)ver1.0:基本形 ・(0-2-1)ver1.1:基本形+ばねの描画 ・(0-2-2)ver1.2:基本形+ばねの描画
- 【0-3】グラフ作成ソフト gnuplot のインストールと使い方
- 【A-1】参考文献
・(A-1-1)OpenGL について
・(A-1-2)VisualC++ について
・(A-1-3)物理シミュレーション
・(A-1-4)数値計算
第1章 様々な運動
- 【1-1】等速度直線運動
- ・(1-1-1)物理量について
- ・(1-1-2)等速度直線運動のアルゴリズムの導出
- ・(1-1-3)等速度直線運動のシミュレーション
- ・(1-1-4)等速度直線運動のグラフ化
- ・(1-1-5)等差数列を用いた等速度直線運動の解析解の導出
- ・(1-1-6)等速度直線運動のシミュレーション結果と解析解との比較
- 【1-2】等加速度直線運動
- ・(1-2-1)等加速度直線運動のアルゴリズムの導出
- ・(1-2-2)等加速度直線運動のシミュレーション
- ・(1-2-3)等加速度直線運動のグラフ化
- ・(1-2-4)階差数列を用いた等加速度直線運動の解析解の導出
- ・(1-2-5)等加速度直線運動のシミュレーション結果と解析解との比較
- ・(1-2-6)差分と微分
- ・(1-2-7)和分と積分
- ・(1-2-8)微分・積分を利用した等加速度直線運動の解析解の導出
- 【1-3】等加加速度直線運動
- ・(1-3-1)等加加速度直線運動のアルゴリズムの導出
- ・(1-3-2)等加加速度直線運動のシミュレーション
- ・(1-3-3)等加加速度直線運動のグラフ化
- ・(1-3-4)等加加速度直線運動の解析解の導出
- ・(1-3-5)等速度直線運動のシミュレーション結果と解析解との比較
- 【1-4】等速度円運動
- ・(1-4-1)等速度円運動のアルゴリズムの導出
- ・(1-4-2)等速度円運動のシミュレーション
- ・(1-4-3)等速度円運動のグラフ化
- ・(1-4-4)等速度円運動の解析解の導出
- ・(1-4-5)等速度円運動のシミュレーション結果と解析解との比較
- 【1-5】等加速度円運動
- ・(1-5-1)等加速度円運動のアルゴリズムの導出
- ・(1-5-2)等加速度円運動のシミュレーション
- ・(1-5-3)等加速度円運動のグラフ化
- ・(1-5-4)等加速度円運動の解析解の導出
- ・(1-5-5)等加速度円運動のシミュレーション結果と解析解との比較
第2章 ニュートンの運動方程式
- 【2-1】重力による運動:自由落下運動
- ・(2-1-1)加速度・力・質量の関係:ニュートンの運動方程式
- ・(2-1-2)重力による運動のアルゴリズムの導出
- ・(2-1-3)重力による自由落下運動のシミュレーション
- ・(2-1-4)重力による自由落下運動のグラフ化
- ・(2-1-5)重力による自由落下運動の解析解の導出
- ・(2-1-6)重力による自由落下運動のシミュレーション結果と解析解との比較
- 【2-2】重力による運動:鉛直投射運動
- ・(2-2-1)重力による鉛直投射運動のシミュレーション
- ・(2-2-2)重力による鉛直投射運動のグラフ化
- ・(2-2-3)重力による鉛直投射運動の解析解
- ・(2-2-4)重力による鉛直投射運動のシミュレーション結果と解析解との比較
- 【2-3】重力による運動:水平投射運動
- ・(2-3-1)重力による水平投射運動のシミュレーション
- ・(2-3-2)重力による水平投射運動のグラフ化
- ・(2-3-3)重力による水平投射運動の解析解
- ・(2-3-4)重力による水平投射運動のシミュレーション結果と解析解との比較
- 【2-4】重力による運動:斜方投射運動
- ・(2-4-1)重力による斜方投射運動のシミュレーション
- ・(2-4-2)重力による斜方投射運動のグラフ化
- ・(2-4-2)重力による斜方投射運動の解析解
- ・(2-4-3)重力による斜方投射運動のシミュレーション結果と解析解との比較
- 【2-5】重力による運動:斜方投射運動2
- ・(2-5-1)重力による斜方投射運動における初速度の分解
- ・(2-5-2)重力による斜方投射運動2のシミュレーション
- ・(2-5-3)重力による斜方投射運動2のグラフ化
- ・(2-5-4)重力による斜方投射運動2の投射角度と飛距離の関係
- ・(2-5-5)重力による斜方投射運動2の解析解
- ・(2-5-6)重力による斜方投射運動2のシミュレーション結果と解析解との比較
- 【2-6】重力+空気抵抗力による運動:自由落下運動
- ・(2-6-1)重力+空気抵抗力による運動の運動方程式とアルゴリズム
- ・(2-6-2)重力+空気抵抗力による自由落下運動のシミュレーション
- ・(2-6-3)重力+空気抵抗力による自由落下運動の解析解の導出
- ・(2-6-4)重力+空気抵抗力による自由落下運動のシミュレーション結果と解析解との比較
- 【2-7】重力+空気抵抗力による運動:斜方投射運動
- ・(2-7-1)重力+空気抵抗力による斜方投射運動のシミュレーション
- ・(2-7-2)重力+空気抵抗力による斜方投射運動の解析解の導出
- ・(2-7-3)重力+空気抵抗力による斜方投射運動のシミュレーション結果と解析解との比較
- 【2-8】ばね弾性力による運動:単振動運動(1次元)
- ・(2-8-1)ばね弾性力による運動の運動方程式(1次元)とアルゴリズム
- ・(2-8-2)ばね弾性力による単振動運動(1次元)のシミュレーション
- ・(2-8-3)ばね弾性力による単振動運動(1次元)の解析解の導出
- ・(2-8-4)ばね弾性力による単振動運動(1次元)のシミュレーション結果と解析解との比較
- 【2-9】ばね弾性力による運動:単振動運動(2次元)
- ・(2-9-1)ばね弾性力による運動の運動方程式と(2次元)アルゴリズム
- ・(2-9-2)ばね弾性力による単振動運動(2次元)のシミュレーション
- ・(2-9-3)ばね弾性力による単振動運動(2次元)の解析解:微分・積分を利用した解析解の導出
- ・(2-9-4)ばね弾性力による単振動運動(2次元)のシミュレーション結果と解析解との比較
- 【2-10】ばね弾性力+重力がある場合の運動:単振動運動(3次元)
- ・(2-10-1)ばね弾性力+重力による運動の運動方程式(3次元)とアルゴリズム
- ・(2-10-2)ばね弾性力+重力による単振動運動(3次元)のシミュレーション
- ・(2-10-3)ばね弾性力+重力による単振動運動(3次元)の解析解の導出
- ・(2-10-4)ばね弾性力+重力による単振動運動(3次元)のシミュレーション結果と解析解との比較
- 【2-11】ばね弾性力+重力がある場合の運動:多重連結ばねの運動
- ・(2-11-1)ばね弾性力の一般化
- ・(2-11-2)ばね弾性力+重力による多重連結ばねの運動の運動方程式
- ・(2-11-3)ばね弾性力+重力による多重連結ばねの運動のアルゴリズムの導出
- ・(2-11-4)ばね弾性力+重力による多重連結ばねの運動のシミュレーション
- 【2-12】ばね弾性力+遠心力がある場合の運動:円運動
- ・(2-12-1)遠心力の導出
- ・(2-12-2)遠心力+ばね弾性力による運動の運動方程式
- ・(2-12-3)遠心力+ばね弾性力による運動のアルゴリズムの導出
- ・(2-12-4)遠心力+ばね弾性力による運動のシミュレーション
- 【2-13】張力+重力がある場合の運動:単振子運動.html
- ・(2-13-1)張力の導出
- ・(2-13-2)張力+重力による多重連結ばねの運動の運動方程式
- ・(2-13-3)張力による多重連結ばねの運動のアルゴリズムの導出
- ・(2-13-4)張力による多重連結ばねの運動のシミュレーション
- 【2-14】万有引力による運動1:楕円運動
- ・(2-14-1)万有引力の導出
- ・(2-14-2)万有引力による運動の運動方程式
- ・(2-14-3)万有引力による運動のアルゴリズムの導出
- ・(2-14-4)初速度と軌道の関係
- ・(2-14-5)万有引力による運動のシミュレーション
- 【2-15】万有引力による運動2:円運動
- ・(2-15-1)万有引力による運動が円運動になるための条件
- ・(2-15-2)初速度の与え方
- ・(2-15-3)万有引力による円運動のシミュレーション
- 【2-16】万有引力による運動3:3体運動
- ・(2-16-1)3体運動とは
- ・(2-16-2)万有引力による3体運動の運動方程式
- ・(2-16-3)万有引力による3体運動のアルゴリズムの導出
- ・(2-16-3)万有引力による3体運動のシミュレーション
第3章 剛体の運動(エネルギー保存則と運動量保存則)
- 【3-1】床と球の衝突
- 【3-2】壁と球の衝突
- 【3-3】円柱と球の衝突(1次元)
- 【3-4】円柱と球の衝突(2次元 x-y平面)
- 【3-4-2】円柱と球の衝突(2次元 x-z平面)
- 【3-5】球と球の衝突(1次元)
- 【3-6】球と球の衝突(2次元)
- 【3-7】球と球の衝突(3次元)
付録
- 【A-1】参考文献
・(A-1-1)OpenGL について
・(A-1-2)VisualC++ について
・(A-1-3)物理シミュレーション
・(A-1-4)数値計算
未分類
力学
- ・ 2体問題のシミュレーション(C言語+ルンゲクッタ法)
- ・ ラグランジュ未定乗数法を用いた2重振子のシミュレーション
- ・ ラグランジュ未定乗数法を用いた球面振子のシミュレーション
- ・ ラグランジュ運動方程式2:極座標を用いた球面振子
- ・ ラグランジュ運動方程式1:極座標を用いた単振子
量子力学
- ・ 1次元量子力学の調和振動子における単一エネルギーの時間発展
- ・ 1次元量子力学の調和振動子における任意の初期状態に対する時間発展
- ・ 1次元量子力学の調和振動子におけるコヒーレント状態の空間分布
- ・ 1次元量子力学の調和振動子におけるコヒーレント状態の時間発展
- ・ 1次元量子力学の調和振動子における n励起状態の運動量表示
- ・ 1次元量子力学の調和振動子における任意の初期運動量分布に対する時間発展
- ・ 1次元量子力学の調和振動子における任意の初期運動量分布+任意の中心座標に対する時間発展
- ・ 1次元量子力学の調和振動子における任意の初期空間分布+任意の中心運動量に対する時間発展(最も簡単な方法)
- ・ 無限に深い2次元井戸型ポテンシャル
- ・ 無限に深い2次元井戸型ポテンシャル2
- ・ 無限に深い2次元井戸型ポテンシャル3
- ・ 無限に深い2次元井戸型ポテンシャル
- ・ 2次元自由粒子
- ・ 無限に深い井戸型ポテンシャルの時間発展2
- ・ シュレディンガー方程式に従う粒子の時間発展2:無限に深い井戸型ポテンシャルの時間発展
- ・ シュレディンガー方程式に従う粒子の時間発展:自由粒子