FORTRAN Library Interface¶
The next enumerations and functions are declared in primme_svds_f77.h
.
sprimme_svds_f77¶
-
sprimme_svds_f77
(svals, svecs, resNorms, primme_svds)¶ Solve a real singular value problem using single precision.
Parameters: - svals(*) (real) -- (output) array at least of size
numSvals
to store the computed singular values; all processes in a parallel run return this local array with the same values. - resNorms(*) (real) -- array at least of size
numSvals
to store the residual norms of the computed triplets; all processes in parallel run return this local array with the same values. - svecs(*) (real) -- array at least of size (
mLocal
+nLocal
) timesnumSvals
to store columnwise the (local part of the) computed left singular vectors and the right singular vectors. - primme_svds (ptr) -- parameters structure.
Returns: error indicator; see Error Codes.
- svals(*) (real) -- (output) array at least of size
cprimme_svds_f77¶
-
cprimme_svds_f77
(svals, svecs, resNorms, primme_svds)¶ Solve a complex singular value problem using single precision.
Parameters: - svals(*) (real) -- (output) array at least of size
numSvals
to store the computed singular values; all processes in a parallel run return this local array with the same values. - resNorms(*) (real) -- array at least of size
numSvals
to store the residual norms of the computed triplets; all processes in parallel run return this local array with the same values. - svecs(*) (complex) -- array at least of size (
mLocal
+nLocal
) timesnumSvals
to store columnwise the (local part of the) computed left singular vectors and the right singular vectors. - primme_svds (ptr) -- parameters structure.
Returns: error indicator; see Error Codes.
- svals(*) (real) -- (output) array at least of size
dprimme_svds_f77¶
-
dprimme_svds_f77
(svals, svecs, resNorms, primme_svds)¶ Solve a real singular value problem using double precision.
Parameters: - svals(*) (double precision) -- (output) array at least of size
numSvals
to store the computed singular values; all processes in a parallel run return this local array with the same values. - resNorms(*) (double precision) -- array at least of size
numSvals
to store the residual norms of the computed triplets; all processes in parallel run return this local array with the same values. - svecs(*) (double precision) -- array at least of size (
mLocal
+nLocal
) timesnumSvals
to store columnwise the (local part of the) computed left singular vectors and the right singular vectors. - primme_svds (ptr) -- parameters structure.
Returns: error indicator; see Error Codes.
- svals(*) (double precision) -- (output) array at least of size
zprimme_svds_f77¶
-
zprimme_svds_f77
(svals, svecs, resNorms, primme_svds)¶ Solve a complex singular value problem using double precision.
Parameters: - svals(*) (double precision) -- (output) array at least of size
numSvals
to store the computed singular values; all processes in a parallel run return this local array with the same values. - resNorms(*) (double precision) -- array at least of size
numSvals
to store the residual norms of the computed triplets; all processes in parallel run return this local array with the same values. - svecs(*) (complex*16) -- array at least of size (
mLocal
+nLocal
) timesnumSvals
to store columnwise the (local part of the) computed left singular vectors and the right singular vectors. - primme_svds (ptr) -- parameters structure.
Returns: error indicator; see Error Codes.
- svals(*) (double precision) -- (output) array at least of size
primme_svds_initialize_f77¶
primme_svds_set_method_f77¶
-
primme_svds_set_method_f77
(method, methodStage1, methodStage2, primme_svds, ierr)¶ Set PRIMME SVDS parameters to one of the preset configurations.
Parameters: - method (integer) --
(input) preset configuration to compute the singular triplets; one of
PRIMME_SVDS_default
, currently set asPRIMME_SVDS_hybrid
.PRIMME_SVDS_normalequations
, compute the eigenvectors of \(A^*A\) or \(A A^*\).PRIMME_SVDS_augmented
, compute the eigenvectors of the augmented matrix, \(\left(\begin{array}{cc} 0 & A^* \\ A & 0 \end{array}\right)\).PRIMME_SVDS_hybrid
, start withPRIMME_SVDS_normalequations
; use the resulting approximate singular vectors as initial vectors forPRIMME_SVDS_augmented
if the required accuracy was not achieved.
- methodStage1 (primme_preset_method) -- (input) preset method to compute the eigenpairs at the first stage; see available values at
primme_set_method_f77()
. - methodStage2 (primme_preset_method) -- (input) preset method to compute the eigenpairs with
the second stage of
PRIMME_SVDS_hybrid
; see available values atprimme_set_method_f77()
. - primme_svds (ptr) -- (input/output) parameters structure.
- ierr (integer) -- (output) if 0, successful; if negative, something went wrong.
- method (integer) --
primme_svds_display_params_f77¶
-
primme_svds_display_params_f77
(primme_svds)¶ Display all printable settings of
primme_svds
into the file descriptoroutputFile
.Parameters: - primme_svds (ptr) -- (input) parameters structure.
primme_svds_free_f77¶
primme_svds_set_member_f77¶
-
primme_svds_set_member_f77
(primme_svds, label, value)¶ Set a value in some field of the parameter structure.
Parameters: - primme_svds (ptr) -- (input) parameters structure.
- label (integer) --
field where to set value. One of:
PRIMME_SVDS_stats_numOuterIterations
PRIMME_SVDS_stats_numRestarts
PRIMME_SVDS_stats_numMatvecs
PRIMME_SVDS_stats_numPreconds
PRIMME_SVDS_stats_elapsedTime
- value -- (input) value to set.
Note
Don't use this function inside PRIMME SVDS's callback functions, e.g.,
matrixMatvec
orapplyPreconditioner
, or in functions called by these functions.
primme_svdstop_get_member_f77¶
-
primme_svdstop_get_member_f77
(primme_svds, label, value)¶ Get the value in some field of the parameter structure.
Parameters: - primme_svds (ptr) -- (input) parameters structure.
- label (integer) -- (input) field where to get value. One of
the detailed in function
primmesvds_top_set_member_f77()
. - value -- (output) value of the field.
Note
Don't use this function inside PRIMME SVDS's callback functions, e.g.,
matrixMatvec
orapplyPreconditioner
, or in functions called by these functions. In those cases useprimme_svds_get_member_f77()
.Note
When
label
is one ofPRIMME_SVDS_matrixMatvec
,PRIMME_SVDS_applyPreconditioner
,PRIMME_SVDS_commInfo
,PRIMME_SVDS_intWork
,PRIMME_SVDS_realWork
,PRIMME_SVDS_matrix
andPRIMME_SVDS_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_svds_set_member_f77(primme_svds, PRIMME_SVDS_commInfo, comm) ... subroutine par_GlobalSumDouble(x,y,k,primme_svds) use iso_c_binding implicit none ... MPI_Comm, pointer :: comm type(c_ptr) :: pcomm call primme_svds_get_member_f77(primme_svds, PRIMME_SVDS_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_svds_get_member_f77¶
-
primme_svds_get_member_f77
(primme_svds, label, value)¶ Get the value in some field of the parameter structure.
Parameters: - primme_svds (ptr) -- (input) parameters structure.
- label (integer) -- (input) field where to get value. One of
the detailed in function
primme_svdstop_set_member_f77()
. - value -- (output) value of the field.
Note
Use this function exclusively inside PRIMME SVDS's callback functions, e.g.,
matrixMatvec
orapplyPreconditioner
, or in functions called by these functions. Otherwise, e.g., from the main program, use the functionprimme_svdstop_get_member_f77()
.Note
When
label
is one ofPRIMME_SVDS_matrixMatvec
,PRIMME_SVDS_applyPreconditioner
,PRIMME_SVDS_commInfo
,PRIMME_SVDS_intWork
,PRIMME_SVDS_realWork
,PRIMME_SVDS_matrix
andPRIMME_SVDS_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_svds_set_member_f77(primme_svds, PRIMME_SVDS_commInfo, comm) ... subroutine par_GlobalSumDouble(x,y,k,primme_svds) use iso_c_binding implicit none ... MPI_Comm, pointer :: comm type(c_ptr) :: pcomm call primme_svds_get_member_f77(primme_svds, PRIMME_SVDS_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.