Ném xiên có lực cản không khí — Mô phỏng và khám phá¶
Mục tiêu học tập¶
Sau bài này, em sẽ:
- Mô tả được sự khác biệt giữa ném xiên lý tưởng (parabola) và ném xiên thực tế khi có lực cản không khí.
- Sử dụng
scipy.integrate.solve_ivpđể giải phương trình vi phân chuyển động khi không có công thức giải tích. - Phân tích ảnh hưởng của hệ số cản $k$ và góc ném $\theta$ lên tầm xa, độ cao và hình dạng quỹ đạo.
- Khám phá rằng góc ném tối ưu $\theta_{opt}$ KHÔNG còn là $45°$ khi có lực cản.
Yêu cầu trước¶
- Đã học cơ học chất điểm lớp 10 (định luật II Newton, phương trình chuyển động).
- Biết khái niệm đạo hàm và tích phân ở mức cơ bản (lớp 11 hoặc bài tự học).
- Quen thuộc với Python mức sơ cấp: định nghĩa hàm, list, vòng lặp
for.
Bối cảnh thực tế¶
Em hãy nghĩ về ba tình huống quen thuộc:
- Bóng rổ ném 3 điểm — bóng vạch một đường cong tới rổ.
- Đạn pháo trong các trận đánh lịch sử — pháo binh phải ngắm góc rất chính xác.
- Hạt mưa rơi từ mây cao 2 km — vì sao không đập chết người ta khi tới mặt đất?
Sách giáo khoa lớp 10 mô tả ném xiên là parabola — quỹ đạo đối xứng đẹp đẽ với góc tối ưu $45°$. Nhưng quan sát thực tế cho thấy:
- Bóng rổ rơi xuống dốc hơn so với khi đi lên.
- Pháo binh thực tế nhắm góc dưới $45°$ để đạt tầm xa lớn nhất.
- Hạt mưa không tăng tốc mãi — đạt một vận tốc tới hạn rồi rơi đều.
Tại sao quỹ đạo thực không phải parabola? Câu trả lời: lực cản không khí. Bài này sẽ mô phỏng và đo lường ảnh hưởng đó.
"""Phần 0 — Chuẩn bị: import thư viện và cấu hình mặc định cho biểu đồ."""
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
import pint
# Cấu hình mặc định: figure rộng hơn, font dễ đọc trên màn hình lớp
plt.rcParams['figure.figsize'] = (10, 6)
plt.rcParams['font.size'] = 11
ureg = pint.UnitRegistry() # dùng cho kiểm tra đơn vị ở phần mở rộng
# In phiên bản thư viện để kết quả tái lập được
import scipy, matplotlib
print(f'numpy {np.__version__}')
print(f'scipy {scipy.__version__}')
print(f'matplotlib {matplotlib.__version__}')
print(f'pint {pint.__version__}')
numpy 2.4.4 scipy 1.17.1 matplotlib 3.10.8 pint 0.25.2
Phần 1 — Ôn tập: ném xiên không có lực cản¶
Khi bỏ qua lực cản không khí, vật ném xiên chỉ chịu trọng lực $\vec{P} = m\vec{g}$. Phương trình chuyển động đơn giản:
$$ \begin{aligned} v_x(t) &= v_0 \cos\theta \\ v_y(t) &= v_0 \sin\theta - g t \\ x(t) &= v_0 \cos\theta \cdot t \\ y(t) &= v_0 \sin\theta \cdot t - \tfrac{1}{2} g t^2 \end{aligned} $$
Từ đó dẫn ra hai đại lượng quan trọng:
- Tầm xa (khi đáp đất, $y = 0$): $$ R = \frac{v_0^2 \sin(2\theta)}{g} $$
- Độ cao cực đại (khi $v_y = 0$): $$ h_{max} = \frac{v_0^2 \sin^2\theta}{2g} $$
⇒ Tầm xa lớn nhất khi $\sin(2\theta) = 1$, tức $\theta = 45°$. Đây là kết quả "quen mặt" trong sách giáo khoa.
def projectile_no_drag(v0, theta_deg, g=9.81, n_pts=200):
"""Tính quỹ đạo ném xiên KHÔNG có lực cản — công thức giải tích.
Tham số:
v0 — vận tốc ban đầu (m/s)
theta_deg — góc ném so với mặt đất (độ)
g — gia tốc trọng trường (m/s²), mặc định 9.81
n_pts — số điểm lấy mẫu trên quỹ đạo
Trả về:
(t, x, y, range, h_max, t_flight)
"""
theta = np.deg2rad(theta_deg)
vx0, vy0 = v0 * np.cos(theta), v0 * np.sin(theta)
t_flight = 2 * vy0 / g # thời gian chạm đất
t = np.linspace(0, t_flight, n_pts)
x = vx0 * t
y = vy0 * t - 0.5 * g * t**2
range_ = vx0 * t_flight
h_max = vy0**2 / (2 * g)
return t, x, y, range_, h_max, t_flight
# Test với tham số chuẩn: v0=20 m/s, θ=45°
t0, x0, y0, R0, hmax0, tf0 = projectile_no_drag(v0=20, theta_deg=45)
plt.plot(x0, y0, 'b-', linewidth=2)
plt.xlabel('x (m)'); plt.ylabel('y (m)')
plt.title('Ném xiên không lực cản — v₀=20 m/s, θ=45°')
plt.grid(alpha=0.3); plt.axis('equal')
plt.show()
print(f'Tầm xa R = {R0:.2f} m')
print(f'Độ cao max h = {hmax0:.2f} m')
print(f'Thời gian bay = {tf0:.2f} s')
Tầm xa R = 40.77 m Độ cao max h = 10.19 m Thời gian bay = 2.88 s
Phần 2 — Lực cản không khí: tuyến tính vs bậc hai¶
Lực cản không khí phụ thuộc vào vận tốc. Có hai mô hình thường dùng:
| Mô hình | Công thức | Khi nào dùng |
|---|---|---|
| Tuyến tính | $\vec{F}_d = -b\vec{v}$ | Vật rất nhỏ, vận tốc thấp (vd: hạt bụi, giọt sương) |
| Bậc hai | $\vec{F}_d = -\tfrac{1}{2}\rho C_d A \lvert\vec{v}\rvert \vec{v}$ | Vật cỡ macro ở tốc độ thường (bóng, đạn, người) |
Bài này dùng mô hình bậc hai vì sát với thực tế ném vật cỡ macro:
$$ \vec{F}_d = -\tfrac{1}{2} \rho C_d A \, |\vec{v}| \, \vec{v} $$
Trong đó:
- $\rho$ — khối lượng riêng không khí ($\approx 1.225$ kg/m³ ở mực nước biển)
- $C_d$ — hệ số cản (không thứ nguyên, phụ thuộc hình dạng vật: bóng cầu $\approx 0.47$, đĩa phẳng $\approx 1.28$)
- $A$ — diện tích mặt cắt ngang vật (m²)
Để gọn, ta gom thành hệ số cản hiệu dụng $k$:
$$ k = \frac{\rho C_d A}{2m} $$
Sao cho gia tốc do cản:
$$ \vec{a}_{drag} = -k \, |\vec{v}| \, \vec{v} \quad \text{(đơn vị: 1/m)} $$
Giá trị $k$ điển hình:
- Bóng tennis: $k \approx 0.005$–$0.01$
- Bóng rổ: $k \approx 0.005$
- Trái cầu lông: $k \approx 0.05$ (cản rất mạnh!)
Phần 3 — Phương trình vi phân: vì sao phải giải số?¶
Áp dụng định luật II Newton cho cả hai trục:
$$ m\frac{d\vec{v}}{dt} = m\vec{g} + \vec{F}_d $$
Tách thành phần $x$ và $y$, chia hai vế cho $m$:
$$ \begin{aligned} \frac{dv_x}{dt} &= -k \, |v| \, v_x \\ \frac{dv_y}{dt} &= -g - k \, |v| \, v_y \end{aligned} \quad\text{với}\quad |v| = \sqrt{v_x^2 + v_y^2} $$
Cộng với hai phương trình động học:
$$ \frac{dx}{dt} = v_x, \qquad \frac{dy}{dt} = v_y $$
Tổng cộng 4 phương trình vi phân bậc 1 với 4 ẩn $[x, y, v_x, v_y]$.
Vì sao không có nghiệm giải tích? Lực cản $-k|v|v$ chứa $|v| = \sqrt{v_x^2 + v_y^2}$ — đây là hàm phi tuyến ghép $v_x$ và $v_y$. Không thể tách biến để tích phân được như trường hợp không cản.
⇒ Phải dùng giải số: chia thời gian thành các bước nhỏ $\Delta t$, dùng máy tính cập nhật $[x, y, v_x, v_y]$ từng bước. Thư viện scipy.integrate.solve_ivp lo việc này, dùng thuật toán Runge-Kutta thích nghi.
def projectile_with_drag(t, state, k, g=9.81):
"""Hệ ODE bậc 1 cho ném xiên có lực cản bậc hai.
state = [x, y, vx, vy]
Trả về [dx/dt, dy/dt, dvx/dt, dvy/dt].
"""
x, y, vx, vy = state
speed = np.sqrt(vx**2 + vy**2)
return [vx, vy, -k * speed * vx, -g - k * speed * vy]
def hit_ground(t, state, k, g=9.81):
"""Sự kiện: y = 0 khi vật chạm đất (chỉ kích hoạt khi đang đi xuống)."""
return state[1]
hit_ground.terminal = True # dừng tích phân khi event xảy ra
hit_ground.direction = -1 # chỉ kích hoạt khi y giảm (đi xuống)
def solve_projectile(v0, theta_deg, k, g=9.81, t_max=20):
"""Wrapper: giải ODE cho v0, góc, hệ số cản — trả quỹ đạo + tầm xa + h_max.
Tham số:
v0 — vận tốc ban đầu (m/s)
theta_deg — góc ném (độ)
k — hệ số cản (1/m)
g — gia tốc trọng trường (m/s²)
t_max — thời gian tối đa cho ODE (giây), bảo hiểm khi chưa rơi
Trả về dict: {t, x, y, range, h_max, t_flight}
"""
theta = np.deg2rad(theta_deg)
state0 = [0, 0, v0 * np.cos(theta), v0 * np.sin(theta)]
sol = solve_ivp(
projectile_with_drag, (0, t_max), state0,
args=(k, g), events=hit_ground,
dense_output=True, max_step=0.01, rtol=1e-8,
)
# Lấy quỹ đạo dày để vẽ mượt
t_dense = np.linspace(0, sol.t[-1], 500)
x_dense, y_dense = sol.sol(t_dense)[0], sol.sol(t_dense)[1]
return {
't': t_dense, 'x': x_dense, 'y': y_dense,
'range': float(x_dense[-1]), 'h_max': float(y_dense.max()),
't_flight': float(sol.t[-1]),
}
# Test với k=0.01 (gần bóng rổ)
res = solve_projectile(v0=20, theta_deg=45, k=0.01)
print(f"Tầm xa = {res['range']:.2f} m")
print(f"Độ cao max = {res['h_max']:.2f} m")
print(f"Thời gian = {res['t_flight']:.2f} s")
Tầm xa = 31.32 m Độ cao max = 8.78 m Thời gian = 2.67 s
# So sánh 4 quỹ đạo: không cản + 3 mức cản tăng dần
v0, theta = 20, 45
_, x_nd, y_nd, *_ = projectile_no_drag(v0, theta)
ks = [0.005, 0.02, 0.05]
colors = ['#2ecc71', '#f39c12', '#e74c3c']
plt.plot(x_nd, y_nd, 'b-', linewidth=2.5, label='Không lực cản')
for k, c in zip(ks, colors):
res = solve_projectile(v0=v0, theta_deg=theta, k=k)
plt.plot(res['x'], res['y'], color=c, linewidth=2,
label=f'k = {k} (R = {res["range"]:.1f} m)')
plt.xlabel('x (m)'); plt.ylabel('y (m)')
plt.title(f'Ảnh hưởng của lực cản — v₀={v0} m/s, θ={theta}°')
plt.legend(loc='upper right')
plt.grid(alpha=0.3); plt.axis('equal')
plt.show()
Phần 4 — Khám phá: lực cản làm gì?¶
Quan sát biểu đồ trên, em có thấy:
| Hệ số cản $k$ | Loại vật | Tầm xa giảm | Hình dạng quỹ đạo |
|---|---|---|---|
| 0 (không cản) | lý tưởng | 0% | parabola đối xứng |
| 0.005 | bóng tennis nhẹ | nhỏ (~10%) | gần parabola, đỉnh hơi lệch trái |
| 0.02 | bóng cứng nặng | trung bình (~30%) | rõ bất đối xứng — đi lên thoải, rơi dốc |
| 0.05 | trái cầu lông | rất lớn (>50%) | đi lên cong, rơi gần như thẳng đứng |
Quan sát: bất đối xứng "đi lên ≠ rơi xuống"¶
Khi không cản, quỹ đạo parabola đối xứng: thời gian đi lên = thời gian rơi xuống, đỉnh đúng giữa quỹ đạo.
Khi có cản, vật đi lên thoải hơn rồi rơi xuống dốc hơn:
- Khi đi lên: trọng lực + cản đều cản chuyển động ⇒ vật chậm rất nhanh.
- Khi rơi xuống: trọng lực kéo xuống, cản kéo lên ⇒ vật đạt vận tốc tới hạn $v_t$ rồi rơi đều (giống hạt mưa).
- Vận tốc nhỏ hơn ⇒ vật bay ngang chậm hơn khi rơi.
Câu hỏi cần trả lời¶
Trong sách lớp 10, góc tối ưu để vật bay xa nhất là $45°$. Khi có lực cản, góc tối ưu có còn là $45°$ không? Hãy đoán trước khi xem phần dưới — em nghĩ $\theta_{opt}$ sẽ lớn hơn hay nhỏ hơn $45°$?
# Quét nhiều góc với cùng v0 và k cố định
v0_test, k_test = 20, 0.02
thetas = [30, 40, 45, 50, 60]
plt.figure(figsize=(10, 6))
for th in thetas:
res = solve_projectile(v0=v0_test, theta_deg=th, k=k_test)
plt.plot(res['x'], res['y'], linewidth=2,
label=f'θ = {th}° (R = {res["range"]:.1f} m)')
plt.xlabel('x (m)'); plt.ylabel('y (m)')
plt.title(f'Quét góc ném — v₀={v0_test} m/s, k={k_test}')
plt.legend(loc='upper right')
plt.grid(alpha=0.3); plt.axis('equal')
plt.show()
def find_optimal_angle(v0, k, theta_min=20, theta_max=70, n=51):
"""Quét góc trong khoảng [theta_min, theta_max] để tìm góc cho tầm xa lớn nhất.
Trả về (thetas, ranges, theta_opt, range_opt).
"""
thetas = np.linspace(theta_min, theta_max, n)
ranges = np.array([
solve_projectile(v0, th, k)['range'] for th in thetas
])
i_opt = int(np.argmax(ranges))
return thetas, ranges, thetas[i_opt], ranges[i_opt]
# Quét 5 mức k để xem θ_opt thay đổi thế nào
v0_scan = 20
k_values = [0, 0.005, 0.01, 0.02, 0.05]
colors = plt.cm.viridis(np.linspace(0, 0.85, len(k_values)))
plt.figure(figsize=(10, 6))
for k, c in zip(k_values, colors):
if k == 0:
# Công thức giải tích cho k=0 — nhanh hơn solve_ivp
thetas = np.linspace(20, 70, 51)
ranges = v0_scan**2 * np.sin(2 * np.deg2rad(thetas)) / 9.81
th_opt, R_opt = 45.0, ranges.max()
else:
thetas, ranges, th_opt, R_opt = find_optimal_angle(v0_scan, k)
plt.plot(thetas, ranges, color=c, linewidth=2,
label=f'k = {k} → θ_opt ≈ {th_opt:.0f}° (R = {R_opt:.1f} m)')
plt.scatter([th_opt], [R_opt], color=c, s=80, zorder=5, edgecolor='black')
plt.xlabel('Góc ném θ (độ)')
plt.ylabel('Tầm xa R (m)')
plt.title(f'Tầm xa theo góc ném — v₀={v0_scan} m/s, đa mức cản')
plt.legend(loc='lower center', fontsize=9)
plt.grid(alpha=0.3)
plt.show()
print('\nKết luận: lực cản càng lớn (k tăng), góc tối ưu θ_opt càng GIẢM dưới 45°.')
print('Lý do: thời gian bay càng dài thì cản tích lũy càng nhiều.')
print('⇒ Góc thấp hơn → đường đi nhanh hơn → cản tích lũy ít hơn.')
Kết luận: lực cản càng lớn (k tăng), góc tối ưu θ_opt càng GIẢM dưới 45°. Lý do: thời gian bay càng dài thì cản tích lũy càng nhiều. ⇒ Góc thấp hơn → đường đi nhanh hơn → cản tích lũy ít hơn.
Bài tập¶
Bài 1 — Bóng rổ ngoài trời¶
Cho bóng rổ với:
- Khối lượng $m = 0.6$ kg
- Bán kính $r = 0.12$ m (⇒ diện tích mặt cắt $A = \pi r^2$)
- Hệ số cản $C_d = 0.5$ (cầu nhẵn)
- Không khí mực nước biển $\rho = 1.225$ kg/m³
Yêu cầu:
- Tính hệ số cản hiệu dụng $k = \dfrac{\rho C_d A}{2m}$.
- Một cầu thủ ném bóng với $v_0 = 9$ m/s, góc $\theta = 50°$. So sánh tầm xa có và không có lực cản. Sai số bao nhiêu phần trăm?
Gợi ý: dùng solve_projectile() ở Phần 3.
Bài 2 — Tìm góc tối ưu cho bóng rổ¶
Vẫn với bóng rổ ở BT1, $v_0 = 9$ m/s. Tìm góc ném tối ưu $\theta_{opt}$ sao cho tầm xa lớn nhất.
Gợi ý: gọi find_optimal_angle() với $k$ vừa tính ở BT1.
So với góc $45°$ "sách giáo khoa", $\theta_{opt}$ cao hơn hay thấp hơn? Bao nhiêu độ?
Bài 3 ★ — Vận tốc tới hạn của hạt mưa¶
Một hạt mưa hình cầu rơi tự do từ mây cao 2 km xuống đất:
- Bán kính $r = 1$ mm $= 0.001$ m
- Khối lượng riêng nước $\rho_{nước} = 1000$ kg/m³ (⇒ tính khối lượng từ thể tích cầu)
- Không khí $\rho = 1.225$ kg/m³, $C_d = 0.47$
Yêu cầu:
- Khi vận tốc hạt mưa không đổi (rơi đều), cản cân bằng trọng lực: $$ m g = \tfrac{1}{2} \rho C_d A v_t^2 $$ Suy ra vận tốc tới hạn: $v_t = \sqrt{\dfrac{2 m g}{\rho C_d A}}$. Tính $v_t$ (đáp án dạng m/s).
- Để xác nhận, dùng
solve_projectile()với $\theta = -90°$ (rơi thẳng đứng), $v_0 = 0$, $k$ tính theo công thức trên. Vẽ đồ thị $v_y(t)$ và quan sát giá trị bão hòa.
Gợi ý: bài này khó hơn — em sẽ cần lắp ghép kiến thức Phần 2 và viết code mới một chút.
Lời giải KHÔNG được cung cấp ở đây. Em hãy thử trước, sau đó so sánh với bạn cùng lớp hoặc giáo viên.
🎯 Tổng kết¶
Sau bài này, em đã:
- Phân biệt ném xiên lý tưởng (parabola, có công thức giải tích) với ném xiên thực tế (cần giải số).
- Lập hệ 4 phương trình vi phân từ định luật II Newton + lực cản bậc hai.
- Sử dụng
scipy.integrate.solve_ivpvới event để dừng khi vật chạm đất. - Quan sát sự bất đối xứng "đi lên ≠ rơi xuống" và bản chất vận tốc tới hạn.
- Khám phá rằng góc ném tối ưu $\theta_{opt}$ giảm dưới $45°$ khi lực cản tăng.
🚀 Liên hệ thực tế¶
- Pháo binh: tầm xa thực tế tối ưu ở góc $\sim 35°$–$42°$, không phải $45°$ như sách dạy.
- Thể thao: cầu thủ bóng rổ ném 3 điểm với góc $\sim 50°$–$55°$ vì cần độ cao đỉnh để bóng rơi vào rổ — đây là tối ưu khác (không phải tầm xa).
- Khí động học: hệ số cản $C_d$ là biến số chính trong thiết kế xe đua, máy bay, tên lửa.
📚 Thử thách mở rộng¶
Em muốn đi xa hơn? Thử các phần mở rộng sau:
- Gió ngang: thêm vector gió $\vec{w}$ vào vận tốc tương đối $\vec{v}_{rel} = \vec{v} - \vec{w}$, lực cản tính theo $\vec{v}_{rel}$.
- Hiệu ứng Magnus: vật quay tròn (như bóng đá xoáy) sinh lực vuông góc với vận tốc — thử mô phỏng "cú đá knuckleball".
- Khí quyển không đồng nhất: $\rho$ giảm theo độ cao $h$ (mô hình $\rho(h) = \rho_0 e^{-h/H}$ với $H \approx 8400$ m). Cần thiết khi mô phỏng đạn pháo bắn cao.
- Pint check dimensional: dùng
pintđể xác minh đơn vị của $k$ thực sự là 1/m, $v_t$ thực sự là m/s.
Vật lý số không chỉ là "máy tính giải bài thay em" — nó là công cụ để khám phá những điều mà công thức giải tích không thể tới được.