ファイルから次のリンクAで取得できる次の入力 numpy 2d-array が与えられた場合、 scipy.sparse.linalg.eigsのような反復ソルバーを使用して固有値のサブセットのみを計算できれば素晴らしいでしょう。hill_mat.npy
まず、コンテキストについて少し説明します。この行列は、 double サイズの同等の固有値問題で線形化されたAサイズの 2 次固有値問題から得られます。次の構造を持っています (青色はゼロです):N2*NA
plt.imshow(np.where(A > 1e-15,1.,0), interpolation='None')
と次の機能:
A shape = (748, 748)
A dtype = float64
A sparsity ratio = 77.64841716949297 %
の実際の寸法はA、この再現可能な小さな例よりもはるかに大きいです。この場合、実際の希薄率と形状は近いと予想され95%ます(5508, 5508)。
結果として得られる の固有値Aは複素数(複素共役ペアになります) であり、モジュラスの虚数部が最小のものに関心があります。
問題: 直接ソルバーを使用する場合:
w_dense = np.linalg.eigvals(A)
idx = np.argsort(abs(w_dense.imag))
w_dense = w_dense[idx]
計算時間は急速に法外なものになります。したがって、スパースアルゴリズムの使用を検討しています:
from scipy.sparse import csc_matrix, linalg as sla
w_sparse = sla.eigs(A, k=100, sigma=0+0j, which='SI', return_eigenvectors=False)
しかし、この方法では ARPACK は固有値を見つけられないようです。scipy/arpack チュートリアルから、 のような小さな固有値を探す場合、kwargwhich = 'SI'を指定して、いわゆるシフト反転モードを使用する必要があります。つまり、アルゴリズムがこれらの固有値を見つけることができる場所を知ることができるようにするためです。それにもかかわらず、私の試みはすべて結果をもたらしませんでした...sigma
この機能の経験が豊富な人が、これを機能させるために手を貸してくれませんか?
コード スニペット全体を以下に示します。
import numpy as np
from matplotlib import pyplot as plt
from scipy.sparse import csc_matrix, linalg as sla
A = np.load('hill_mat.npy')
print('A shape =', A.shape)
print('A dtype =', A.dtype)
print('A sparsity ratio =',(np.product(A.shape) - np.count_nonzero(A)) / np.product(A.shape) *100, '%')
# quick look at the structure of A
plt.imshow(np.where(A > 1e-15,1.,0), interpolation='None')
# direct
w_dense = np.linalg.eigvals(A)
idx = np.argsort(abs(w_dense.imag))
w_dense = w_dense[idx]
# sparse
w_sparse = sla.eigs(csc_matrix(A), k=100, sigma=0+0j, which='SI', return_eigenvectors=False)