# 相疊的木塊
> 作者：王一哲
> 日期：2018/4/22
<br />
這是會出現在基礎物理2B下的例題，出現的章節甚至橫跨動量守恆、功與能量、碰撞等三章，題目敘述如下：
> 木塊A質量為 $2m$、木塊B質量為 $m$，兩個木塊相疊，桌面與木塊A之間沒有摩擦力，A的初速度為 $v$。兩個木塊之間的動摩擦係數為 $\mu_k$，若A的長度夠長，A、B最後能以相同的速度前進，請問：(1) A、B的末速為何？ (2) 經過多久之後兩者速度相同。
<img style="display: block; margin-left: auto; margin-right: auto" height="50%" width="50%" src="https://i.imgur.com/AenXbEl.png">
<div style="text-align:center">兩木塊相疊示意圖</div>
<br />
題目有兩種版本，第一種是像上圖的版本，木塊B的初位置在木塊A的最右側且原為靜止，木塊A的初速度向右；第二種則是木塊B的初位置在木塊A的最左側且初速度向右，木塊A原為靜止。這次的目標是將運動過程中木塊A、B的速度 - 時間關係圖、能量 - 時間關係圖畫出來。
<img style="display: block; margin-left: auto; margin-right: auto" height="100%" width="100%" src="https://i.imgur.com/CP8ftNO.png">
<div style="text-align:center">模擬程式畫面截圖</div>
<br />
## 程式 13-1.相疊的木塊完全非彈性碰撞 （[取得程式碼](https://github.com/YiZheWangTw/VPythonTutorial/blob/master/13.%E7%9B%B8%E7%96%8A%E6%9C%A8%E5%A1%8A/13-1_blocks_collision.py)） （[GlowScript 網站動畫連結](http://www.glowscript.org/#/user/yizhe/folder/Public/program/13-1blockscollision)）
```python=
"""
VPython教學: 13-1.相疊的木塊完全非彈性碰撞
Ver. 1: 2017/6/16
Ver. 2: 2018/1/4 修改為 Python 3.X 版程式碼
Ver. 3: 2018/2/26 增加 v-t 圖
Ver. 4: 2019/9/8
作者: 王一哲
"""
from vpython import *
"""
1. 參數設定, 設定變數及初始值
(1) m1 = 0.2, v1 = 0.0, m2 = 0.1, v2 = 2.0
(2) m1 = 0.2, v1 = 1.0, m2 = 0.1, v2 = 0.0
(3) m1 = 0.2, v1 = 0.0, m2 = 0.1, v2 = 2.5 => b2 會飛出去
(4) m1 = 0.2, v1 = 2.0, m2 = 0.1, v2 = 0.0 => b1 到畫面右側時還沒有達到等速度
"""
d1 , h1, w1 = 1.8, 0.2, 0.2 # 下方木塊1的長度 = 1.8m, 高度 = 0.2m, 寬度 = 0.2 m
d2 , h2, w2 = 0.2, 0.2, 0.2 # 上方木塊2的長度 = 0.2m, 高度 = 0.2m, 寬度 = 0.2 m
m1, v1, c1 = 0.2, 0.0, color.red # 下方木塊1的質量 = 0.2 kg, 初速 = 0.0 m/s, 紅色
m2, v2, c2 = 0.1, 2.0, color.green # 上方木塊2的質量 = 0.1 kg, 初速 = 2.0 m/s, 綠色
xmax, xmin = 2.0, -2.0 # x 軸範圍
g = 9.8 # 重力加速度 = 9.8 m/s^2
mu = 0.1 # 動摩擦係數
dt = 0.0005 # 畫面更新的時間間隔，單位為s, 原為0.001但不夠準確, 故改為0.0005
t = 0 # 模擬所經過的時間 ，單位為s，初始值為0
bx = 0 # 計算 b2 初位置用的變數
i, te = 0, -1 # 記錄 b1、b2 達到等速度經過時間用的變數
"""
2. 畫面設定
"""
# 產生動畫視窗
scene = canvas(title="Two Blocks", width=800, height=300, center=vec(0, 0.4, 0), background=vec(0, 0.6, 0.6))
# 產生地板
floor = box(pos=vec(0, -0.5*h2, 0), size=vec(xmax - xmin, 0.05, 0.8), color=color.blue)
# 產生下方木塊 b1, 位於畫面最左側 xmin + d1/2, 初速度 v1
b1 = box(pos=vec(xmin + d1/2, 0, 0), size=vec(d1, h1, w1), color = c1, v=vec(v1, 0, 0))
# 產生上方木塊 b2, 若 v2 >= v1 則位於畫面最左側 xmin + d2/2; 若 v2 < v1, 則位於 b1 最右側 xmin + d1 - d2/2, 初速度 v2
if(v2 >= v1): bx = xmin + 0.5*d2
else: bx = xmin + d1 - 0.5*d2
b2 = box(pos=vec(bx, h1, 0), size=vec(d2, h2, w2), color=c2, v=vec(v2, 0, 0))
# 繪圖部分
gd = graph(title="<i>E</i> - <i>t</i> plot", x=0, y=300, width=600, height=450, xtitle="<i>t</i> (s)",
ytitle="red: <i>K</i><sub>1</sub>, green: <i>K</i><sub>2</sub>, blue: <i>E</i> (J)")
kt1 = gcurve(graph=gd, color=c1)
kt2 = gcurve(graph=gd, color=c2)
et = gcurve(graph=gd, color=color.blue)
gd2 = graph(title="<i>v</i> - <i>t</i> plot", x=0, y=750, width=600, height=450, xtitle="<i>t</i> (s)",
ytitle="red: <i>v</i><sub>1</sub>, green: <i>v</i><sub>2</sub> (m/s)")
vt1 = gcurve(graph=gd2, color=c1)
vt2 = gcurve(graph=gd2, color=c2)
"""
3. 物體運動部分, 重複執行直到: (1) b1 抵達邊緣時, (2) b2 抵達 b1 邊緣時
"""
while((b1.pos.x <= xmax - d1/2) and (b2.pos.x + d2/2 <= b1.pos.x + d1/2 + 0.001)):
#動畫執行頻率, 每秒500次
rate(500)
#由 b1, b2 速度判斷動摩擦力量值及加速度方向,
if(b2.v.x > b1.v.x):
force = mu * m2 * g
b1.a = vec(force / m1, 0, 0)
b2.a = vec(-force / m2, 0, 0)
elif(b2.v.x < b1.v.x):
force = mu * m2 * g
b1.a = vec(-force / m1, 0, 0)
b2.a = vec(force / m2, 0, 0)
else:
force = 0
b1.a = vec(0, 0, 0)
b2.a = vec(0, 0, 0)
# 記錄 b1、b2 達到等速度需要的時間
if(abs(b2.v.x - b1.v.x) < 0.0005 and i == 0):
te = t
i = 1
#計算 b1, b2 的速度, 位置
b1.v += b1.a * dt
b2.v += b2.a * dt
b1.pos += b1.v * dt
b2.pos += b2.v * dt
# 計算 b1, b2 動能及系統力學能, 畫能量 - 時間關係圖
k1 = 0.5 * m1 * mag2(b1.v)
k2 = 0.5 * m2 * mag2(b2.v)
e = k1 + k2
kt1.plot(pos = (t, k1))
kt2.plot(pos = (t, k2))
et.plot(pos = (t, e))
# 畫速度 - 時間關係圖
vt1.plot(pos = (t, b1.v.x))
vt2.plot(pos = (t, b2.v.x))
#更新時間
t += dt
#結束 while 迴圈, 印出末速及 b1、b2 達到等速度需要的時間
print("v1 = ", b1.v.x)
print("v2 = ", b2.v.x)
print("te = ", te)
print("end")
```
<br />
### 參數設定
在此設定變數為木塊的長、寬、高、質量、初速度、顏色、x軸的範圍、重力加速度、動摩擦係數、時間、時間間隔，其中時間間隔 dt 設定為 0.0005，這是因為設定為 0.001 時的計算誤差較大，故選擇較小的數值。
<br />
### 畫面設定
大部分的功能在之前的動畫中都已經使用過了，應該只需要解釋比較特別的部分：**由 v1、v2 決定 b2 的初位置**。當 v2 >= v1 時，b2 的初位置要在 b1 的左側，否則在執行動畫之後 b2 會立刻飛出 b1；反之，當 v2 < v1 時，b2 的初位置要在 b1 的右側，因此需要以下兩行程式碼計算 b2 的初位置：
```python
if(v2 >= v1): bx = xmin + d2/2
else: bx = xmin + d1 - d2/2
```
<br />
### 物體運動
1. 為了使動畫重覆執行，直到以下 2 種情況發生為止：
1. b1 抵達畫面邊緣時
2. b2 抵達 b1 邊緣時
因此 while 迴圈中的條件設定為
```python
((b1.pos.x <= xmax - d1/2) and (b2.pos.x + d2/2 <= b1.pos.x + d1/2 + 0.001)
```
第二個條件中加上 0.001 是為了避免當 v1 > v2 時，動畫一開始執行就停止運作。
2. b1、b2 之間的動摩擦力方向與當時的速度有關，共有 3 種狀況：
1. 當 b2.v.x > b1.v.x 時，動摩擦力對 b1 的作用方向向右、對 b2 的作用方向向左。
2. 當 b2.v.x < b1.v.x 時，動摩擦力對 b1 的作用方向向左、對 b2 的作用方向向右。
3. 當 b2.v.x = b1.v.x 時，動摩擦力為 0。
因此需要第 60 ~ 71 行的程式碼，利用 if … elif … else 決定動摩擦力的量值及方向。其實 else 的部分可以省略，因為在程式中 b1.v.x 和 b2.v.x 都是浮點數，幾乎不可能相等。
3. 為了記錄 b1、b2 達到等速度經過的時間，需要第 73 ~ 75 行程式碼。由於 b1.v.x 和 b2.v.x 幾乎不可能相等，因此條件設定為 **abs(b2.v.x - b1.v.x) < 0.0005**。由於只需要記錄 b1.v.x 和 b2.v.x 幾乎相等的第一個時刻，因此條件中再加上 i == 0，並在執行 if 中的程式碼時設定 i = 1，就不會再執行第二次。
4. 計算木塊的加速度，更新速度、位置，計算木塊的動能，畫速度 - 時間關係圖及能量 - 時間關係圖，更新時間。
5. 最後印出木塊末速及 b1、b2 達到等速度需要的時間。
<br />
### 模擬結果
b1、b2 之間的動摩擦力與加速度量值分別為
$$ f_k = \mu_k N = \mu_k m_2 g~~~~~a_1 = \frac{\mu_k m_2 g}{m_1}~~~~~a_2 = \mu_k g $$
可以由上式及動量守恆計算木塊的末速及達成等速度需要的時間。我測試了 4 種不同的條件：
1. m1 = 0.2, v1 = 0.0, m2 = 0.1, v2 = 2.0
理論值：v1’ = v2’ = 2/3, te ≈ 1.360544
模擬結果：v1’ = 0.6668900000000263, v2’ = 0.6662199999997995, te = 1.3604999999999057
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/RFP7XEG.png">
<div style="text-align:center">條件1的模擬結果：速度 - 時間關係圖</div>
<br />
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/h3yc84S.png">
<div style="text-align:center">條件1的模擬結果：能量 - 時間關係圖</div>
<br />
2. m1 = 0.2, v1 = 1.0, m2 = 0.1, v2 = 0.0
理論值：v1’ = v2’ = 2/3, te ≈ 0.680272
模擬結果：v1’ = 0.666554999999931, v2’ = 0.6668899999999899, te = 0.6799999999999805
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/G9Rmt5P.png">
<div style="text-align:center">條件2的模擬結果：速度 - 時間關係圖</div>
<br />
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/08Nmo0r.png">
<div style="text-align:center">條件2的模擬結果：能量 - 時間關係圖</div>
<br />
3. m1 = 0.2, v1 = 0.0, m2 = 0.1, v2 = 2.5
理論值：b2 會飛出去，不會達到等速度
模擬結果：b2 會飛出 b1 時 v1’ = 0.41943999999999326, v2’ = 1.6611199999998263
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/AOb1wmM.png">
<div style="text-align:center">條件3的模擬結果：速度 - 時間關係圖</div>
<br />
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/hBoDeWv.png">
<div style="text-align:center">條件3的模擬結果：能量 - 時間關係圖</div>
<br />
4. m1 = 0.2, v1 = 2.0, m2 = 0.1, v2 = 0.0
理論值：b1 到達畫面右側時兩木塊還未達到等速度
模擬結果：b1 到達畫面右側時 v1’ = 1.357854999999867, v2’ = 1.2842900000000423
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/phdvIdQ.png">
<div style="text-align:center">條件4的模擬結果：速度 - 時間關係圖</div>
<br />
<img style="display: block; margin-left: auto; margin-right: auto" height="80%" width="80%" src="https://imgur.com/OCERphs.png">
<div style="text-align:center">條件4的模擬結果：能量 - 時間關係圖</div>
<br />
## VPython官方說明書
1. **canvas**: http://www.glowscript.org/docs/VPythonDocs/canvas.html
2. **box**: http://www.glowscript.org/docs/VPythonDocs/box.html
3. **graph**: http://www.glowscript.org/docs/VPythonDocs/graph.html
---
###### tags:`VPython`
Sign in with Wallet
Connect another wallet