Parameter Description¶
primme_svds_params¶
-
primme_svds_params
¶ Structure to set the problem matrix and the solver options.
-
PRIMME_INT
m
¶ Number of rows of the matrix.
Input/output:
primme_initialize()
sets this field to 0;this field is read bydprimme()
.
-
PRIMME_INT
n
¶ Number of columns of the matrix.
Input/output:
primme_initialize()
sets this field to 0;this field is read bydprimme()
.
-
void
(*matrixMatvec)
(void *x, PRIMME_INT ldx, void *y, PRIMME_INT ldy, int *blockSize, int *transpose, primme_svds_params *primme_svds, int *ierr)¶ Block matrix-multivector multiplication, \(y = A x\) if
transpose
is zero, and \(y = A^*x\) otherwise.Parameters: - x -- input array.
- ldx -- leading dimension of
x
. - y -- output array.
- ldy -- leading dimension of
y
. - blockSize -- number of columns in
x
andy
. - transpose -- if non-zero, the transpose A should be applied.
- primme_svds -- parameters structure.
- ierr -- output error code; if it is set to non-zero, the current call to PRIMME will stop.
If
transpose
is zero, thenx
andy
are arrays of dimensionsnLocal
xblockSize
andmLocal
xblockSize
respectively. Elsewhere they have dimensionsmLocal
xblockSize
andnLocal
xblockSize
. Both arrays are in column-major order (elements in the same column with consecutive row indices are consecutive in memory).The actual type of
x
andy
depends on which function is being calling. Fordprimme_svds()
, it isdouble
, forzprimme_svds()
it isPRIMME_COMPLEX_DOUBLE
, forsprimme_svds()
it isfloat
and forcprimme_svds()
it isPRIMME_COMPLEX_FLOAT
.Input/output:
primme_initialize()
sets this field to NULL;this field is read bydprimme_svds()
andzprimme_svds()
.Note
Integer arguments are passed by reference to make easier the interface to other languages (like Fortran).
-
void
(*applyPreconditioner)
(void *x, PRIMME_INT ldx, void *y, PRIMME_INT ldy, int *blockSize, int *mode, primme_svds_params *primme_svds, int *ierr)¶ Block preconditioner-multivector application, \(y = M^{-1}x\) for finding singular values close to \(\sigma\). Depending on
mode
, \(M\) is expected to be an approximation of the following operators:primme_svds_op_AtA
: \(M \approx A^*Ax - \sigma^2 I\),primme_svds_op_AAt
: \(M \approx AA^*x - \sigma^2 I\),primme_svds_op_augmented
: \(M \approx \left(\begin{array}{cc} 0 & A^* \\ A & 0 \end{array}\right) - \sigma I\).
Parameters: - x -- input array.
- ldx -- leading dimension of
x
. - y -- output array.
- ldy -- leading dimension of
y
. - blockSize -- number of columns in
x
andy
. - mode -- one of
primme_svds_op_AtA
,primme_svds_op_AAt
orprimme_svds_op_augmented
. - primme_svds -- parameters structure.
- ierr -- output error code; if it is set to non-zero, the current call to PRIMME will stop.
If
mode
isprimme_svds_op_AtA
, thenx
andy
are arrays of dimensionsnLocal
xblockSize
; if mode isprimme_svds_op_AAt
, they aremLocal
xblockSize
; and otherwise they are (mLocal
+nLocal
) xblockSize
. Both arrays are in column-major order (elements in the same column with consecutive row indices are consecutive in memory).The actual type of
x
andy
depends on which function is being calling. Fordprimme_svds()
, it isdouble
, forzprimme_svds()
it isPRIMME_COMPLEX_DOUBLE
, forsprimme_svds()
it isfloat
and forcprimme_svds()
it isPRIMME_COMPLEX_FLOAT
.Input/output:
primme_initialize()
sets this field to NULL;this field is read bydprimme_svds()
andzprimme_svds()
.
-
int
numProcs
¶ Number of processes calling
dprimme_svds()
orzprimme_svds()
in parallel.Input/output:
primme_initialize()
sets this field to 1;this field is read bydprimme()
andzprimme_svds()
.
-
int
procID
¶ The identity of the local process within a parallel execution calling
dprimme_svds()
orzprimme_svds()
. Only the process with id 0 prints information.Input/output:
primme_svds_initialize()
sets this field to 0;dprimme_svds()
sets this field to 0 ifnumProcs
is 1;this field is read bydprimme_svds()
andzprimme_svds()
.
-
PRIMME_INT
mLocal
¶ Number of local rows on this process. The value depends on how the matrix and preconditioner is distributed along the processes.
Input/output:
primme_svds_initialize()
sets this field to -1;this field is read bydprimme_svds()
andzprimme_svds()
.See also:
matrixMatvec
andapplyPreconditioner
.
-
PRIMME_INT
nLocal
¶ Number of local columns on this process. The value depends on how the matrix and preconditioner is distributed along the processes.
Input/output:
primme_svds_initialize()
sets this field to -1;this field is read bydprimme_svds()
andzprimme_svds()
.
-
void *
commInfo
¶ A pointer to whatever parallel environment structures needed. For example, with MPI, it could be a pointer to the MPI communicator. PRIMME does not use this. It is available for possible use in user functions defined in
matrixMatvec
,applyPreconditioner
andglobalSumReal
.Input/output:
primme_svds_initialize()
sets this field to NULL;
-
void
(*globalSumReal)
(double *sendBuf, double *recvBuf, int *count, primme_svds_params *primme_svds, int *ierr)¶ Global sum reduction function. No need to set for sequential programs.
Parameters: - sendBuf -- array of size
count
with the local input values. - recvBuf -- array of size
count
with the global output values so that the i-th element of recvBuf is the sum over all processes of the i-th element ofsendBuf
. - count -- array size of
sendBuf
andrecvBuf
. - primme_svds -- parameters structure.
- ierr -- output error code; if it is set to non-zero, the current call to PRIMME will stop.
The actual type of
sendBuf
andrecvBuf
depends on which function is being calling. Fordprimme_svds()
andzprimme_svds()
it isdouble
, and forsprimme_svds()
andcprimme_svds()
it isfloat
. Note thatcount
is the number of values of the actual type.Input/output:
primme_svds_initialize()
sets this field to an internal function;this field is read bydprimme_svds()
andzprimme_svds()
.When MPI is used, this can be a simply wrapper to MPI_Allreduce() as shown below:
void par_GlobalSumForDouble(void *sendBuf, void *recvBuf, int *count, primme_svds_params *primme_svds, int *ierr) { MPI_Comm communicator = *(MPI_Comm *) primme_svds->commInfo; if (MPI_Allreduce(sendBuf, recvBuf, *count, MPI_DOUBLE, MPI_SUM, communicator) == MPI_SUCCESS) { *ierr = 0; } else { *ierr = 1; } }
When calling
sprimme_svds()
andcprimme_svds()
replaceMPI_DOUBLE
by`MPI_FLOAT
.- sendBuf -- array of size
-
int
numSvals
¶ Number of singular triplets wanted.
Input/output:
primme_svds_initialize()
sets this field to 1;
-
primme_svds_target
target
¶ Which singular values to find:
primme_svds_smallest
- Smallest singular values;
targetShifts
is ignored. primme_svds_largest
- Largest singular values;
targetShifts
is ignored. primme_svds_closest_abs
- Closest in absolute value to the shifts in
targetShifts
.
Input/output:
primme_svds_initialize()
sets this field toprimme_svds_smallest
;this field is read bydprimme_svds()
andzprimme_svds()
.
-
int
numTargetShifts
¶ Size of the array
targetShifts
. Used only whentarget
isprimme_svds_closest_abs
. The default values is 0.Input/output:
primme_svds_initialize()
sets this field to 0;this field is read bydprimme_svds()
andzprimme_svds()
.
-
double *
targetShifts
¶ Array of shifts, at least of size
numTargetShifts
. Used only whentarget
isprimme_svds_closest_abs
.Singular values are computed in order so that the i-th singular value is the closest to the i-th shift. If
numTargetShifts
<numSvals
, the last shift given is used for all the remaining i's.Input/output:
primme_svds_initialize()
sets this field to NULL;this field is read bydprimme_svds()
andzprimme_svds()
.Note
Eventually this is used by
dprimme()
andzprimme()
. Please see considerations oftargetShifts
.
-
int
printLevel
¶ The level of message reporting from the code. All output is written in
outputFile
.One of:
0: silent.
1: print some error messages when these occur.
2: as in 1, and info about targeted singular triplets when they are marked as converged:
#Converged $1 sval[ $2 ]= $3 norm $4 Mvecs $5 Time $7 stage $10
or locked:
#Lock striplet[ $1 ]= $3 norm $4 Mvecs $5 Time $7 stage $10
3: as in 2, and info about targeted singular triplets every outer iteration:
OUT $6 conv $1 blk $8 MV $5 Sec $7 SV $3 |r| $4 stage $10
Also, if using
PRIMME_DYNAMIC
, show JDQMR/GD+k performance ratio and the current method in use.4: as in 3, and info about targeted singular triplets every inner iteration:
INN MV $5 Sec $7 Sval $3 Lin|r| $9 SV|r| $4 stage $10
5: as in 4, and verbose info about certain choices of the algorithm.
Output key:
$1: Number of converged triplets up to now.$2: The index of the triplet currently converged.$3: The singular value.$4: Its residual norm.$5: The current number of matrix-vector products.$6: The current number of outer iterations.$7: The current elapsed time.$8: Index within the block of the targeted triplet.$9: QMR norm of the linear system residual.$10: stage (1 or 2)In parallel programs, when
printLevel
is 0 to 4 onlyprocID
0 produces output. ForprintLevel
5 output can be produced in any of the parallel calls.Input/output:
primme_svds_initialize()
sets this field to 1;this field is read bydprimme_svds()
andzprimme_svds()
.Note
Convergence history for plotting may be produced simply by:
grep OUT outpufile | awk '{print $8" "$14}' > out grep INN outpufile | awk '{print $3" "$11}' > inn
Or in gnuplot:
plot 'out' w lp, 'inn' w lp
-
double
aNorm
¶ An estimate of the 2-norm of \(A\), which is used in the default convergence criterion (see
eps
).If
aNorm
is less than or equal to 0, the code uses the largest absolute Ritz value seen. On return,aNorm
is then replaced with that value.Input/output:
primme_svds_initialize()
sets this field to 0.0;this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
double
eps
¶ If
convTestFun
is NULL, a triplet \((u,\sigma,v)\) is marked as converged when \(\sqrt{\|A v - \sigma u\|^2 + \|A^* u - \sigma v\|^2}\) is less thaneps
*aNorm
, or close to the minimum tolerance that the selected method can achieve in the given machine precision. See Preset Methods.The default value is machine precision times \(10^4\).
Input/output:
primme_svds_initialize()
sets this field to 0.0;this field is read bydprimme_svds()
andzprimme_svds()
.
-
FILE *
outputFile
¶ Opened file to write down the output.
Input/output:
primme_svds_initialize()
sets this field to the standard output;
-
int
locking
¶ If set to 1, the underneath eigensolvers will use hard locking. See
locking
.Input/output:
primme_svds_initialize()
sets this field to -1;written byprimme_svds_set_method()
(see Preset Methods);this field is read bydprimme_svds()
andzprimme_svds()
.
-
int
initSize
¶ On input, the number of initial vector guesses provided in
svecs
argument indprimme_svds()
andzprimme_svds()
.On output,
initSize
holds the number of converged triplets. Withoutlocking
allnumSvals
approximations are insvecs
but only the firstinitSize
are converged.During execution, it holds the current number of converged triplets.
Input/output:
primme_svds_initialize()
sets this field to 0;this field is read and written bydprimme_svds()
andzprimme_svds()
.-
int
numOrthoConst
¶ Number of vectors to be used as external orthogonalization constraints. The left and the right vector constraints are provided as input of the
svecs
argument insprimme_svds()
or other variant, and must be orthonormal.PRIMME SVDS finds new triplets orthogonal to these constraints (equivalent to solving the problem \((I-UU^*)A(I-VV^*)\) where \(U\) and \(V\) are the given left and right constraint vectors). This is a handy feature if some singular triplets are already known, or for finding more triplets after a call to
dprimme_svds()
orzprimme_svds()
, possibly with different parameters (see an example inTEST/exsvd_zseq.c
).Input/output:
primme_svds_initialize()
sets this field to 0;this field is read bydprimme_svds()
andzprimme_svds()
.
-
int
-
int
maxBasisSize
¶ The maximum basis size allowed in the main iteration. This has memory implications.
Input/output:
primme_svds_initialize()
sets this field to 0;this field is read and written byprimme_svds_set_method()
(see Preset Methods);this field is read bydprimme_svds()
andzprimme_svds()
.
-
int
maxBlockSize
¶ The maximum block size the code will try to use.
The user should set this based on the architecture specifics of the target computer, as well as any a priori knowledge of multiplicities. The code does not require that
maxBlockSize
> 1 to find multiple triplets. For some methods, keeping to 1 yields the best overall performance.Input/output:
primme_svds_initialize()
sets this field to 1;this field is read and written byprimme_svds_set_method()
(see Preset Methods);this field is read bydprimme_svds()
andzprimme_svds()
.
-
PRIMME_INT
maxMatvecs
¶ Maximum number of matrix vector multiplications (approximately half the number of preconditioning operations) that the code is allowed to perform before it exits.
Input/output:
primme_svds_initialize()
sets this field toINT_MAX
;this field is read bydprimme_svds()
andzprimme_svds()
.
-
int
intWorkSize
¶ If
dprimme_svds()
orzprimme_svds()
is called with all arguments as NULL except forprimme_svds_params
then it returns immediately withintWorkSize
containing the size in bytes of the integer workspace that will be required by the parameters set.Otherwise if
intWorkSize
is not 0, it should be the size of the integer work array in bytes that the user provides inintWork
. IfintWorkSize
is 0, the code will allocate the required space, which can be freed later by callingprimme_svds_free()
.Input/output:
primme_svds_initialize()
sets this field to 0;this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
size_t
realWorkSize
¶ If
dprimme_svds()
orzprimme_svds()
is called with all arguments as NULL except forprimme_svds_params
then it returns immediately withrealWorkSize
containing the size in bytes of the real workspace that will be required by the parameters set.Otherwise if
realWorkSize
is not 0, it should be the size of the real work array in bytes that the user provides inrealWork
. IfrealWorkSize
is 0, the code will allocate the required space, which can be freed later by callingprimme_svds_free()
.Input/output:
primme_svds_initialize()
sets this field to 0;this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
int *
intWork
¶ Integer work array.
If NULL, the code will allocate its own workspace. If the provided space is not enough, the code will return the error code
-21
.Input/output:
primme_svds_initialize()
sets this field to NULL;this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
void *
realWork
¶ Real work array.
If NULL, the code will allocate its own workspace. If the provided space is not enough, the code will return the error code
-20
.Input/output:
primme_svds_initialize()
sets this field to NULL;this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
PRIMME_INT
iseed
¶ The
PRIMME_INT iseed[4]
is an array with the seeds needed by the LAPACK dlarnv and zlarnv.The default value is an array with values -1, -1, -1 and -1. In that case,
iseed
is set based on the value ofprocID
to avoid every parallel process generating the same sequence of pseudorandom numbers.Input/output:
primme_svds_initialize()
sets this field to[-1, -1, -1, -1]
;this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
void *
matrix
¶ This field may be used to pass any required information in the matrix-vector product
matrixMatvec
.Input/output:
primme_svds_initialize()
sets this field to NULL;
-
void *
preconditioner
¶ This field may be used to pass any required information in the preconditioner function
applyPreconditioner
.Input/output:
primme_svds_initialize()
sets this field to NULL;
-
int
precondition
¶ Set to 1 to use preconditioning. Make sure
applyPreconditioner
is not NULL then!Input/output:
primme_svds_initialize()
sets this field to 0;this field is read and written byprimme_svds_set_method()
(see Preset Methods);this field is read bydprimme_svds()
andzprimme_svds()
.
-
primme_svds_op_operator
method
¶ Select the equivalent eigenvalue problem that will be solved:
primme_svds_op_AtA
: \(A^*Ax = \sigma^2 x\),primme_svds_op_AAt
: \(AA^*x = \sigma^2 x\),primme_svds_op_augmented
: \(\left(\begin{array}{cc} 0 & A^* \\ A & 0 \end{array}\right) x = \sigma x\).
The options for this solver are stored in
primme
.Input/output:
primme_svds_initialize()
sets this field toprimme_svds_op_none
;this field is read and written byprimme_svds_set_method()
(see Preset Methods);this field is read bydprimme_svds()
andzprimme_svds()
.
-
primme_svds_op_operator
methodStage2
¶ Select the equivalent eigenvalue problem that will be solved to refine the solution. The allowed options are
primme_svds_op_none
to not refine the solution andprimme_svds_op_augmented
to refine the solution by solving the augmented problem with the current solution as the initial vectors. Seemethod
.The options for this solver are stored in
primmeStage2
.Input/output:
primme_svds_initialize()
sets this field toprimme_svds_op_none
;this field is read and written byprimme_svds_set_method()
(see Preset Methods);this field is read bydprimme_svds()
andzprimme_svds()
.
-
primme_params
primme
¶ Parameter structure storing the options for underneath eigensolver that will be called at the first stage. See
method
.Input/output:
primme_svds_initialize()
initialize this structure;this field is read and written byprimme_svds_set_method()
(see Preset Methods);this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
primme_params
primmeStage2
¶ Parameter structure storing the options for underneath eigensolver that will be called at the second stage. See
methodStage2
.Input/output:
primme_svds_initialize()
initialize this structure;this field is read and written byprimme_svds_set_method()
(see Preset Methods);this field is read and written bydprimme_svds()
andzprimme_svds()
.
-
void
(*convTestFun)
(double *sval, void *leftsvec, void *rightsvec, double *rNorm, int *isconv, primme_svds_params *primme_svds, int *ierr)¶ Function that evaluates if the approximate triplet has converged. If NULL, it is used the default convergence criteria (see
eps
).Parameters: - sval -- the approximate singular value to evaluate.
- leftsvec -- one dimensional array of size
mLocal
containing the approximate left singular vector; it can be NULL. - rightsvec -- one dimensional array of size
nLocal
containing the approximate right singular vector; it can be NULL. - resNorm -- the norm of the residual vector.
- isconv -- (output) the function sets zero if the pair is not converged and non zero otherwise.
- primme_svds -- parameters structure.
- ierr -- output error code; if it is set to non-zero, the current call to PRIMME will stop.
The actual type of
leftsvec
andrightsvec
depends on which function is being calling. Fordprimme_svds()
, it isdouble
, forzprimme_svds()
it isPRIMME_COMPLEX_DOUBLE
, forsprimme_svds()
it isfloat
and forcprimme_svds()
it isPRIMME_COMPLEX_FLOAT
.Warning
When solving the augmented problem (for the method
primme_svds_augmented
and at the second stage in the methodprimme_svds_hybrid
), the given residual vector normresNorm
is an approximation of the actual residual. Alsoleftsvec
andrightsvec
may not have length 1.Input/output:
svds_primme_initialize()
sets this field to NULL;this field is read and written bydprimme_svds()
.
-
void *
convtest
¶ This field may be used to pass any required information to the function
convTestFun
.Input/output:
primme_svds_initialize()
sets this field to NULL;
-
void
(*monitorFun)
(void *basisSvals, int *basisSize, int *basisFlags, int *iblock, int *blockSize, void *basisNorms, int *numConverged, void *lockedSvals, int *numLocked, int *lockedFlags, void *lockedNorms, int *inner_its, void *LSRes, primme_event *event, int *stage, primme_svds_params *primme_svds, int *ierr)¶ Convergence monitor. Used to customize how to report solver information during execution (stage, iteration number, matvecs, time, residual norms, targets, etc).
Parameters: - basisSvals -- array with approximate singular values of the basis.
- basisSize -- size of the arrays
basisSvals
,basisFlags
andbasisNorms
. - basisFlags -- state of every approximate triplet in the basis.
- iblock -- indices of the approximate triplet in the block.
- blockSize -- size of array
iblock
. - basisNorms -- array with residual norms of the triplets in the basis.
- numConverged -- number of triplets converged in the basis plus the number of the locked triplets (note that this value isn't monotonic).
- lockedSvals -- array with the locked triplets.
- numLocked -- size of the arrays
lockedSvals
,lockedFlags
andlockedNorms
. - lockedFlags -- state of each locked triplets.
- lockedNorms -- array with residual norms of the locked triplets.
- inner_its -- number of performed QMR iterations in the current correction equation.
- LSRes -- residual norm of the linear system at the current QMR iteration.
- event -- event reported.
- stage --
0
for first stage,1
for second stage. - primme_svds -- parameters structure; the counter in
stats
are updated with the current number of matrix-vector products, iterations, elapsed time, etc., since start. - ierr -- output error code; if it is set to non-zero, the current call to PRIMME will stop.
This function is called at the next events:
*event == primme_event_outer_iteration
: every outer iterations.It is provided
basisSvals
,basisSize
,basisFlags
,iblock
andblockSize
.basisNorms[iblock[i]]
has the residual norms for the selected triplets in the block. PRIMME avoids computing the residual of soft-locked triplets,basisNorms[i]
fori<iblock[0]
. So those values may correspond to previous iterations. The valuesbasisNorms[i]
fori>iblock[blockSize-1]
are not valid.If
locking
is enabled,lockedSvals
,numLocked
,lockedFlags
andlockedNorms
are also provided.inner_its
andLSRes
are not provided.*event == primme_event_inner_iteration
: every QMR iteration.basisSvals[0]
andbasisNorms[0]
provides the approximate singular value and the residual norm of the triplet which is improved in the current correction equation. IfconvTest
isprimme_adaptive
orprimme_adaptive_ETolerance
,basisSvals[0]
, andbasisNorms[0]
are updated every QMR iteration.inner_its
andLSRes
are also provided.lockedSvals
,numLocked
,lockedFlags
, andlockedNorms
may not be provided.*event == primme_event_converged
: a new triplet in the basis passed the convergence criterioniblock[0]
is the index of the newly converged triplet in the basis which will be locked or soft locked. The following are provided:basisSvals
,basisSize
,basisFlags
andblockSize[0]==1
.lockedSvals
,numLocked
,lockedFlags
andlockedNorms
may not be provided.inner_its
andLSRes
are not provided.*event == primme_event_locked
: a new triplet added to the locked singular vectors.lockedSvals
,numLocked
,lockedFlags
andlockedNorms
are provided. The last element oflockedSvals
,lockedFlags
andlockedNorms
corresponds to the recent locked triplet.basisSvals
,numConverged
,basisFlags
andbasisNorms
may not be provided.inner_its
andLSRes
are not be provided.
The values of
basisFlags
andlockedFlags
are:0
: unconverged.1
: internal use; only inbasisFlags
.2
: passed convergence test (seeeps
).3
: converged because the solver may not be able to reduce the residual norm further.
Input/output:
primme_initialize()
sets this field to NULL;dprimme_svds()
sets this field to an internal function if it is NULL;this field is read bydprimme_svds()
andzprimme_svds()
.
-
void *
monitor
¶ This field may be used to pass any required information to the function
monitorFun
.Input/output:
primme_svds_initialize()
sets this field to NULL;
-
PRIMME_INT
stats.numOuterIterations
¶ Hold the number of outer iterations.
Input/output:
primme_svds_initialize()
sets this field to 0;written bydprimme_svds()
andzprimme_svds()
.
-
PRIMME_INT
stats.numRestarts
¶ Hold the number of restarts.
Input/output:
primme_svds_initialize()
sets this field to 0;written bydprimme_svds()
andzprimme_svds()
.
-
PRIMME_INT
stats.numMatvecs
¶ Hold how many vectors the operator in
matrixMatvec
has been applied on.Input/output:
primme_svds_initialize()
sets this field to 0;written bydprimme_svds()
andzprimme_svds()
.
-
PRIMME_INT
stats.numPreconds
¶ Hold how many vectors the operator in
applyPreconditioner
has been applied on.Input/output:
primme_svds_initialize()
sets this field to 0;written bydprimme_svds()
andzprimme_svds()
.
-
double
stats.elapsedTime
¶ Hold the wall clock time spent by the call to
dprimme_svds()
orzprimme_svds()
.Input/output:
primme_svds_initialize()
sets this field to 0;written bydprimme_svds()
andzprimme_svds()
.
-
PRIMME_INT
Preset Methods¶
-
primme_svds_preset_method
¶ -
primme_svds_default
¶ Set as
primme_svds_hybrid
.
-
primme_svds_normalequations
¶ Solve the equivalent eigenvalue problem \(A^*A V = \Sigma^2 V\) and computes \(U\) by normalizing the vectors \(AV\). If
m
is smaller thann
, \(AA^*\) is solved instead.With
primme_svds_normalequations
primme_svds_set_method()
setsmethod
toprimme_svds_op_AtA
ifm
is larger or equal thann
, and toprimme_svds_op_AAt
otherwise; andmethodStage2
is set toprimme_svds_op_none
.The minimum residual norm that this method can achieve is \(\|A\|\epsilon\sigma^{-1}\), where \(\epsilon\) is the machine precision and \(\sigma\) the required singular value.
-
primme_svds_augmented
¶ Solve the equivalent eigenvalue problem \(\left(\begin{array}{cc} 0 & A^* \\ A & 0 \end{array}\right) X = \sigma X\) with \(X = \left(\begin{array}{cc}V\\U\end{array}\right)\).
With
primme_svds_augmented
primme_svds_set_method()
setsmethod
toprimme_svds_op_augmented
andmethodStage2
toprimme_svds_op_none
.The minimum residual norm that this method can achieve is \(\|A\|\epsilon\), where \(\epsilon\) is the machine precision. However it may not return triplets with singular values smaller than \(\|A\|\epsilon\).
-
primme_svds_hybrid
¶ First solve the equivalent normal equations (see
primme_svds_normalequations
) and then refine the solution solving the augmented problem (seeprimme_svds_augmented
).With
primme_svds_normalequations
primme_svds_set_method()
setsmethod
toprimme_svds_op_AtA
ifm
is larger or equal thann
, and toprimme_svds_op_AAt
otherwise; andmethodStage2
is set toprimme_svds_op_augmented
.The minimum residual norm that this method can achieve is \(\|A\|\epsilon\), where \(\epsilon\) is the machine precision. However it may not return triplets with singular values smaller than \(\|A\|\epsilon\) if
eps
is smaller than \(\|A\|\epsilon\sigma^{-1}\).
-
Error Codes¶
The functions dprimme_svds()
and zprimme_svds()
return one of the next values:
- 0: success,
- 1: reported only amount of required memory,
- -1: failed in allocating int or real workspace,
- -2: malloc failed in allocating a permutation integer array,
- -3: main_iter() encountered problem; the calling stack of the functions where the error occurred was printed in 'stderr',
- -4:
primme_svds
is NULL, - -5: Wrong value for
m
orn
ormLocal
ornLocal
, - -6: Wrong value for
numProcs
, - -7:
matrixMatvec
is not set, - -8:
applyPreconditioner
is not set butprecondition
== 1 , - -9:
numProcs
>1 butglobalSumReal
is not set, - -10: Wrong value for
numSvals
, it's larger than min(m
,n
), - -11: Wrong value for
numSvals
, it's smaller than 1, - -13: Wrong value for
target
, - -14: Wrong value for
method
, - -15: Not supported combination of method and
methodStage2
, - -16: Wrong value for
printLevel
, - -17:
svals
is not set, - -18:
svecs
is not set, - -19:
resNorms
is not set - -20: not enough memory for
realWork
- -21: not enough memory for
intWork
- -100 up to -199: eigensolver error from first stage; see the value plus 100 in Error Codes.
- -200 up to -299: eigensolver error from second stage; see the value plus 200 in Error Codes.