• 24小时客服在线
  • info@qimumu.top

标签:python

用Python来计算

位置导航:首页/ 技术杂谈(投稿)


最近做毕业设计的时候需要计算机器人运动学和动力学,但是这玩意在计算过程中需要计算很多矩阵,这个就很头疼。一般这种计算直接给MATLAB,但是它太太太大了,界面也太太太丑了,而且一运行就会很卡,再加上MATLAB也很贵,因此用了python的sympy库。用了之后发现这个是一个宝库,遇到问题,别问,问就python因此写一下。

Sympy是啥

这是一个 可以免费使用 的用python写的库(相当于打包好的一堆函数),然后你可以调用这些函数来完成诸如多项式求值、求极限、解方程、求积分、微分方程、级数展开、矩阵运算等等计算问题。sympy(symbol-python)主要进行的是符号计算,就是化简表达式之类的。

sympy有以下功能:

  • 表达式化简,求值,变形
  • 微积分
  • 求解线性或非线性方程
  • 求解微分方程或差分方程
  • 矩阵运算
  • 数学公式的 TeX 或 LaTeX 显示

矩阵运算、微积分和表达式化简一般来说比较常用,LaTeX输出和AxMath结合可以直接把结果输出到Word里面,或者Markdown。

符号计算应用

矩阵应用

这种用例子看最快了。首先看个矩阵运算

假设一个点p[px py pz]^T,它先绕X轴旋转α度,再绕Z轴旋转β度,求坐标,很容易得到新点的表达式为

其中Rot(axis, angle)表示旋转矩阵。

要求解这个首先需要导入sympy这个库,为了它的输出更加美观使用init_printing函数

import sympy as sym
sym.init_printing(use_latex="mathjax")

然后表示出你设置的符号,像上面的px, py, pz, α , β 就是的

px, py, pz, alpha, beta = sym.symbols(r"p_x, p_y, p_z, \alpha, \beta")

其中sympy用symbols函数申明符号,用逗号隔开,支持LaTeX语法,这里在字符串前面加r是为了防止一些转义字符引起bug

符号定义好了就可以定义矩阵的,使用sympy的Matrix类创建矩阵,为了方便使用,这里用定义三个函数给出(其实只用写一个然后坐标轮换就可以了)

def RX(theta):
    return sym.Matrix([
        [1, 0, 0],
        [0, sym.cos(theta), -sym.sin(theta)],
        [0, sym.sin(theta), sym.cos(theta)],
    ])
def RY(theta):
    return sym.Matrix([
        [sym.cos(theta), 0,sym.sin(theta)],
        [0, 1, 0],
        [-sym.sin(theta), 0,sym.cos(theta)],
    ])
def RZ(theta):
    return sym.Matrix([
        [sym.cos(theta), -sym.sin(theta), 0],
        [sym.sin(theta), sym.cos(theta), 0],
        [0, 0, 1],
    ])

符号和矩阵定义好了,下面就是进行计算

p_after = RZ(beta)*RX(alpha)*sym.Matrix([px, py, pz])
p_after #这一行用于输出

这样运行一下就可以了,输出结果如下

还可以通过latex函数输出latex代码如下,用【\】替换【\\】就可以粘贴到AxMath里了。

矩阵+积分

这里用一个惯性矩矩阵的例子来说明,比如说要求下面这个方块的惯性矩矩阵

惯性矩矩阵的定义如下

这时候就要用到积分了,sympy中计算积分的函数是 integrate (function)可以计算不定积分,加入定义域后可以计算定积分 integrate (function, (x, x_low, x_up))

由于这里是三次积分,就需要进行三次套娃,例如Ixx可以这么写

x, y, z, rho = sym.symbols(r"x, y, z, \rho")
w, l, h = sym.symbols("w, l, h") 
Ixx = sym.integrate(sym.integrate(sym.integrate((y**2 + z**2)*rho, (x, 0, w)), (y, 0, l)), (z, 0, h))

最终写成代码并输出是这样的

和课本上是一样的(注意m=ρwlh)

Fu_Qingchen, LearningBySharing2020