干货警告:用python解方程的那些事

图:作者上台笔记本电脑的主板及本文涉及的部分代码

时光飞逝,转眼间已是2018年的初秋。304TechStudio公众号创办以来,一直稀稀拉拉没有发太多的文章,也是略显尴尬。所以从这个夏天开始,我们打算恢复在这里的更新,并推出周期性更新的栏目:周闻。

这个栏目主要有两大块含义:一是技术性内容,包括几位同学学习了解计算机相关技术过程中的收获和总结,以及几位同学所学的化学及相关专业的内容;二是随感,内容主要是日常生活中、读书看帖时的所见所感。这个栏目的名称有两重含义:一则取“周详”、“周全”的含义,以示栏目内容涉猎广泛;二则取“每周”的含义,表明此栏目以周为单位更新。

周闻 技术 | 干货警告:用python解方程的那些事

最近在和同学合作用Python开发过程中,有个函数需要一个解五元一次方程组的过程,这里也成为了目前整个算法的决速步骤,因此在最近集中解决了一下用Python已有的库解方程组的问题。

虽然最后还是把方程直接解出来了,不过整个调试的过程还是有不少所得是值得记录的。

关于Numpy、Scipy和Sympy库中解方程或方程组的方法

(1) 如何调用三个库中解方程的方法

(2) 官方文档中对scipy.linalg和numpy.linalg的比较

补充阅读

(1) 关于scipy、numpy和sympy中解方程或方程组相关方法的文章和教程

(2) 利用Python Constraint 模块解决简单的约束性编程问题

(3) Python环境下的8种简单线性回归算法

(4) 线性方程组解法概观

参考链接

关于Numpy、Scipy和Sympy库中解方程或方程组的方法

(1)如何调用三个库中解方程的方法

Numpy、Scipy和Sympy都是Python环境中常用的科学计算库,它们都可以通过安装Anaconda来安装。

首先,Numpy和Scipy的线性代数部分scipy.linalgnumpy.linalg提供了相同的方法solve()解线性方程或方程组。

参数中第一个矩阵是方程的系数或方程组的系数矩阵,第二个矩阵是方程或方程组中每个方程等号右侧的参数项,函数返回的解则是一个numpy数组,即ndarray。

由此可知,如果系数矩阵不是方阵、系数矩阵与常数项矩阵不匹配,以及方程组无解或存在任意解的情况都会引发报错,因此在调用此方法前应该对输入的矩阵进行检查。

使用举例(注意scipy的调用方法,以及结果的精确程度):

importnumpy

importscipy

fromscipyimportlinalg

print("a : ",a)

print("a[0] :",a[])

print("b : ",b)

print("c : ",c)

a : [2. 1.]

另外scipy.optimize部分提供了fsolve()root()方法用于解非线性方程或方程组。输入参数既可以直接输入待解的函数和初值,也可以输入待解函数、初值和可变参数,返回的结果也是numpy数组ndarray。举例如下:

fromscipy.optimizeimportfsolve

importmath

importnumpyasnp

deff0(x):

returnmath.cos(x) - x**2

result_0 = fsolve(f0,1)

print("result_0 :", result_0)

deff1(x,*arg):

f1_1 = x**2*arg[] + x*arg[1] + arg[2]

f1_2 =2*x +9

return(f1_1-f1_2)

arg_f1 = (1,1,2)

result_1 = fsolve(f1,-0.5,arg_f1)

print("result_1 :", result_1)

deff2(x):

x0 =float(x[])

x1 =float(x[1])

x2 =float(x[2])

return[

5*x1+3,

4*x0*x0 -2*math.sin(x1*x2),

x1*x2-1.5

]

x0 = [1,1,1]

result_2 = fsolve(f2,x0)

print("result_2 :",result_2)

result_1: [-2.1925824]

Python代码缩进应为四个空格,但在不同设备上浏览本文中代码文字部分时,实际显示的空格数量可能会有变化,可见在微信推送中用空格维持格式并不是个明智的选择。

需要注意的是,由于该方法使用迭代法求解,初值的设置对于能否顺利求解还是比较重要的,同时还要在调用前确定方程有解。

scipy.optimize. root()同scipy.optimize. fsolve()不同的是,这个方法支持调用不同的解方程方法以进行比较。这在官方文档中有详细的说明:

https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.root.html

fromscipy.optimizeimportleastsq

importmath

importnumpyasnp

deff3(x,*arg):

f3_1 = x**2+x+2

f3_2 =2*x +9

return(f3_1-f3_2)

result_3 = leastsq(f3,)

print("result_3 :", result_3)

result_3 :(array([-2.1925824]), 2)

Sympy是Python另一个常用的数值计算库,由其文档可知其中有solveset()、linsolve()、nonlineslove()、roots()、dsolve()和solve()等函数,可满足解不同方程(包括微分方程)或方程组的需求。其中sympy.solve()函数最为常用。

文档页面:

https://docs.sympy.org/latest/tutorial/solvers.html#solving-equations-algebraically

sympy.solve()在调用时要先输入的数组内容是齐次方程或方程组,即方程等号右侧为0,第二个数组是前面方程或方程组的未知数,用Symbol()方法初始化。

fromsympyimport*

x = Symbol('x')

y = Symbol('y')

z = Symbol('z')

a = solve([x**2- y**2/exp(x)], [x,y])

print(a)

xx = Symbol('xx')

yy = Symbol('yy')

b =solve([xx+yy-3,2*xx+5*yy-12],[xx,yy])

print(b)

xxx = Symbol('xxx')

yyy = Symbol('yyy')

c =solve([xxx+yyy-3,2*xxx+2*yyy-12],[xxx,yyy])

print(c)

xxxx = Symbol('xxxx')

yyyy = Symbol('yyyy')

d = solve([xxxx-3,*xxxx+*yyyy-],[xxxx,yyyy])

print(d)

[]

[]

注意方程解的形式,注意当方程组无和存在在任意解时解的形式。

使用Sympy的其他方法求解不同函数的情况作者精力所限,在此不再赘述,只是通过这个例子见微知著,提醒读者调用时注意。

importscipy

fromscipyimportlinalg

print(a)

print(b)

[-4.4408921e-16 1.0000000e-01]

[-0.17 0.22]

(2)官方文档中对scipy.linalg和numpy.linalg的比较

在Scipy v1.1.0的官方文档中,Linear Algebra (scipy.linalg)页面认为scipy.linalg比numpy.linalg模块有更多的方法,并更多地使用了BLAS/LAPACK运算因此速度更快;而在同一版本文档下另一页面中, numpy.linalg比scipy.linalg有更多的线性代数函数。所以二者真正孰优孰劣,还有待于进一步关注。

Linear Algebra (scipy.linalg)

https://docs.scipy.org/doc/scipy/reference/tutorial/linalg.html

Linear algebra (scipy.linalg)

https://docs.scipy.org/doc/scipy/reference/linalg.html

Linear algebra (numpy.linalg)

https://docs.scipy.org/doc/numpy/reference/routines.linalg.html

补充阅读

(1)关于scipy、numpy和sympy中解方程或方程组相关方法的文章和教程

用Python做科学计算

http://old.sebug.net/paper/books/scipydoc/index.html

这是一个简单的numpy、scipy、sympy等python科学计算库的教程,可以作为入门参考。

python/scipy求解非线性方程(fsolve/root)

使用 Python 解数学方程

https://www.cnblogs.com/baby123/p/6296629.html

这是一篇关于sympy的简单介绍

Welcome to SymPy’s documentation!

https://docs.sympy.org/latest/index.html

使用Python+SymPy 求解线性方程组

https://www.jianshu.com/p/f7bec09515fd

这篇文章提到了用Sympy求解存在自由变量的方程组,但作者尝试似乎并不能运行

什么是Condition Number(条件数)?

https://www.cnblogs.com/2008nmj/p/8820401.html

(2)利用PythonConstraint 模块解决简单的约束性编程问题

这是在准备资料过程中搜索到的一个有趣的问题,作者精力和水平所限,在此不做过多探讨。

用Python怎么解这种类型的问题?

https://www.zhihu.com/question/35520855/answer/63116891

python constraint 模块

python-constraint Constraint SolvingProblem resolver for Python

https://labix.org/python-constraint

(3) Python环境下的8种简单线性回归算法

(4)线性方程组解法概观

《用Python学习数值分析--解方程》

《用Python学习数值分析--解方程组》

线性方程组的直接解法(python)

参考链接

python用fsolve、leastsq对非线性方程组进行求解

使用scipy来解非线性方程

python求非线性方程的解/非线性方程组的解

https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.leastsq.html

https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.fsolve.html

Quick examples

https://github.com/sympy/sympy/wiki/Quick-examples

Optimization and root finding(scipy.optimize)

https://docs.scipy.org/doc/scipy/reference/optimize.html

Optimization (scipy.optimize)

https://docs.scipy.org/doc/scipy/reference/tutorial/optimize.html

https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.fsolve.html

Nonlinear solvers

https://docs.scipy.org/doc/scipy/reference/optimize.nonlin.html

https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.solve.html

https://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.solve.html

  • 发表于:
  • 原文链接https://kuaibao.qq.com/s/20180831G1N6X800?refer=cp_1026
  • 腾讯「云+社区」是腾讯内容开放平台帐号(企鹅号)传播渠道之一,根据《腾讯内容开放平台服务协议》转载发布内容。
  • 如有侵权,请联系 yunjia_community@tencent.com 删除。

扫码关注云+社区

领取腾讯云代金券