Generate matrices A, with random Gaussian entries, B, a Toeplitz matrix, where A ∈Rn×m and B ∈Rm×m, for n = 200, m = 500.
import random
from scipy.linalg import toeplitz
import numpy as np
n = 200
m = 500
#A n*m
#B m*m
A = np.random.normal(size=(n,m))
B = toeplitz(np.random.random((m,)), np.random.random((m,)))
Exercise 9.1:
Matrix operations Calculate A + A, AA>,A>A and AB. Write a function that computes A(B−λI) for any λ.
print("Exercise 9.1")
print("A+A:", A + A)
print("AAT:", np.dot(A, A.T))
print("ATA:", np.dot(A.T, A))
print("AB:", np.dot(A, B))
Exercise 9.2:
Solving a linear system Generate a vector b with m entries and solve Bx = b.
print("Exercise 9.2")
b = np.random.random(size=(m,1))
x = np.linalg.solve(B,b)
print("Bx=b:")
print("b:",b)
print("x:",x)
Exercise 9.3:
Norms Compute the Frobenius norm of A: kAkF and the infinity norm of B: kBk∞. Also find the largest and smallest singular values of B.
print("Exercise 9.3")
print("AF:",np.linalg.norm(A,ord='fro'))
print("BInfinity:", np.linalg.norm(B, np.inf))
e, v = np.linalg.eig(B)
print("max eigen", e.max())
print("min eigen", e.min())
Exercise 9.4:
Power iteration Generate a matrix Z, n × n, with Gaussian entries, and use the power iteration to find the largest eigenvalue and corresponding eigenvector of Z. How many iterations are needed till convergence?
Optional: use the time.clock() method to compare computation time when varying n.
print("Exercise 9.4")
Z = np.random.normal(size=(n,n))
v = np.ones((n,))
def lpower(A, v, eps=1e-6, N=1e8):
count = 1
y = v/v[np.argmax(np.abs(v))]
x = np.dot(A, v)
b = x[np.argmax(np.abs(x))]
if (np.abs(x/b - y) < eps).all():
t = b
return (t, y, count)
while (np.abs(x/b - y) > eps).any() and count < N:
count = count+1
m = x[np.argmax(np.abs(x))]
y = x/m
x = A.dot(y)
b = x[np.argmax(np.abs(x))]
return (b, y, count)
e, vec, count= lpower(Z, v)
print("特征值:", e)
print("特征向量:",vec)
print("迭代次数为:", count)
关于迭代次数与n的关系,理论上迭代次数很大程度上依赖于特征值的特征。但是当n变大时,程序就很难运行出结果了。
Exercise 9.5:
Singular values Generate an n×n matrix, denoted by C, where each entry is 1 with probability p and 0 otherwise. Use the linear algebra library of Scipy to compute the singular values of C. What can you say about the relationship between n, p and the largest singular value?
print("Exercise 9.5")
import scipy
p = 0.5
C = np.where(np.random.random((n,n)) < p, 1, 0)
U, s, V = scipy.linalg.svd(C)
print(n)
print("n*p=",p*n)
print("s.max()=",s.max())
通过观察,n*p近似于最大的奇异值。
Exercise 9.6:
Nearest neighbor Write a function that takes a value z and an array A and finds the element in A that is closest to z. The function should return the closest value, not index.
Hint: Use the built-in functionality of Numpy rather than writing code to find this value manually. In particular, use brackets and argmin.
print("Exercise 9.6")
def findNearestEle(A, z):
n, n = A.shape
index = np.argmin(np.abs(A-z))
return A[index//n, index%n]