Inverse power method for solving the eigenvalue problem in Python,inverseeigenvalue,''' lam,x =


''' lam,x = inversePower(a,s,tol=1.0e-6).    Inverse power method for solving the eigenvalue problem    [a]{x} = lam{x}. Returns 'lam' closest to 's' and the    corresponding eigenvector {x}.'''from numpy import zeros,dot,identityfrom LUdecomp import *from math import sqrtfrom random import randomdef inversePower(a,s,tol=1.0e-6):    n = len(a)    aStar = a - identity(n)*s   # Form [a*] = [a] - s[I]    aStar = LUdecomp(aStar)     # Decompose [a*]    x = zeros(n)    for i in range(n):          # Seed [x] with random numbers        x[i] = random()    xMag = sqrt(dot(x,x))       # Normalize [x]    x =x/xMag    for i in range(50):         # Begin iterations              xOld = x.copy()         # Save current [x]        x = LUsolve(aStar,x)    # Solve [a*][x] = [xOld]        xMag = sqrt(dot(x,x))   # Normalize [x]        x = x/xMag        if dot(xOld,x) < 0.0:   # Detect change in sign of [x]            sign = -1.0            x = -x        else: sign = 1.0        if sqrt(dot(xOld - x,xOld - x)) < tol:            return s + sign/xMag,x    print 'Inverse power method did not converge'

评论关闭