Sympy.solve求解三重摆拉格朗日方程时持续运行不停止

核心原因 #

你的拉格朗日方程推导逻辑没有问题，双摆场景能正常运行已经验证了这一点。长时间跑不出结果的核心原因是**sympy.solve是通用符号求解器，不适合直接求解这类结构明确的线性符号方程组**：

• 你构造的LE1/LE2/LE3确实是关于the1_dd、the2_dd、the3_dd的线性方程组，但solve不会自动识别这个线性结构，默认会按通用非线性方程的逻辑遍历所有消元路径，还会执行大量你无法通过参数完全关闭的隐式化简
• 双摆场景下方程组是2阶的，展开后符号项少，哪怕求解器走冗余路径也能很快算出结果；三重摆的3阶方程组展开后有上百个三角函数交叉项，通用消元的计算量随阶数指数级上涨，很容易出现几小时跑不出结果的情况。

修复方案 #

针对线性方程组的结构做定向求解，速度可以提升几个数量级，通常几秒就能得到结果：

1. 用SymPy内置的linear_eq_to_matrix方法，直接从三个拉格朗日方程中提取未知数的系数矩阵M和常数项向量b，这个方法专门针对线性方程组做了优化，不会做无意义的非线性求解尝试
2. 通过线性代数求解M * acc = b得到三个角加速度的解析解，替代原来的通用solve调用

替换原有求解部分的代码即可，后续接lambdify转数值函数做仿真的逻辑和双摆场景完全兼容：

# 替换原代码中的 sols = smp.solve(...) 行
from sympy import linear_eq_to_matrix
# 待求解的未知数：三个角加速度
acc_vars = [the1_dd, the2_dd, the3_dd]
# 提取线性方程组的系数矩阵和常数项
M, b = linear_eq_to_matrix([LE1, LE2, LE3], acc_vars)
# 求解线性方程组得到角加速度解析解
sols = M.solve(b)

额外优化提示 #

如果后续需要做数值仿真，不需要对解析解做全量simplify化简，直接用lambdify把sols转成numpy兼容的数值函数即可，能省掉大量符号化简的时间。

内容的提问来源于stack exchange，提问作者Oriol