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) times numSvals 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.

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) times numSvals 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.

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) times numSvals 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.

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) times numSvals 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.

primme_svds_initialize_f77

primme_svds_initialize_f77(primme_svds)

Set PRIMME SVDS parameters structure to the default values.

Parameters:
  • primme_svds (ptr) – (output) parameters structure.

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:

primme_svds_display_params_f77

primme_svds_display_params_f77(primme_svds)

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

Parameters:
  • primme_svds (ptr) – (input) parameters structure.

primme_svds_free_f77

primme_svds_free_f77(primme_svds)

Free memory allocated by PRIMME SVDS and delete all values set.

Parameters:
  • primme_svds (ptr) – (input/output) parameters structure.

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:

Note

Don’t use this function inside PRIMME SVDS’s callback functions, e.g., matrixMatvec or applyPreconditioner, 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 or applyPreconditioner, or in functions called by these functions. In those cases use primme_svds_get_member_f77().

Note

When label is one of PRIMME_SVDS_matrixMatvec, PRIMME_SVDS_applyPreconditioner, PRIMME_SVDS_commInfo, PRIMME_SVDS_intWork, PRIMME_SVDS_realWork, PRIMME_SVDS_matrix and PRIMME_SVDS_preconditioner, the returned value 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 or applyPreconditioner, or in functions called by these functions. Otherwise, e.g., from the main program, use the function primme_svdstop_get_member_f77().

Note

When label is one of PRIMME_SVDS_matrixMatvec, PRIMME_SVDS_applyPreconditioner, PRIMME_SVDS_commInfo, PRIMME_SVDS_intWork, PRIMME_SVDS_realWork, PRIMME_SVDS_matrix and PRIMME_SVDS_preconditioner, the returned value 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.