简介《理论力学哈工大版》PPT课件面向高校工科学生及相关自学者聚焦质点动力学核心章节梳理静力学、运动学与动力学的关系系统讲解牛顿三定律、国际单位制、质点运动微分方程及动力学的两类基本问题有助于建立理论力学分析与解题框架。压缩包共1个PPT文件大小约1010KB单文件课件适合课堂讲解与自主复习便于直接投屏浏览或打印成讲义。内容以第十章质点动力学为主线配有曲柄连杆机构等典型例题逐步演示如何由运动方程求力、由力求运动并说明直角坐标、自然坐标等坐标系的选用思路可帮助读者理解抽象概念并应用于工程实际问题。资源已有5626人学习及下载适合正在学习理论力学或准备考研的工科生用于概念巩固与例题研习。1. 理论力学哈工大版为何仍是仿真工程师的底层接口物理引擎跑飞了、机械臂动力学模型输出 NaN、卫星姿态仿真发散——这些问题表面出在数值积分步长上深挖下去多半是力学建模那一步就偏了。哈尔滨工业大学理论力学教研组编的《理论力学》被工科院校使用了几十年教材从静力学、运动学到动力学始终在做同一件事把一个机械问题逐句化成不依赖直觉的方程。反过来看IT 工程师调试仿真代码时需要的正是同一能力把受力关系、运动约束转换成计算机能迭代的数值模型。这本书适合两类人一类是工作中需要补刚体动力学基础写物理引擎、做多体仿真的另一类是学生时代学过高等数学现在要靠自学把力学捡起来目的不是考试而是能读懂开源动力学库里的偏导与雅可比。这篇文章的重点不在复述教材例题而是给出一条可执行路径——顺着哈工大版的章节顺序每一步配上能跑的 Python 代码让「教材推导」与「数值仿真」在同一个工作流里对上账。2. 哈工大版静力学的数值化力系简化与姿态矩阵计算2.1 力系简化在哈工大版中的本质是向量投影哈工大版静力学篇对平面任意力系的标准处理流程是先把所有力沿直角坐标轴分解再分别对某点取矩得到主矢和主矩。对一个写代码的人来说这套流程天然就是矩阵与数组运算。设平面上有 n 个力每个力可以用分量描述作用点坐标、Fx、Fy以及该力对参考点产生的矩 Mz。整组力系的主矢就是 Fx 列与 Fy 列的求和主矩就是 Mz 列的求和。教材里反复强调的「力线平移」在数值上等价于给每个力的矩分量补上一项 ±F×d。这个视角最直接的好处是把「力的多边形法则」这种几何操作变成可复现的数组运算。调参、增力、删力都只改一行数据。图解法再熟练遇到三十个力的复杂桁架也没法手工画。2.2 用 NumPy 实现平面任意力系简化下面这段代码接受一组力返回该力系对指定简化中心的主矢和主矩import numpy as np def reduce_force_system(forces, ref_point(0.0, 0.0)): 平面任意力系简化 forces: list of (x, y, Fx, Fy) 每个元素表示力的作用点坐标以及力在 x/y 方向的分量 ref_point: 简化中心坐标 返回: (主矢, 主矩) R np.zeros(2) M 0.0 x0, y0 ref_point for (x, y, fx, fy) in forces: R np.array([fx, fy]) # 力对简化中心取矩M_z rx * Fy - ry * Fx M (x - x0) * fy - (y - y0) * fx return R, M # 例三力作用于同一构件求对点 O 的简化结果 forces [ (0.0, 0.0, 100.0, 0.0), # 沿 x 方向的 100N 力作用在原点 (2.0, 0.0, 0.0, -50.0), # 构件右端向下的 50N 力 (1.0, 1.0, 30.0, 40.0), # 斜向力分量 30N / 40N ] R, M reduce_force_system(forces) print(主矢 , R, 主矩 , M)这段代码里矩的计算用的是教材上的简化中心定义先求力作用点到简化中心的矢径 r再算 r × F 的 z 分量。注意正负号由叉乘顺序决定哈工大版默认逆时针为正代码沿用这个约定。实际工程中如果力的数据来自 CAD 或传感器坐标系务必先确认坐标轴方向和单位制否则符号问题会在后续动力学求解中被成倍放大。输出结果中主矢是二维数组、主矩是标量。判断力系平衡只需要同时检查两者的模是否为零。浮点场景下不要用 0.0判断一般取 1e-6 以下的容差即可。2.3 运动学描述的两种路线与选型哈工大版运动学篇把刚体运动分成平动、定轴转动和平面一般运动。描述平面运动刚体上一点的速度教材给了两条路线矢量法和解析法。两种方法在代码里呈现的风格完全不同对比项矢量法解析法数学工具矢量叉乘、方向余弦链式求导、雅可比矩阵适用场景四连杆、曲柄滑块机构机械臂、多体树状结构方程形式显式矢量方程递推的位置/速度关系数值稳定性依赖几何构造精度依赖导数计算完备性教材典型工具速度瞬心法、基点法直角坐标法、自然轴系法写仿真代码时我会把两种方法交叉验证机构简单时先用手推的矢量方程做真值再让计算机用解析法求雅可比两者误差不超过 1e-9 才继续向下做动力学。这个方法能很快发现坐标系正方向定义不一致的问题。2.4 姿态矩阵与欧拉角的数值表达进入动力学之前姿态表示是绕不开的一关。哈工大版在空间运动部分给出了方向余弦矩阵、欧拉角与角速度矢量之间的关系。平面运动里一个角度 θ 就够到了空间运动方向余弦矩阵 R 的三列分别是动系三轴在定系中的投影。绕 z 轴旋转 θ 的旋转矩阵是import numpy as np def rot_z(theta): c, s np.cos(theta), np.sin(theta) return np.array([[c, -s, 0], [s, c, 0], [0, 0, 1]])三个方向余弦矩阵连乘得到完整姿态。连乘顺序对应旋转顺序不同教材的约定不一样自己推导矩阵乘法时很容易写反。实际写代码我建议统一用 scipy.spatial.transform.Rotation 的from_euler直接指定旋转顺序例如xyz表示依次绕定系 x、y、z 轴旋转把次序交给库去管避免手写矩阵乘出差错。欧拉角的万向锁问题在仿真里很常见。判断条件是第二个旋转角接近 ±90° 时第一次与第三次旋转不可区分此时继续用欧拉角做四元数变换会输出 NaN。处理方式通常是切换到四元数表示或改变旋转顺序。3. 用 SymPy 复现哈工大推导平面运动刚体的瞬心与角速度合成3.1 符号计算如何配合哈工大版矢量推导哈工大版教材的推导有一个特点坐标系的选取非常讲究选得好公式极简选不好代数展开能占满三页纸。手工推导时经常因为丢一项白算一个小时。用 SymPy 做符号推导不是要替代手工能力而是把「选坐标系」这个物理判断保留给人把「展开化简」这个机械劳动交给机器。对平面机构我习惯先定义所有几何变量为符号比如杆长 l1、l2角度 θ1、θ2直接写出位置约束方程再对时间求导得到速度与加速度关系。这正是哈工大版运动学解析法的标准路径区别只在把纸笔换成符号引擎推导过程每一步都可审计。3.2 速度瞬心法的 SymPy 验证先看教材里最经典的速度瞬心法。平面运动中每一瞬时都存在速度为零的点——速度瞬心刚体上任意一点的速度等于角速度矢量叉乘该点相对瞬心的矢径。用 SymPy 验证曲柄滑块机构的角速度关系import sympy as sp from sympy.physics.mechanics import dynamicsymbols t dynamicsymbols._t theta1, theta2 dynamicsymbols(theta1 theta2) l1, l2 sp.symbols(l1 l2, positiveTrue) # A 为曲柄固定铰链B 为曲柄与连杆连接点C 为滑块 A sp.Matrix([0, 0]) B A sp.Matrix([l1 * sp.cos(theta1), l1 * sp.sin(theta1)]) C B sp.Matrix([l2 * sp.cos(theta2), l2 * sp.sin(theta2)]) # 滑块在水平导轨上运动C 的纵向位移恒为 0 constraint sp.Eq(C[1], 0) # 对时间求导得到速度约束 vel_constraint sp.diff(constraint.lhs - constraint.rhs, t) # 解出连杆角速度 theta2 theta2_dot sp.solve(vel_constraint, sp.diff(theta2, t)) print(theta2_dot)这段代码先写出位置闭链方程再对时间求导得到的符号解就是连杆角速度与曲柄角速度的关系。l1、l2用positiveTrue声明能帮助 SymPy 在化简时去除绝对值等冗余分支theta1、theta2必须用dynamicsymbols声明为时间函数否则求导链会断。如果解出来的表达式里残留Derivative(theta1(t), t)说明你把时间函数写成了普通符号需要回头改声明方式。3.3 角速度合成定理的参考系转换角速度合成定理的教材表述是绝对角速度 牵连角速度 相对角速度。用 SymPy 验证符号关系容易但工程里真正的坑在于这个矢量加法在三维空间必须统一参考系不能在分量上直接相加。import numpy as np def cross_matrix(w): 把角速度向量转成叉乘矩阵用于与向量做矢量积 return np.array([[0, -w[2], w[1]], [w[2], 0, -w[0]], [-w[1], w[0], 0]]) # 牵连角速度载体自身的旋转 omega_e np.array([0.0, 0.0, 1.2]) # 相对角速度动系中测得的关节角速度 omega_r_body np.array([0.5, 0.0, 0.0]) # 动系相对定系绕 z 轴旋转 30 度方向余弦矩阵 theta np.deg2rad(30.0) R np.array([[np.cos(theta), -np.sin(theta), 0], [np.sin(theta), np.cos(theta), 0], [0, 0, 1]]) # 相对角速度转换到定系后再合成 omega_r_world R omega_r_body omega_a_world omega_e omega_r_world print(绝对角速度(定系) , omega_a_world)注意这里的合成顺序如果相对角速度是在动系里表示的必须先左乘方向余弦矩阵 R 转到定系再与牵连角速度相加。哈工大版教材在推导过程中默认各角速度都在同一参考系下表示新手照搬公式时容易漏掉参考系转换这一步导致机械臂末端速度场错得离谱。cross_matrix函数在后续处理陀螺力矩时也会用到提前封装好。3.4 瞬心法在数值仿真中的失效边界速度瞬心法在机构学里非常好用但它只在「瞬时」成立。进入动力学仿真后瞬心位置随时间变化而且当刚体接近纯平动时瞬心会跑到无穷远处数值上表现为一个巨大的坐标值。此时若继续用瞬心法计算速度浮点误差会剧烈放大。我一般会在程序里加保护当角速度模长小于某个阈值比如 1e-4 rad/s时自动切换到平动假设把瞬心法退化为基点法计算。这个判断放在状态更新函数之前防止积分器在接近奇异的位置浪费时间。4. 把哈工大版分析力学落地为仿真循环拉格朗日方程与 RK4 数值积分4.1 为什么拉格朗日方程更适合写进代码牛顿欧拉法对每个刚体列写力与力矩平衡方程铰链内部的约束力也作为未知量出现在方程组中拉格朗日方程把系统看作一个整体用广义坐标和动能/势能写出运动方程理想约束力自动消去。对机械臂路径规划、卫星姿控这类场景约束力本身常常是需要求解的值但建模阶段先消去约束状态变量更少代码也更不容易因为约束方程写错而崩溃。哈工大版把拉格朗日方程放在动力学与分析力学部分标准形式是d/dt(∂L/∂q̇) − ∂L/∂q Q其中 L T − VQ 是非保守广义力。用 SymPy 推导的好处是动能和势能写完后所有偏导自动生成不会出现漏项。4.2 用 SymPy 推导二连杆广义坐标方程以平面二连杆为例两杆质量 m1、m2长度 l1、l2广义坐标取关节角 θ1、θ2。先写出动能与势能import sympy as sp from sympy.physics.mechanics import dynamicsymbols t dynamicsymbols._t theta1, theta2 dynamicsymbols(theta1 theta2) l1, l2, m1, m2, g sp.symbols(l1 l2 m1 m2 g, positiveTrue) # 杆1质心位置近似取杆中点 x1 (l1 / 2) * sp.cos(theta1) y1 (l1 / 2) * sp.sin(theta1) # 杆2质心位置受两个关节角共同影响 x2 l1 * sp.cos(theta1) (l2 / 2) * sp.cos(theta1 theta2) y2 l1 * sp.sin(theta1) (l2 / 2) * sp.sin(theta1 theta2) # 速度平方 v1_sq sp.diff(x1, t) ** 2 sp.diff(y1, t) ** 2 v2_sq sp.diff(x2, t) ** 2 sp.diff(y2, t) ** 2 # 动能细杆相对质心的转动惯量取 m*l^2/12先加平动项 T sp.Rational(1, 2) * m1 * v1_sq sp.Rational(1, 2) * m2 * v2_sq V m1 * g * y1 m2 * g * y2 L T - V print(sp.simplify(L))这里用dynamicsymbols让 θ1、θ2 自动成为时间函数求导不会断链。动能里暂时只写了平动项若考虑细杆绕质心的转动需要为每根杆补充Rational(1, 2) * I * omega**2。表达式中用sp.Rational(1, 2)而不是1 / 2因为 Python 的浮点除法会污染符号表达式导致后续化简变慢或精度下降。下一步由拉格朗日方程得到运动方程再用sp.solve解出二阶导数# 对两个广义坐标分别列拉格朗日方程 lag_eq1 sp.diff(sp.diff(L, theta1.diff(t)), t) - sp.diff(L, theta1) lag_eq2 sp.diff(sp.diff(L, theta2.diff(t)), t) - sp.diff(L, theta2) sols sp.solve([lag_eq1, lag_eq2], [theta1.diff(t, 2), theta2.diff(t, 2)])sols是一个字典键是 θ1、θ2 的二阶导符号值是对应的符号表达式。这一步输出的是显式二阶常微分方程组可以直接喂给数值积分器。注意求解得到的表达式里可能包含大量 cos/sin 组合不需要手工化简交给下一步的 lambdify 即可。4.3 把符号方程 lambdify 成数值函数并跑 RK4符号表达式要变成可调用的数值函数用sp.lambdifyf_func sp.lambdify( (theta1, theta2, theta1.diff(t), theta2.diff(t), l1, l2, m1, m2, g), [sols[theta1.diff(t, 2)], sols[theta2.diff(t, 2)]], numpy )lambdify的第一个参数是自变量元组必须包含所有符号变量第二个参数是输出表达式列表第三个参数numpy指定后端让输出函数能直接接受 numpy 数组做批量求值。这一步把符号层的成果编译成数值层可用的函数后续所有求解都在数值层完成不再碰符号。然后写标准 RK4 积分器。状态向量取[theta1, theta2, omega1, omega2]import numpy as np def state_deriv(state, t, params): th1, th2, w1, w2 state l1, l2, m1, m2, g params a1, a2 f_func(th1, th2, w1, w2, l1, l2, m1, m2, g) return np.array([w1, w2, a1, a2]) def rk4_step(f, state, t, dt, params): k1 f(state, t, params) k2 f(state 0.5 * dt * k1, t 0.5 * dt, params) k3 f(state 0.5 * dt * k2, t 0.5 * dt, params) k4 f(state dt * k3, t dt, params) return state dt / 6.0 * (k1 2 * k2 2 * k3 k4)rk4_step里四个 k 值对应龙格库塔四阶法的四次斜率估计。dt 的选择要同时满足数值稳定性与精度对杆长 1 米、质量 1 千克量级的系统dt 取 1e-3 秒通常足够。如果能量曲线出现高频抖动先减半 dt 观察结果是否显著变化不要一上来就换更高阶积分器。4.4 教材例题向仿真参数迁移的注意事项哈工大版课后题给的大部分是理想化参数不计摩擦、不计阻尼、杆件视为均质细杆。做仿真时建议先用教材例题的原始参数跑一遍确认自己的数值结果和教材解析解对得上再逐步加入阻尼项和驱动项。教材答案通常保留 2~3 位小数对照时允许 1e-2 量级的误差如果量级明显不对多半是符号推导时丢了某个交叉项而不是数值积分的问题。5. 收尾技巧用能量守恒校验哈工大版例题的数值化边界5.1 能量漂移检测数值积分一定会引入能量漂移RK4 在无阻尼系统中能量误差通常随时间缓慢增长。检验仿真代码正确性的最直接办法是打印总能量 TV 的时间序列def total_energy(state, params): th1, th2, w1, w2 state l1, l2, m1, m2, g params # 建议直接用 4.2 节符号表达式 lambdify 得到 T 和 V 的数值函数 # 这里给出手写形式来展示结构 T_lin 0.5 * m1 * ((l1 / 2) * w1) ** 2 V_lin m1 * g * (l1 / 2) * np.sin(th1) return T_lin_sp V_lin_sp # 实际替换为完整表达式如果能量在 1 秒仿真时间内变化超过初始值的 1%先怀疑积分步长再怀疑约束方程是否被正确消去最后才检查代码本身的浮点问题。教材中的无阻尼系统理论上能量守恒所以这个校验在哈工大版例题上特别有效。5.2 教材例题当成回归测试的三种做法第一种抽 3 个典型例题——一个静力学力系平衡、一个运动学速度瞬心求速度、一个动力学拉格朗日方程或动量矩定理把教材的数值答案写成断言放进 CI每次改完代码跑一遍。第二种保留 Symbolic 推导脚本任何模型修改都先让符号层输出新方程再重新 lambdify绝不手改数值函数避免两边不一致。第三种记录每次仿真的能量漂移率作为指标如果某次改动导致漂移率翻倍说明这次改动大概率破坏了动力学模型的一致性。5.3 固定验证顺序把前面所有步骤串成一个固定顺序先跑静力学例对比主矢主矩再跑平面运动学例对比瞬心速度最后跑完整动力学仿真对比能量曲线。三道全过哈工大版的核心内容就真正变成了你自己的调试能力。同样一个四连杆机构用瞬心法算出的速度分布和完整的拉格朗日方程仿真结果如果相差超过 1e-6不要急着在数值积分里找原因先回头检查第 3.2 节的约束在符号层是否被正确求导。这个检查顺序节省的调试时间往往比重读整章教材更有效。本文还有配套的精品资源点击获取
