FNLMath : Eigensystem Analysis : Eigenvalues and Eigenvectors Computed with ARPACK
Eigenvalues and Eigenvectors Computed with ARPACK
First see Using ARPACK for Ordinary and Generalized Eigenvalue Problems in the Usage Notes section of this chapter. We describe here the Fortran 2003 usage of four basic problem types. There must be compiler support for the object-oriented features of Fortran 2003 to use these routines.
The generalized eigenvalue problem Ax = λBx requires that some eigenvalues and eigenvectors be computed. The organization of the user-written function for matrix-vector products depends on the part of the eigenvalue spectrum that is desired.
For an ordinary problem with A symmetric and B = I, the eigenvalues of largest or smallest magnitude can be computed by providing the operator products ω = Ax. Here x is an input vector and ω is the result of applying the linear operator A to x. This process is repeated several times within the Arnoldi algorithm, and the net result is a few eigenvalues of A and the corresponding eigenvectors.
For a generalized problem, it is useful and efficient to consider a shift value σ and the ordinary eigenvalue problem Cx  (A - σB)-1 Bx = vx. The matrix pencil A - σB is non-singular. The purpose of the user-written function is to provide results for the individual operator products ω = Bx, ω =  (A - σB)-1x, and ω = Ax. Usually the inverse matrix product will be computed by solving linear systems, where the matrix pencil is the coefficient matrix. The desired eigenvalues of this ordinary problem satisfy λj = σ + vj-1.
In the special case that B is positive definite, well-conditioned, and symmetric, one may compute the Cholesky decomposition B =RTR and then solve the ordinary eigenvalue problem Cy  R-T AR-1y = λy. The product operation required by the Arnoldi algorithm, ω = Cx, is performed in steps: Solve Rz = x for z, compute y = Az, and solve RTx = y for ω. The eigenvectors, Y, of C are transformed to those of the generalized problem, X, by solving RX = Y for X.
The operations required by ARPACK codes are returned as array functions. An array of input values, x, will yield an output array, y. These functions are written by the user. They must be written according to an abstract interface, given below. There are two user functions, double precision real and complex, that we support for the eigenvalue problem, and a third for the singular value decomposition, using double precision real data only. This interface, the named or enumerated constants that describe what is needed, and the eigenvalue codes are in the module ARPACK_INT. We use the notation: DKIND=kind(1.D0) to specify two double precision data types: REAL(DKIND) and COMPLEX(DKIND). The interface SVDMV(...) is for the singular value decomposition products only. For that problem the components EXTYPE%MROWS and EXTYPE%NCOLS are switched between the operator sizes M and N to account for computing y = AMxNx or y = ATx.
The Abstract Interfaces for User-Written Array Functions
Abstract Interface

IMPORT DKIND, ARPACKBASE
REAL(DKIND), INTENT(INOUT) :: X(:)
CLASS (ARPACKBASE), INTENT(INOUT) :: EXTYPE
REAL(DKIND) Y(SIZE(X))
END FUNCTION
FUNCTION ZMV (X, TASK, EXTYPE) RESULT(Y)
IMPORT DKIND, ARPACKBASE
CLASS (ARPACKBASE), INTENT(INOUT) :: EXTYPE
COMPLEX (DKIND), INTENT(INOUT) :: X(:)
COMPLEX (DKIND) Y(SIZE(X))
END FUNCTION
FUNCTION SVDMV (X, TASK, EXTYPE) RESULT(Y)
IMPORT DKIND, ARPACKBASE
CLASS (ARPACKBASE), INTENT(INOUT) :: EXTYPE
REAL (DKIND), INTENT(INOUT) :: X(EXTYPE%NCOLS)
REAL (DKIND) Y(EXTYPE%MROWS)
END FUNCTION

End Interface