FORTRAN 77 Library Interface¶
The next enumerations and functions are declared in primme_f77.h
.
-
type
ptr
¶ Fortran datatype with the same size as a pointer. Use
integer*4
when compiling in 32 bits andinteger*8
in 64 bits.
primme_initialize_f77¶
-
subroutine
primme_initialize_f77
(primme)¶ Allocate and initialize a PRIMME parameters structure to the default values.
After calling
dprimme_f77()
(or a variant), callprimme_free_f77()
to release allocated resources by PRIMME.- Parameters
primme [ptr] :: (output) parameters structure.
primme_set_method_f77¶
-
subroutine
primme_set_method_f77
(method, primme, ierr)¶ Set PRIMME parameters to one of the preset configurations.
- Parameters
method [integer] ::
(input) 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 [ptr] :: (input) parameters structure.
ierr [integer] :: (output) if 0, successful; if negative, something went wrong.
primme_free_f77¶
sprimme_f77¶
-
subroutine
sprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a real symmetric standard or generalized eigenproblem.
All arrays should be hosted on CPU. The computations are performed on CPU (see
magma_sprimme_f77()
for using GPUs).- Parameters
evals (*) [real] :: (output) array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [real] :: (input/output) 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] :: (output) array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 2.0.
cprimme_f77¶
-
subroutine
cprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a Hermitian standard or generalized eigenproblem.
All arrays should be hosted on CPU. The computations are performed on CPU (see
magma_cprimme_f77()
for using GPUs).- Parameters
evals (*) [real] :: (output) array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex real] :: (input/output) 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] :: (output) array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: (input) parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 2.0.
dprimme_f77¶
-
subroutine
dprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a real symmetric standard or generalized eigenproblem.
All arrays should be hosted on CPU. The computations are performed on CPU (see
magma_dprimme_f77()
for using GPUs).- Parameters
evals (*) [double precision] :: (output) array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [double precision] :: (input/output) 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 (*) [double precision] :: (output) array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
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
.
zprimme_f77¶
-
subroutine
zprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a Hermitian standard or generalized eigenproblem.
All arrays should be hosted on CPU. The computations are performed on CPU (see
magma_zprimme_f77()
for using GPUs).- Parameters
evals (*) [double precision] :: (output) array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex double precision] :: (input/output) 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 (*) [double precision] :: (output) array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.
magma_sprimme_f77¶
-
subroutine
magma_sprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a real symmetric standard or generalized eigenproblem.
Most of the computations are performed on GPU (see
sprimme_f77()
for using only the CPU).- Parameters
evals (*) [real] :: (output) CPU array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [real] :: (input/output) 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] :: (output) CPU array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
magma_cprimme_f77¶
-
subroutine
magma_cprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a Hermitian standard or generalized eigenproblem.
Most of the computations are performed on GPU (see
cprimme_f77()
for using only the CPU).- Parameters
evals (*) [real] :: (output) CPU array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex real] :: (input/output) 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] :: (output) CPU array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: (input) parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
magma_dprimme_f77¶
-
subroutine
magma_dprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a real symmetric standard or generalized eigenproblem.
Most of the computations are performed on GPU (see
dprimme_f77()
for using only the CPU).- Parameters
evals (*) [double precision] :: (output) CPU array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [double precision] :: (input/output) 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 (*) [double precision] :: (output) CPU array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
magma_zprimme_f77¶
-
subroutine
magma_zprimme_f77
(evals, evecs, resNorms, primme, ierr)¶ Solve a Hermitian standard or generalized eigenproblem.
Most of the computations are performed on GPU (see
zprimme_f77()
for using only the CPU).- Parameters
evals (*) [double precision] :: (output) CPU array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex double precision] :: (input/output) 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 (*) [double precision] :: (output) CPU array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: (input) parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
cprimme_normal_f77¶
-
subroutine
cprimme_normal_f77
(evals, evecs, resNorms, primme, ierr)¶ 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_cprimme_normal_f77()
for using GPUs).- Parameters
evals (*) [real] :: (output) array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex real] :: (input/output) 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] :: (output) array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: (input) parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
zprimme_normal_f77¶
-
subroutine
zprimme_normal_f77
(evals, evecs, resNorms, primme, ierr)¶ 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_f77()
for using GPUs).- Parameters
evals (*) [double precision] :: (output) array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex double precision] :: (input/output) 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 (*) [double precision] :: (output) array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
magma_cprimme_normal_f77¶
-
subroutine
magma_cprimme_normal_f77
(evals, evecs, resNorms, primme, ierr)¶ 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
cprimme_normal_f77()
for using only the CPU).- Parameters
evals (*) [real] :: (output) CPU array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex real] :: (input/output) 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] :: (output) CPU array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: (input) parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
magma_zprimme_normal_f77¶
-
subroutine
magma_zprimme_normal_f77
(evals, evecs, resNorms, primme, ierr)¶ 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
zprimme_normal_f77()
for using only the CPU).- Parameters
evals (*) [double precision] :: (output) CPU array at least of size
numEvals
to store the computed eigenvalues; all parallel calls return the same value in this array.evecs (*) [complex double precision] :: (input/output) 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 (*) [double precision] :: (output) CPU array at least of size
numEvals
to store the residual norms of the computed eigenpairs; all parallel calls return the same value in this array.primme [ptr] :: (input) parameters structure.
ierr [integer] :: (output) error indicator; see Error Codes.
Further descriptions of evals, evecs, and resNorms on notes in subroutine
dprimme_f77()
.New in version 3.0.
primme_set_member_f77¶
-
subroutine
primme_set_member_f77
(primme, label, value)¶ Set a value in some field of the parameter structure.
- Parameters
primme [ptr] :: (input) parameters structure.
label [integer] ::
field where to set value. One of:
PRIMME_stats_flopsDense
PRIMME_stats_timeDense
PRIMME_stats_estimateMaxEVal
PRIMME_stats_estimateBNorm
PRIMME_stats_estimateInvBNorm
value ::
(input) value to set.
If the type of the option is integer (
int
,PRIMME_INT
,size_t
), the type ofvalue
should be as long asPRIMME_INT
, which isinteger*8
by default.
Note
Don’t use this subroutine inside PRIMME’s callback functions, e.g.,
matrixMatvec
orapplyPreconditioner
, or in functions called by these functions.
primmetop_get_member_f77¶
-
subroutine
primmetop_get_member_f77
(primme, label, value)¶ Get the value in some field of the parameter structure.
- Parameters
primme [ptr] :: (input) parameters structure.
label [integer] :: (input) field where to get value. One of the detailed in subroutine
primmetop_set_member_f77()
.value ::
(output) value of the field.
If the type of the option is integer (
int
,PRIMME_INT
,size_t
), the type ofvalue
should be as long asPRIMME_INT
, which isinteger*8
by default.
Note
Don’t use this subroutine inside PRIMME’s callback functions, e.g.,
matrixMatvec
orapplyPreconditioner
, or in functions called by these functions. In those cases useprimme_get_member_f77()
.Note
When
label
is one ofPRIMME_matrixMatvec
,PRIMME_applyPreconditioner
,PRIMME_commInfo
,PRIMME_matrix
andPRIMME_preconditioner
, the returnedvalue
is a C pointer (void*
). Use Fortran pointer or other extensions to deal with it. For instance:use iso_c_binding MPI_Comm comm comm = MPI_COMM_WORLD call primme_set_member_f77(primme, PRIMME_commInfo, comm) ... subroutine par_GlobalSumDouble(x,y,k,primme) use iso_c_binding implicit none ... MPI_Comm, pointer :: comm type(c_ptr) :: pcomm call primme_get_member_f77(primme, PRIMME_commInfo, pcomm) call c_f_pointer(pcomm, comm) call MPI_Allreduce(x,y,k,MPI_DOUBLE,MPI_SUM,comm,ierr)
Most users would not need to retrieve these pointers in their programs.
primmetop_get_prec_shift_f77¶
-
subroutine
primmetop_get_prec_shift_f77
(primme, index, value)¶ Get the value in some position of the array
ShiftsForPreconditioner
.- Parameters
primme [ptr] :: (input) parameters structure.
index [integer] :: (input) position of the array; the first position is 1.
value :: (output) value of the array at that position.
primme_get_member_f77¶
-
subroutine
primme_get_member_f77
(primme, label, value)¶ Get the value in some field of the parameter structure.
- Parameters
primme [ptr] :: (input) parameters structure.
label [integer] :: (input) field where to get value. One of the detailed in subroutine
primmetop_set_member_f77()
.value ::
(output) value of the field.
If the type of the option is integer (
int
,PRIMME_INT
,size_t
), the type ofvalue
should be as long asPRIMME_INT
, which isinteger*8
by default.
Note
Use this subroutine exclusively inside PRIMME’s callback functions, e.g.,
matrixMatvec
orapplyPreconditioner
, or in functions called by these functions. Otherwise, e.g., from the main program, use the subroutineprimmetop_get_member_f77()
.Note
When
label
is one ofPRIMME_matrixMatvec
,PRIMME_applyPreconditioner
,PRIMME_commInfo
,PRIMME_matrix
andPRIMME_preconditioner
, the returnedvalue
is a C pointer (void*
). Use Fortran pointer or other extensions to deal with it. For instance:use iso_c_binding MPI_Comm comm comm = MPI_COMM_WORLD call primme_set_member_f77(primme, PRIMME_commInfo, comm) ... subroutine par_GlobalSumDouble(x,y,k,primme) use iso_c_binding implicit none ... MPI_Comm, pointer :: comm type(c_ptr) :: pcomm call primme_get_member_f77(primme, PRIMME_commInfo, pcomm) call c_f_pointer(pcomm, comm) call MPI_Allreduce(x,y,k,MPI_DOUBLE,MPI_SUM,comm,ierr)
Most users would not need to retrieve these pointers in their programs.
primme_get_prec_shift_f77¶
-
subroutine
primme_get_prec_shift_f77
(primme, index, value)¶ Get the value in some position of the array
ShiftsForPreconditioner
.- Parameters
primme [ptr] :: (input) parameters structure.
index [integer] :: (input) position of the array; the first position is 1.
value :: (output) value of the array at that position.
Note
Use this subroutine exclusively inside the subroutine
matrixMatvec
,massMatrixMatvec
, orapplyPreconditioner
. Otherwise use the subroutineprimmetop_get_prec_shift_f77()
.