下邊代碼段是關(guān)于N lowest eigenvalues of the tridiagonal matrix in python的代碼袄简,希望能對(duì)碼農(nóng)們有較大用。
''' r = lamRange(d,c,N).
? ? Returns the sequence {r[0],r[1],...,r[N]} that
? ? separates the N lowest eigenvalues of the tridiagonal
? ? matrix [A] = [cdc]; that is, r[i] < lam[i] < r[i+1].
'''
from numpy import ones
def lamRange(d,c,N):
? ? lamMin,lamMax = gerschgorin(d,c)
? ? r = ones(N+1)
? ? r[0] = lamMin
? # Search for eigenvalues in descending order?
? ? for k in range(N,0,-1):
? ? ? # First bisection of interval(lamMin,lamMax)
? ? ? ? lam = (lamMax + lamMin)/2.0
? ? ? ? h = (lamMax - lamMin)/2.0
? ? ? ? for i in range(1000):
? ? ? ? ? # Find number of eigenvalues less than lam
? ? ? ? ? ? p = sturmSeq(d,c,lam)
? ? ? ? ? ? numLam = numLambdas(p)
? ? ? ? ? # Bisect again & find the half containing lam
? ? ? ? ? ? h = h/2.0
? ? ? ? ? ? if numLam < k: lam = lam + h
? ? ? ? ? ? elif numLam > k: lam = lam - h
? ? ? ? ? ? else: break
? ? ? # If eigenvalue located, change the upper limit
? ? ? # of search and record it in [r]
? ? ? ? lamMax = lam
? ? ? ? r[k] = lam
? ? return r