别再死记硬背！用Python模拟刚体转动，轻松理解角速度、转动惯量（附Jupyter代码）

用Python模拟刚体转动从代码中理解物理本质刚体转动是经典力学中令人着迷又颇具挑战性的主题。传统教学中学生们往往被要求死记硬背各种转动惯量公式在抽象的例题中挣扎。但物理本质上是一门实验科学而今天我们可以用Python代码搭建数字实验室通过动画和交互式模拟直观感受角速度、力矩与转动惯量的动态关系。1. 为什么用代码学习刚体转动翻开任何一本大学物理教材刚体转动章节通常充斥着积分符号和代数推导。虽然数学严谨但这种呈现方式容易让学习者陷入符号操作而忽略物理图像。当我们用Python将转动过程可视化时抽象概念突然变得触手可及角速度不再只是ω符号而是可以看到圆盘转速变化的动画转动惯量不再只是Jmr²而是能观察到质量分布改变时转动状态的实时响应力矩的作用效果可以直接对比不同施力方式的影响import numpy as np import matplotlib.pyplot as plt from matplotlib.animation import FuncAnimation # 初始化圆盘参数 radius 1.0 mass 2.0 angular_velocity 0.0 torque 0.5 # 恒定力矩 moment_of_inertia 0.5 * mass * radius**2 # 实心圆盘转动惯量 dt 0.05 # 时间步长 # 创建图形和轴 fig, ax plt.subplots(figsize(8, 8)) ax.set_xlim(-1.5, 1.5) ax.set_ylim(-1.5, 1.5) ax.set_aspect(equal) ax.grid(True) # 绘制静止圆盘 circle plt.Circle((0, 0), radius, fillFalse, linewidth2) ax.add_patch(circle) # 标记参考点 point, ax.plot([0, radius], [0, 0], r-) angle_text ax.text(0.05, 0.95, , transformax.transAxes) velocity_text ax.text(0.05, 0.90, , transformax.transAxes) def update(frame): global angular_velocity # 计算角加速度 (τ Iα) angular_acceleration torque / moment_of_inertia # 更新角速度 (ω ω₀ αt) angular_velocity angular_acceleration * dt # 计算角位移 (θ θ₀ ωt) angle angular_velocity * frame * dt # 更新圆盘显示 point.set_data([0, radius*np.cos(angle)], [0, radius*np.sin(angle)]) angle_text.set_text(f角度: {np.degrees(angle):.1f}°) velocity_text.set_text(f角速度: {angular_velocity:.2f} rad/s) return point, angle_text, velocity_text ani FuncAnimation(fig, update, frames200, interval50, blitTrue) plt.show()这段基础代码模拟了恒定力矩作用下实心圆盘的加速转动过程。运行后会看到红色标记线随圆盘转动而移动实时显示当前转角和角速度数值角速度随时间线性增加匀角加速转动2. 转动惯量的可视化理解转动惯量是刚体转动中的质量它决定了物体抵抗转动状态改变的难易程度。传统教学中学生需要记忆各种几何形状的转动惯量公式但通过模拟可以直观理解其物理意义。2.1 质量分布的影响让我们修改前面的代码比较不同质量分布的转动惯量def calculate_moment_of_inertia(shape, mass, radius): 计算不同形状的转动惯量 if shape disk: return 0.5 * mass * radius**2 elif shape ring: return mass * radius**2 elif shape sphere: return 0.4 * mass * radius**2 else: raise ValueError(未知形状) # 比较相同质量不同形状的角加速度 mass 2.0 radius 1.0 torque 0.5 shapes [disk, ring, sphere] for shape in shapes: I calculate_moment_of_inertia(shape, mass, radius) alpha torque / I print(f{shape}的转动惯量: {I:.2f}, 角加速度: {alpha:.2f} rad/s²)输出结果disk的转动惯量: 1.00, 角加速度: 0.50 rad/s² ring的转动惯量: 2.00, 角加速度: 0.25 rad/s² sphere的转动惯量: 0.80, 角加速度: 0.62 rad/s²这个简单实验展示了质量分布离轴越远转动惯量越大圆环圆盘球体相同力矩下转动惯量越大角加速度越小2.2 平行轴定理的验证平行轴定理指出刚体绕任意轴的转动惯量等于绕通过质心的平行轴的转动惯量加上质量与两轴距离平方的乘积$$ I I_{\text{cm}} md^2 $$我们可以用数值计算验证这一定理def rod_moment_of_inertia(length, mass, distance_from_end): 计算细棒绕不同轴的转动惯量 # 绕端点 if distance_from_end 0: return (1/3) * mass * length**2 # 绕中心 elif distance_from_end length/2: return (1/12) * mass * length**2 # 绕任意平行轴 else: I_cm (1/12) * mass * length**2 d abs(length/2 - distance_from_end) return I_cm mass * d**2 # 参数设置 length 2.0 # 细棒长度 mass 1.0 # 细棒质量 positions np.linspace(0, length, 20) # 从一端到另一端的转轴位置 # 计算各位置的转动惯量 I_values [rod_moment_of_inertia(length, mass, p) for p in positions] # 绘制结果 plt.figure(figsize(10, 6)) plt.plot(positions, I_values, b-, linewidth2) plt.axvline(xlength/2, colorr, linestyle--, label质心位置) plt.xlabel(转轴位置 (m), fontsize12) plt.ylabel(转动惯量 (kg·m²), fontsize12) plt.title(细棒转动惯量与转轴位置的关系, fontsize14) plt.grid(True) plt.legend(fontsize12) plt.show()这张图清晰展示了当转轴位于细棒端点时转动惯量最大转动惯量随转轴向质心移动而减小在质心处达到最小值$I_{cm}$转轴偏离质心时转动惯量对称增加3. 动态变化的转动惯量经典例题小虫落在旋转细杆上展示了转动惯量动态变化的场景。让我们用Python模拟这一过程# 模拟小虫在细杆上的运动 length 2.0 # 细杆长度 mass_rod 0.5 # 细杆质量 mass_bug 0.1 # 小虫质量 omega0 2.0 # 初始角速度 # 小虫运动参数 bug_speed 0.05 # 小虫爬行速度 (m/s) total_time 20.0 # 总模拟时间 dt 0.1 # 时间步长 times np.arange(0, total_time, dt) # 初始化数组存储结果 positions [] omegas [] # 初始条件 bug_pos length/2 # 小虫起始位置细杆中点 omega omega0 for t in times: # 计算当前转动惯量 (细杆 小虫) I_rod (1/12) * mass_rod * length**2 I_bug mass_bug * bug_pos**2 I_total I_rod I_bug # 角动量守恒: L Iω 常数 L I_total * omega # 小虫移动 bug_pos - bug_speed * dt if bug_pos 0: bug_pos 0 # 小虫到达端点 # 更新转动惯量 I_rod (1/12) * mass_rod * length**2 I_bug mass_bug * bug_pos**2 I_total_new I_rod I_bug # 根据角动量守恒计算新角速度 omega L / I_total_new # 存储结果 positions.append(bug_pos) omegas.append(omega) # 绘制结果 fig, (ax1, ax2) plt.subplots(2, 1, figsize(10, 8), sharexTrue) ax1.plot(times, positions, b-) ax1.set_ylabel(小虫位置 (m), fontsize12) ax1.grid(True) ax2.plot(times, omegas, r-) ax2.set_xlabel(时间 (s), fontsize12) ax2.set_ylabel(角速度 (rad/s), fontsize12) ax2.grid(True) plt.suptitle(小虫爬行对细杆转动的影响, fontsize14) plt.show()模拟结果展示了两个关键现象随着小虫向转轴移动系统转动惯量减小根据角动量守恒角速度相应增加这种数值模拟比代数推导更直观地展示了物理过程特别是转动惯量与角速度的动态平衡关系。4. 复杂场景变力矩与非匀质刚体现实中的转动问题往往比教科书例题复杂得多。Python模拟可以轻松处理这些情况4.1 随时间变化的力矩# 变力矩下的转动模拟 def variable_torque(t): 随时间变化的力矩函数 if t 5: return 0.2 # 初始小力矩 elif 5 t 10: return 0.5 # 增大力矩 else: return -0.3 # 反向力矩 # 模拟参数 I 1.0 # 转动惯量 total_time 15.0 dt 0.05 times np.arange(0, total_time, dt) # 初始化 omega 0.0 theta 0.0 omegas [] thetas [] torques [variable_torque(t) for t in times] # 模拟循环 for t, torque in zip(times, torques): alpha torque / I omega alpha * dt theta omega * dt omegas.append(omega) thetas.append(theta) # 绘制结果 fig, (ax1, ax2, ax3) plt.subplots(3, 1, figsize(10, 10), sharexTrue) ax1.plot(times, torques, g-) ax1.set_ylabel(力矩 (N·m), fontsize12) ax1.grid(True) ax2.plot(times, omegas, b-) ax2.set_ylabel(角速度 (rad/s), fontsize12) ax2.grid(True) ax3.plot(times, thetas, r-) ax3.set_xlabel(时间 (s), fontsize12) ax3.set_ylabel(角位移 (rad), fontsize12) ax3.grid(True) plt.suptitle(变力矩下的转动响应, fontsize14) plt.show()这个模拟展示了力矩变化时角加速度的即时响应角速度是角加速度的时间积分角位移是角速度的时间积分反向力矩导致减速运动4.2 非匀质刚体的转动现实中刚体的质量分布往往不均匀。我们可以用离散化方法计算这类刚体的转动惯量# 非匀质细棒的转动惯量计算 length 2.0 # 细棒长度 n_points 1000 # 离散点数 positions np.linspace(0, length, n_points) # 非均匀质量密度函数 def density(x): 线性变化的质量密度 return 1.0 0.5 * x # 密度沿细棒线性增加 # 计算总质量 dx length / n_points dm_values density(positions) * dx total_mass np.sum(dm_values) # 计算绕端点的转动惯量 I_end np.sum(dm_values * positions**2) # 计算绕质心的转动惯量 # 首先计算质心位置 center_of_mass np.sum(dm_values * positions) / total_mass # 然后计算绕质心的转动惯量 I_cm np.sum(dm_values * (positions - center_of_mass)**2) print(f总质量: {total_mass:.3f} kg) print(f质心位置: {center_of_mass:.3f} m) print(f绕端点转动惯量: {I_end:.3f} kg·m²) print(f绕质心转动惯量: {I_cm:.3f} kg·m²) print(f验证平行轴定理: {I_cm total_mass*center_of_mass**2:.3f} ≈ {I_end:.3f})输出示例总质量: 2.000 kg 质心位置: 1.167 m 绕端点转动惯量: 2.667 kg·m² 绕质心转动惯量: 0.556 kg·m² 验证平行轴定理: 2.667 ≈ 2.667这种方法可以推广到任意质量分布的刚体只需修改density函数即可。对于更复杂的二维或三维刚体可以采用类似的离散化方法。