This example illustrates how to solve the standard symmetric-definite eigenvalue problem for a symmetric matrix
Given an
A solution for this problem is given by a pair
-
$X$ is an$n \times n$ orthogonal matrix containing (as columns) the eigenvectors$x_i$ for$i = 0, \dots, n-1$ and -
$\Lambda$ is an$n \times n$ diagonal matrix containing the eigenvalues$\lambda_i$ for$i = 0, \dots, n-1$
such that
- Declare and initialize a number of constants for the input matrix.
- Allocate and initialize the host matrix
$A$ . - Allocate device memory and copy input data from host to device.
- Create a hipSOLVER compatibility API handle and query the working space size.
- Invoke
hipsolverDnDsyevdxto calculate the eigenvalues and eigenvectors of matrix$A$ . - Copy the resulting devInfo value
d_syevdx_infovalue to the host and check if the operation was successful. - Copy the resulting eigenvalues and eigenvectors to the host and print the eigenvalues.
- Validate the solution by checking if
$X \cdot \Lambda - A \cdot X$ is the zero matrix using the hipBLAS API. - Free device memory and the handles.
- Print validation result.
-
hipSOLVER is initialized by calling
hipsolverDnCreate(hipsolverHandle_t*)and it is terminated by callinghipsolverDnDestroy(hipsolverHandle_t). -
In this example
hipsolverDnHandle_tis used instead ofhipsolverHandle_t.hipsolverDnHandle_tis actually a typedef ofhipsolverHandle_t, so they can be used equivalently. -
hipsolverDn[SD]syevdxcomputes the eigenvalues of an$n \times n$ symmetric matrix$A$ . The correct function signature should be chosen based on the datatype of the input matrix:-
S(single-precision real:float) -
D(double-precision real:double)
For single- and double-precision complex values, the function
hipsolverDn[CZ]heevdx(...)is available in hipSOLVER.In this example, a double-precision real input matrix is used, in which case the function accepts the following parameters:
hipsolverHandle_t handle-
hipsolverEigMode_t jobz: Specifies whether the eigenvectors should also be calculated besides the eigenvalues. The following values are accepted:-
HIPSOLVER_EIG_MODE_NOVECTOR: Calculate the eigenvalues only. -
HIPSOLVER_EIG_MODE_VECTOR: Calculate both the eigenvalues and the eigenvectors. The eigenvectors are calculated by a divide and conquer algorithm and are written to the memory location specified by*A.
-
-
hipsolverEigRange_t range: Specifies a range of eigenvalues to be returned. The following values are accepted:-
HIPSOLVER_EIG_RANGE_ALL: The whole spectrum is returned. -
HIPSOLVER_EIG_RANGE_V: Only the eigenvalues in the interval(vl, vu]are returned.vl$>$ vumust be satisfied. -
HIPSOLVER_EIG_RANGE_I: Only the eigenvalues from theil-th to theiu-th are returned.$1$ $\leq$ il$\leq$ iu$\leq$ $n$ must be satisfied.
-
-
hipSolverFillMode_t uplo: Specifies whether the upper or lower triangle of the symmetric matrix is stored. The following values are accepted:-
HIPSOLVER_FILL_MODE_UPPER: The provided*Apointer points to the upper triangle matrix data. -
HIPSOLVER_FILL_MODE_LOWER: The provided*Apointer points to the lower triangle matrix data.
-
-
int n: Number of rows and columns of$A$ . -
double *A: Pointer to matrix$A$ in device memory. -
int lda: Leading dimension of matrix$A$ . -
double vl: Lower bound of the interval to be searched for eigenvalues ifrange=HIPSOLVER_EIG_RANGE_V. -
double vu: Upper bound of the interval to be searched for eigenvalues ifrange=HIPSOLVER_EIG_RANGE_V. -
int il: Smallest index of the eigenvalues to be returned ifrange=HIPSOLVER_EIG_RANGE_I. -
int iu: Largest index of the eigenvalues to be returned ifrange=HIPSOLVER_EIG_RANGE_I. -
int *nev: Number of eigenvalues returned. -
double *W: Pointer to array$W$ in device memory, where the resulting eigenvalues are written. -
double *work: Pointer to working space in device memory. -
int lwork: Size of the working space. -
int *devInfo: Convergence result of the function in device memory.- If 0, the algorithm converged.
- If greater than 0 (
devInfo = ifor$0 < i \leq n$ ), thendevInfoeigenvectors did not converge. - For CUDA backend, if lesser than 0 (
devInfo = -ifor$0 < i \leq n$ ) then the the$i^{th}$ parameter is wrong (not counting the handle).
Return type:
hipsolverStatus_t. -
-
hipsolverDn[SD]syevdxinternally calls tocusolverDn[SD]syevdxfor CUDA backend and to a rocSOLVER's internalsyevxfunction (not the one from the public API) for ROCm backend, as nohipsolver[SD]syevdxfunction exists yet in hipSOLVER regular API. -
hipsolverDn[SD]syevdx_bufferSizeallows to obtain the size (in bytes) needed for the working space for thehipsolverDn[SD]syevdxfunction. The character matched in[SD]coincides with the one inhipsolverDn[SD]syevdx.This function accepts the following input parameters:
hipsolverHandle_t handle-
hipsolverEigMode_t jobz: Specifies whether the eigenvectors should also be calculated besides the eigenvalues. -
hipsolverEigRange_t range: Specifies a range of eigenvalues to be returned. -
hipSolverFillMode_t uplo: Specifies whether the upper or lower triangle of the symmetric matrix is stored. -
int n: Number of rows and columns of$A$ . -
double *A: Pointer to matrix$A$ in device memory. -
int lda: Leading dimension of matrix$A$ . -
double vl: Lower bound of the interval to be searched for eigenvalues ifrange=HIPSOLVER_EIG_RANGE_V. -
double vu: Upper bound of the interval to be searched for eigenvalues ifrange=HIPSOLVER_EIG_RANGE_V. -
int il: Smallest index of the eigenvalues to be returned ifrange=HIPSOLVER_EIG_RANGE_I. -
int iu: Largest index of the eigenvalues to be returned ifrange=HIPSOLVER_EIG_RANGE_I. -
int *nev: Number of eigenvalues returned. -
double *W: Pointer to array$W$ in device memory, where the resulting eigenvalues are written. -
int *lwork: The required buffer size is written to this location.
Return type:
hipsolverStatus_t.
- For validating the solution we have used the hipBLAS functions
hipblasDdgmmandhipblasDgemm.hipblasDgemmcomputes a general scaled matrix multiplication$\left(C = \alpha \cdot A \cdot B + \beta \cdot C\right)$ and is showcased (strided-batched and with single-precision real type) in the gemm_strided_batched example.-
The
hipblas[SDCZ]dgmmfunction performs a matrix--matrix operation between a diagonal matrix and a general$m \times n$ matrix. The order of the multiplication can be determined using ahipblasSideMode_ttype parameter:-
HIPBLAS_SIDE_RIGHT: the operation performed is$C = A \cdot diag(x)$ . This is the one used in the example for computing$X \cdot \Lambda$ . -
HIPBLAS_SIDE_LEFT: the operation performed is$C = diag(x) \cdot A$ .
The correct function signature should be chosen based on the datatype of the input matrices:
-
S(single-precision real:float) -
D(double-precision real:double) -
C(single-precision complex:hipblasComplex) -
Z(double-precision complex:hipblasDoubleComplex).
Return type:
hipblasStatus_t. -
-
- The
hipblasPointerMode_ttype controls whether scalar parameters must be allocated on the host (HIPBLAS_POINTER_MODE_HOST) or on the device (HIPBLAS_POINTER_MODE_DEVICE). It is set by usinghipblasSetPointerMode.
HIPSOLVER_EIG_MODE_VECTORHIPSOLVER_FILL_MODE_UPPER
HIPSOLVER_EIG_RANGE_IhipsolverDnCreatehipsolverDnDestroyhipsolverDnDsyevdxhipsolverDnDsyevdx_bufferSizehipsolverDnHandle_t
HIPBLAS_OP_NHIPBLAS_POINTER_MODE_HOSTHIPBLAS_SIDE_RIGHThipblasCreatehipblasDdgmmhipblasDestroyhipblasDgemmhipblasHandle_thipblasSetPointerMode
hipFreehipMallochipMemcpyhipMemcpyDeviceToHosthipMemcpyHostToDevice