熱門文章

  • 安裝及測試 VPython
    作者：王一哲 日期：2018/3/5 下載並安裝Python Python 是一種直譯式的高階程式語言，功能強大且程式碼較為簡潔，個人認為比 C 語言容易上手。在 Google 上搜尋 Python，第1個搜尋結果通常會是 Python 官方網站的首頁 。目前最新的...
  • 使用 SciPy 計算快速傅立葉變換
    作者：王一哲 日期：2021/3/28 快速傅立葉變換 快速傅立葉變換 (fast Fourier transform, FFT) 是一種用來計算離散傅立葉變換 (discrete Fourier transform, DFT) 及其逆變換的計算方法，目前常用的是庫利－...
  • Word圖表標號及自動產生目錄
    Word圖表標號及自動產生目錄 日期：2017/12/25   我會寫這篇文章是因為某人寫論文時遇到困難，我只好又幫忙測試。其實我很討厭 Microsoft Word，平常只用 LibreOffice Writer。以下是我找出來的方法，使用的版本為 Word 20...

2018年5月5日 星期六

自由落下兩球碰撞

作者：王一哲
日期：2018/5/5




若將兩個彈性極佳的球上下疊放，且下方的球質量較大，將兩球從某個高度釋放，兩球受到重力的作用由靜止開始向下加速，當下方的球撞到地面反彈後會發生什麼事？在 YouTube 上有相當多實驗影片，例如 "Stacked Ball Drop" 。我們可以利用之前寫好的〈自由落下〉及〈一維彈性碰撞〉模擬程式來處理這個問題。

以下共有兩個程式，16-1 有兩個球，16-2 改為三個球。其實可以讓畫只執行到最上方的球飛出去之後就停止，因為在真實的情境下很難讓最上方的球被撞擊後往正上方飛出，掉落下後發生第二次碰撞的機率極低。


程式16-1畫面截圖




程式16-2畫面截圖



程式 16-1.自由落下兩球碰撞


取得程式碼
GlowScript 網站動畫連結

"""
 VPython教學: 16-1.自由落下兩球碰撞
 Ver. 1: 2018/2/28
 Ver. 2: 2019/9/14
 作者: 王一哲
"""
from vpython import *

"""
 1. 參數設定, 設定變數及初始值
    (1) r1 = 2, m1 = 2, r2 = 1, m2 = 1
    (2) r1 = 1, m1 = 1, r2 = 2, m2 = 2
"""
r1, m1, c1 = 2, 2, color.red     # 小球1半徑, 質量, 顏色
r2, m2, c2 = 1, 1, color.green   # 小球2半徑, 質量, 顏色
h, L = 15, 40     # 小球離地高度, 地板邊長
g = 9.8           # 重力加速度 9.8 m/s^2
t, dt = 0, 0.001  # 時間, 時間間隔

"""
 2. 畫面設定
"""
# 產生動畫視窗
scene = canvas(title="Free Fall and Collision", width=600, height=600, x=0, y=0,
               center=vec(0, h/2, 0), background=vec(0, 0.6, 0.6), range=L)
# 產生地板
floor = box(pos=vec(0, -r1, 0), size=vec(L, 0.01, L), texture=textures.wood)
# 產生小球並設定初速度、加速度
b1 = sphere(pos=vec(0, h, 0), radius=r1, color=c1, m=m1, v=vec(0, 0, 0), a=vec(0, -g, 0))
b2 = sphere(pos=vec(0, h+r1+r2, 0), radius=r2, color=c2, m=m2, v=vec(0, 0, 0), a=vec(0, -g, 0))
# y-t plot
gd = graph(title="<i>y</i>-<i>t</i> plot", width=600, height=450, x=0, y=600,
           xtitle="<i>t</i> (s)", ytitle="red: <i>y</i><sub>1</sub>, green: <i>y</i><sub>2</sub> (m)")
yt1 = gcurve(graph=gd, color=c1)
yt2 = gcurve(graph=gd, color=c2)
# v-t plot
gd2 = graph(title="<i>v</i>-<i>t</i> plot", width=600, height=450, x=0, y=1050,
            xtitle="<i>t</i> (s)", ytitle="red: <i>v</i><sub>2</sub>, green: <i>v</i><sub>2</sub> (m/s)")
vt1 = gcurve(graph=gd2, color=c1)
vt2 = gcurve(graph=gd2, color=c2)
# 一維彈性碰撞速度公式
def af_col_v(m1, m2, v1, v2):
    v1_prime = (m1-m2)/(m1+m2)*v1 + (2*m2)/(m1+m2)*v2
    v2_prime = (2*m1)/(m1+m2)*v1 + (m2-m1)/(m1+m2)*v2
    return (v1_prime, v2_prime)

"""
 3. 物體運動部分, 小球觸地時反彈
"""
while(t < 20):
    rate(1000)
# 更新小球速度、位置 
    b1.v += b1.a*dt
    b1.pos += b1.v*dt
    b2.v += b2.a*dt
    b2.pos += b2.v*dt
# 繪製小球 y-t、v-t 圖
    yt1.plot(pos=(t, b1.pos.y))
    yt2.plot(pos=(t, b2.pos.y))
    vt1.plot(pos=(t, b1.v.y))
    vt2.plot(pos=(t, b2.v.y))
# 若 b1 撞到地板則反彈
    if(b1.pos.y <= 0 and b1.v.y < 0): b1.v.y = -b1.v.y
# 若 b1、b2 相撞則計算撞後速度並重新指定給 vy1, vy2
    if(b2.pos.y - b1.pos.y <= r1 + r2):
        b1.v.y, b2.v.y = af_col_v(b1.m, b2.m, b1.v.y, b2.v.y)
# 更新時間
    t += dt




參數設定


在此設定變數為小球的半徑、質量、顏色，小球離地高度、地板邊長、重力加速度、時間、時間間隔，對應的變數名稱請參考程式碼。




畫面設定


產生動畫視窗、地板、小球、繪圖視窗的程式碼在之前動畫當中已經出現很多次，這裡就不再贅述。比較特別之處在於動畫視窗中多加了一個選項range=L，這是為了限制動畫視窗顯示範圍，使畫面不會因為上方的小球反彈到太高的地方而自動縮小。




物體運動部分


  1. 更新小球速度、位置。
  2. 繪製小球 y-t、v-t 圖。
  3. 若b1.pos.y <= 0 代表 b1 撞到地板，將速度加上負號使 b1 反彈。
  4. 若b2.pos.y - b1.pos.y <= r1 + r2 代表 b1、b2 相撞，將 m1、b1.v.y、m2、b2.v.y 代入自訂函式 af_col_v 中計算撞後速度，再重新指定給 b1.v.y、b2.v.y。



模擬結果


我測試了 2 種不同的條件：

  1. r1 = 2, m1 = 2, r2 = 1, m2 = 1
  2. r1 = 1, m1 = 1, r2 = 2, m2 = 2

當下方的小球 b1 較重時，b2 反彈的高度會遠超過原本的離地高度。當下方的小球 b1 較輕時，b1 會不斷地在地板和 b2 之間來回反彈。




條件1的模擬結果：高度 - 時間關係圖




條件1的模擬結果：速度 - 時間關係圖




條件2的模擬結果：高度 - 時間關係圖




條件2的模擬結果：速度 - 時間關係圖



程式 16-2.自由落下三球碰撞


取得程式碼
GlowScript 網站動畫連結

"""
 VPython教學: 16-2.自由落下三球碰撞
 Ver. 1: 2018/2/28
 Ver. 2: 2019/9/14
 作者: 王一哲
"""
from vpython import *

"""
 1. 參數設定, 設定變數及初始值
    (1) r1 = 3, m1 = 3, r2 = 2, m2 = 2, r3 = 1, m3 = 1
    (2) r1 = 1, m1 = 1, r2 = 2, m2 = 2, r3 = 3, m3 = 3
"""
r1, m1, c1 = 3, 3, color.red     # 小球1半徑, 質量, 顏色
r2, m2, c2 = 2, 2, color.green   # 小球2半徑, 質量, 顏色
r3, m3, c3 = 1, 1, color.blue    # 小球3半徑, 質量, 顏色
h, L = 15, 50      # 小球離地高度, 地板邊長
g = 9.8            # 重力加速度 9.8 m/s^2
t, dt = 0, 0.001   # 時間, 時間間隔

"""
 2. 畫面設定
"""
# 產生動畫視窗
scene = canvas(title="Free Fall and Collision", width=600, height=600, x=0, y=0,
               center=vec(0, h/2, 0), background=vec(0, 0.6, 0.6), range=L)
# 產生地板
floor = box(pos=vec(0, -r1, 0), size=vec(L, 0.01, L), texture=textures.wood)
# 產生小球並設定初速度、加速度
b1 = sphere(pos=vec(0, h, 0), radius=r1, color=c1, m=m1, v=vec(0, 0, 0), a=vec(0, -g, 0))
b2 = sphere(pos=vec(0, h+r1+r2, 0), radius=r2, color=c2, m=m2, v=vec(0, 0, 0), a=vec(0, -g, 0))
b3 = sphere(pos=vec(0, h+r1+2*r2+r3, 0), radius=r3, color=c3, m=m3, v=vec(0, 0, 0), a=vec(0, -g, 0))
# y-t plot
gd = graph(title="<i>y</i>-<i>t</i> plot", width=600, height=450, x=0, y=600,
           xtitle="<i>t</i> (s)", ytitle="red: <i>y</i><sub>1</sub>, green: <i>y</i><sub>2</sub>, blue: <i>y</i><sub>3</sub> (m)")
yt1 = gcurve(graph=gd, color=c1)
yt2 = gcurve(graph=gd, color=c2)
yt3 = gcurve(graph=gd, color=c3)
# v-t plot
gd2 = graph(title="<i>v</i>-<i>t</i> plot", width=600, height=450, x=0, y=1050,
            xtitle="<i>t</i> (s)", ytitle="red: <i>v</i><sub>2</sub>, green: <i>v</i><sub>2</sub>, blue: <i>v</i><sub>3</sub> (m/s)")
vt1 = gcurve(graph=gd2, color=c1)
vt2 = gcurve(graph=gd2, color=c2)
vt3 = gcurve(graph=gd2, color=c3)
# 一維彈性碰撞速度公式
def af_col_v(m1, m2, v1, v2):
    v1_prime = (m1-m2)/(m1+m2)*v1 + (2*m2)/(m1+m2)*v2
    v2_prime = (2*m1)/(m1+m2)*v1 + (m2-m1)/(m1+m2)*v2
    return (v1_prime, v2_prime)

"""
 3. 物體運動部分, 小球觸地時反彈
"""
while(t < 20):
    rate(1000)
# 更新小球速度、位置
    b1.v += b1.a*dt
    b1.pos += b1.v*dt
    b2.v += b2.a*dt
    b2.pos += b2.v*dt
    b3.v += b3.a*dt
    b3.pos += b3.v*dt
# 繪製小球 y-t、v-t 圖
    yt1.plot(pos=(t, b1.pos.y))
    yt2.plot(pos=(t, b2.pos.y))
    yt3.plot(pos=(t, b3.pos.y))
    vt1.plot(pos=(t, b1.v.y))
    vt2.plot(pos=(t, b2.v.y))
    vt3.plot(pos=(t, b3.v.y))
# 若 b1 撞到地板則反彈
    if(b1.pos.y <= 0 and b1.v.y < 0): b1.v.y = -b1.v.y
# 若 b1、b2 相撞則計算撞後速度並重新指定給 vy1, vy2
    if(b2.pos.y - b1.pos.y <= r1 + r2):
        b1.v.y, b2.v.y = af_col_v(b1.m, b2.m, b1.v.y, b2.v.y)
# 若 b2、b3 相撞則計算撞後速度並重新指定給 vy2, vy3
    if(b3.pos.y - b2.pos.y <= r2 + r3):
        b2.v.y, b3.v.y = af_col_v(b2.m, b3.m, b2.v.y, b3.v.y)
# 更新時間
    t += dt




程式設計部分


程式 16-2 和 16-1 幾乎相同，只是多了一個小球 b3。




模擬結果


我測試了 2 種不同的條件：

  1. r1 = 3, m1 = 3, r2 = 2, m2 = 2, r3 = 1, m3 = 1
  2. r1 = 1, m1 = 1, r2 = 2, m2 = 2, r3 = 3, m3 = 3

當下方的小球 b1、b2 比 b3 重時，b3 反彈的高度會遠超過原本的離地高度。當下方的小球 b1 較輕時，b1 會不斷地在地板和 b2 之間來回反彈。





條件1的模擬結果：高度 - 時間關係圖




條件1的模擬結果：速度 - 時間關係圖




條件2的模擬結果：高度 - 時間關係圖




條件2的模擬結果：速度 - 時間關係圖



結語


當我們寫過較多的程式之後，寫程式的速度會開始變快，熟練度較高是原因之一，同時也更懂得利用之前寫好的程式碼，例如這次就是利用了自由落下和一維彈性碰撞的程式。接下來可以試著發揮創意，將之前寫過的程式碼組合成有趣的東西。




VPython官方說明書


  1. canvas: http://www.glowscript.org/docs/VPythonDocs/canvas.html
  2. box: http://www.glowscript.org/docs/VPythonDocs/box.html
  3. sphere: http://www.glowscript.org/docs/VPythonDocs/sphere.html
  4. graph: http://www.glowscript.org/docs/VPythonDocs/graph.html




HackMD 版本連結：https://hackmd.io/@yizhewang/B1NnE5Hz7
張貼者： 於
標籤：

沒有留言:

張貼留言

訂閱： 張貼留言 (Atom)