ARPACK_COMPLEX
Compute some eigenvalues and eigenvectors of the generalized eigenvalue problem Ax = λBx. This can be used for the case B = I. The values for A,B are real or complex. When the values are complex the eigenvalues may be complex and are not expected to occur in complex conjugate pairs.
Required Arguments
N — The dimension of the problem. (Input)
F — User-supplied
FUNCTION to return matrix-vector operations or linear solutions. This user function is written corresponding to the
abstract interface for the function
ZMV(…). The usage is
F (
X,
TASK,
EXTYPE), where
Function Return Value
F — An array of length N containing matrix-vector operations or linear equations solutions. Operations provided as code in the function F will be made depending upon the value of argument TASK.
Required Arguments
X — An array of length N containing the vector to which the operator will be applied. (Input)
TASK — An enumerated type which specifies the operation to be performed. (Input)
TASK is an enumerated integer value, use-associated from the module ARPACK_INT. It will be one of the following:
Value | Description |
|---|
ARPACK_Prepare | Take initial steps to prepare for the operations that follow. These steps can include defining the data for the matrices, factorizations for upcoming linear system solves, or recording the vectors used in the operations. |
ARPACK_A_x | y = Ax |
ARPACK_B_x | y = Bx |
EXTYPE — A derived type of the extensible class ARPACKBASE, which may be used to pass additional information to/from the user-supplied function. (Input/Output)
The user must include a USE ARPACK_INT statement in the calling program to define this derived type. If EXTYPE is not included as an argument to ARPACK_COMPLEX it should be ignored in the user-function, F.
The function F must be written according to the abstract interface for ZMV. If F is not use-associated nor contained in the calling program, declare it with PROCEDURE(ZMV) F.
ZVALUES — A complex array of eigenvalues. (Output)
The value NEV=size(ZVALUES) defines the number of eigenvalues to be computed. The calling program declares or allocates the array ZVALUES(1:NEV). The number of eigenvalues computed accurately is optionally available as the component EXTYPE%NACC of the base class EXTYPE.
Optional Arguments
PLACE — Defines the output content of VALUES. (Input)
PLACE specifies the placement within the spectrum for the required eigenvalues. PLACE can be one of the following enumerated integers as defined in ARPACK_INT:
Value |
|---|
ARPACK_Largest_Magnitude |
ARPACK_Smallest_Magnitude |
ARPACK_Largest_Real_Parts |
ARPACK_Smallest_Real_Parts |
ARPACK_Largest_Imag_Parts |
ARPACK_Smallest_Imag_Parts |
Default: PLACE = ARPACK_Largest_Magnitude.
TYPE — Defines the eigenvalue problem as either a standard or generalized eigenvalue problem. (Input)
TYPE can be one of the following enumerated integers as defined in ARPACK_INT:
Value | Description |
|---|
ARPACK_Standard | Ax = λx |
ARPACK_Generalized | Ax = λBx |
Default: TYPE = ARPACK_Standard.
CATEGORY — CATEGORY and TYPE define the operation sequence provided in the user-written function. (Input)
CATEGORY can be one of the following enumerated integers as defined in ARPACK_INT:
Value | Description |
|---|
ARPACK_Regular | y = Ax |
ARPACK_Regular_Inverse | y = Ax, y = Bx, y = B-1x |
ARPACK_Shift_Invert | y = Ax, y = (A - σ I)-1x, y = (A -σ B)-1x |
ARPACK_Complex_Part_Shift_Invert | y = Ax, y = Bx, y = Re{(A - σ B)-1x}
y = Ax, y = Bx, y = Im{(A - σ B)-1x} |
Default: CATEGORY = ARPACK_Regular.
EXTYPE — A derived type of the extensible class ARPACKBASE, which may be used to pass additional information to/from the user-supplied function. (Input/Output)
The user must include a USE ARPACK_INT statement in the calling program to define this derived type. If EXTYPE is not included as an argument to ARPACK_COMPLEX it must still be supplied as an argument to user-function, F, but is not used.
VECTORS — An allocatable array of approximate eigenvectors. (Output)
It is not necessary to allocate
VECTORS(
:,:). If this argument is used the allocation occurs within the routine
ARPACK_NONSYMMETRIC. The output sizes are
VECTORS(
1:N,1:NCV). The second dimension value is
NCV=min(N, max(
FACTOR_MAXNCV*NEV,NEV+1)), where the value
FACTOR_MAXNCV is a component of the base class,
ARPACKBASE. The first
NEV columns of
VECTORS(
:,:) represent the eigenvectors
wj (see
Comments).
FORTRAN 2003 Interface
Generic: ARPACK_COMPLEX (N, F,ZVALUES [,…])
Specific: The specific interface name is Z_ARPACK_COMPLEX.
Fortran 90 Interface
A Fortran 90 compiler version is not available for this routine.
Description
Routine ARPACK_COMPLEX calls ARPACK subroutines to compute partial eigenvalue-eigenvector decompositions for complex matrices. The ARPACK routines are dzaupd and dzeupd (see ARPACK Users’ Guide, SIAM Publications, (1998)), which use “reverse communication” to obtain the required matrix-vector operations for this approximation. Responses to these requests are made by calling the user-written function F. By including the class object EXTYPE as an argument to this function, user data and procedure pointers are available for the evaluations. A user code must extend the base class EXTYPE to include the extra data and procedure pointers.
Comments
For purposes of checking results the complex residual rj = Awj - λjwj should be small in norm relative to the norm of A. For the generalized eigenvalue problem Ax = λBx the eigenvalues are optionally computed based on the Raleigh Quotient. Because of the shifts used, only the eigenvectors may be computed. The eigenvalues are returned based on solving Awj = λBwj for λj,
where
The eigenvalue λj is finite and valid if the denominator is non-zero. The Raleigh Quotient for eigenvalues of generalized problems is used only when vectors are requested and the user has requested it be used with the base class component EXTYPE%RALEIGH_QUOTIENT = .TRUE. This is the component’s default value.
Example
This example is a variation of the first
example for
ARPACK_SYMMETRIC. We approximate eigenvalues and eigenfunctions of the Laplacian operator
on the unit square, [0,1] × [0,1]. But now the parameter ρ is complex. Thus the eigenvalues and eigenfunctions are complex.
(Example arpack_complex_ex1.f90)
MODULE ARPACK_COMPLEX_EX1_INT
USE ARPACK_INT
IMPLICIT NONE
TYPE, EXTENDS(ARPACKBASE), PUBLIC :: ARPACKBASE_EXT
REAL(DKIND) :: H =0._DKIND
REAL(DKIND) :: HSQ=0._DKIND
COMPLEX(DKIND) :: RHO=(0._DKIND,0._DKIND)
COMPLEX(DKIND) :: DL=(0._DKIND,0._DKIND)
COMPLEX(DKIND) :: DD=(0._DKIND,0._DKIND)
COMPLEX(DKIND) :: DU=(0._DKIND,0._DKIND)
INTEGER :: NX=0
INTEGER :: NV=0
END TYPE ARPACKBASE_EXT
CONTAINS
FUNCTION FZ1(X, TASK, EXTYPE) RESULT(Y)
USE UMACH_INT
CLASS (ARPACKBASE), INTENT(INOUT) :: EXTYPE
COMPLEX (DKIND), INTENT(INOUT) :: X(:)
INTEGER, INTENT(IN) :: TASK
COMPLEX (DKIND) Y(size(X))
COMPLEX (DKIND) DT(3)
REAL(DKIND) :: ONE=1._DKIND
INTEGER J, NOUT
CALL UMACH(2, NOUT)
SELECT TYPE(EXTYPE)
TYPE IS(ARPACKBASE_EXT)
ASSOCIATE(NX => EXTYPE % NX,&
H => EXTYPE % H,&
HSQ => EXTYPE % HSQ,&
RHO => EXTYPE % RHO,&
DL => EXTYPE % DL,&
DD => EXTYPE % DD,&
DU => EXTYPE % DU,&
NV => EXTYPE % NV)
SELECT CASE(TASK)
CASE(ARPACK_A_x)
! Computes y <-- A*x, where A is the N**2 by N**2 block
! tridiagonal matrix deriving from (Laplacian u) + rho*(du/dx).
DT=(/DL,DD,DU/)
Y(1:NX)=T(NX,X(1:NX),DT) - X(NX+1:2*NX)/HSQ
DO J=NX+1,NX**2-NX,NX
Y(J:J+NX-1)=T(NX,X(J:J+NX-1),DT) &
- (X(J-NX:J-1)+ X(J+NX:J+2*NX-1))/HSQ
END DO
Y((NX-1)*NX+1:NX**2)= - X((NX-1)*NX-NX+1:(NX-1)*NX) &
/ HSQ + T(NX,X((NX-1)*NX+1:NX**2),DT)
! Total the number of matrix-vector products.
NV=NV+1
CASE(ARPACK_Prepare)
! Define 1/H**2, etc. so they are available in the evaluator.
H = ONE/REAL(NX+1,DKIND)
HSQ = H**2
DD = (4.0D+0, 0.0D+0) / HSQ
DL = -ONE/HSQ - (5.0D-1, 0.0D+0) *RHO/H
DU = -ONE/HSQ + (5.0D-1, 0.0D+0) *RHO/H
NV = 0
CASE DEFAULT
WRITE(nout,*) TASK, ': INVALID TASK REQUESTED'
STOP 'IMSL_ERROR_WRONG_OPERATION'
END SELECT
END ASSOCIATE
END SELECT
END FUNCTION
FUNCTION T(NX, X, DT)RESULT(V)
INTEGER, INTENT(IN) :: NX
COMPLEX(DKIND), INTENT(IN) :: X(:), DT(3)
COMPLEX(DKIND) :: V(NX)
INTEGER J
ASSOCIATE(DL => DT(1),&
DD => DT(2),&
DU => DT(3))
V(1) = DD*X(1) + DU*X(2)
DO J = 2,NX-1
V(J) = DL*X(J-1) + DD*X(J) + DU*X(J+1)
END DO
V(NX) = DL*X(NX-1) + DD*X(NX)
END ASSOCIATE
END FUNCTION
END MODULE
! Compute the largest magnitude eigenvalues of a discrete Laplacian,
! based on second order divided differences.
! The matrix used is obtained from the standard central difference
! discretization of the convection-diffusion operator
! (Laplacian u) + rho*(du / dx)
! on the unit squre 0,1x0,1 with zero Dirichlet boundary
! conditions.
USE ARPACK_COMPLEX_EX1_INT
USE UMACH_INT
USE WRCRN_INT
INTEGER, PARAMETER :: NEV=6
INTEGER :: J, N, NOUT
COMPLEX(DKIND) :: VALUES(NEV)
COMPLEX(DKIND), ALLOCATABLE :: RES(:), EF(:,:)
COMPLEX(DKIND), ALLOCATABLE :: VECTORS(:,:)
REAL(DKIND) NORM
LOGICAL SMALL, SOLVED
TYPE(ARPACKBASE_EXT) EX
ASSOCIATE(NX => EX % NX, &
NV => EX % NV, &
RHO => EX % RHO,&
NACC => EX % NACC)
CALL UMACH(2, NOUT)
NX=10
RHO=(100._DKIND,1._DKIND)
! Define size of matrix problem.
N=NX**2
! Note that VECTORS(:,:) does not need to be allocated
! in the calling program. That happens within the
! routine ARPACK_COMPLEX(). It is OK to do this but
! the sizes (N,NCV) are determined in ARPACK_COMPLEX.
CALL ARPACK_COMPLEX(N, FZ1, VALUES, EXTYPE=EX, VECTORS=VECTORS)
WRITE(NOUT, *) 'Number of eigenvalues requested, and accurate'
WRITE(NOUT, *) '---------------------------------------------'
WRITE(NOUT, '(5X, I4, 5X, I4)') NEV, NACC
WRITE(NOUT, *) 'Number of Matrix-Vector Products Required, ZEX-1'
WRITE(NOUT, *) '------------------------------------------------'
WRITE(NOUT, '(5X, I4)') NV
CALL WRCRN ('Largest Magnitude Operator Eigenvalues', VALUES)
! Check residuals, A*vectors = values*vectors:
ALLOCATE(RES(N))
DO J=1,NACC
RES=FZ1(VECTORS(:,J),ARPACK_A_x,EX)-VALUES(J)*VECTORS(:,J)
NORM=maxval(abs(RES))
SMALL=(NORM <= ABS(VALUES(J))*SQRT(EPSILON(NORM)))
IF(J==1) SOLVED=SMALL
SOLVED=SOLVED .and. SMALL
END DO
IF(SOLVED) THEN
WRITE(NOUT,'(A///)') &
'All Ritz Values and Vectors have small residuals.'
ELSE
WRITE(NOUT,'(A///)') &
'Some Ritz Values and Vectors have large residuals.'
END IF
END ASSOCIATE
END
Output
Number of eigenvalues requested, and accurate
---------------------------------------------
6 6
Number of Matrix-Vector Products Required, ZEX-1
------------------------------------------------
475
Largest Magnitude Operator Eigenvalues
1 ( 727.0,-1029.6)
2 ( 705.4, 1029.6)
3 ( 698.4,-1029.6)
4 ( 676.8, 1029.6)
5 ( 653.3,-1029.6)
6 ( 631.7, 1029.6)
All Ritz Values and Vectors have small residuals.
