[ Python | 物理シミュレーション ]
路面をはねる質点の運動をUnrealで表示

 

質点の運動というと物理(力学)の分野で初めに勉強する内容なので,

高校で物理を履修した方には馴染みのある内容だと思いますが,

今回ご紹介する内容はより発展的な内容です.

 

高校物理の範囲ではまず,質点を静止状態から落下させる自由落下や,

ある初速度で質点を投げる斜方投射について勉強すると思います

 

 

自由落下や斜方投射の計算は手計算で簡単にできるので,Pythonを使って計算しても面白くないため今回はしません.

 

ここでは,ある初速度で投げられた質点が床に当たってバウンドしながら移動する挙動を,Pythonを使って計算してみたいと思います.

 

 

バウンドする質点の挙動は高校物理の範囲でも扱われますが,
反発係数という衝突した質点がどのくらいの勢いで跳ね返されるかという概念を使って計算します.

 

反発係数を使って計算できれば良いのですが,現実世界の現象を計算したい場合に,反発係数では表現できない状況も出てきます.

 

そこでここでは,質点のバウンドの挙動の計算に反発係数は使わずに,物体間の接触を表現する要素を使って,より柔軟に様々な状況を計算可能なプログラムをご紹介します.

 

 

通常こういった物理現象を計算する事をシミュレーションと呼び,シミュレーションに使う計算プログラムをシミュレーションモデルと呼びます.覚えておくと良いでしょう.

 

 

1. 計算プログラム

 

まずは計算プログラムを示します.

実行すると床の上を跳ねる質点の運動が計算できます.

 


Pythonで計算した結果をUnrealで表示させています

 

Unrealの利用方法については今後別の記事でご紹介します

 

2. 計算プログラムの解説

 

2.1. 必要なモジュールをインポート | 1行目~6行名

 

まずは必要なモジュールをインポートします.

 

このプログラムで利用するモジュールと目的は以下です.

 

・[モジュール] pyplot (パッケージ:matplotlib) / [目的] グラフの作成

・[モジュール] csv / [目的] csvファイルの作成

・[モジュール] numpy / [目的] csvファイルの作成のための配列の作成

 

2.2. モデルのパラメータを設定 | 7行目~13行名

 

重力加速度,質点の質量,質点と床の接触を表す要素であるばね定数と粘性係数を設定します.

 

ここで,質点と床の接触を表す要素について説明します.

 

質点と床の接触を表す要素

 

接触を表す最も基本的で一般的な方法は,接触する物体どうしの間にばねとダンパを想定する方法です.

 

 

質点が床よりも上の位置にある場合はばねとダンパは考慮せず,床よりも下にある場合に考慮することで,接触を表現することができます.

 

物体と物体が接触する場合,接触後にお互いの物体は弾かれますが,弾かれる挙動や,どの程度の勢いで弾かれるか等を,この方法により表現することができます.

 

ばねにはばね定数,ダンパには粘性係数を設定します.

 

ばねにより生じる力とダンパにより生じる力は,質点と床の間の相対距離をz,相対速度をvz,ばね定数と粘性係数をそれぞれk,cとすると,以下のようにして計算することができます.

 

(ばねで生じる力)

Fspring = k * z

 

(ダンパで生じる力)

Fdamper = c * vz

 

(トータルで生じる力)

Ftotal = Fspring + Fdamper

 

質点が床よりも下にある時だけ,Ftotalを計算します.

 

簡単ですね.

 

場合によってはもっと複雑なことを考えないといけない場合がありますが,まずは今回紹介する最も簡単な方法を覚えておくとよいでしょう.

 

今回はkzを100000N/m,cを50Ns/mと設定しています.

 

この値には特に根拠はなく,適当に設定しています.

 

質点がどの程度の勢いでバウンドするかによって,計算を流す人が決めるべき値です.

 

粘性係数を0Ns/mとすると,質点が床に衝突する直前の速度で,衝突後に床から跳ね返ります.

 

粘性係数を0Ns/mとすると,衝突時のエネルギの損失がないことになるため,同じ勢いで跳ね返ります.

 

高校物理で習う,反発係数が1という状況は,ここでいう粘性係数が0の状態を表します.

 

2.3. 初期設定 | 14行目~25行名

 

種々の初期値を設定します.

 

各行のそれぞれの意味は,プログラムに書いてあります.

 

このシミュレーションでは,時間を少しずつ増やしていき,微小時間経過する毎の質点の位置を計算しています.

18行目のdtは,計算の時間幅を表していて,この値が小さいほどより厳密な計算をすることになりますが,計算により多くの時間を要します.計算が終わるまで沢山待たなければなりません.

 

また逆にdtを大きな値とすると,計算時間は短くなりますが,計算精度は落ちます.

 

正しい計算をしたい場合はdtを小さくしないといけませんが,計算時間とのバランスを見て,程よい値を決める必要があります.適した値は状況により異なるので,答えはありません.

 

21行目から25行目で設定している初期位置や初速度の値は,想定したい状況によって変える値です.好きな値を自由に設定します.

 

2.4. メインパート | 26行目~63行名

 

計算の主要な部分です.

 

初期設定で設定した計算終了時間まで,設定した計算の時間幅dtずつ時間を足していって,各時間ごとの質点の位置を繰り返し計算します.

 

(32行目)

繰り返し計算をするために,関数whileを使います.

 

(33行目~37行目)

運動方程式を立てています.

 

質点のz座標(高さ)が0よりも小さい場合は,重力とばね,ダンパで生じる力を足した値を系に作用する正味の力としています.

 

質点のz座標(高さ)が0よりも大きい場合は,重力のみを系に作用する正味の力としています.

 

[point]
・ ニュートンの運動方程式を使う方法は物体の運動を計算する上で基本的な方法
・ 加速度は速度の時間微分
・ 速度は距離の時間微分
・ 加速度を時間で積分すると速度
・ 速度を時間で積分すると距離
・ 加速度を時間で積分すると速度
このあたりは物理の基本法則なのであえて説明しません

(39行目~41行目)

 

x軸(前後方向)とz軸(上下方向)の,各軸方向の加速度を計算しています.

 

x軸方向の力は発生していないと想定しているため,加速度を0としています.

 

(43行目~48行目)

 

微小時間dtの間の,各軸方向の速度の変化量dvxとdvzをそれぞれ計算しています.

 

速度の初期値(微小時間経過する直前の値)の値に,変化量dvxとdvzをそれぞれ足すことで,微小時間後の各軸方向の速度を計算しています.

 

(50行目~58行目)

 

速度の場合と同様に,まずは微小時間dtの間の,各軸方向の位置の変化量dxとdzをそれぞれ計算します.

 

位置の初期値(微小時間経過する直前の値)の値に,変化量dxとdzをそれぞれ足すことで,微小時間後の各軸方向の位置を計算しています.

 

これである時刻の質点の各軸方向の加速度と,速度,位置を計算することができました.

 

58行目で,次の時刻の計算をするために,時間の初期値timeに微小時間dtを足します.

 

(60行目~63行目)

 

計算した座標x,zと時間timeを,結果の配列の末尾に加えます.

 

結果の配列を作ることで,計算開始から終了までの,一連の質点の速度と位置を記録することができます.

 

63行目が実行された後,シミュレーション内部の時間timeが,終了時間stoptimeと一致するまで,再度32行目まで戻って同様の計算を繰り返します.

 

2.5. グラフの作成 | 64行目~72行名

 

結果のグラフを作成します.

 

グラフの作成に,61行目と62行目で作成した質点の位置のリストを使います.

 

2.6. csvファイルの作成 | 73行目~最終行

 

結果をcsvファイルに保存します.

 

pythonのカレントフォルダに,data.csvというファイル名で,結果を出力します.

 

3. まとめ

 

床の上をバウンドする質点の挙動を,pythonを使ってシミュレーションする方法を紹介しました.

 

とても簡単なプログラムなので,物理シミュレーションの初心者の方には,分かりやすい内容だったと思います.

 

もし何かわからないことがあれば,以下のコメント欄でご質問ください.

 

Have a good Python life 🙂

コメントを残す

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です