Scipy是一个用于数学、科学、工程领域的常用软件包,可以处理插值、积分、优化、图像处理、常微分方程数值解的求解、信号处理等问题。它用于有效计算Numpy矩阵,使Numpy和Scipy协同工作,高效解决问题。
一、文件输入/输出
1、以matlab文件形式保存和加载数组
import numpy as np
from scipy import io as spio
data = np.ones((3, 3))
spio.savemat('file.mat', {'a':data})
data2 = spio.loadmat('file.mat', struct_as_record=True)
print(data2['a'])
2、以图像形式保存和加载数组
import numpy as np
from scipy import misc
data = np.ones((3, 3))
misc.imsave('file.jpg', data)
data2 = misc.imread('file.jpg')
print(data2)
二、物理和数学常数
scipy的constants模块包含了大量用于科学计算的常数
from scipy import constants as C
#圆周率
print(C.pi)
#黄金比例
print(C.golden)
#真空中的光速
print(C.c)
#普朗克常数
print(C.h)
#一英里等于多少米
print(C.mile)
#一英寸等于多少米
print(C.inch)
#一度等于多少弧度
print(C.degree)
#一分钟等于多少秒
print(C.minute)
#标准重力加速度
print(C.g)
三、线性代数
1、求解线性方程组
import numpy as np
from scipy import linalg
# 2x - y = 0
# x + 3y = 7
a = np.array([[2, -1],
[1, 3]])
b = np.array([0, 7])
# 求解线性方程组
x = linalg.solve(a, b)
print(x)
2、计算方阵A的行列式 |A|
import numpy as np
from scipy import linalg
# |A| = det(A) = 1*4-2*3
A = np.array([[1, 2],
[3, 4]])
# 计算方阵A的行列式
x = linalg.det(A)
print(x)
3、计算方阵A的逆
import numpy as np
from scipy import linalg
A = np.array([[1, 2],
[3, 4]])
# 计算方阵A的逆
iA = linalg.inv(A)
print(iA)
4、计算方阵A的特征值和特征向量
A(方阵A)* v(特征向量) = λ(特征值)* v(特征向量)
import numpy as np
from scipy import linalg
A = np.array([[1, 2],
[3, 4]])
# 计算方阵A的特征值和特征向量
λ, V = linalg.eig(A)
print(λ) # 特征值
print(v) # 特征向量
5、奇异值分解
奇异值分解(SVD)可以被认为是特征值问题扩展到非方阵的矩阵。
A = U(左奇异矩阵) * sigma(奇异值)* VT(右奇异矩阵的转置)
import numpy as np
from scipy import linalg
A = np.array([[1, 2],
[3, 4],
[5, 6]])
# 计算矩阵A的左奇异矩阵、奇异值、右奇异矩阵转置
U, sigma, VT = linalg.svd(A)
print(sigma) # 奇异值
print(U) # 左奇异矩阵
print(VT.T) # 右奇异矩阵
四、数值积分
1、单积分
# 园的面积 x^2 + y^2 = 1
import numpy as np
from scipy import integrate
def half_circle(x):
return (1-x**2)**0.5
pi_half, err = integrate.quad(half_circle, -1, 1)
print(2*pi_half)
# 除了自定义的函数外,quad函数还可以积分单个变量的标准预定义NumPy函数,如exp,sin和cos
i, err = integrate.quad(np.sin, 0, 1)
print(i)
2、双重积分
# 球的体积 x^2 + y^2 + z^2 = 1
import numpy as np
from scipy import integrate
def half_circle(x):
return (1-x**2)**0.5
def half_sphere(x, y):
return (1-x**2-y**2)**0.5
half_sphere, err = integrate.dblquad(half_sphere, -1, 1,
lambda x: -half_circle(x),
lambda x: half_circle(x)
)
print(2*half_sphere)
五、数据插值
1、一维插值
import numpy as np
from scipy import interpolate
import matplotlib.pyplot as plt
x = np.linspace(0, 4, 12)
y = np.cos(x**2/3+4)
f1 = interpolate.interp1d(x, y,kind = 'linear')
f2 = interpolate.interp1d(x, y, kind = 'cubic')
xnew = np.linspace(0, 4, 30)
plt.plot(x, y, 'o' , xnew, f1(xnew), '-', xnew, f2(xnew), '--')
plt.legend(['data', 'linear', 'cubic','nearest'], loc = 'best')
plt.show()
2、二维插值
import numpy as np
from scipy import interpolate
import matplotlib.pyplot as plt
import matplotlib.cm as cm
def func(x, y):
return (x + y) * np.exp(-5.0 * (x ** 2 + y ** 2))
x = np.linspace(-1,1,15)
y = np.linspace(-1,1,15)
x, y= np.meshgrid(x, y) # 15*15的网格
fvals = func(x, y) # 计算每个网格点上的函数值
plt.subplot(121)
im1=plt.imshow(fvals, extent=[-1,1,-1,1], cmap=cm.hot, interpolation='nearest', origin="lower")#pl.cm.jet
plt.colorbar(im1)
########################################################
# 二维插值
newfunc = interpolate.interp2d(x, y, fvals, kind='cubic')
xnew = np.linspace(-1,1,100)
ynew = np.linspace(-1,1,100)
fnew = newfunc(xnew, ynew) # 计算100*100的网格上的插值
plt.subplot(122)
im2=plt.imshow(fnew, extent=[-1,1,-1,1], cmap=cm.hot, interpolation='nearest', origin="lower")
plt.colorbar(im2)
plt.show()