FORTRAN 90 Library Interface¶
New in version 3.0.
The next enumerations and functions are declared in primme_f90.inc
.
-
type
c_ptr
¶ Fortran datatype for C pointers in module iso_c_binding.
-
subroutine
primme_eigs_matvec
(x, ldx, y, ldy, blockSize, primme, ierr)¶ Abstract interface for the callbacks
matrixMatvec
,massMatrixMatvec
, andapplyPreconditioner
.- Parameters
x (ldx,*) [type(*),in] :: matrix of size
nLocal
xblockSize
in column-major order with leading dimensionldx
.ldx [c_int64_t] :: the leading dimension of the array
x
.y (ldy,*) [type(*),out] :: matrix of size
nLocal
xblockSize
in column-major order with leading dimensionldy
.ldy [c_int64_t] :: the leading dimension of the array
y
.blockSize [c_int,in] :: number of columns in
x
andy
.primme [c_ptr,in] :: parameters structure created by
primme_params_create()
.ierr [c_int,out] :: output error code; if it is set to non-zero, the current call to PRIMME will stop.
See more details about the precision and type for x and y in the documentation of the callbacks.
primme_params_create¶
-
function
primme_params_create
()¶ Allocate and initialize a parameters structure to the default values.
After calling
xprimme()
(or a variant), callprimme_params_destroy()
to release allocated resources by PRIMME.- Return
primme_params_create [c_ptr] :: pointer to a parameters structure.
primme_set_method¶
-
function
primme_set_method
(method, primme)¶ Set PRIMME parameters to one of the preset configurations.
- Parameters
method [c_int,in] ::
preset configuration. One of:
PRIMME_DYNAMIC
PRIMME_DEFAULT_MIN_TIME
PRIMME_DEFAULT_MIN_MATVECS
PRIMME_Arnoldi
PRIMME_GD
PRIMME_GD_plusK
PRIMME_GD_Olsen_plusK
PRIMME_JD_Olsen_plusK
PRIMME_RQI
PRIMME_JDQR
PRIMME_JDQMR
PRIMME_JDQMR_ETol
PRIMME_STEEPEST_DESCENT
PRIMME_LOBPCG_OrthoBasis
PRIMME_LOBPCG_OrthoBasis_Window
See
primme_preset_method
.primme [c_ptr,in] :: parameters structure created by
primme_params_create()
.
- Return
primme_set_method [c_int] :: nonzero value if the call is not successful.
primme_params_destroy¶
-
function
primme_params_destroy
(primme)¶ Free memory allocated by PRIMME associated to a parameters structure created with
primme_params_create()
.- Parameters
primme [c_ptr] :: parameters structure.
- Return
primme_params_destroy :: nonzero value if the call is not successful.
xprimme¶
-
function
xprimme
(evals, evecs, resNorms, 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_xprimme()
for using GPUs).- Parameters
evals (*) [real(kind),out] :: 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 (*) [real(kind) or complex(kind)] :: array at least of size
nLocal
times (numOrthoConst
+numEvals
) with leading dimensionldevecs
to store column-wise the (local part for this process of the) computed eigenvectors.resNorms (*) [real(kind),out] :: 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 [c_ptr,in] :: parameters structure created by
primme_params_create()
.
- Return
xprimme [c_int] :: error indicator; see Error Codes.
The arrays
evals
,evecs
, andresNorms
should have the same kind.On input,
evecs
should start with the content of thenumOrthoConst
vectors, followed by theinitSize
vectors.On return, the i-th eigenvector starts at evecs((
numOrthoConst
+ i - 1)*ldevecs
+ 1), with value evals(i) and associated residual 2-norm resNorms(i). The first vector has index i=1. The number of eigenpairs marked as converged (seeeps
) is returned oninitSize
. 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
+ 1 and going up to eithernumEvals
or the last resNorms(i) with non-negative value.All internal operations are performed at the same precision than
evecs
unless the user setsinternalPrecision
otherwise.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
, andconvTestFun
.
magma_xprimme¶
-
function
magma_xprimme
(evals, evecs, resNorms, primme)¶ Solve a real symmetric/Hermitian standard or generalized eigenproblem.
Most of the computations are performed on GPU (see
xprimme()
for using only the CPU).- Parameters
evals (*) [real(kind),out] :: 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 dimensionldevecs
to store column-wise the (local part for this process of the) computed eigenvectors.resNorms (*) [real(kind),out] :: 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 [c_ptr,in] :: parameters structure created by
primme_params_create()
.
- Return
magma_xprimme [c_int] :: error indicator; see Error Codes.
The arrays
evals
,evecs
, andresNorms
should have the same kind.Further descriptions of evals, evecs, and resNorms on notes in function
xprimme()
.
xprimme_normal¶
-
function
xprimme_normal
(evals, evecs, resNorms, 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_xprimme_normal()
for using GPUs).- Parameters
evals (*) [complex(kind),out] :: 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 dimensionldevecs
to store column-wise the (local part for this process of the) computed eigenvectors.resNorms (*) [real(kind),out] :: 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 [c_ptr,in] :: parameters structure created by
primme_params_create()
.
- Return
xprimme_normal [c_int] :: error indicator; see Error Codes.
The arrays
evals
,evecs
, andresNorms
should have the same kind.Further descriptions of evals, evecs, and resNorms on notes in function
xprimme()
.
magma_xprimme_normal¶
-
function
magma_xprimme_normal
(evals, evecs, resNorms, primme)¶ Solve a normal standard eigenproblem, which may not be Hermitian.
Most of the arrays are stored on GPU, and also most of the computations are done on GPU (see
xprimme()
for using only the CPU).- Parameters
evals (*) [complex(kind),out] :: 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 dimensionldevecs
to store column-wise the (local part for this process of the) computed eigenvectors.resNorms (*) [real(kind),out] :: 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 [c_ptr,in] :: parameters structure created by
primme_params_create()
.
- Return
magma_xprimme_normal [c_int] :: error indicator; see Error Codes.
The arrays
evals
,evecs
, andresNorms
should have the same kind.Further descriptions of evals, evecs, and resNorms on notes in function
xprimme()
.
primme_set_member¶
-
function
primme_set_member
(primme, label, value)¶ Set a value in some field of the parameter structure.
- Parameters
primme [c_ptr,in] :: parameters structure created by
primme_params_create()
.label [c_int,in] ::
field where to set value. One of:
PRIMME_stats_flopsDense
PRIMME_stats_timeDense
PRIMME_stats_estimateMaxEVal
PRIMME_stats_estimateBNorm
PRIMME_stats_estimateInvBNorm
value [in] :: value to set. The allowed types are c_int64, c_double, c_ptr, c_funptr and
procedure(primme_eigs_matvec)
- Return
primme_set_member [c_int] ::
nonzero value if the call is not successful.
Examples:
type(c_ptr) :: primme integer :: ierr ... integer(c_int64_t) :: n = 100 ierr = primme_set_member(primme, PRIMME_n, n) ierr = primme_set_member(primme, PRIMME_correctionParams_precondition, 1_c_int64_t) real(c_double) :: tol = 1.0D-12 ierr = primme_set_member(primme, PRIMME_eps, tol) integer(c_int64_t), parameter :: numTargetShifts = 2 real(c_double) :: TargetShifts(numTargetShifts) = (/3.0D0, 5.1D0/) ierr = primme_set_member(primme, PRIMME_numTargetShifts, numTargetShifts) ierr = primme_set_member(primme, PRIMME_targetShifts, TargetShifts) ierr = primme_set_member(primme, PRIMME_target, primme_closest_abs) procedure(primme_eigs_matvec) :: MV, ApplyPrecon ierr = primme_set_member(primme, PRIMME_matrixMatvec, MV) ierr = primme_set_member(primme, PRIMME_applyPreconditioner, c_funloc(ApplyPrecon))
primme_get_member¶
-
function
primme_get_member
(primme, label, value)¶ Get the value in some field of the parameter structure.
- Parameters
primme [c_ptr,in] :: parameters structure created by
primme_params_create()
.label [integer,in] :: field where to get value. One of the detailed in function
primme_set_member()
.value [out] :: value of the field. The allowed types are c_int64, c_double, and c_ptr.
- Return
primme_get_member [c_int] :: nonzero value if the call is not successful.
Examples:
type(c_ptr) :: primme integer :: ierr ... integer(c_int64_t) :: n ierr = primme_get_member(primme, PRIMME_n, n) real(c_double) :: aNorm ierr = primme_get_member(primme, PRIMME_aNorm, aNorm) real(c_double), pointer :: shifts(:) type(c_ptr) :: pshifts ierr = primme_get_member(primme, PRIMME_ShiftsForPreconditioner, pshifts) call c_f_pointer(pshifts, shifts, shape=[k])