1. 程式人生 > >高階程式設計技術_課後作業(十四)

高階程式設計技術_課後作業(十四)

Numpy

下面為作業文件中Numpy這一章的內容

題目

解答

Exercise 9.1

高斯隨機分佈矩陣的生成:

numpy.random.normal(loc=0.0, scale=1.0, size=None)  ,其中loc為均值,scale為方差,size為輸出形狀乘法:numpy中*是矩陣元素逐個計算,np.dot(a,b)是按照矩陣乘法的運演算法則來運算toeplitz--生成託普利茲矩陣
T=toeplitz(c,r)生成非對稱託普利茲矩陣,將c作為第一列,r作為第一行,若c(1)與r(1)不相等,則使用c(1)作為矩陣的第一個元素,同時列印一條警告資訊。

程式碼:

import numpy as np
from scipy.linalg import toeplitz

n = 200
m = 500

A = np.random.normal(size = (n,m))

print("A+A\n")
print(A+A)

print("AA'\n")
print(np.dot(A, A.T))

print("A'A\n")
print(np.dot(A.T, A))

c = np.random.normal(size = (1,500))
r = np.random.normal(size = (1,500))
B = toeplitz(c,r)

print("AB\n")
print(np.dot(A,B))

def A_B_I(r):
    A = np.random.normal(size = (n,m))
    c = np.random.normal(size = (1,m))
    r = np.random.normal(size = (1,m))
    B = toeplitz(c,r)
    I = np.identity(m)
    return (np.dot(A, (B - r*I)))
print("A(B-rI)\n")
print(A_B_I(5))

結果截圖:


Exercise 9.2

程式碼:

import numpy as np
from scipy.linalg import toeplitz

m = 500

c = np.random.normal(size = (1,m))
r = np.random.normal(size = (1,m))
B = toeplitz(c,r)

b = np.random.permutation(m)

ans = np.linalg.solve(B, b)

print("Ans:\n")
print(ans)

結果截圖:


Exercise 9.3

程式碼:

import numpy as np
from scipy.linalg import toeplitz

n = 200
m = 500

A = np.random.normal(size = (n,m))

c = np.random.normal(size = (1,500))
r = np.random.normal(size = (1,500))
B = toeplitz(c,r)

frobenius = np.linalg.norm(A, "fro")
infinity = np.linalg.norm(B, np.inf)
e = np.linalg.eigvals(B)

print("the Frobenius norm of A: ")
print(frobenius)

print("\nthe infinity norm of B")
print(infinity)

print("\nthe largest singular values of B: ")
print(max(e))

print("\nthe smallest singular values of B: ")
print(min(e))

結果截圖:


Exercise 9.4

程式碼:

import time
import numpy as np
from scipy.linalg import toeplitz

for n in range(2, 6):
    start_CPU = time.clock()
    i = 1
    eig = 0

    Z = np.random.normal(1,1,size = (n,n))
    u = np.random.randn(n)
    print("##Test of " + str(n))

    while True:
        v = np.dot(Z, u)
        last = eig
        #eig = np.linalg.norm(v, np.inf)
        eig = v[np.argmax(np.abs(v))]
        u = v / eig
        if (i != 1) & (abs(last - eig) < (1/2)*(10**(-5))):
            break
        i = i + 1

    print("迭代次數:")
    print(i)
    end_CPU = time.clock()
    print("\nThe time of " + str(n) + " is " + str(end_CPU - start_CPU))
    print("\n")
    print("The lagest eigenvalue of Z")
    print(eig)
    print("\n")
    print("Its corresponding eigenvector")
    print(u)
    print("\n")

結果截圖:


Exercise 9.5

程式碼:

import time
import numpy as np
from scipy.linalg import svdvals

n = 4
p = .5
C = np.random.binomial(1, p, size = (n, n))
print(C)
lsv = svdvals(C)
print(lsv.max())

結果截圖:


Maybe the largest singular value = n * p

Exercise 9.6

程式碼:

import numpy as np


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

n = 5
A =  np.random.normal(size = (n,1))
z = np.random.randint(1)
print(z)
ans = find_Closest_Element(A, z)
#ans = find_Closest_Element(ans, z)
print(A)
print(ans)

結果截圖: