NumPy练习

最新推荐文章于 2026-05-01 06:28:53 发布

原创最新推荐文章于 2026-05-01 06:28:53 发布 · 415 阅读

0 ·

本内容遵循CC 4.0 BY-SA版权协议

Generate matrices A, with random Gaussian entries, B, a Toeplitz matrix, where A ∈ R n×m and B ∈ R m×m ,

for n = 200, m = 500.

生成n×m矩阵A，其中各元素是服从高斯分布的随机数。生成常对角矩阵B。

代码：

import numpy as np
import scipy.linalg as la

n = 200
m = 500
A = np.mat(np.random.randn(n, m))
col_1 = np.random.random(size=(m, 1))
row_1 = np.random.random(size=(1, m))
B = np.mat(la.toeplitz(col_1, row_1))

Exercise 9.1: Matrix operations

Calculate A + A, AA > ,A > A and AB. Write a function that computes A(B − λI) for any λ.

代码：

def cal(lambdaa):
    return A*(B-lambdaa*I)

print(
    'A = '+str(A),
    'B = '+str(B),
    'A+A = '+str(A+A),
    'A*A.T = '+str(A*A.T),
    'A.T*A = '+str(A.T*A),
    'A*B = '+str(A*B),
    'I = '+str(I),
    'A(B − λI) =' + str(cal(5)),
    sep='\n\n'
)

输出：

A = [[ 0.04195854 -0.38710441  0.02414068 ...  0.12998791 -1.09805947
  -0.49228629]
 [ 0.28992152 -0.25311269  1.51203024 ...  0.87529587  0.42869605
  -1.01772641]
 [-0.01180153 -0.52115087 -1.06669412 ...  0.09279688 -0.18493147
   0.47616983]
 ...
 [-1.59083133  2.15117569  0.91646871 ... -0.07697636  0.58731201
   0.35030513]
 [-0.47906178  0.50324141 -0.14811297 ... -2.20791176 -1.22901378
  -0.64334274]
 [ 0.91530925 -1.58723489 -0.03392294 ... -0.91057087 -1.15192596
   0.89204535]]


B = [[0.82404747 0.34030171 0.04510667 ... 0.97886067 0.96370754 0.22601283]
 [0.77107993 0.82404747 0.34030171 ... 0.52657226 0.97886067 0.96370754]
 [0.01290139 0.77107993 0.82404747 ... 0.38751267 0.52657226 0.97886067]
 ...
 [0.89640885 0.48097235 0.85224666 ... 0.82404747 0.34030171 0.04510667]
 [0.99096391 0.89640885 0.48097235 ... 0.77107993 0.82404747 0.34030171]
 [0.42959832 0.99096391 0.89640885 ... 0.01290139 0.77107993 0.82404747]]


A+A = [[ 0.08391708 -0.77420882  0.04828136 ...  0.25997582 -2.19611895
  -0.98457258]
 [ 0.57984304 -0.50622538  3.02406047 ...  1.75059174  0.85739211
  -2.03545282]
 [-0.02360306 -1.04230174 -2.13338824 ...  0.18559375 -0.36986295
   0.95233966]
 ...
 [-3.18166265  4.30235138  1.83293741 ... -0.15395271  1.17462402
   0.70061025]
 [-0.95812355  1.00648281 -0.29622593 ... -4.41582353 -2.45802756
  -1.28668548]
 [ 1.8306185  -3.17446977 -0.06784589 ... -1.82114174 -2.30385192
   1.7840907 ]]


A*A.T = [[ 5.27537748e+02  1.34696479e-01 -1.24809925e+01 ... -3.93112193e+01
   2.15343154e+01  4.26919125e+01]
 [ 1.34696479e-01  5.41951934e+02  1.90018025e+01 ... -1.85704375e+01
   3.20290916e+01 -1.35193393e+01]
 [-1.24809925e+01  1.90018025e+01  4.60985958e+02 ...  3.89682737e+01
  -1.37186793e+01  3.72866375e+01]
 ...
 [-3.93112193e+01 -1.85704375e+01  3.89682737e+01 ...  4.88881188e+02
  -1.58172861e+01  3.50264454e+00]
 [ 2.15343154e+01  3.20290916e+01 -1.37186793e+01 ... -1.58172861e+01
   4.71088145e+02 -1.91876202e+01]
 [ 4.26919125e+01 -1.35193393e+01  3.72866375e+01 ...  3.50264454e+00
  -1.91876202e+01  5.07528739e+02]]


A.T*A = [[181.59339518  -3.4769892   13.47451816 ...  -2.85963471   4.77597727
   29.82138842]
 [ -3.4769892  205.4581638   -3.30189924 ...   0.36806718  20.29013161
  -20.0666659 ]
 [ 13.47451816  -3.30189924 181.18403999 ...  11.05256883 -16.15320599
    4.64771307]
 ...
 [ -2.85963471   0.36806718  11.05256883 ... 204.89356816  -1.18547402
   17.16772537]
 [  4.77597727  20.29013161 -16.15320599 ...  -1.18547402 198.06877941
  -13.69891785]
 [ 29.82138842 -20.0666659    4.64771307 ...  17.16772537 -13.69891785
  181.26927499]]


A*B = [[  4.25658264   6.68511215  14.31291597 ...  15.29729322  22.04987366
    7.30773622]
 [ -4.38481245   8.37493151  12.30946257 ...  12.13157489   7.13387606
    6.31249827]
 [ -6.37549539   4.61405134   2.29595035 ...   7.74256197  -8.35348835
   -9.42119234]
 ...
 [  8.90184728   9.78111981   2.04445731 ...   8.66266684  -3.54664847
    5.99144587]
 [  0.78702146  -6.53679439 -12.00846036 ...   4.2966042   -5.94174445
    9.52998304]
 [-20.97034986 -11.18665034  -3.79174184 ... -21.56821008  -5.90990218
  -14.0203237 ]]


I = [[1. 0. 0. ... 0. 0. 0.]
 [0. 1. 0. ... 0. 0. 0.]
 [0. 0. 1. ... 0. 0. 0.]
 ...
 [0. 0. 0. ... 1. 0. 0.]
 [0. 0. 0. ... 0. 1. 0.]
 [0. 0. 0. ... 0. 0. 1.]]


A(B − λI) =[[  4.04678994   8.62063419  14.19221256 ...  14.64735367  27.54017104
    9.76916766]
 [ -5.83442004   9.64049495   4.74931139 ...   7.75509554   4.99039579
   11.40113032]
 [ -6.31648775   7.21980568   7.62942095 ...   7.27857758  -7.42883099
  -11.8020415 ]
 ...
 [ 16.85600391  -0.97475862  -2.53788623 ...   9.04754861  -6.48320852
    4.23992024]
 [  3.18233034  -9.05300142 -11.26789553 ...  15.33616302   0.20332445
   12.74669674]
 [-25.54689611  -3.2504759   -3.62212712 ... -17.01535573  -0.15027239
  -18.48055046]]

Exercise 9.2: Solving a linear system

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

代码：

b = np.random.random(size=(m, 1))

x = np.mat(np.linalg.solve(B, b))
print('x = '+str(x))

输出：

x = [[ 4.38225617e-01]
 [-2.62604279e-02]
 [ 5.77628412e-02]
 ...
 [-3.41521726e-01]
 [ 9.07959804e-02]
 [-5.83334792e-01]]

经验证Bx = b

Exercise 9.3: Norms
Compute the Frobenius norm of A: kAk F and the infinity norm of B: kBk ∞ . Also find the largest and

smallest singular values of B.

代码：

AF = np.linalg.norm(A)
Binf = np.linalg.norm(B, np.inf)
s = la.svdvals(B)
print(
    'Frobenius norm of A: '+str(AF),
    'Infinity norm of B: '+str(Binf),
    'The largest singular values of B: '+str(s.max()),
    'The smallest singular values of B: ' + str(s.min()),
    sep = '\n'
)

输出：

Frobenius norm of A: 315.95558444399796
Infinity norm of B: 259.47901617509274
The largest singular values of B: 253.9171417278695
The smallest singular values of B: 0.04916465790673367

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 time

def power_iteration(A):
    n, d = A.shape
    count = 0

    v = np.ones(d)
    ev = 0

    begin_time = time.clock()

    while True:
        count += 1
        Av = A.dot(v)
        ev_new = np.linalg.norm(Av)
        v_new = Av / ev_new

        if np.abs(ev - ev_new) < 0.01:
            break

        v = v_new
        ev = ev_new

    end_time = time.clock()
    print('Time: ' + str(end_time-begin_time))

    return ev_new, v_new, count

n = 200
Z = np.random.randn(n, n)

ev, v, count = power_iteration(Z)
print(
    'eigenvalue: ' + str(ev),
    'eigenvector: ' + str(v),
    'number of iteration: ' + str(count),
    sep='\n'
)

输出：

Time: 0.0018198326227696195
eigenvalue: 15.81618091989844
eigenvector: [-0.03700825  0.0029305   0.01227961  0.02983836  0.01829982  0.02772856
  0.02107871  0.09893484  0.0350257   0.04419609  0.04727178  0.056736
 -0.04377146  0.09618859 -0.01862322  0.01417478 -0.02309586 -0.00194714
 -0.03120478 -0.02046295 -0.0968875  -0.03675232  0.03318669  0.10242163
  0.05797829  0.07278333  0.02818215 -0.01875538  0.12128094 -0.07273369
 -0.01526538 -0.17105178 -0.11733502 -0.01788494 -0.0601756  -0.06721134
  0.09800221  0.01506051  0.02491852 -0.12336477  0.02428129 -0.02547432
 -0.03982852  0.09438698 -0.09507363 -0.10752629 -0.03039604  0.02625876
  0.02946467 -0.03132434 -0.00859429  0.02732496 -0.02247782  0.06659113
 -0.01402789  0.04534168  0.03565373  0.00755469 -0.0681837  -0.06179226
 -0.01214626 -0.03232888 -0.00067058 -0.06898048  0.02596048 -0.04140793
 -0.13078973  0.05529428 -0.05294517  0.00409866  0.01130606  0.08446122
  0.16050908  0.01138752  0.03961415 -0.1426926   0.05989797 -0.13517472
  0.00196125 -0.06839017  0.01876062  0.00799628  0.07576671  0.13489427
 -0.08642153  0.01101552 -0.01891506 -0.03754352 -0.11884894  0.05107088
  0.00951774 -0.02215805 -0.04574779 -0.08041138  0.08242144 -0.19971271
 -0.03669991 -0.0268403   0.07627959 -0.09874723 -0.03435932 -0.00241873
 -0.06670095 -0.01555391  0.05338037 -0.21958027 -0.05404186  0.03915557
 -0.10399761  0.08424245  0.07007432  0.02019906 -0.00034468 -0.01125532
 -0.04812383 -0.07710452  0.09902331  0.00081299 -0.08204022 -0.05218078
  0.01894426 -0.0911309  -0.03287158  0.02815841  0.02207783 -0.05711911
  0.00441251 -0.21093863  0.0020175  -0.02082109 -0.11669596 -0.06667732
  0.07087201  0.03347427  0.02888406 -0.02505937 -0.09300211  0.07533929
  0.05001491 -0.02080881  0.09205245 -0.07718376 -0.01268754  0.0559997
  0.08342342  0.08266444 -0.05132904 -0.04049663 -0.11326853 -0.14141689
  0.1126093  -0.01487698  0.09328524 -0.01179341  0.09013115  0.01899532
  0.04128708  0.04779245  0.01369805 -0.05740939  0.04079265 -0.00880687
  0.02980492 -0.12020824  0.02103489 -0.04582447  0.05758544 -0.00679323
  0.04331108  0.03973211 -0.01544934  0.09382855 -0.06072975 -0.0369994
  0.08413063  0.09077527  0.08635241  0.04281156  0.04731971 -0.04836914
  0.00405065  0.02900649  0.08715658  0.09654352  0.05667637  0.00851858
  0.02996912 -0.12203476 -0.06564825  0.0525448   0.03626988  0.10106428
  0.09465614 -0.08401987 -0.20462801 -0.09707087 -0.07242656 -0.0545679
 -0.14721702  0.02890988]
number of iteration: 28

所需时间如下图，矩阵大小从100×100到2000×2000，步长为100，对每个大小运行10次取均值。

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 matplotlib.pyplot as plt
from matplotlib import cm
from mpl_toolkits.mplot3d import Axes3D

def svd(n, p):
    C = np.random.binomial(1, p, n)
    for _ in range(n-1):
        row = np.random.binomial(1, p, n)
        C = np.row_stack((C, row))
    return la.svdvals(C)

X = np.arange(100, 1100, 100)
Y = np.arange(0.1, 1.1, 0.1)
X, Y = np.meshgrid(X, Y)
Z = [svd(x,y).max() for x, y in zip(X.flat, Y.flat)]
Z = np.reshape(Z, (10,10))

fig = plt.figure()
ax = Axes3D(fig)
ax.plot_surface(X, Y, Z, rstride=1, cstride=1, cmap=cm.viridis)
plt.xlabel('n')
plt.ylabel('p')

plt.show()

输出：

可知最大奇异值虽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.

代码：

z = np.random.random()

a = np.reshape(A, (1, 100000))
index = np.argmin(np.abs((a-z)))
print(
    'z = ' + str(z),
    'the element in A that is closest to z: ' + str(a[0,index]),
    sep = '\n'
)

输出：

z = 0.7783356000465671
the element in A that is closest to z: 0.7782910522477339