行列式の値を計算するプログラムをPythonで作ってみた。(2)

 前回のプログラムをちょっと改良してみました。(^_^;
 1次と4次のときは展開結果を利用して、それ以外はLaplace展開を用いて計算するようにしました。(^_^;
 ここで、4次の行列式の展開式は次式を用いました。(^_^;

\left|\begin{array}{cc}a&b&s&t\\c&d&u&v\\e&f&w&x\\g&h&y&z\end{array}\right|\\=+\left|\begin{array}{cc}a&b\\c&d\end{array}\right|\left|\begin{array}{cc}w&x\\y&z\end{array}\right|-\left|\begin{array}{cc}a&b\\e&f\end{array}\right|\left|\begin{array}{cc}u&v\\y&z\end{array}\right|+\left|\begin{array}{cc}a&b\\g&h\end{array}\right|\left|\begin{array}{cc}u&v\\w&x\end{array}\right|+\left|\begin{array}{cc}c&d\\e&f\end{array}\right|\left|\begin{array}{cc}s&t\\y&z\end{array}\right|-\left|\begin{array}{cc}c&d\\g&h\end{array}\right|\left|\begin{array}{cc}s&t\\w&x\end{array}\right|+\left|\begin{array}{cc}e&f\\g&h\end{array}\right|\left|\begin{array}{cc}s&t\\u&v\end{array}\right|\\=+\left|\begin{array}{cc}a&b\\c&d\end{array}\right|\left|\begin{array}{cc}w&x\\y&z\end{array}\right|+\left|\begin{array}{cc}c&d\\e&f\end{array}\right|\left|\begin{array}{cc}s&t\\y&z\end{array}\right|+\left|\begin{array}{cc}e&f\\g&h\end{array}\right|\left|\begin{array}{cc}s&t\\u&v\end{array}\right|-\left|\begin{array}{cc}a&b\\e&f\end{array}\right|\left|\begin{array}{cc}u&v\\y&z\end{array}\right|-\left|\begin{array}{cc}c&d\\g&h\end{array}\right|\left|\begin{array}{cc}s&t\\w&x\end{array}\right|+\left|\begin{array}{cc}a&b\\g&h\end{array}\right|\left|\begin{array}{cc}u&v\\w&x\end{array}\right|

 実行結果は、倍速ぐらいになりました。再帰から早く抜けるほど速くなりそうですが、展開結果を書くのが大変そうなのでこの辺にしておきます。(^_^;
 ちなみに、RingsInTheSquare.py等で、detLap()を使いたい場合は、次のようにして使うとよいです。

from DetLaplace2 import detLap

● DetLaplace2.py

# coding: UTF-8
# DetLaplace2.py

def det4(m):
    r = (m[0][0]*m[1][1]-m[0][1]*m[1][0])*(m[2][2]*m[3][3]-m[2][3]*m[3][2])
    r-= (m[0][0]*m[1][2]-m[0][2]*m[1][0])*(m[2][1]*m[3][3]-m[2][3]*m[3][1])
    r+= (m[0][0]*m[1][3]-m[0][3]*m[1][0])*(m[2][1]*m[3][2]-m[2][2]*m[3][1])
    r+= (m[0][1]*m[1][2]-m[0][2]*m[1][1])*(m[2][0]*m[3][3]-m[2][3]*m[3][0])
    r-= (m[0][1]*m[1][3]-m[0][3]*m[1][1])*(m[2][0]*m[3][2]-m[2][2]*m[3][0])
    r+= (m[0][2]*m[1][3]-m[0][3]*m[1][2])*(m[2][0]*m[3][1]-m[2][1]*m[3][0])
    return r

def detLap(m):
    n = len(m)
    if n!=len(m[0]): return float('nan')
    if n==1: return m[0][0]
    if n==4: return det4(m)
    # 1次と4次のときは展開結果を利用、それ以外はLaplace展開
    sum = 0
    for i in range(n):
        mnr = [[m[k][j] for j in range(n) if j!=i] for k in range(1,n)]
        sum+= (-1)**i*m[0][i]*detLap(mnr)   # Laplace展開
    return sum

def main():
    import numpy as np
    import random
    from time import time
    tm=time()  # Timer Start
    for n in range(1,10+1):
        m = [[random.randint(0,99) for x in range(n)] for y in range(n)]
        a = np.array(m)
        print("%2d : %s | %s"%(n,detLap(m),np.linalg.det(a)))

    tm=time()-tm  # Timer Stop
    print("Runtime : %.3f [sec]"%tm)

if __name__ == '__main__':
    main()

●実行結果

 1 : 12 | 12.0
 2 : 2446 | 2446.0
 3 : -21276 | -21276.0
 4 : -74043 | -74043.0
 5 : 722025219 | 722025219.0
 6 : 2459179794 | 2459179794.0
 7 : -9473195486322 | -9.47319548632e+12
 8 : 284054910623389 | 2.84054910623e+14
 9 : -110893831223072229 | -1.10893831223e+17
10 : -570726085044361050 | -5.70726085044e+17
Runtime : 3.522 [sec]

※参考URL
行列式の値を計算するプログラムをPythonで作ってみた。
行列式の展開式を文字列で求めるプログラムをPythonで作ってみた。
行列式の値を計算するプログラムをPythonで作ってみた。(3)

Pythonスタートブック

Pythonスタートブック