7. 固有値と固有ベクトル

Pythonで学ぶ線形代数学

問7.4.

code: prob1.py

from sympy import *

from numpy.random import choice

N = -3, -2, -1, 1, 2, 3

def f(truth):

while True:

A = Matrix(choice(N, (2, 2)))

eigenvals = A.eigenvals()

if len(eigenvals) == 2 and 0 not in eigenvals:

if all(x.is_real for x in eigenvals) == truth:

print(eigenvals)

return A

ソースコード: prob1.py

実行例:

code: python

>> A = f(); A

{1 - sqrt(6): 1, 1 + sqrt(6): 1}

Matrix([

-1, 1,

2, 3])

>> B = f(False); B

{-1/2 - sqrt(3)*I/2: 1, -1/2 + sqrt(3)*I/2: 1}

Matrix([

-2, -1,

3, 1])

$ \begin{bmatrix}1&1\\0&1\\\end{bmatrix}の固有値、固有ベクトル

code: eig2.py

import sympy as sp

import numpy as np

A = 1, 1], [0, 1

a = sp.Matrix(A).eigenvects()

b = np.linalg.eig(A)

print(f'''SymPy

eigen value: {a00}

multiplicity: {a01}

eigen vector:

{a020}

NumPy

eigen values: {b00}, {b01}

eigen vectors:

{b1:, 0}

{b1:, 1}''')

ソースコード: eig2.py

実行結果:

code: python

SymPy

eigen value: 1

multiplicity: 2

eigen vector:

Matrix(1], [0)

NumPy

eigen values: 1.0, 1.0

eigen vectors:

1. 0.

-1.00000000e+00 2.22044605e-16

問7.5.

code: prob2.py

from sympy import *

from numpy.random import choice

D = -5, -4, -3, -2, -1, 1, 2, 3, 4, 5

def f():

while True:

A = Matrix(choice(D, (3, 3)))

cp = A.charpoly(Symbol('lmd'))

F = factor_list(cp)

if len(F1) == 3:

print(f'det(A - lmd*I) = {factor(cp)}\nA = {A}')

return A

ソースコード: prob2.py

実行例:

code: python

>> A = f(); A

det(A - lmd*I) = (lmd - 3)*(lmd - 1)*(lmd + 6)

A = Matrix(-4, -2, -3], -5, -1, -3, [-4, 4, 3)

Matrix([

-4, -2, -3,

-5, -1, -3,

-4, 4, 3])

問7.10.

code: prob3.py

from sympy import Matrix

from numpy.random import choice

N = -3, -2, -1, 1, 2, 3

def g(symmetric=True):

if symmetric:

a, b, d = choice(N, 3)

return Matrix(a, b], [b, d)

else:

a, b = choice(N, 2)

return Matrix(a, b], [-b, a)

ソーシコード: prob3.py

実行例:

code: python

>> A = g(); A

Matrix([

-1, -2,

-2, -3])

>> B = g(False); B

Matrix([

-1, -2,

2, -1])

行列ノルム

code: unitcircle.py

from numpy import array, arange, pi, sin, cos, isreal

from numpy.linalg import eig, norm

import matplotlib.pyplot as plt

def arrow(p, v, c=(0, 0, 0), w=0.02):

plt.quiver(p0, p1, v0, v1, units='xy', scale=1,

color=c, width=w)

n = 0

A = [array(1, -2], [2, 2),

array(3, 1], [1, 3),

array(2, 1], [0, 2),

array(2, 1], [0, 3)]

T = arange(0, 2 * pi, pi / 500)

U = array((cos(t), sin(t)) for t in T)

plt.plot(U:, 0, U:, 1)

V = array([An.dot(u) for u in U])

plt.plot(V:, 0, V:, 1)

for u, v in zip(U::20, V::20):

arrow(u, v - u)

o = array(0, 0)

Lmd, Vec = eig(An)

if isreal(Lmd0):

arrow(o, Lmd0 * Vec:, 0, c=(1, 0, 0), w=0.1)

if isreal(Lmd1):

arrow(o, Lmd1 * Vec:, 1, c=(0, 1, 0), w=0.1)

plt.axis('scaled'), plt.xlim(-4, 4), plt.ylim(-4, 4), plt.show()

ソースコード: unitcircle.py

https://gyazo.com/0c81b661e1c26f2b255fd9d0f08038d0https://gyazo.com/330c0e62a0add6e3d8733f87c24fa714https://gyazo.com/8434a5584ec8c539282fecb72cec7105https://gyazo.com/dead3cd14b738c65e197a3012edf1a04

スペクトル半径、数域半径

code: matrixnorm.py

from numpy import array, arange, pi, sin, cos

from numpy.linalg import eig, norm

M = [array(1, 2], [2, 1),

array(1, 2], [-2, 1),

array(1, 2], [3, 4)]

T = arange(0, 2 * pi, pi / 500)

U = array((cos(t), sin(t)) for t in T)

for A in M:

r1 = max(abs((A.dot(u)).dot(u)) for u in U)

r2 = max([abs(e) for e in eig(A)0])

r3 = max(norm(A.dot(u)) for u in U)

print(f'{A}: num={r1:.2f}, spec={r2:.2f}, norm={r3:.2f}')

ソースコード: matrixnorm.py

実行結果:

code: python

[1 2

2 1]: num=3.00, spec=3.00, norm=3.00

[ 1 2

-2 1]: num=1.00, spec=2.24, norm=2.24

[1 2

3 4]: num=5.42, spec=5.37, norm=5.46

行列の指数関数

code: exp_np.py

from numpy import matrix, e, exp, diag

from numpy.linalg import eigh

A = matrix(1, 2], [2, 1)

m, B = 1, 0

for n in range(10):

B += A ** n / m

m *= n + 1

print(B)

a = eigh(A)

S, V = diag(e**a0), a1

print(V * S * V.H)

print(exp(A))

ソースコード: exp_np.py

実行結果:

code: python

[10.21563602 9.84775683

9.84775683 10.21563602]

[10.22670818 9.85882874

9.85882874 10.22670818]

[2.71828183 7.3890561

7.3890561 2.71828183]