C Library Interface

The PRIMME interface is composed of the following functions. To solve real symmetric and complex Hermitian problems, standard \(A x = \lambda x\) and generalized \(A x = \lambda B x\), call:

int dprimme (double *evals, double *evecs, double *resNorms,
                        primme_params *primme)
int zprimme (double *evals, PRIMME_COMPLEX_DOUBLE *evecs,
                 double *resNorms, primme_params *primme)

There are more versions for matrix problems working in other precisions:

Precision

Real

Complex

half

hprimme() hsprimme()

kprimme() ksprimme()

single

sprimme()

cprimme()

double

dprimme()

zprimme()

To solve standard eigenproblems with normal but not necessarily Hermitian matrices call:

int zprimme_normal (PRIMME_COMPLEX_DOUBLE *evals,
                 PRIMME_COMPLEX_DOUBLE *evecs,
                 double *resNorms, primme_params *primme)

There are more versions for matrix problems working in other precisions:

Precision

Complex

half

kprimme_normal() kcprimme_normal()

single

cprimme_normal()

double

zprimme_normal()

Other useful functions:

void primme_initialize (primme_params *primme)
int primme_set_method (primme_preset_method method,
                                                     primme_params *params)
void primme_display_params (primme_params primme)
void primme_free (primme_params *primme)

PRIMME stores its data on the structure primme_params. See Parameters Guide for an introduction about its fields.

Running

To use PRIMME, follow these basic steps.

  1. Include:

    #include "primme.h"   /* header file is required to run primme */
    
  2. Initialize a PRIMME parameters structure for default settings:

    primme_params primme;
    primme_initialize (&primme);
  3. Set problem parameters (see also Parameters Guide), and, optionally, set one of the preset methods:

    primme.matrixMatvec = LaplacianMatrixMatvec; /* MV product */
    primme.n = 100;                   /* set problem dimension */
    primme.numEvals = 10;       /* Number of wanted eigenpairs */
    ret = primme_set_method (method, &primme);
    ...
  4. Then call the solver:

    ret = dprimme (evals, evecs, resNorms, &primme);

    The call arguments are:

    • evals, array to return the found eigenvalues;

    • evecs, array to return the found eigenvectors;

    • resNorms, array to return the residual norms of the found eigenpairs; and

    • ret, returned error code.

  5. To free the work arrays in PRIMME:

    primme_free (&primme);

See usage examples at the directory examples.

Parameters Guide

PRIMME stores the data on the structure primme_params, which has the next fields:

Basic
PRIMME_INT n, matrix dimension.
void (* matrixMatvec )(...), matrix-vector product.
void (* massMatrixMatvec )(...), mass matrix-vector product (null for standard problems).
int numEvals, how many eigenpairs to find.
primme_target target, which eigenvalues to find.
int numTargetShifts, for targeting interior eigenpairs.
double * targetShifts
double eps, tolerance of the residual norm of converged eigenpairs.

For parallel programs
int numProcs, number of processes
int procID, rank of this process
PRIMME_INT nLocal, number of rows stored in this process
void (* globalSumReal )(...), sum reduction among processes
void (* broadcastReal )(...), broadcast array among processes

Accelerate the convergence
void (* applyPreconditioner )(...), preconditioner-vector product.
int initSize, initial vectors as approximate solutions.

User data
void * commInfo
void * matrix
void * massMatrix
void * convtest
void * monitor

Advanced options
PRIMME_INT ldevecs, leading dimension of the evecs.
int numOrthoConst, orthogonal constrains to the eigenvectors.
PRIMME_INT maxMatvecs
PRIMME_INT iseed [4]
double aNorm
double BNorm
double invBNorm
FILE * outputFile
primme_init initBasisMode
struct projection_params projectionParams
struct restarting_params restartingParams
struct correction_params correctionParams
struct primme_stats stats
void (* convTestFun )(...), custom convergence criterion.
PRIMME_INT ldOPs, leading dimension to use in matrixMatvec.
void (* monitorFun )(...), custom convergence history.
primme_op_datatype matrixMatvec_type
primme_op_datatype massMatrixMatvec_type
primme_op_datatype applyPreconditioner_type
primme_op_datatype globalSumReal_type
primme_op_datatype broadcastReal_type
primme_op_datatype internalPrecision
primme_orth orth

PRIMME requires the user to set at least the dimension of the matrix (n) and the matrix-vector product (matrixMatvec), as they define the problem to be solved. For parallel programs, nLocal, procID and globalSumReal are also required.

In addition, most users would want to specify how many eigenpairs to find, numEvals, and provide a preconditioner applyPreconditioner (if available).

It is useful to have set all these before calling primme_set_method(). Also, if users have a preference on maxBasisSize, maxBlockSize, etc, they should also provide them into primme_params prior to the primme_set_method() call. This helps primme_set_method() make the right choice on other parameters. It is sometimes useful to check the actual parameters that PRIMME is going to use (before calling it) or used (on return) by printing them with primme_display_params().

Interface Description

The next enumerations and functions are declared in primme.h.

?primme

int hprimme(PRIMME_HALF *evals, PRIMME_HALF *evecs, PRIMME_HALF *resNorms, primme_params *primme)
int hsprimme(float *evals, PRIMME_HALF *evecs, float *resNorms, primme_params *primme)
int kprimme(PRIMME_HALF *evals, PRIMME_COMPLEX_HALF *evecs, PRIMME_HALF *resNorms, primme_params *primme)
int ksprimme(float *evals, PRIMME_COMPLEX_HALF *evecs, float *resNorms, primme_params *primme)

New in version 3.0.

int sprimme(float *evals, float *evecs, float *resNorms, primme_params *primme)
int cprimme(float *evals, PRIMME_COMPLEX_FLOAT *evecs, float *resNorms, primme_params *primme)

New in version 2.0.

int dprimme(double *evals, double *evecs, double *resNorms, primme_params *primme)
int zprimme(double *evals, PRIMME_COMPLEX_DOUBLE *evecs, double *resNorms, primme_params *primme)

Solve a real symmetric/Hermitian standard or generalized eigenproblem.

All arrays should be hosted on CPU. The computations are performed on CPU (see magma_dprimme() for using GPUs).

Parameters
  • evals – array at least of size numEvals to store the computed eigenvalues; all processes in a parallel run return this local array with the same values.

  • evecs – array at least of size nLocal times (numOrthoConst + numEvals) with leading dimension ldevecs to store column-wise the (local part for this process of the) computed eigenvectors.

  • resNorms – array at least of size numEvals to store the residual norms of the computed eigenpairs; all processes in parallel run return this local array with the same values.

  • primme – parameters structure.

Returns

error indicator; see Error Codes.

On input, evecs should start with the content of the numOrthoConst vectors, followed by the initSize vectors.

On return, the i-th eigenvector starts at evecs[( numOrthoConst + i)* ldevecs ], with value evals[i] and associated residual 2-norm resNorms[i]. The first vector has index i=0. The number of eigenpairs marked as converged (see eps) is returned on initSize. Since version 4.0, if the returned error code is PRIMME_MAIN_ITER_FAILURE, PRIMME may return also unconverged eigenpairs and its residual norms in evecs, evals, and resNorms starting at i=|initSize| and going up to either numEvals-1 or the last resNorms[i] with non-negative value.

All internal operations are performed at the same precision than evecs unless the user sets internalPrecision otherwise. The functions hsprimme() and ksprimme() perform all computations in half precision by default and report the eigenvalues and the residual norms in single precision. These functions may help in applications that may be not built with a compiler supporting half precision.

The type and precision of the callbacks is also the same as evecs. Although this can be changed. See details for matrixMatvec, massMatrixMatvec, applyPreconditioner, globalSumReal, broadcastReal, and convTestFun.

magma_?primme

int magma_hprimme(PRIMME_HALF *evals, PRIMME_HALF *evecs, PRIMME_HALF *resNorms, primme_params *primme)
int magma_hsprimme(float *evals, PRIMME_HALF *evecs, float *resNorms, primme_params *primme)
int magma_kprimme(PRIMME_HALF *evals, PRIMME_COMPLEX_HALF *evecs, PRIMME_HALF *resNorms, primme_params *primme)
int magma_sprimme(float *evals, float *evecs, float *resNorms, primme_params *primme)
int magma_ksprimme(float *evals, PRIMME_COMPLEX_HALF *evecs, float *resNorms, primme_params *primme)
int magma_cprimme(float *evals, PRIMME_COMPLEX_FLOAT *evecs, float *resNorms, primme_params *primme)
int magma_dprimme(double *evals, double *evecs, double *resNorms, primme_params *primme)
int magma_zprimme(double *evals, PRIMME_COMPLEX_DOUBLE *evecs, double *resNorms, primme_params *primme)

Solve a real symmetric/Hermitian standard or generalized eigenproblem.

Most of the computations are performed on GPU (see dprimme() for using only the CPU).

Parameters
  • evals – CPU array at least of size numEvals to store the computed eigenvalues; all processes in a parallel run return this local array with the same values.

  • evecs – GPU array at least of size nLocal times (numOrthoConst + numEvals) with leading dimension ldevecs to store column-wise the (local part for this process of the) computed eigenvectors.

  • resNorms – CPU array at least of size numEvals to store the residual norms of the computed eigenpairs; all processes in parallel run return this local array with the same values.

  • primme – parameters structure.

Returns

error indicator; see Error Codes.

On input, evecs should start with the content of the numOrthoConst vectors, followed by the initSize vectors.

On return, the i-th eigenvector starts at evecs[( numOrthoConst + i)* ldevecs ], with value evals[i] and associated residual 2-norm resNorms[i]. The first vector has index i=0. The number of eigenpairs marked as converged (see eps) is returned on initSize. Since version 4.0, if the returned error code is PRIMME_MAIN_ITER_FAILURE, PRIMME may return also unconverged eigenpairs and its residual norms in evecs, evals, and resNorms starting at i=|initSize| and going up to either numEvals-1 or the last resNorms[i] with non-negative value.

All internal operations are performed at the same precision than evecs unless the user sets internalPrecision otherwise. The functions hsprimme() and ksprimme() perform all computations in half precision by default and report the eigenvalues and the residual norms in single precision. These functions may help in applications that may be not built with a compiler supporting half precision.

The type and precision of the callbacks is also the same as evecs. Although this can be changed. See details for matrixMatvec, massMatrixMatvec, applyPreconditioner, globalSumReal, broadcastReal, and convTestFun.

New in version 3.0.

?primme_normal

int kprimme_normal(PRIMME_COMPLEX_HALF *evals, PRIMME_COMPLEX_HALF *evecs, PRIMME_HALF *resNorms, primme_params *primme)
int kcprimme_normal(PRIMME_COMPLEX_FLOAT *evals, PRIMME_COMPLEX_HALF *evecs, float *resNorms, primme_params *primme)
int cprimme_normal(PRIMME_COMPLEX_FLOAT *evals, PRIMME_COMPLEX_FLOAT *evecs, float *resNorms, primme_params *primme)
int zprimme_normal(PRIMME_COMPLEX_DOUBLE *evals, PRIMME_COMPLEX_DOUBLE *evecs, double *resNorms, primme_params *primme)

Solve a normal standard eigenproblem, which may not be Hermitian.

All arrays should be hosted on CPU. The computations are performed on CPU (see magma_zprimme_normal() for using GPUs).

Parameters
  • evals – array at least of size numEvals to store the computed eigenvalues; all processes in a parallel run return this local array with the same values.

  • evecs – array at least of size nLocal times (numOrthoConst + numEvals) with leading dimension ldevecs to store column-wise the (local part for this process of the) computed eigenvectors.

  • resNorms – array at least of size numEvals to store the residual norms of the computed eigenpairs; all processes in parallel run return this local array with the same values.

  • primme – parameters structure.

Returns

error indicator; see Error Codes.

On input, evecs should start with the content of the numOrthoConst vectors, followed by the initSize vectors.

On return, the i-th eigenvector starts at evecs[( numOrthoConst + i)* ldevecs ], with value evals[i] and associated residual 2-norm resNorms[i]. The first vector has index i=0. The number of eigenpairs marked as converged (see eps) is returned on initSize. Since version 4.0, if the returned error code is PRIMME_MAIN_ITER_FAILURE, PRIMME may return also unconverged eigenpairs and its residual norms in evecs, evals, and resNorms starting at i=|initSize| and going up to either numEvals-1 or the last resNorms[i] with non-negative value.

All internal operations are performed at the same precision than evecs unless the user sets internalPrecision otherwise. The functions hsprimme() and ksprimme() perform all computations in half precision by default and report the eigenvalues and the residual norms in single precision. These functions may help in applications that may be not built with a compiler supporting half precision.

The type and precision of the callbacks is also the same as evecs. Although this can be changed. See details for matrixMatvec, applyPreconditioner, globalSumReal, broadcastReal, and convTestFun.

New in version 3.0.

magma_?primme_normal

int magma_kprimme_normal(PRIMME_COMPLEX_HALF *evals, PRIMME_COMPLEX_HALF *evecs, PRIMME_HALF *resNorms, primme_params *primme)
int magma_kcprimme_normal(PRIMME_COMPLEX_FLOAT *evals, PRIMME_COMPLEX_HALF *evecs, float *resNorms, primme_params *primme)
int magma_cprimme_normal(PRIMME_COMPLEX_FLOAT *evals, PRIMME_COMPLEX_FLOAT *evecs, float *resNorms, primme_params *primme)
int magma_zprimme_normal(PRIMME_COMPLEX_DOUBLE *evals, PRIMME_COMPLEX_DOUBLE *evecs, double *resNorms, primme_params *primme)

Solve a normal standard eigenproblem, which may not be Hermitian.

Most of the computations are performed on GPU (see zprimme_normal() for using only the CPU).

Parameters
  • evals – CPU array at least of size numEvals to store the computed eigenvalues; all processes in a parallel run return this local array with the same values.

  • evecs – GPU array at least of size nLocal times (numOrthoConst + numEvals) with leading dimension ldevecs to store column-wise the (local part for this process of the) computed eigenvectors.

  • resNorms – CPU array at least of size numEvals to store the residual norms of the computed eigenpairs; all processes in parallel run return this local array with the same values.

  • primme – parameters structure.

Returns

error indicator; see Error Codes.

On input, evecs should start with the content of the numOrthoConst vectors, followed by the initSize vectors.

On return, the i-th eigenvector starts at evecs[( numOrthoConst + i)* ldevecs ], with value evals[i] and associated residual 2-norm resNorms[i]. The first vector has index i=0. The number of eigenpairs marked as converged (see eps) is returned on initSize. Since version 4.0, if the returned error code is PRIMME_MAIN_ITER_FAILURE, PRIMME may return also unconverged eigenpairs and its residual norms in evecs, evals, and resNorms starting at i=|initSize| and going up to either numEvals-1 or the last resNorms[i] with non-negative value.

All internal operations are performed at the same precision than evecs unless the user sets internalPrecision otherwise. The functions hsprimme() and ksprimme() perform all computations in half precision by default and report the eigenvalues and the residual norms in single precision. These functions may help in applications that may be not built with a compiler supporting half precision.

The type and precision of the callbacks is also the same as evecs. Although this can be changed. See details for matrixMatvec, applyPreconditioner, globalSumReal, broadcastReal, and convTestFun.

New in version 3.0.

primme_initialize

void primme_initialize(primme_params *primme)

Initialize a PRIMME parameters structure to the default values.

After calling dprimme() (or a variant), call primme_free() to release allocated resources by PRIMME.

Parameters
  • primme – parameters structure.

Example:

primme_params primme;
primme_initialize(&primme);

primme.n = 100;
...
dprimme(evals, evecs, rnorms, &primme);
...

primme_free(&primme);

See the alternative function primme_params_create() that also allocates the primme_params structure.

primme_params_create

primme_params *primme_params_create(void)

Allocate and initialize a parameters structure to the default values.

After calling dprimme() (or a variant), call primme_params_destroy() to release allocated resources by PRIMME.

Returns

pointer to a parameters structure.

Example:

primme_params *primme = primme_params_create();

primme->n = 100;
...
dprimme(evals, evecs, rnorms, primme);
...

primme_params_destroy(primme);

See the alternative function primme_initialize() that only initializes the structure.

New in version 3.0.

primme_set_method

int primme_set_method(primme_preset_method method, primme_params *primme)

Set PRIMME parameters to one of the preset configurations.

Parameters

See also Preset Methods.

primme_display_params

void primme_display_params(primme_params primme)

Display all printable settings of primme into the file descriptor outputFile.

Parameters
  • primme – parameters structure.

primme_free

void primme_free(primme_params *primme)

Free memory allocated by PRIMME.

Parameters
  • primme – parameters structure.

primme_params_destroy

int primme_params_destroy(primme_params *primme)

Free memory allocated by PRIMME associated to a parameters structure created with primme_params_create().

Parameters
  • primme – parameters structure.

Returns

nonzero value if the call is not successful.

New in version 3.0.