/linalg/5. 行列の基本変形と不変量

generated at 2/17/2025, 12:24:57 PM
5. 行列の基本変形と不変量
Pythonで学ぶ線形代数学

基本行列による単位立方体の変化
 elementary_vp.pyfrom vpython import *
import numpy as np

o = vec(0, 0, 0)
x, y, z = vec(1, 0, 0), vec(0, 1, 0), vec(0, 0, 1)
X, Y, Z = [o, y, z, y+z], [o, z, x, z+x], [o, x, y, x+y]
box(pos=(x+y+z)/2)
for v in [x, y, z]:
    curve(pos=[-v, 3*v], color=v)


def T(A, u):
    return vec(*np.dot(A, (u.x, u.y, u.z)))


E = [[1, 0, 0], [0, 1, 0], [0, 0, 1]]
E1 = [[1, 2, 0], [0, 1, 0], [0, 0, 1]]
E2 = [[0, 1, 0], [1, 0, 0], [0, 0, 1]]
E3 = [[2, 0, 0], [0, 1, 0], [0, 0, 1]]

for u, U in [(x, X), (y, Y), (z, Z)]:
    for v in U:
        A = E
        curve(pos=[T(A, v), T(A, u+v)], color=u)
ソースコード: elementary_vp.py



基本行列
 elementary_sp.pyfrom sympy import Matrix, var

var('x y a11 a12 a13 a21 a22 a23 a31 a32 a33')
E1 = Matrix([[1, x, 0], [0, 1, 0], [0, 0, 1]])
E2 = Matrix([[0, 1, 0], [1, 0, 0], [0, 0, 1]])
E3 = Matrix([[1, 0, 0], [0, y, 0], [0, 0, 1]])
A = Matrix([[a11, a12, a13], [a21, a22, a23], [a31, a32, a33]])
ソースコード: elementary_sp.py

対話モードでの実行例（基本行列を左から掛けた場合））
 pyhon>>> E1 * A
Matrix([
[a11 + a21*x, a12 + a22*x, a13 + a23*x],
[        a21,         a22,         a23],
[        a31,         a32,         a33]])
>>> E2 * A
Matrix([
[a21, a22, a23],
[a11, a12, a13],
[a31, a32, a33]])
>>> E3 * A
Matrix([
[  a11,   a12,   a13],
[a21*y, a22*y, a23*y],
[  a31,   a32,   a33]])

対話モードでの実行例（基本行列を右から掛けた場合））
 python>>> A * E1
Matrix([
[a11, a11*x + a12, a13],
[a21, a21*x + a22, a23],
[a31, a31*x + a32, a33]])
>>> A * E2
Matrix([
[a12, a11, a13],
[a22, a21, a23],
[a32, a31, a33]])
>>> A * E3
Matrix([
[a11, a12*y, a13],
[a21, a22*y, a23],
[a31, a32*y, a33]])


階数の問題生成プログラム
 prob_rank.pyfrom numpy.random import seed, choice, permutation
from sympy import Matrix


def f(P, m1, m2, n):
    if n > min(m1, m2):
        return Matrix(choice(P, (m1, m2)))
    else:
        while True:
            X, Y = choice(P, (m1, n)), choice(P, (n, m2))
            A = Matrix(X.dot(Y))
            if A.rank() == n:
                return A


m1, m2 = 3, 4
seed(2020)
for i in permutation(max(m1, m2)):
    print(f([-3, -2, -1, 1, 2, 3], m1, m2, i+1))
ソースコード: prob_rank.py

実行例
 pythonMatrix([[4, 5, -5, 4], [8, 1, -17, 16], [0, -3, -3, 0]])
Matrix([[-10, 4, -8, -8], [7, -2, 8, 4], [0, -1, -3, 2]])
Matrix([[-2, 2, -2, 1], [-1, -3, -2, -2], [-1, -1, -2, -3]])
Matrix([[3, -1, -3, -1], [9, -3, -9, -3], [-3, 1, 3, 1]])


定義に従って行列式を計算する
 determinant.pyfrom functools import reduce


def P(n):
    if n == 1:
        return [([0], 1)]
    else:
        Q = []
        for p, s in P(n-1):
            Q.append((p + [n-1], s))
            for i in range(n-1):
                q = p + [n-1]
                q[i], q[-1] = q[-1], q[i]
                Q.append((q, -1*s))
        return Q


def prod(L):
    return reduce(lambda x, y: x*y, L)


def det(A):
    n = len(A)
    a = sum([s*prod([A[i][p[i]] for i in range(n)]) for p, s in P(n)])
    return a


if __name__ == '__main__':
    A = [[1, 2], [2, 3]]
    B = [[1, 2], [2, 4]]
    C = [[1, 2, 3], [2, 3, 4], [3, 4, 5]]
    D = [[1, 2, 3], [2, 3, 1], [3, 1, 2]]
    print(det(A), det(B), det(C), det(D))
ソースコード: determinant.py

実行結果
 python-1 0 0 -18


問5.5の解答例(n=5)の場合
 random_matrix.pyfrom numpy.random import seed, randint
from numpy.linalg import det


def f(n):
    S = 0
    for i in range(10000):
        D = det(randint(0, 10, (n, n)))
        if abs(D) < 1.0e-10:
            print(D)


f(5)
ソースコード: random_matrix.py

実行例
 python0.0
7.087663789206989e-13
0.0
0.0
-7.505107646466066e-13
-1.3429257705865887e-12
-1.0924594562311528e-12
-1.3522516439934362e-12
6.03961325396085e-13
非正則行列でも、行列式の値が正確に0にならないことがあることに注意


行列式（逆行列）の問題生成プログラム
 prob_det.pyfrom numpy.random import seed, choice, permutation
from sympy import Matrix


def f(P, m, p):
    while True:
        A = Matrix(choice(P, (m, m)))
        if p == 0:
            if A.det() == 0:
                return A
        elif A.det() != 0:
            return A


m = 3
seed(2020)
for p in permutation(2):
    print(f([-3, -2, -1, 1, 2, 3], m, p))
ソースコード: prob_det.py

実行例
 pythonMatrix([[-3, 1, 1], [1, 3, 1], [-3, 3, -3]])
Matrix([[-2, 1, -1], [-3, 2, 3], [3, -2, -3]])


連立方程式の問題生成プログラム
 prob_eqn.pyfrom numpy.random import seed, choice, shuffle
from sympy import Matrix, latex, solve, zeros
from sympy.abc import x, y, z


def f(P, m, n):
    while True:
        A = Matrix(choice(P, (3, 4)))
        if A[:, :3].rank() == m and A.rank() == n:
            break
    A, b = A[:, :3], A[:, 3]
    u = Matrix([[x], [y], [z]])
    print(f'{latex(A)}{latex(u)}={latex(b)}')
    print(solve(A*u - b, [x, y, z]))


seed(2020)
m, n = 2, 2
f(range(2, 10), m, n)
ソースコード: prob_eqn.py

実行例
 python\left[\begin{matrix}7 & 8 & 8\\5 & 8 & 6\\8 & 8 & 9\end{matrix}\right]\left[\begin{matrix}x\\y\\z\end{matrix}\right]=\left[\begin{matrix}3\\5\\2\end{matrix}\right]
{y: 5/4 - z/8, x: -z - 1}


余因子行列を用いた逆行列の計算
 inv.pyfrom numpy import array, linalg, random

n = 5
A = random.randint(0, 10, (n, n))
K = [[j for j in range(n) if j != i] for i in range(n)]
B = array([[(-1) ** (i+j) * linalg.det(A[K[i], :][:, K[j]])
            for i in range(n)] for j in range(n)])

print(A.dot(B/linalg.det(A)))
ソースコード: inv.py

実行結果
 python[[ 1.00000000e+00 -2.66453526e-15 -2.66453526e-15  1.59872116e-14
  -1.15463195e-14]
 [-5.55111512e-16  1.00000000e+00 -1.33226763e-15  4.44089210e-15
  -2.22044605e-15]
 [-4.44089210e-16 -6.21724894e-15  1.00000000e+00  2.13162821e-14
  -1.95399252e-14]
 [-1.77635684e-15 -5.32907052e-15 -1.66533454e-15  1.00000000e+00
  -1.77635684e-14]
 [-1.77635684e-15 -6.21724894e-15 -1.99840144e-15  1.42108547e-14
   1.00000000e+00]]
単位行列になっている（非対角要素にわずかな誤差を含む）