matlab lanczos算法用来计算大型稀疏矩阵的最大最小本征值及相应的本征矢量

function varargout = irbleigs(varargin) % IRBLEIGS: Finds a few eigenvalues and eigenvectors of a Hermitian matrix. % % IRBLEIGS will find a few eigenvalues and eigenvectors for either the % standard eigenvalue problem A*x = lambda*x or the generalized eigenvalue % problem A*x = lambda*M*x, where A is a sparse Hermitian matrix and the % matrix M is positive definite. % % [V,D,PRGINF] = IRBLEIGS(A,OPTIONS) % [V,D,PRGINF] = IRBLEIGS(A,M,OPTIONS) % [V,D,PRGINF] = IRBLEIGS('Afunc',N,OPTIONS) % [V,D,PRGINF] = IRBLEIGS('Afunc',N,M,OPTIONS) % % The first input argument into IRBLEIGS can be a numeric matrix A or an M-file % ('Afunc') that computes the product A*X, where X is a (N x blsz) matrix. If A is % passed as an M-file then the (N x blsz) matrix X is the first input argument, the % second input argument is N, (the size of the matrix A), and the third input argument % is blsz, i.e. Afunc(X,N,blsz). For the generalized eigenvalue problem the matrix M is % positive definite and passed only as a numeric matrix. IRBLEIGS will compute the % Cholesky factorization of the matrix M by calling the internal MATLAB function CHOL. % M may be a dense or sparse matrix. In all the implementations IRBLEIGS(A,...) can % be replaced with IRBLEIGS('Afunc',N,...). % % NOTE: If the M-file 'Afunc' requires additional input parameters, e.g. a data structure, % use the option FUNPAR to pass any additional parameters to your function (X,N,blsz,FUNPAR). % % OUTPUT OPTIONS: % --------------- % % I.) IRBLEIGS(A) or IRBLEIGS(A,M) % Displays the desired eigenvalues. % % II.) D = IRBLEIGS(A) or D = IRBLEIGS(A,M) % Returns the desired eigenvalues in the vector D. % % III.) [V,D] = IRBLEIGS(A) or [V,D] = IRBLEIGS(A,M) % D is a diagonal matrix that contains the desired eigenvalues along the % diagonal and the matrix V contains the corresponding eigenvectors, such % that A*V = V*D or A*V = M*V*D. If IRBLEIGS reaches the maximum number of % iterations before convergence then V = [] and D = []. Use output option % IV.) to get approximations to the Ritz pairs. % % IV.) [V,D,PRGINF] = IRBLEIGS(A) or [V,D,PRGINF] = IRBLEIGS(A,M) % This option returns the same as (III) plus a two dimensional array PRGINF % that reports if the algorithm converges and the number of matrix vector % products. If PRGINF(1) = 0 then this implies normal return: all eigenvalues have % converged. If PRGINF(1) = 1 then the maximum number of iterations have been % reached before all desired eigenvalues have converged. PRGINF(2) contains the % number of matrix vector products used by the code. If the maximum number of % iterations are reached (PRGINF(1) = 1), then the matrices V and D contain any % eigenpairs that have converged plus the last Ritz pair approximation for the % eigenpairs that have not converged. % % INPUT OPTIONS: % -------------- % % ... = IRBLEIGS(A,OPTS) or ... = IRBLEIGS(A,M,OPTS) % OPTS is a structure containing input parameters. The input parameters can % be given in any order. The structure OPTS may contain some or all of the % following input parameters. The string for the input parameters can contain % upper or lower case characters. % % INPUT PARAMETER DESCRIPTION % % OPTS.BLSZ Block size of the Lanczos tridiagonal matrix. % DEFAULT VALUE BLSZ = 3 % % OPTS.CHOLM Indicates if the Cholesky factorization of the matrix M is available. If % the Cholesky factorization matrix R is available then set CHOLM = 1 and % replace the input matrix M with R where M = R'*R. % DEFAULT VALUE CHOLM = 0 % % OPTS.PERMM Permutation vector for the Cholesky factorization of M(PERMM,PERMM). % When the input matrix M is replaced with R where M(PERMM,PERMM)=R'*R % then the vector PERMM is the permutation vector. % DEFAULT VALUE PERMM=1:N % % OPTS.DISPR Indicates if K Ritz values and residuals are to be displayed on each % iteration. Set positive to display the Ritz values and residuals on % each iteration. % DEFAULT VALUE DISPR = 0 % % OPTS.EIGVEC A matrix of converged eigenvectors. % DEFAULT VALUE EIGVEC = [] % % OPTS.ENDPT Three letter string specifying the location of the interior end- % points for the dampening interval(s). % 'FLT' - Let the interior endpoints float. % 'MON' - Interior endpoints are chosen so that the size of the % dampening interval is increasing. This creates a nested % sequence of intervals. The interior endpoint will approach % the closest Ritz value in the undampened part of the spectrum % to the dampening interval. % DEFAULT VALUE ENDPT = 'MON' (If SIGMA = 'LE' or 'SE'.) % DEFAULT VALUE ENDPT = 'FLT' (If SIGMA = a numeric value NVAL.) % % OPTS.FUNPAR If A is passed as a M-file then FUNPAR contains any additional parameters % that the M-file requires in order to compute the matrix vector product. % FUNPAR can be passed as any type, numeric, character, data structure, etc. The % M-file must take the input parameters in the following order (X,n,blsz,FUNPAR). % DEFAULT VALUE FUNPAR =[] % % OPTS.K Number of desired eigenvalues. % DEFAULT VALUE K = 3 % % OPTS.MAXIT Maximum number of iterations, i.e. maximum number of block Lanczos restarts. % DEFAULT VALUE MAXIT = 100 % % OPTS.MAXDPOL Numeric value indicating the maximum degree of the dampening % polynomial allowed. % DEFAULT VALUE MAXDPOL = 200 (If SIGMA = 'LE' or 'SE'.) % DEFAULT VALUE MAXDPOL = N (If SIGMA = a numeric value NVAL.) % % OPTS.NBLS Number of blocks in the Lanczos tridiagonal matrix. The program may increase % NBLS to ensure certain requirements in [1] are satisfied. A warning message % will be displayed if NBLS increases. % DEFAULT VALUE NBLS = 3 % % OPTS.SIGMA Two letter string or numeric value specifying the location % of the desired eigenvalues. % 'SE' Smallest Real eigenvalues. % 'LE' Largest Real eigenvalues. % NVAL A numeric value. The program searches for the K closest % eigenvalues to the numeric value NVAL. % DEFAULT VALUE SIGMA = 'LE' % % OPTS.SIZINT Size of the dampening interval. Value of 1 indicates consecutive % Ritz values are used to determine the endpoints of the dampening % interval. Value of 2 indicates endpoints are chosen from Ritz % values that are seprated by a single Ritz value. A value of 3 % indicates endpoints are chosen from Ritz values that are seprated % by two Ritz values. Etc. The minimum value is 1 and the maximum % value is (