使用SymPy Mechanics推导3D摆运动方程遇NaN问题求助

问题出在SphericalJoint的child_point参数定义上：你直接用了bob.x，但此时摆锤（bob）的本体坐标系还没通过关节与天花板（ceiling）的参考系建立关联，这种未定义的坐标引用会导致符号计算冲突，最终生成含NaN的质量矩阵。

修正方案如下：

• 改用天花板的参考系来定义摆锤的悬挂点位置（因为关节是连接天花板和摆锤，悬挂点相对于天花板的位置是固定的）
• 不需要手动定义q和u的导数，JointsMethod会自动处理这些关系

修正后的完整代码：

from sympy import symbols, simplify
from sympy.physics.mechanics import Body, SphericalJoint, JointsMethod, dynamicsymbols, mechanics_printing, mprint

mechanics_printing()

# 定义常量和动力学符号
l, m, g = symbols('l m g')
q1, q2, q3 = dynamicsymbols('q1 q2 q3')
u1, u2, u3 = dynamicsymbols('u1 u2 u3')

# 创建参考体和摆锤体
ceiling = Body('C')
bob = Body('U', mass=m)

# 关键修正：用天花板的参考系定义child_point，而非摆锤自身的未关联坐标系
joint = SphericalJoint(
    'P1', 
    ceiling, 
    bob, 
    child_point=-l * ceiling.x,  # 悬挂点相对于天花板沿x轴负方向l处
    coordinates=[q1, q2, q3], 
    speeds=[u1, u2, u3]
)

# 施加重力（用天花板的y轴方向）
bob.apply_force(-m * g * ceiling.y)

# 生成运动方程
method = JointsMethod(ceiling, joint)
method.form_eoms()

# 计算并输出结果
mm = method.mass_matrix_full
fo = method.forcing_full
mmi = simplify(mm).inv()
qudots = mmi * fo
mprint(simplify(qudots))

运行这段代码后，质量矩阵不再包含NaN，能正常推导得到3D摆的运动微分方程。

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