(十五）Numpy

本次作业的主角毫无疑问就是numpy了。所以首先要做的就是将numpy加载进来，所以简单起见，也将此模块加载进来，还有之后产生随机数和记录时间的功能，所以还需要添加random和time:

import numpy as np

import time

import random

 

接下来就是对含有正态矩阵A的构建以及托普利兹矩阵B的构建：

矩阵A利用np.random.normal来构建，第一个参数是平均值，第二个是标准差，第三个则是大小（200行500列）

 

得到结果：

 

 

 

托普利兹矩阵则利用自己写的函数实现（scipy模块在本地无法使用lingal，原因未知）：

 

得到B矩阵：

 

 

9.1:矩阵相加相减都可以用运算符直接运算，但是乘法需要使用numpy.dot(a,b)

temp=A+A
print("A+A:")
print(temp)

temp=np.dot(A,A.T)
print("A*A.T:")
print(temp)

temp=np.dot(A.T,A)
print("A.T*A")
print(temp)

temp=np.dot(A,B)
print("A*B:")
print(temp)

 

 

9.2:利用函数numpy.linalg.inv得到B的逆，等式两边同时左乘B-1就可以得到X，X是一个n个元素的向量（此时n为200）

b=np.array(range(1000,1500))
print("the x of Bx=b is:")
print(np.dot(np.linalg.inv(B),b))

 

 

 

9.3：

print("|A|F:")
print(np.linalg.norm(A,ord='fro'))
print("")
print("||B||:")
print(np.linalg.norm(B,ord=np.inf))

U,sigma,VT=np.linalg.svd(B)

print("the max singular of B is;")
print(np.max(sigma))
print("the min singular of B is;")
print(np.min(sigma))

 

结果如下：

 

9.4:首先设置一个足够小的数字作为判别是否收敛的条件（这里取了0.5*10^-4），接着利用求特征值的算法（参考文献：https://en.wikipedia.org/wiki/Power_iteration

）即可得到最大特征根

Z=np.random.normal(100,100,n*n).reshape((n,n))
print("Z:")
print(Z)
y1=np.ones(n)
y2=np.zeros(n)
x=y1
t1=time.clock()
while abs(np.max(y1-y2))>1/2*(0.0001):
    y1=y2
    y2=np.dot(Z,x)
    x=y2/np.linalg.norm(y2)
t2=time.clock()
print("time cost:"+str(t2-t1))
print("x=")
print(x)
print("max eigenvalue:")
print(max(x))

 

结果如下：

 

 

9.5:利用9.3的求奇异值的方法得到奇异值。

p=float(input("The p :"))
temp=[1 if random.random()<p else 0 for i in range(0,n*n)]
C=np.array(temp).reshape(n,n)
print("C:")
print(C)
U,sigma,VT=np.linalg.svd(C)
e=np.max(sigma)
print("the eigenvalue is:")
print(e)

 

其中一次的结果：

 

经过反复测试，得到最大奇异值约等于n*p这一结论。

9.6:利用数组可以直接运用运算符“-”去减一个常数的特点，配合argmin可以很简洁的完成！

def argMin(A,z):
    return A[0][np.argmin(np.abs(A - z))]

K=np.random.randn(1,1000)
z=input("input the num:")
z=float(z)
a=argMin(K,z)
print(a)

其中一次输入测试：