SymPy Physics/Mechanics 中的 Kane 法KanesMethod 构建多体系统运动方程的完整指南【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympy本文基于 SymPy 官方文档 Kanes Method in Physics/Mechanics 撰写完整梳理sympy.physics.mechanics模块中 Kane 法的方程体系5 类基本方程组、KanesMethod对象的构造流程与参数、非完整约束与辅助速度的处理方式并结合 kane.py 源码实现解释质量矩阵与力矢量forcing vector的生成原理、linearize线性化接口与线性求解器选择。读完本文你将能够针对含坐标/速度依赖关系与运动约束的多体系统用 SymPy 符号推导出形如M(q,t)·udot f(q, qdot, u, t)的运动方程并转换为数值积分器可用的状态空间形式。一、Kane 法在 SymPy 中的定位sympy.physics.mechanics 模块的目标是为复杂多体系统创建符号运动方程。根据 经典力学文档索引前向动力学forward dynamics问题通常分为 5 步描述系统的几何与构型指定系统如何运动包括运动约束描述系统所受外力与力矩按牛顿第二定律组合上述信息整理方程使其可以被积分求解轨迹。SymPy 承担的是第 4、5 步用户需要提供完整的自由度描述即自由体图的符号化模块负责按牛顿第二定律完成繁琐的代数推导输出适合通用 ODE 求解器积分的方程组。模块内实现了 Kane 法与 Lagrange 法两种方程推导方法Kane 法的核心实现类是 KanesMethod。二、方程体系五类基本方程文档首先指出sympy.physics.mechanics假设描述一个系统需要5 类基本方程组完整约束holonomic constraints$\mathbf{f_h}(q, t) 0$非完整约束non-holonomic constraints$\mathbf{k_{nh}}(q, t), u \mathbf{f_{nh}}(q, t) 0$运动学微分方程kinematic differential equations$\mathbf{k_{k\dot{q}}}(q, t), \dot{q} \mathbf{k_{ku}}(q, t), u \mathbf{f_k}(q, t) 0$动力学方程dynamic equations$\mathbf{k_d}(q, t), \dot{u} \mathbf{f_d}(q, \dot{q}, u, t) 0$微分后的非完整方程differentiated non-holonomic equations$\mathbf{k_{dnh}}(q, t), \dot{u} \mathbf{f_{dnh}}(q, \dot{q}, u, t) 0$几点使用要点原文档的明确说明完整约束只用于线性化过程。模块假设完整约束过于复杂、无法解出依赖坐标若能轻松解出完整约束应当重新定义问题、改用更小的坐标集合或者直接提供对时间微分后的约束。Kane 法形成两个表达式 $F_r$ 与 $F_r^*$二者之和为零$F_r F_r^* 0$其中 $F_r$ 为广义主动力、$F_r^*$ 为广义惯性力参考 Kane Levinson,Dynamics Theory and Applications, 1985。模块将表达式整理为线性质量矩阵形式$$\mathbf{M}(q, t), \dot{u} \mathbf{f}(q, \dot{q}, u, t)$$对含有 $o$ 个总速度与 $m$ 个运动约束的非完整系统最终得到o − m 个方程。质量矩阵与力矢量按下式增广$$\mathbf{M}(q, t) \begin{bmatrix} \mathbf{k_d}(q, t) \ \mathbf{k_{dnh}}(q, t) \end{bmatrix}, \qquad \mathbf{f}{(forcing)}(q, \dot{q}, u, t) \begin{bmatrix} -\mathbf{f_d}(q, \dot{q}, u, t) \ -\mathbf{f{dnh}}(q, \dot{q}, u, t) \end{bmatrix}$$这一增广方式在源码中有直接对应mass_matrix 属性返回Matrix([self._k_d, self._k_dnh])forcing 返回-Matrix([self._f_d, self._f_dnh])与文档给出的矩阵定义完全一致。三、构造 KanesMethod 对象运动方程的推导从创建KanesMethod对象开始。初始化时必须提供惯性参考系以及坐标、速度等系统基本信息 from sympy.physics.mechanics import * N ReferenceFrame(N) q1, q2, u1, u2 dynamicsymbols(q1 q2 u1 u2) q1d, q2d, u1d, u2d dynamicsymbols(q1 q2 u1 u2, 1) KM KanesMethod(N, [q1, q2], [u1, u2])注意源码中的强制校验__init__会检查frame是否为ReferenceFrame实例不是则抛出TypeError(An inertial ReferenceFrame must be supplied)见 kane.py#L287-L289。dynamicsymbols(q1 q2 u1 u2, 1)中的第二个参数表示生成其一阶时间导数符号即q1d表示 $q_1$。3.1 依赖坐标与依赖速度的顺序要求当系统存在依赖坐标和依赖速度时顺序必须正确独立量在前依赖量在后或者显式使用关键字参数。文档示例中假设q2依赖、u2和u3依赖 q1, q2, q3, q4 dynamicsymbols(q1 q2 q3 q4) u1, u2, u3, u4 dynamicsymbols(u1 u2 u3 u4) # 假设 q2 依赖u2 和 u3 依赖 # 需要先提供约束方程 KM KanesMethod(N, [q1, q3, q4], [u1, u4])从源码结构看KanesMethod.q 与 KanesMethod.u 属性的 docstring 明确说明坐标矩阵是独立广义坐标在前、依赖坐标在后的列向量速度矩阵同理。后续的约束矩阵分割也依赖这个顺序——_initialize_constraint_matrices中按num_dof num_speeds - num_dep_speeds将速度约束的系数矩阵B切成独立列与依赖列kane.py#L413-L424。3.2 辅助速度auxiliary speeds若存在辅助速度需在构造时显式声明。示例中u4为辅助速度 KM KanesMethod(N, [q1, q3, q4], [u1, u2, u3], u_auxiliary[u4])辅助速度用于处理不贡献主动力的约束力问题。从源码看当存在辅助速度时kanes_equations 会内部构造一个新的KanesMethod实例以辅助速度作为 u 向量分别计算辅助的 $F_r$ 与 $F_r^*$并将二者与主方程拼接结果可通过auxiliary_eqs属性访问。若存在依赖速度新实例还会继承原实例的速度约束与加速度约束self._k_nh * self.u self._f_nh等。四、运动学微分方程运动学微分方程必须提供形式为一组等于零的表达式列表它把广义速度与广义坐标的时间导数线性关联起来。最简单的例子 kd [q1d - u2, q2d - u2]文档建议打开mechanics_printing()让表达式显著变短也可以改用mprint和mpprint命令。源码层面这些方程在 _initialize_kindiffeq_matrices 中被解析通过kdeqs.jacobian(u)与kdeqs.jacobian(qdot)分别提取系数矩阵 $\mathbf{k_{ku}}$ 和 $\mathbf{k_{k\dot{q}}}$将 $u0$、$\dot{q}0$ 代入得到 $\mathbf{f_k}$会校验方程关于 $u$ 和 $\dot{q}$ 的线性性若发现非线性项直接抛出ValueError随后用线性求解器解出 $\dot{q}$得到显式形式 $\mathbf{k_{k\dot{q}}^{-1}\mathbf{k_{ku}}, u I, \dot{q} \mathbf{k_{k\dot{q}}^{-1}}\mathbf{f_k} 0$并构建内部映射_qdot_u_map即 $\dot{q} \mapsto u$ 的字典供后续所有表达式替换使用。由此得到便捷方法 kindiffdict mechanics_printing(pretty_printFalse) KM.kindiffdict() {q1: u1, q2: u2, q3: u3, q4: u4}返回的字典把每个 $\dot{q}$ 映射到对应速度的表达式这正是 ODE 求解器状态向量中 $\dot{q}$ 分量的来源。关于显式/隐式两种运动学形式KanesMethod提供explicit_kinematics参数默认True显式形式下运动学质量矩阵 $\mathbf{k_{k\dot{q}}} \mathbf{I}$设为False则保留隐式形式$\mathbf{k_{k\dot{q}}}$ 不必是单位阵对非平凡运动学可得到更紧凑的方程。对应属性 mass_matrix_kin 与 forcing_kin 会依据该标志返回不同矩阵。五、非完整约束与依赖速度的完整示例存在非完整约束时需要指定依赖速度依赖坐标只在系统线性化时起作用。约束以等于零的表达式列表提供。文档给出了一个平凡运动约束 构型约束的完整示例 N ReferenceFrame(N) q1, q2, q3, q4 dynamicsymbols(q1 q2 q3 q4) q1d, q2d, q3d, q4d dynamicsymbols(q1 q2 q3 q4, 1) u1, u2, u3, u4 dynamicsymbols(u1 u2 u3 u4) # 假设 q2 依赖u2 和 u3 依赖 speed_cons [u2 - u1, u3 - u1 - u4] coord_cons [q2 - q1] q_ind [q1, q3, q4] q_dep [q2] u_ind [u1, u4] u_dep [u2, u3] kd [q1d - u1, q2d - u2, q3d - u3, q4d - u4] KM KanesMethod(N, q_ind, u_ind, kd, ... q_dependentq_dep, ... configuration_constraintscoord_cons, ... u_dependentu_dep, ... velocity_constraintsspeed_cons)这里涉及的构造参数均可在 KanesMethod 的 docstring 中查到完整说明参数含义校验规则q_ind/u_ind独立广义坐标 / 速度必须是可迭代的dynamicsymbol非可迭代对象抛TypeErrorkd_eqs运动学微分方程关于 $u$、$\dot{q}$ 线性数量须与坐标数相等非线性抛ValueErrorq_dependent依赖广义坐标数量须与configuration_constraints相等configuration_constraints完整构型约束仅用于线性化数量须与依赖坐标相等u_dependent依赖广义速度数量须与速度约束数相等velocity_constraints速度约束非完整约束 微分后完整约束的组合与nonholonomic_constraints二选一acceleration_constraints加速度约束默认取速度约束的时间导数数量须与依赖速度相等u_auxiliary辅助广义速度见 3.2 节kd_eqs_solver/constraint_solver线性方程组求解方法默认LU见 7.1 节从源码的_initialize_constraint_matriceskane.py#L336-L433可以看到几处值得注意的实现细节速度约束的自动组合如果同时提供configuration_constraints完整约束与nonholonomic_constraints严格非完整约束模块会自动将完整约束对时间微分通过_qdot_u_map替换 $\dot{q}$后与非完整约束垂直拼接col_join形成速度约束若只提供其中之一则直接使用。加速度约束的自动推导未显式提供时按_f_dnh ∂_t(∂k_nh/∂t)·u ∂_t(f_nh)从速度约束推导kane.py#L398-L404。依赖速度消去矩阵Ars将速度约束写成 $B u C 0$ 后按顺序把 $B$ 分成独立列 $B_{ind}$ 与依赖列 $B_{dep}$令 $A_{rs} -B_{dep}^{-1} B_{ind}$使依赖速度可由独立速度表示$u_{dep} A_{rs} u_{ind}$忽略 $C$ 项。若 $B_{dep}$ 奇异抛出带有明确提示的NonInvertibleMatrixErrorArs is singular, maybe due to a wrong set of dependent speeds. Try a different set of dependent speeds.——这正是要求依赖量必须排在独立量之后的底层原因。六、生成运动方程kane_equations 与质量矩阵最后一步是提供刚体/质点列表以及力/力矩列表。力列表由(Point, Vector)或(ReferenceFrame, Vector)二元组构成前者表示某点处的力后者表示作用在坐标系上的力矩 N ReferenceFrame(N) q, u dynamicsymbols(q u) qd, ud dynamicsymbols(q u, 1) P Point(P) P.set_vel(N, u * N.x) Pa Particle(Pa, P, 5) BL [Pa] FL [(P, 7 * N.x)] KM KanesMethod(N, [q], [u], [qd - u]) (fr, frstar) KM.kanes_equations(BL, FL) KM.mass_matrix Matrix([[5]]) KM.forcing Matrix([[7]])解读这个结果一个质量为 5 的质点在惯性系 $N$ 中沿 $x$ 方向以速度 $u$ 运动受常力 $7 N.x$。Kane 方程整理后为 $5,\dot{u} 7$因此mass_matrix为 $[5]$、forcing为 $[7]$。6.1 Fr 与 Fr* 的构造原理kanes_equations 内部依次调用_form_fr(loads)与_form_frstar(bodies)_form_fr广义主动力用 partial_velocity 对力作用点/系的速度列表计算各速度方向上的偏速度再与力做点积求和$F_r \sum_j \frac{\partial v_j}{\partial u_r} \cdot F_j$。若存在依赖速度则用 $A_{rs}^T$ 把依赖分量折叠进独立分量。_form_frstar广义惯性力先组装 $MM$即 $k_d$与nonMM不含 $\dot{u}$ 的项质点的惯性力为 $M, v M, a$刚体则额外包含质心平动项与欧拉力矩项 $\dot{I}\cdot\omega I\cdot\alpha \omega \times (I\omega)$见 kane.py#L603-L657。最终fr_star -(MM·udot nonMM)即文档中 $\mathbf{k_d}\dot{u} \mathbf{f_d} 0$ 的具体体现其中self._k_d MM、self._f_d -(self._fr - nonMM)。存在运动约束时如文档所述质量矩阵被 $k_{dnh}(q, t)$ 矩阵增广力矢量被 $f_{dnh}(q, \dot{q}, u, t)$ 增广——即第二节给出的分块结构。6.2 完整的full质量矩阵文档还介绍了 full 质量矩阵与 full 力矢量它们把运动学微分方程并入规模是(n o) × (n o) 的方阵覆盖全部坐标与速度 KM.mass_matrix_full Matrix([ [1, 0], [0, 5]]) KM.forcing_full Matrix([ [u], [7]])对应实现见 mass_matrix_full按[k_kqdot, 0; 0, M]的分块结构拼接本例中 $n o 1$$\dot{q} - u 0$ 给出第一行 $[1, 0]$。forcing_full则是运动学力 $[u]$来自 $\dot{q} u$与动力学力 $[7]$ 的堆叠。七、从方程到可积分形式线性化与求解器7.1 线性化接口文档最后鼓励读者探索KanesMethod的示例以加深理解结合源码线性化是这套流程的自然终点。KanesMethod提供两个入口linearize直接返回线性化结果。若A_and_BFalse默认返回 $M, A, B, r$对应 $M [q, u]^T A [q_{ind}, u_{ind}]^T B r$若A_and_BTrue则返回状态空间形式 $dx A x B r$$x [q_{ind}, u_{ind}]^T$。to_linearizer返回 Linearizer 实例便于围绕不同工作点op_point更高效地重复线性化。完整约束正是这里发挥作用的地方——这与文档第二节完整约束只用于线性化的说明呼应。文档中 docstring 自带的弹簧-质量-阻尼器示例演示了完整链路 KM KanesMethod(N, q_ind[q], u_ind[u], kd_eqskd) (fr, frstar) KM.kanes_equations(BL, FL) MM KM.mass_matrix forcing KM.forcing rhs MM.inv() * forcing rhs Matrix([[(-c*u(t) - k*q(t))/m]]) KM.linearize(A_and_BTrue)[0] Matrix([ [ 0, 1], [-k/m, -c/m]])对 $x [q, u]^T$、$r$ 为空的情形得到标准状态空间矩阵 $\begin{bmatrix} 0 1 \ -k/m -c/m \end{bmatrix}$。另外 rhs 方法直接输出一阶形式 $x f(q, u, r, p, t)$ 的右端注释明确说明这是大多数数值 ODE 积分器所需的形态其中 $\dot{u}$ 分量默认用mass_matrix.LUsolve(forcing)求解。7.2 线性求解器的选择构造KanesMethod时可为两组内部线性方程组分别指定求解器kd_eqs_solver解运动学微分方程与constraint_solver解速度约束两者默认都是LU。根据 docstring 的 Notes 部分默认LU运算量相对较少但可能出现除零错误导致nanCRAMER当 LU 出现除零时可改用 Cramer 法则求解运算更多更慢但默认每个解分量只产生一次除法自定义 callable可传入形如f(A, rhs)的函数例如lambda A, b: simplify(Matrix.LUsolve(A, b))用于化简结果或lambda A, b: tuple(linsolve((A, b)))[0]利用linsolve求唯一解。字符串形式必须是 MatrixBase.solve 支持的有效方法名。八、验证与延伸阅读Kane 法的行为由一组分层测试文件覆盖test_kane.py、test_kane2.py、test_kane3.py、test_kane4.py、test_kane5.py分别对应不同规模与约束类型的系统速度约束线性性的专门测试见 test_linearity_of_velocity_constraints.py。这些测试文件是验证构造顺序、约束数量匹配、矩阵增广等行为的最佳参照。配套文档均在仓库中从 经典力学文档目录 出发Kane 法说明本文所依据的原文档Lagrange 法说明线性化说明质量体质点/刚体说明高级话题含辅助速度、不贡献力等参考页模块中的其他实现文件也值得按需查阅particle.pyParticle含linear_momentum等动量接口、rigidbody.py刚体与中心惯性、linearize.pyLinearizer类、functions.pymsubs、find_dynamicsymbols、_parse_linear_solver等内部工具。九、小结标准工作流把文档与源码串起来用KanesMethod推导多体系统运动方程的标准工作流是建系与建模创建惯性ReferenceFrame用dynamicsymbols生成坐标/速度及其导数用Point.set_vel完成运动学描述构造Particle/RigidBody列表组装 KanesMethod按独立在前、依赖在后的顺序给出q_ind、u_ind提供kd_eqs必要时补齐q_dependent、configuration_constraints、u_dependent、velocity_constraints与u_auxiliary生成方程调用kane_equations(BL, FL)得到 $F_r$、$F_r^*$通过mass_matrix/forcing或含运动学的mass_matrix_full/forcing_full读取增广后的 $\mathbf{M}\dot{u} \mathbf{f}$转为可积分形式用kindiffdict()得到 $\dot{q}$ 分量、rhs()得到一阶右端或linearize(A_and_BTrue)/to_linearizer()得到线性化状态空间形式交给数值 ODE 求解器积分。整个流程中所有矩阵提取、约束消元与线性求解都在符号层面完成这正是该模块被设计为前向动力学第 4、5 步的实现方式。【免费下载链接】sympyA computer algebra system written in pure Python项目地址: https://gitcode.com/GitHub_Trending/sy/sympy创作声明：本文部分内容由AI辅助生成（AIGC），仅供参考