Skip to content

The MRRR algorithm for multi-core and SMP systems

License

Notifications You must be signed in to change notification settings

SebastianAchilles/mr3smp

 
 

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 

History

45 Commits
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

Repository files navigation

These files provide the routine 'mrrr'. Its purpose is to compute all 
or a subset of eigenvalues and optionally eigenvectors of a symmetric 
tridiagonal matrix. 'mrrr' is based on the algorithm of Multiple 
Relatively Robust Representations (MRRR). The routine thereby targets 
multi-core processors and SMP systems made out of them. The 
implementation is based on LAPACK's routine 'dstemr'.


Building the archive libmrrr.a:
--------------------------------
1) In the folder INSTALL, depending on the compiler used, 
   replace the file 'make.inc' (which assumes GNU gcc) 
   with the corresponding file. 
   For example if you want to use Intel's compilers: 
   $cp ./INSTALL/make.inc.intel ./make.inc
   If you do not find a file matching your compiler, then use any of 
   the files as a basis to edit.
2) Edit the file 'make.inc' according to your system. Especially, set 
   the variable INCLAPACK to zero if you plan to link to LAPACK and 
   BLAS, instead of compiling the necessary routines and including 
   them in the archive. (LAPACK 3.3 or equivalent is recommended.)
3) Run make.

Known issues:
-------------
- PROBLEM: Spinlocks are not supported
  SOLUTION: Compile the library with SPINLOCK_SUPPORT=0 in 'make.inc'
- PROBLEM: The C99 feature of complex numbers is not supported
  SOLUTION: Use a newer compiler or compile the library with
  COMPLEX_SUPPORT=0 (default) in 'make.inc'. In the latter case the 
  dense Hermitian routines will not be available in the archive

Using libmrrr.a:
----------------
The folder EXAMPLES contains simple examples of how to use the routine 
'mrrr' in C and Fortran code. Edit the 'Makefile' in these folders if 
you do not use the GNU compilers and run 'make' to compile the 
examples.

In general, the code that calls 'mrrr' needs to be linked to the 
library 'libmrrr.a', which is created in the LIB folder. 

Below are given the C and Fortran prototypes of the function 'mrrr'.
For more information please see 'INCLUDE/mrrr.h'.


######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
# int mrrr(char *jobz, char *range, int *n, double  *D,              #
#          double *E, double *vl, double *vu, int *il, int *iu,      #
#          int *tryrac, int *m, double *W, double *Z, int *ldz,      #
#	   int *Zsupp);                                              #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# jobz              "N" or "n" - compute only eigenvalues            #
#                   "V" or "v" - compute also eigenvectors           #
#                   "C" or "c" - count the maximal number of         #
#                                locally computed eigenvectors       #
# range             "A" or "a" - all                                 #
#                   "V" or "v" - by interval: (VL,VU]                #
#                   "I" or "i" - by index:     IL-IU                 #
# n                 Matrix size                                      #
# ldz               Leading dimension of eigenvector matrix Z;       #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# D (double[n])     Diagonal elements of tridiagonal T.              #
#                   (On output the array will be overwritten).       #
# E (double[n])     Off-diagonal elements of tridiagonal T.          #
#                   First n-1 elements contain off-diagonals,        #
#                   the last element can have an abitrary value.     #
#                   (On output the array will be overwritten.)       #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index 'il' to 'iu'#
#                   are computed.                                    #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index 'il' to 'iu'# 
#                   are computed.                                    #
# tryrac            0 - do not try to achieve high relative accuracy.#
#                   1 - relative accuracy will be attempted;         #
#                       on output it is set to zero if high relative #
#                       accuracy is not achieved.                    #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed. #
#                   If jobz="C", 'm' will be set to the number of    #
#                   eigenvalues/eigenvectors that will be computed.  #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues.   #
# Z                 Eigenvectors                                     #
# (double[n*m])     Enough space must be provided to store the       #
#                   vectors. 'm' should be bigger or equal           #
#                   to 'n' for range="A" and 'iu-il+1' for range="I".#
#                   For range="V" the minimum 'm' can be queried     #
#                   using jobz="C".                                  #
# Zsupp             Support of eigenvectors, which is given by       #
# (double[2*n])     Z[2*i] to Z[2*i+1] for the i-th eigenvector      #
#                   stored locally (1-based indexing).               #
#                                                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# MRRR(JOBZ, RANGE, N, D, E, VL, VU, IL, IU, TRYRAC, M, W,           #
#      Z, LDZ, ZSUPP, INFO)                                          #
#                                                                    #
# CHARACTER        JOBZ, RANGE                                       #
# INTEGER          N, IL, IU, TRYRAC, M, LDZ, ZSUPP(*), INFO         #
# DOUBLE PRECISION D(*), D(*), VL, VU, W(*), Z(*,*)                  #
######################################################################

The number of threads created can be specified via the environment 
variable PMR_NUM_THREADS. In case this variable is not 
defined, then the routine will create as many threads as specified by 
the variable DEFAULT_NUM_THREADS (which is set in 'mrrr.h'). 




ADDITIONAL ROUTINES FOR DENSE STANDARD AND GENERALIZED EIGENPROBLEMS
-------------------------------------------------------------------- 
For convinience there are also routines for the standard and 
generalized symmetric and Hermitian eigenproblem added that make use 
of the multi-threaded MRRR. The functionality is similar to LAPACK's 
routines for this operations. The routines using packed storage should 
be avoided whenever possible, since their performance is much worse.

The routines are (using a LAPACK like notation):
----------------------------------------------------------------------
Standard eigenvector problem A*x = lambda*x:
Computing all or a subset of eigenvalues, and optionally 
eigenvectors
----------------------------------------------------------------------
dsyeig: A is symmetric
zheeig: A is Hermitian
dspeig: A is symmetric and stored in packed storage
zhpeig: A is Hermitian and stored in packed storage
----------------------------------------------------------------------
----------------------------------------------------------------------
Standard eigenvector problem A*x = lambda*B*x or A*B*x = lambda*x 
or B*A*x = lambda*x:
Computing all or a subset of eigenvalues, and optionally 
eigenvectors
----------------------------------------------------------------------
dsygeig: A and B are symmetric and B is definite
zhegeig: A and B are Hermitian and B is definite
dspgeig: A and B are symmetric, stored in packed storage, 
         and B is definite
zhpgeig: A and B are Hermitian, stored in packed storage,
         and B is definite
----------------------------------------------------------------------




The dense standard eigenvalue problem:
--------------------------------------
For the reduction and back-transformation the routines call the 
corresponding LAPACK routines. Therefore LAPACK is required to run 
the routines. Simply code snippets of how to call the functions are 
included in the examples folder.
The prototypes are similar to LAPACK's "dsyevx" and are as follows.

The symmetric case:

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
# int dsyeig(char *jobz, char *range, char *uplo, int *n, double *A, #
#            int *lda, double *vl, double *vu, int *il, int *iu,     #
#	     int *m, double *W, double *Z, int *ldz);                #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A stored                 #
#                   "U" - Lower triangle of A stored                 #
# n                 Order of the matrix A                            #
# lda               Leading dimension of matrix A;                   #
#                   often equal to matrix size n                     #
# ldz               Leading dimension of eigenvector matrix Z;       #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# A                 On entry symmetric input matrix, stored          #
# (double[lda*n])   in column major ordering. Depending on the       #
#                   value of 'uplo' only the upper or lower          #
#                   triangular part is referenced.                   #
#                   (On output the array will be overwritten).       #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed.                                    #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed.                                    #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed. #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues.   #
# Z                 Eigenvectors.                                    #
# (double[n*m])     Enough space must be provided to store the       #
#                   vectors. 'm' should be bigger or equal           #
#                   to 'n' for range="A" or "V" and 'iu-il+1' for    #
#                   range="I".                                       #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                 0 - Success                                        #
#                 1 - Wrong input parameter                          #
#                 2 - Misc errors                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# DSYEIG(JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, M, W,         #
#        Z, LDZ, INFO)                                               #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDA, LDZ, INFO                      #
# DOUBLE PRECISION A(*,*), VL, VU, W(*), Z(*,*)                      #
#                                                                    #
######################################################################



The symmetric case (packed storage):

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
#  int dspeig(char *jobz, char *range, char *uplo, int *n,           #
#             double *AP, double *vl, double *vu, int *il, int *iu,  #
#	      int *m, double *W, double *Z, int *ldz);               #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A stored                 #
#                   "U" - Lower triangle of A stored                 #
# n                 Order of the matrix A                            #
# ldz               Leading dimension of eigenvector matrix Z;       #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# AP (double[s])    On entry symmetric input matrix, stored          #
# s = (n*(n+1))/2   in column major ordering. Depending on the       #
#                   value of 'uplo' only the upper or lower          #
#                   triangular part is stored.                       #
#                   (On output the array will be overwritten).       #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed.                                    #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #    
#                   are computed.                                    #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed  #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues    #
# Z                 Eigenvectors                                     #
# (double[n*m])     Enough space must be provided to store the       #
#                   vectors. 'm' should be bigger or equal           #
#                   to 'n' for range="A" or "V" and 'iu-il+1' for    #
#                   range="I"                                        #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                   0 - Success                                      #
#                   1 - Wrong input parameter                        #
#                   2 - Misc errors                                  #
#                                                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# DSPEIG(JOBZ, RANGE, UPLO, N, AP, VL, VU, IL, IU, M, W,             #
#        Z, LDZ, INFO)                                               #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDZ, INFO                           #
# DOUBLE PRECISION AP(*), VL, VU, W(*), Z(*,*)                       #
#                                                                    #
######################################################################



The Hermitian case (only added to 'libmrrr.a' if specified in 
'make.inc'):

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
# int zheeig(char *jobz, char *range, char *uplo, int *n,            #
#	   double complex *A, int *lda, double *vl, double *vu,      #
#	   int *il, int *iu, int *m, double *W, double complex *Z,   #
#	   int *ldz);                                                #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A stored                 #
#                   "U" - Lower triangle of A stored                 #
# n                 Order of the matrix A                            #
# lda               Leading dimension of matrix A;                   #
#                   often equal to matrix size n                     #
# ldz               Leading dimension of eigenvector matrix Z;       #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# A                 On entry Hermitian input matrix, stored          #
# (double complex   in column major ordering. Depending on the       #
#  [lda*n])         value of 'uplo' only the upper or lower          #
#                   triangular part is referenced.                   #
#                   (On output the array will be overwritten).       #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed.                                    #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed.                                    #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed. #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues.   #
# Z                 Eigenvectors.                                    #
# (double complex   Enough space must be provided to store the       #
#  [n*m])           vectors. 'm' should be bigger or equal           #
#                   to 'n' for range="A" or "V" and 'iu-il+1' for    #
#                   range="I".                                       #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                 0 - Success                                        #
#                 1 - Wrong input parameter                          #
#                 2 - Misc errors                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# DSYEIG(JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, M, W,         #
#        Z, LDZ, INFO)                                               #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDA, LDZ, INFO                      #
# DOUBLE PRECISION VL, VU, W(*)                                      #
# COMPLEX*16       A(*,*), Z(*,*)                                    #
#                                                                    #
######################################################################



The Hermitian case (packed storage; only added to if specified in 
'make.inc'):

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
#  int zhpeig(char *jobz, char *range, char *uplo, int *n,           #
#	   double complex *AP, double *vl, double *vu,               #
#	   int *il, int *iu, int *m, double *W,                      #
#	   double complex *Z, int *ldz);                             #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A stored                 #
#                   "U" - Lower triangle of A stored                 #
# n                 Order of the matrix A                            #
# ldz               Leading dimension of eigenvector matrix Z;       #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# AP (double        On entry Hermitian input matrix, stored          #
#     complex[s])   in packed storage by columns. Depending on the   #
# s = (n*(n+1))/2   value of 'uplo' only the upper or lower          #
#                   triangular part is stored.                       #
#                   (On output the array will be overwritten).       #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined.                      #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues.                        #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    # 
#                   are computed.                                    #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'.                         #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed.                                    #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed. #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues.   #
# Z                 Eigenvectors                                     #
# (double           Enough space must be provided to store the       #
#  complex[n*m])    vectors. 'm' should be bigger or equal           #
#                   to 'n' for range="A" or "V" and 'iu-il+1' for    #
#                   range="I".                                       #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                   0 - Success                                      #
#                   1 - Wrong input parameter                        #
#                   2 - Misc errors                                  #
#                                                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# ZHPEIG(JOBZ, RANGE, UPLO, N, AP, VL, VU, IL, IU, M, W,             #
#             Z, LDZ, INFO)                                          #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDZ, INFO                           #
# DOUBLE PRECISION VL, VU, W(*)                                      #
# COMPLEX*16       AP(*), Z(*,*)                                     #
#                                                                    #
######################################################################




The dense generalized definite eigenvalue problem:
--------------------------------------------------
Using LAPACK notation the routines are called DSYGEIG and ZHEGEIG.
Simply code snippets of how to call the functions are included in the 
examples folder.
The problems that can be solved have the following form:
A*x = lambda*B*x or A*B*x = lambda*x or B*A*x = lambda*x,
where A and B are symmetric/Hermitian and B is definite.

Symmetric-definite case:

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
# int dsygeig(int *itype, char *jobz, char *range, char *uplo,       #
#             int *n, double *A, int *lda, double *B, int *ldb,      #
#	      double *vl, double *vu, int *il, int *iu, int *m,      #
#	      double *W);                                            #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# itype              1  - A*x = lambda*B*x                           #
#                    2  - A*B*x = lambda*x                           #
#                    3  - B*A*x = lambda*x                           #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A and B stored           #
#                   "U" - Lower triangle of A and B stored           #
# n                 Order of the matrix A and B                      #
# lda               Leading dimension of matrix A;                   #
#                   often equal to matrix size n                     #
# ldb               Leading dimension of eigenvector matrix B;       #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# A                 On entry symmetric input matrix, stored          #
# (double[lda*n])   in column major ordering. Depending on the       #
#                   value of 'uplo' only the upper or lower          #
#                   triangular part is referenced                    #
#                   On output the array will contain the             #
#                   'm' computed eigenvectors                        #
# B                 On entry symmetric definite input matrix,        #
# (double[ldb*n])   stored in column major ordering.                 #
#                   Depending on the value of 'uplo' only the upper  #
#                   or lower triangular part is referenced           #
#                   On output overwritten                            #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input      #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input      #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed  #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues    #
#                                                                    #
# NOTICE:           The routine will allocate work space of size     #
#                   double[n*n] for range="A" or "V" and double[m*n] #
#                   for range="I"                                    #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                 0 - Success                                        #
#                 1 - Wrong input parameter                          #
#                 2 - Misc errors                                    #
#                                                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# DSYGEIG(ITYPE, JOBZ, RANGE, UPLO, N, A, LDA, B, LDB,               #
#         VL, VU, IL, IU, M, W, INFO);                               #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDA, LDB, INFO                      #
# DOUBLE PRECISION A(*,*), B(*,*), VL, VU, W(*)                      #
#                                                                    #
######################################################################



Symmetric-definite case (packed storage):

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
# int dspgeig(int *itype, char *jobz, char *range, char *uplo,       #
#             int *n, double *AP, double *BP, double *vl,            #
#             double *vu, int *il, int *iu, int *m, double *W,       #
#             double *Z, int *ldz);                                  #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# itype              1  - A*x = lambda*B*x                           #
#                    2  - A*B*x = lambda*x                           #
#                    3  - B*A*x = lambda*x                           #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A and B stored           #
#                   "U" - Lower triangle of A and B stored           #
# n                 Order of the matrix A and B                      #
# ldz               Leading dimension of matrix Z;                   #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# AP                On entry symmetric input matrix, stored in       #
# (double[s])       packed format by columns. Depending on the       #
# s = (n*(n+1))/2   value of 'uplo' only the upper or lower          #
#                   triangular part is stored                        #
#                   On output the array will contain the             #
#                   'm' computed eigenvectors                        #
# BP                On entry symmetric definite input matrix, stored #
# (double[s])       in packed format by columns.                     #
# s = (n*(n+1))/2   Depending on the value of 'uplo' only the upper  #
#                   or lower triangular part is stored               #
#                   On output overwritten                            #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input      #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input      #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed  #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues    #
# Z (double[m*n])   Eigenvectors                                     #
#                   Enough space must be provided to store the       #
#                   vectors. 'm' should be bigger or equal           #
#                   to 'n' for range="A" or "V" and 'iu-il+1' for    #
#                   range="I".                                       #
#                                                                    #
# NOTICE:           The routine will allocate work space of size     #
#                   double[n*n] for range="A" or "V" and double[m*n] #
#                   for range="I"                                    #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                 0 - Success                                        #
#                 1 - Wrong input parameter                          #
#                 2 - Misc errors                                    #
#                                                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# DSPGEIG(ITYPE, JOBZ, RANGE, UPLO, N, AP, BP,                       #
#         VL, VU, IL, IU, M, W, Z, LDZ, INFO)                        #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDZ, INFO                           #
# DOUBLE PRECISION AP(*), BP(*), VL, VU, W(*), Z(*,*)                #
#                                                                    #
######################################################################



Hermitian-definite case:

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
# int zhegeig(int *itype, char *jobz, char *range, char *uplo,       #
#             int *n, double complex *A, int *lda,                   #
#	      double complex *B, int *ldb, double *vl, double *vu,   #
#	      int *il, int *iu, int *m, double *W);                  #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# itype              1  - A*x = lambda*B*x                           #
#                    2  - A*B*x = lambda*x                           #
#                    3  - B*A*x = lambda*x                           #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A and B stored           #
#                   "U" - Lower triangle of A and B stored           #
# n                 Order of the matrix A and B                      #
# lda               Leading dimension of matrix A;                   #
#                   often equal to matrix size n                     #
# ldb               Leading dimension of eigenvector matrix B;       #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# A                 On entry symmetric input matrix, stored          #
# (double           in column major ordering. Depending on the       #
#  complex[lda*n])  value of 'uplo' only the upper or lower          #
#                   triangular part is referenced                    #
#                   On output the array will contain the             #
#                   'm' computed eigenvectors                        #
# B                 On entry symmetric definite input matrix,        #
# (double           stored in column major ordering.                 #
#  complex[ldb*n])  Depending on the value of 'uplo' only the upper  #
#                   or lower triangular part is referenced           #
#                   On output overwritten                            #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input      #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input      #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed  #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues    #
#                                                                    #
# NOTICE:           The routine will allocate work space of size     #
#                   double[n*n] for range="A" or "V" and double[m*n] #
#                   for range="I"                                    #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                 0 - Success                                        #
#                 1 - Wrong input parameter                          #
#                 2 - Misc errors                                    #
#                                                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# ZHEGEIG(ITYPE, JOBZ, RANGE, UPLO, N, A, LDA, B, LDB,               #
#         VL, VU, IL, IU, M, W, INFO);                               #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDA, LDB, INFO                      #
# DOUBLE PRECISION VL, VU, W(*)                                      #
# COMPLEX*16       A(*,*), B(*,*)                                    #
#                                                                    #
######################################################################
 


Hermitian-definite case (packed storage):

######################################################################
# C function prototype:                                              #
###################################################################### 
#                                                                    #
# int zhpgeig(int *itype, char *jobz, char *range, char *uplo,       #
#             int *n, double complex *AP, double complex *BP,        #
#	      double *vl, double *vu, int *il, int *iu, int *m,      #
#	      double *W, double complex *Z, int *ldz);               #
#                                                                    #
# Arguments:                                                         #
# ----------                                                         #
#                                                                    #
# INPUTS:                                                            #
# -------                                                            #
# itype              1  - A*x = lambda*B*x                           #
#                    2  - A*B*x = lambda*x                           #
#                    3  - B*A*x = lambda*x                           #
# jobz              "N" - compute only eigenvalues                   #
#                   "V" - compute also eigenvectors                  #
# range             "A" - all                                        #
#                   "V" - by interval: (VL,VU]                       #
#                   "I" - by index:     IL-IU                        #
# uplo              "L" - Upper triangle of A and B stored           #
#                   "U" - Lower triangle of A and B stored           #
# n                 Order of the matrix A and B                      #
# ldz               Leading dimension of matrix Z;                   #
#                   often equal to matrix size n                     #
#                                                                    #
# INPUT + OUTPUT:                                                    #
# ---------------                                                    #
# AP                On entry Hermitian input matrix, stored in       #
# (double           packed format by columns. Depending on the       #
#  complex[s])      value of 'uplo' only the upper or lower          #
# s = (n*(n+1))/2   triangular part is stored                        #
#                   On output the array will contain the             #
#                   'm' computed eigenvectors                        #
# BP                On entry Hermitian definite input matrix, stored #
# (double           in packed format by columns.                     #
#  complex[s])      Depending on the value of 'uplo' only the upper  #
# s = (n*(n+1))/2   or lower triangular part is stored               #
#                   On output overwritten                            #
# vl                If range="V", lower bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input      #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# vu                If range="V", upper bound of interval            #
#                   (vl,vu], on output refined                       #
#                   If range="A" or "I" not referenced as input.     #
#                   On output the interval (vl,vu] contains ALL      #
#                   the computed eigenvalues                         #
# il                If range="I", lower index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input.     #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
# iu                If range="I", upper index (1-based indexing) of  #
#                   the subset 'il' to 'iu'                          #
#                   If range="A" or "V" not referenced as input      #
#                   On output the eigenvalues with index il to iu    #
#                   are computed                                     #
#                                                                    #
# OUTPUT:                                                            #
# -------                                                            #
# m                 Number of eigenvalues and eigenvectors computed  #
# W (double[n])     Eigenvalues                                      #
#                   The first 'm' entries contain the eigenvalues    #
# Z (double         Eigenvectors                                     #
#    complex[m*n])  Enough space must be provided to store the       #
#                   vectors. 'm' should be bigger or equal           #
#                   to 'n' for range="A" or "V" and 'iu-il+1' for    #
#                   range="I".                                       #
#                                                                    #
# NOTICE:           The routine will allocate work space of size     #
#                   double[n*n] for range="A" or "V" and double[m*n] #
#                   for range="I"                                    #
#                                                                    #
# RETURN VALUE:                                                      #
# -------------                                                      #
#                 0 - Success                                        #
#                 1 - Wrong input parameter                          #
#                 2 - Misc errors                                    #
#                                                                    #
######################################################################


######################################################################
# Fortran prototype:                                                 #
######################################################################
#                                                                    #
# ZHPGEIG(ITYPE, JOBZ, RANGE, UPLO, N, AP, BP,                       #
#         VL, VU, IL, IU, M, W, Z, LDZ, INFO)                        #
#                                                                    #
# CHARACTER        JOBZ, RANGE, UPLO                                 #
# INTEGER          N, IL, IU, M, LDZ, INFO                           #
# DOUBLE PRECISION VL, VU, W(*)                                      #
# COMPLEX*16       AP(*), BP(*), Z(*,*)                              #
#                                                                    #
######################################################################



COMMENTS AND BUGS:
petschow@aices.rwth-aachen.de

About

The MRRR algorithm for multi-core and SMP systems

Resources

License

Stars

Watchers

Forks

Releases

No releases published

Packages

No packages published

Languages

  • Fortran 58.6%
  • C 38.5%
  • Makefile 2.2%
  • Other 0.7%