第11周作业 #高级编程技术

numpy

---

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 numpy as np
from scipy.linalg import toeplitz
A = np.random.normal(2, 1, [200, 500])
B = toeplitz(np.random.normal(2, 1, 500), np.random.normal(2, 1, 500))

Exercise 9.1

Matrix operations Calculate A + A, AAT, ATA and AB. Write a function that computes A(B−λI) for any λ.

print("A:", A)
print("B:", B)
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))

def foo(λ):
    return np.dot(A, B - λ)

运行结果：(数据量用了2×5和5×5)

A: [[1.65428653 1.29250801 0.12205231 0.74224551 2.11893814]
 [1.31702596 1.57987211 1.76251052 1.09167503 0.72512259]]
B: [[1.75756406 2.38712337 0.4477556  3.09135636 1.92525955]
 [2.42381625 1.75756406 2.38712337 0.4477556  3.09135636]
 [2.31593935 2.42381625 1.75756406 2.38712337 0.4477556 ]
 [2.61638874 2.31593935 2.42381625 1.75756406 2.38712337]
 [2.29753941 2.61638874 2.31593935 2.42381625 1.75756406]]
A + A: [[3.30857306 2.58501602 0.24410461 1.48449102 4.23787628]
 [2.63405192 3.15974422 3.52502104 2.18335006 1.45024519]]
AAT: [[9.46296488 6.78263494]
 [6.78263494 9.05455374]]
ATA: [[4.47122129 4.21891117 2.52318159 2.6656511  4.4603361 ]
 [4.21891117 4.16657284 2.94229479 2.6840652  3.88434548]
 [2.52318159 2.94229479 3.12134009 2.0146815  1.53665748]
 [2.6656511  2.6840652  2.0146815  1.74268277 2.36437055]
 [4.4603361  3.88434548 1.53665748 2.36437055 5.01570162]]
AB: [[13.13332887 13.77944548 10.74700582 12.42453151 12.73118452]
 [14.74818863 14.61808725 11.78414015 12.66237593 12.08915085]]

---

Exercise 9.2

Solving a linear system Generate a vector b with m entries and solve Bx = b.

print("B:")
print(B)

b = np.random.random(5)
print("b:", b)
x = np.linalg.solve(B, b)
print("x:", x)
print("np.dot(B, x):\n", np.dot(B, x))

运行结果：(数据量用了5×5)

B:
[[ 2.5121108   1.84558016  1.69360915  1.14977797  3.51985693]
 [ 2.89274206  2.5121108   1.84558016  1.69360915  1.14977797]
 [ 3.59327592  2.89274206  2.5121108   1.84558016  1.69360915]
 [ 0.87142587  3.59327592  2.89274206  2.5121108   1.84558016]
 [-0.09840123  0.87142587  3.59327592  2.89274206  2.5121108 ]]
b: [0.63848187 0.88284346 0.00751235 0.26442866 0.93755242]
x: [ 0.23215661 -0.41407927 -2.29703059  2.99870062  0.35851666]
np.dot(B, x):
 [0.63848187 0.88284346 0.00751235 0.26442866 0.93755242]

---

Exercise 9.3

Norms Compute the Frobenius norm of A: ||A||F and the infinity norm of B: ||B||∞. Also find the largest and smallest singular values of B.

print("A:")
print(A)
print("B:")
print(B)

print("the Frobenius norm of A:", np.linalg.norm(A))
print("the infinity norm of B:", np.linalg.norm(B, np.inf))
print("the largest singular value of B:", np.linalg.norm(B, 2))
print("the smallest singular value of B:", np.linalg.norm(B, -2))

运行结果：(数据量用了2×5和5×5)

A:
[[ 0.89763757  3.90870489  2.24489208  1.06493451 -0.38395911]
 [ 1.9064736  -0.60809743  2.44124628  3.29251693  2.95935474]]
B:
[[3.26704861 1.56438763 1.90880029 3.5782664  2.64795173]
 [2.55789806 3.26704861 1.56438763 1.90880029 3.5782664 ]
 [3.60588272 2.55789806 3.26704861 1.56438763 1.90880029]
 [1.30917891 3.60588272 2.55789806 3.26704861 1.56438763]
 [0.96130538 1.30917891 3.60588272 2.55789806 3.26704861]]
the Frobenius norm of A: 7.208837168639553
the infinity norm of B: 12.966454652614956
the largest singular value of B: 12.56650752516994
the smallest singular value of B: 1.8246744060618487

---

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.

import numpy as np
import time
def foo(n):
    Z = np.random.normal(0, 1, [n, n])
    b_k = np.random.rand(Z.shape[0])
    b_k1 = np.random.rand(Z.shape[0])
    while b_k.all() != b_k1.all():
        b_k1 = np.dot(Z, b_k)
        b_k1_norm = np.linalg.norm(b_k1)
        b_k = b_k1 / b_k1_norm
    # print(b_k)
    # print(np.dot(Z, b_k))
    print("n:", n, "time:", time.clock())
t = time.clock()
foo(5)
foo(10)
foo(100)
foo(1000)
foo(10000)

运行结果：

n: 5 time: 0.00014279107049498438
n: 10 time: 0.0002451910413678816
n: 100 time: 0.001004088603281464
n: 1000 time: 0.07183528623351858
n: 10000 time: 5.034103741410491

---

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?

import numpy as np
import scipy as sp
# print(help(sp.linalg.eigvals))
def foo(p, n):
    C = np.random.binomial(1, p, [n, n])
    # print(sp.linalg.eig(C))
    print("p:", p, "n:", n, max(sp.linalg.eigvals(C)))


foo(.5, 10)
foo(.5, 20)
foo(.5, 100)
foo(.3, 10)
foo(.3, 20)
foo(.3, 100)
foo(.9, 10)
foo(.9, 20)
foo(.9, 100)

运行结果：

p: 0.5 n: 10 (5.641844752967631+0j)
p: 0.5 n: 20 (9.242385892004453+0j)
p: 0.5 n: 100 (50.14720097026784+0j)
p: 0.3 n: 10 (3.456889014392734+0j)
p: 0.3 n: 20 (6.498898769869614+0j)
p: 0.3 n: 100 (30.05949038152027+0j)
p: 0.9 n: 10 (9.010646243839297+0j)
p: 0.9 n: 20 (17.458772002141625+0j)
p: 0.9 n: 100 (90.72610340531895+0j)

规律：

$$\lim_{n\to∞}max(eigvals)=np$$

---

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.

import numpy as np

print(help(np.argmin))
def foo(A, z):
    B = np.abs(A - z)
    return A[np.argmin(B)]

a = np.random.randint(0, 100, 6)
print(a)
b = np.random.randint(100)
print("b:", b, "Nearest neighbor:", foo(a,b))

运行结果：

[17 93 87 77 20 26]
b: 9 Nearest neighbor: 17

---