Developer Reference for Intel® oneAPI Math Kernel Library - C

Developer Reference

Contents

Extended Eigensolver RCI Interface Description

The Extended Eigensolver RCI interfaces can be used to solve standard or generalized eigenvalue problems, and are independent of the format of the matrices. As mentioned earlier, the Extended Eigensolver algorithm is based on the contour integration techniques of the matrix resolvent
G
(
σ
 )= (
σ
B
 -
A
)
-1
 over a circle. For solving a generalized eigenvalue problem, Extended Eigensolver has to perform one or more of the following operations at each contour point denoted below by
Z
e
 :
  • Factorize the matrix
    (
    Z
    e
     *
    B
     -
    A
    )
  • Solve the linear system
    (
    Z
    e
     *
    B
     -
    A
    )
    X
     =
    Y
     or
    (
    Z
    e
     *
    B
     -
    A
    )
    H
    X
     =
    Y
     with multiple right hand sides, where H means transpose conjugate
  • Matrix-matrix multiply
    B
    X
     =
    Y
     or
    A
    X
     =
    Y
For solving a standard eigenvalue problem, replace the matrix
B
 with the identity matrix
I
.
The primary aim of RCI interfaces is to isolate these operations: the linear system solver, factorization of the matrix resolvent at each contour point, and matrix-matrix multiplication. This gives universality to RCI interfaces as they are independent of data structures and the specific implementation of the operations like matrix-vector multiplication or inner system solvers. However, this approach requires some additional effort when calling the interface. In particular, operations listed above are performed by routines that you supply on data structures that you find most appropriate for the problem at hand.
To initialize an Extended Eigensolver RCI routine, set the job indicator (
ijob
) parameter to the value -1. When the routine requires the results of an operation, it generates a special value of
ijob
 to indicate the operation that needs to be performed. The routine also returns
ze
, the coordinate along the complex contour, the values of array
work
 or
workc
, and the number of columns to be used. Your subroutine then must perform the operation at the given contour point
ze
, store the results in prescribed array, and return control to the Extended Eigensolver RCI routine.
The following pseudocode shows the general scheme for using the Extended Eigensolver RCI functionality for a real symmetric problem: 
    ijob=-1; // initialization
    do while (ijob!=0) {
       ?feast_srci(&ijob, &N, &Ze, work, workc, Aq,  Bq,
                   fpm, &epsout, &loop, &Emin, &Emax, &M0, E, lambda, &q, res, &info);

       switch(ijob) {
       case 10:  // Factorize the complex matrix (ZeB-A)
                 break;
     
       case 11: // Solve the complex linear system (ZeB-A)x=workc
                // Put result in workc
                 break;
     
       case 30: // Perform multiplication A by Q
i
..Q
j
 columns of Q
NxM0

                // where i = fpm[23] and j = fpm[23]+fpm[24]−1
                // Q
i
..Q
j
 located in q starting from q+N*(i-1)
                 break;
     
       case 40: // Perform multiplication B by Q
i
..Q
j
 columns of Q
NxM0

                // where i = fpm[23] and j = fpm[23]+fpm[24]−1
                // Q
i
..Q
j
 located in q starting from q+N*(i-1)
                // Result is stored in work+N*(i-1)

                 break;
       }
     }
The
?
 option in
?feast
 in the pseudocode given above should be replaced by either
s
 or
d
, depending on the matrix data type of the eigenvalue system.
The next pseudocode shows the general scheme for using the Extended Eigensolver RCI functionality for a complex Hermitian problem: 
    ijob=-1; // initialization
    while (ijob!=0) {
       ?feast_hrci(&ijob, &N, &Ze, work, workc, Aq,  Bq,
                   fpm, &epsout, &loop, &Emin, &Emax, &M0, E, lambda, &q, res, &info);

       switch (ijob) {
       case 10: // Factorize the complex matrix (ZeB-A)
       break;

       case 11: // Solve the linear system (ZeB−A)y=workc
                // Put result in workc
       break;

       case 20: // Factorize (if needed by case 21) the complex matrix     (ZeB−A)ˆH
	               // ATTENTION: This option requires additional memory storage
	               // (i.e . the resulting matrix from case 10 cannot be overwritten)
       break;

       case 21: // Solve the linear system (ZeB−A)ˆHy=workc
                // Put result in  workc
	               // REMARK: case 20 becomes obsolete if this solve can be performed
               	// using the factorization in case 10
       break;

       case 30: // Multiply A by Q
i
..Q
j
 columns of Q
NxM0
,
                // where i = fpm[23] and j = fpm[23]+fpm[24]−1
                // Q
i
..Q
j
 located in q starting from q+N*(i-1)
                // Result is stored in work+N*(i-1) 
       break;

       case 40: // Perform multiplication B by Q
i
..Q
j
 columns of Q
NxM0

                // where i = fpm[23] and j = fpm[23]+fpm[24]−1
                // Q
i
..Q
j
 located in q starting from q+N*(i-1)
                // Result is stored in work+N*(i-1)
       break;

       }
    }
end do
The
?
 option in
?feast
 in the pseudocode given above should be replaced by either
c
 or
z
, depending on the matrix data type of the eigenvalue system.
If
case 20
 can be avoided, performance could be up to twice as fast, and Extended Eigensolver functionality would use half of the memory.
If an iterative solver is used along with a preconditioner, the factorization of the preconditioner could be performed with
ijob
 = 10
 (and
ijob
 = 20
 if applicable) for a given value of
Z
e
, and the associated iterative solve would then be performed with
ijob
 = 11
 (and
ijob
 = 21
 if applicable).
Product and Performance Information
Performance varies by use, configuration and other factors. Learn more at www.Intel.com/PerformanceIndex.
Notice revision #20201201

Product and Performance Information

1

Performance varies by use, configuration and other factors. Learn more at www.Intel.com/PerformanceIndex.