[已收纳]python之科学计算笔记

15195775117

#变量引用易错：
a=[1,2,3]
b=a#这里是引用，不是复制
a.append(4)
print(b)#输出：[1, 2, 3, 4]
#用array比较符合习惯：
a=np.array([1,2,3])
b=a
a=np.append(a,4)
print(b)#输出：[1 2 3]


切片：x=range(10);y=np.arange(1,20)；
print(x[1:-1:3],y[0::3])
结果range(1, 9, 3) [ 1  4  7 10 13 16 19]
逆序：x=np.arange(5)
y=x[-1:0:-1]或x=list(reversed(x))结果[4 3 2 1]
切片与原数组共享内存，改变切片也就改变了原数组：
x=np.arange(8);y=x[0::2];y[0]=999;print(y,x)
结果[999   2   4   6] [999   1   2   3   4   5   6   7]
共享内存与不同享内存的区别：
x=np.arange(5)
y=x[range(3)];y[0]=999;print(x)
y=x[0:3];y[0]=999;print(x)
输出：[0 1 2 3 4]，[999   1   2   3   4]
可见，y=x[0:3]会共享，x[[0,1,2]]不会


#画心形
import matplotlib.pyplot as plt
x,y=np.mgrid[-2:2:500j,-2:2:500j]
z=(x**2+y**2-1)**3-x**2*y**3
plt.contourf(x,y,z,levels=[-1,0],colors=['red'])
plt.show()



像plt可以设置字体一样，数字显示也可以提前预设：
np.set_printoptions(precision=2)#小数只显示2位

python的int类型任意大，float类型是双精度64位。
查数据类型：a=np.arange(12);print(a.dtype)
numpy的ndarray只支持单一数据类型，这种“缺点”类似IDL：
np.array([1,2,'a'])，结果['1' '2' 'a']
获取numpy中所有数据类型：set(np.typeDict.values())
类型转换：x=np.arange(3);y=x.astype(np.str);print(y)
结果：['0' '1' '2']


修改数组大小：
a=np.arange(12)
a.shape=3,4
print(a)
[[ 0  1  2  3]
[ 4  5  6  7]
[ 8  9 10 11]]
也可以用a.shape=3,-1，列数会自动算出。
还可以用reshape函数b=a.reshape(3,4)


#创建时指定数据类型,可以指定为内置类型，也可以是np类型:
#float16,float32,float64,float128...
x=np.array([1,2,3],dtype=np.int32)
x=np.array([1,2,3],dtype=np.float)
x=np.array([1,2,3],dtype=complex)
例子:x=np.array(32767+1,dtype=np.int16)
print(x) #-32768


等比数列，类似linspace,从10^0到10^2，3个数:
x=np.logspace(0,2,3,base=10);print(x)结果[1. 10. 100.]
base是底，默认10.
全1阵：x=np.ones((2,3))


x=10-np.arange(10);p=[1,3,5,7]
print(x[p]);结果[9 7 5 3]
numpy有着类似IDL的数组处理方式，python内置数学库不行
注意：不要只把np.arange当索引序列使用，可以更灵活些：
np.arange(10,5,-1)


切片的同时变维：
x=np.arange(21)*10;a=np.arange(1,21,2)
b=np.arange(2,21,2);p=np.vstack((a,b))
y=x[p];print(y)
# [[ 10 30 50 70 90 110 130 150 170 190]
# [ 20 40 60 80 100 120 140 160 180 200]]

切片变量：x=np.arange(10);p=np.s_[::2];print(x[p])
#结果[0 2 4 6 8]，与x[::2]相同

对值的大小进行过滤：
x=np.random.randint(0,11,10)
p=[x>5];print(p);print(x[p])
# [array([True,False,False,True,True,
# False,True,False,False,False])]
# [ 9  7 10  6]
这里是支持多个检索条件的：
x=np.random.randint(0,11,10);print(x)
p=[(x>5)&(x<8)];print(x[p])
#[ 2  3  6  2 10  0  4  3  2  7]   [6 7]
原来之前错怪np.where了：
x=np.random.randint(0,11,10);print(x)
p=[(x>5)&(x<8)];print(x[p])
p=np.where((x>5)&(x<8));print(x[p])
#[0 0 6 5 0 0 2 6 8 6]
#[6 6 6]
#[6 6 6]

# 对角线：
x=np.arange(100).reshape(10,10)
print(x[np.arange(10),np.arange(10)])

矩阵生成：
x=10*np.arange(3).reshape(-1,1)+np.arange(4)
[[ 0  1  2  3]
[10 11 12 13]
[20 21 22 23]]

numpy比math库快10倍！

#以下每个函数都可用axis分行列计算:
a=np.arange(12).reshape(3,4)
x=np.sum(a)
x=np.mean(a)
x=np.mean(a,axis=0)#每列求均值
x=np.mean(a,axis=1)#每行求均值
x=np.median(a)
x=np.var(a)#方差
x=np.std(a)#标准差
x=np.product(a)#元素乘积，可用于阶乘
x=np.min(a)
x=np.max(a)
p=np.argmin(a)#返回的索引是一维化的索引
p=np.argmax(a)#返回的索引是一维化的索引
p=np.argsort(a,axis=0)#排序

#求集合的另一种方法:
a=np.random.randint(1,11,(10,))
x,p=np.unique(a,return_index=True)
print('排序索引',p);print('排序值',a[p])
# 排序索引 [2 8 5 1 0 7]
# 排序值 [1 3 4 5 8 9]


#直方图:
a=np.random.randn(20)
#density=False表示返回值为各区间频数
x=np.histogram(a,bins=5,density=False)
print('区间元素个数=',x[0])
print('区间界=',x[1])
区间元素个数= [2 6 8 3 1]
区间界= [-1.83 -1.02 -0.21  0.6   1.41  2.22]
区间可以自定义、非均匀:
a=np.random.randn(20)
x=np.histogram(a,bins=[-1,-0.9,0.9,1],density=False)
print('区间元素个数=',x[0])
print('区间界=',x[1])
# 区间元素个数= [ 0 13 0]
# 区间界= [-1. -0.9 0.9 1. ]


数值积分：
from scipy.integrate import quad
def f(x):
    return np.sqrt(1-x**2)
y=quad(f,-1,1)[0]*2 #第二个返回值是啥？
print(y) #3.1415926535897967

#符号计算库SymPy:
from sympy import symbols,integrate,sqrt
x=symbols('x')
y=integrate(sqrt(1-x**2),(x,-1,1))*2
print(y) #输出结果:pi


def g(x):
    return x**2-3*x+2
a=np.array([1,-3,2])#系数降幂
f=np.poly1d(a)#一元多项式函数，此处返回值同g函数
print('f和g一样:',f(3),g(3))#6 6
#符号计算:
df=f.deriv()#微分，返回值也是poly1d型
print('微分函数:',df)#微分，输出的写法有上标
print('微分函数带入值计算:',df(3))
int_f=f.integ()#积分，返回值也是poly1d型
print('积分函数带入值计算:',int_f(3))
root=np.roots(f);print('根:',root)#求根
print('由根反求多项式系数:',np.poly(root))
# f和g一样: 2 2
# 微分函数:  
# 2 x - 3
# 微分函数带入值计算: 3
# 积分函数带入值计算: 1.5
# 根: [2. 1.]
# 由根反求多项式系数: [ 1. -3.  2.]
#matlab中求根的也是roots,反求系数的也是poly


#多项式拟合
x=np.linspace(-np.pi,np.pi,100)
y=np.sin(x)
plt.plot(x,y)
deg=np.arange(3,8)#最高次从3到7
for i in deg:
    coef=np.polyfit(x,y,i)#拟合
    y2=np.polyval(coef,x)#多项式(依据系数)代入自变量序列
    plt.plot(x,y2,label='多项式次数'+str(i))
    plt.legend()
plt.show()

之前还单独发过一个帖子：python之多项式拟合


#poly1d是降序,np.polynomial模块的多项式系数是升序的:
#Polynomial,Chebyshev,Legendre
from numpy.polynomial import Polynomial,Chebyshev
def g(x):
    return 1.0-2*x+x**3
f=Polynomial([1,-2,0,1])
print(f(2),g(2))
df=f.deriv();print('导函数系数',df)


# #三次样条二维插值
# func = interpolate.interp2d(x[px],y[py],zone, kind='cubic')
# x2,y2= np.meshgrid(x[px],y[py])
# # plt.scatter(x2,y2,c=y2,cmap='rainbow')
# rows=x2.shape[0]
# cols=x2.shape[1]
# x3=x2.reshape(rows*cols)
# y3=y2.reshape(rows*cols)
# z3= func(x3,y3)
# z2=z3.reshape(rows,cols)

ゼヅ欢チニ

好东西，谢谢分享