图:作者上台笔记本电脑的主板及本文涉及的部分代码
时光飞逝,转眼间已是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.linalg和numpy.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
领取专属 10元无门槛券
私享最新 技术干货