前回のプログラムをちょっと改良してみました。(^_^;
1次と4次のときは展開結果を利用して、それ以外はLaplace展開を用いて計算するようにしました。(^_^;
ここで、4次の行列式の展開式は次式を用いました。(^_^;
実行結果は、倍速ぐらいになりました。再帰から早く抜けるほど速くなりそうですが、展開結果を書くのが大変そうなのでこの辺にしておきます。(^_^;
ちなみに、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)
- 作者: 辻真吾
- 出版社/メーカー: 技術評論社
- 発売日: 2010/04/24
- メディア: 大型本
- 購入: 19人 クリック: 199回
- この商品を含むブログ (59件) を見る