Package slepc4py :: Module SLEPc :: Class SVD
[hide private]
[frames] | no frames]

Class SVD



SVD

Nested Classes [hide private]
  ConvergedReason
SVD convergence reasons
  ErrorType
SVD error type to assess accuracy of computed solutions
  Type
SVD types
  Which
SVD desired piece of spectrum
Instance Methods [hide private]
a new object with type S, a subtype of T

__new__(S, ...)
 
appendOptionsPrefix(self, prefix)
Appends to the prefix used for searching for all SVD options in the database.
 
cancelMonitor(self)
Clears all monitors for an SVD object.
 
computeError(self, int i, etype=None)
Computes the error (based on the residual norm) associated with the i-th singular triplet.
 
create(self, comm=None)
Creates the SVD object.
 
destroy(self)
Destroys the SVD object.
 
errorView(self, etype=None, Viewer viewer=None)
Displays the errors associated with the computed solution (as well as the eigenvalues).
 
getBV(self)
Obtain the basis vectors objects associated to the SVD object.
 
getConverged(self)
Gets the number of converged singular triplets.
 
getConvergedReason(self)
Gets the reason why the `solve()` iteration was stopped.
 
getCrossEPS(self)
Retrieve the eigensolver object (`EPS`) associated to the singular value solver.
 
getCyclicEPS(self)
Retrieve the eigensolver object (`EPS`) associated to the singular value solver.
 
getCyclicExplicitMatrix(self)
Returns the flag indicating if ``H(A) = [ 0 A ; A^T 0 ]`` is built explicitly.
 
getDimensions(self)
Gets the number of singular values to compute and the dimension of the subspace.
 
getImplicitTranspose(self)
Gets the mode used to handle the transpose of the matrix associated with the singular value problem.
 
getIterationNumber(self)
Gets the current iteration number.
 
getOperator(self)
Gets the matrix associated with the singular value problem.
 
getOptionsPrefix(self)
Gets the prefix used for searching for all SVD options in the database.
 
getSingularTriplet(self, int i, Vec U=None, Vec V=None)
Gets the i-th triplet of the singular value decomposition as computed by `solve()`.
 
getTolerances(self)
Gets the tolerance and maximum iteration count used by the default SVD convergence tests.
 
getType(self)
Gets the SVD type of this object.
 
getValue(self, int i)
Gets the i-th singular value as computed by `solve()`.
 
getVectors(self, int i, Vec U, Vec V)
Gets the i-th left and right singular vectors as computed by `solve()`.
 
getWhichSingularTriplets(self)
Returns which singular triplets are to be sought.
 
reset(self)
Resets the SVD object.
 
setBV(self, BV V, BV U=None)
Associates basis vectors objects to the SVD solver.
 
setCrossEPS(self, EPS eps)
Associate an eigensolver object (`EPS`) to the singular value solver.
 
setCyclicEPS(self, EPS eps)
Associate an eigensolver object (`EPS`) to the singular value solver.
 
setCyclicExplicitMatrix(self, flag=True)
Indicate if the eigensolver operator ``H(A) = [ 0 A ; A^T 0 ]`` must be computed explicitly.
 
setDimensions(self, nsv=None, ncv=None, mpd=None)
Sets the number of singular values to compute and the dimension of the subspace.
 
setFromOptions(self)
Sets SVD options from the options database.
 
setImplicitTranspose(self, mode)
Indicates how to handle the transpose of the matrix associated with the singular value problem.
 
setInitialSpaces(self, spaceright=None, spaceleft=None)
Sets the initial spaces from which the SVD solver starts to iterate.
 
setLanczosOneSide(self, flag=True)
Indicate if the variant of the Lanczos method to be used is one-sided or two-sided.
 
setOperator(self, Mat A)
Sets the matrix associated with the singular value problem.
 
setOptionsPrefix(self, prefix)
Sets the prefix used for searching for all SVD options in the database.
 
setTRLanczosOneSide(self, flag=True)
Indicate if the variant of the thick-restart Lanczos method to be used is one-sided or two-sided.
 
setTolerances(self, tol=None, max_it=None)
Sets the tolerance and maximum iteration count used by the default SVD convergence tests.
 
setType(self, svd_type)
Selects the particular solver to be used in the SVD object.
 
setUp(self)
Sets up all the internal data structures necessary for the execution of the singular value solver.
 
setWhichSingularTriplets(self, which)
Specifies which singular triplets are to be sought.
 
solve(self)
Solves the singular value problem.
 
view(self, Viewer viewer=None)
Prints the SVD data structure.

Inherited from petsc4py.PETSc.Object: __copy__, __deepcopy__, __eq__, __ge__, __gt__, __le__, __lt__, __ne__, __nonzero__, compose, decRef, getAttr, getClassId, getClassName, getComm, getDict, getName, getRefCount, getTabLevel, incRef, incrementTabLevel, query, setAttr, setName, setTabLevel, stateIncrease, viewFromOptions

Inherited from object: __delattr__, __format__, __getattribute__, __hash__, __init__, __reduce__, __reduce_ex__, __repr__, __setattr__, __sizeof__, __str__, __subclasshook__

Properties [hide private]
  bv
  max_it
  tol
  transpose_mode
  which

Inherited from petsc4py.PETSc.Object: classid, comm, fortran, handle, klass, name, prefix, refcount, type

Inherited from object: __class__

Method Details [hide private]

__new__(S, ...)

 


Returns:
a new object with type S, a subtype of T

Overrides: object.__new__

appendOptionsPrefix(self, prefix)

 
Appends to the prefix used for searching for all SVD options
in the database.

Parameters
----------
prefix: string
        The prefix string to prepend to all SVD option requests.

computeError(self, int i, etype=None)

 
Computes the error (based on the residual norm) associated with the i-th
singular triplet.

Parameters
----------
i: int
   Index of the solution to be considered.
etype: `SVD.ErrorType` enumerate
   The error type to compute.

Returns
-------
e: real
   The relative error bound, computed in various ways from the residual norm
   ``sqrt(n1^2+n2^2)`` where ``n1 = ||A*v-sigma*u||_2``,
   ``n2 = ||A^T*u-sigma*v||_2``, ``sigma`` is the singular
   value, ``u`` and ``v`` are the left and right singular
   vectors.

Notes
-----
The index ``i`` should be a value between ``0`` and
``nconv-1`` (see `getConverged()`).

create(self, comm=None)

 
Creates the SVD object.

Parameters
----------
comm: Comm, optional
      MPI communicator; if not provided, it defaults to all
      processes.

destroy(self)

 
Destroys the SVD object.

Overrides: petsc4py.PETSc.Object.destroy

errorView(self, etype=None, Viewer viewer=None)

 
Displays the errors associated with the computed solution
(as well as the eigenvalues).

Parameters
----------
etype: `SVD.ErrorType` enumerate, optional
   The error type to compute.
viewer: Viewer, optional.
        Visualization context; if not provided, the standard
        output is used.

Notes
-----
By default, this function checks the error of all eigenpairs and prints
the eigenvalues if all of them are below the requested tolerance.
If the viewer has format ``ASCII_INFO_DETAIL`` then a table with
eigenvalues and corresponding errors is printed.

getBV(self)

 
Obtain the basis vectors objects associated to the SVD object.

Returns
-------
V: BV
    The basis vectors context for right singular vectors.
U: BV
    The basis vectors context for left singular vectors.

getConverged(self)

 
Gets the number of converged singular triplets.

Returns
-------
nconv: int
       Number of converged singular triplets.

Notes
-----
This function should be called after `solve()` has finished.

getConvergedReason(self)

 
Gets the reason why the `solve()` iteration was stopped.

Returns
-------
reason: `SVD.ConvergedReason` enumerate
        Negative value indicates diverged, positive value
        converged.

getCrossEPS(self)

 
Retrieve the eigensolver object (`EPS`) associated to the
singular value solver.

Returns
-------
eps: EPS
     The eigensolver object.

getCyclicEPS(self)

 
Retrieve the eigensolver object (`EPS`) associated to the
singular value solver.

Returns
-------
eps: EPS
     The eigensolver object.

getCyclicExplicitMatrix(self)

 
Returns the flag indicating if ``H(A) = [ 0 A ; A^T 0 ]`` is
built explicitly.

Returns
-------
flag: boolean
      True if ``H(A)`` is built explicitly.

getDimensions(self)

 
Gets the number of singular values to compute and the
dimension of the subspace.

Returns
-------
nsv: int
     Number of singular values to compute.
ncv: int
     Maximum dimension of the subspace to be used by the
     solver.
mpd: int
     Maximum dimension allowed for the projected problem.

getImplicitTranspose(self)

 
Gets the mode used to handle the transpose of the matrix
associated with the singular value problem.

Returns
-------
impl: boolean
      How to handle the transpose (implicitly or not).

getIterationNumber(self)

 
Gets the current iteration number. If the call to `solve()` is
complete, then it returns the number of iterations carried out
by the solution method.

Returns
-------
its: int
     Iteration number.

getOperator(self)

 
Gets the matrix associated with the singular value problem.

Returns
-------
A: Mat
   The matrix associated with the singular value problem.

getOptionsPrefix(self)

 
Gets the prefix used for searching for all SVD options in the
database.

Returns
-------
prefix: string
        The prefix string set for this SVD object.

Overrides: petsc4py.PETSc.Object.getOptionsPrefix

getSingularTriplet(self, int i, Vec U=None, Vec V=None)

 
Gets the i-th triplet of the singular value decomposition as
computed by `solve()`. The solution consists of the singular
value and its left and right singular vectors.

Parameters
----------
i: int
   Index of the solution to be obtained.
U: Vec
   Placeholder for the returned left singular vector.
V: Vec
   Placeholder for the returned right singular vector.

Returns
-------
s: float
   The computed singular value.

Notes
-----
The index ``i`` should be a value between ``0`` and
``nconv-1`` (see `getConverged()`. Singular triplets are
indexed according to the ordering criterion established with
`setWhichSingularTriplets()`.

getTolerances(self)

 
Gets the tolerance and maximum iteration count used by the
default SVD convergence tests.

Returns
-------
tol: float
     The convergence tolerance.
max_it: int
     The maximum number of iterations

getType(self)

 
Gets the SVD type of this object.

Returns
-------
type: `SVD.Type` enumerate
      The solver currently being used.

Overrides: petsc4py.PETSc.Object.getType

getValue(self, int i)

 
Gets the i-th singular value as computed by `solve()`.

Parameters
----------
i: int
   Index of the solution to be obtained.

Returns
-------
s: float
   The computed singular value.

Notes
-----
The index ``i`` should be a value between ``0`` and
``nconv-1`` (see `getConverged()`. Singular triplets are
indexed according to the ordering criterion established with
`setWhichSingularTriplets()`.

getVectors(self, int i, Vec U, Vec V)

 
Gets the i-th left and right singular vectors as computed by
`solve()`.

Parameters
----------
i: int
   Index of the solution to be obtained.
U: Vec
   Placeholder for the returned left singular vector.
V: Vec
   Placeholder for the returned right singular vector.

Notes
-----
The index ``i`` should be a value between ``0`` and
``nconv-1`` (see `getConverged()`. Singular triplets are
indexed according to the ordering criterion established with
`setWhichSingularTriplets()`.

getWhichSingularTriplets(self)

 
Returns which singular triplets are to be sought.

Returns
-------
which: `SVD.Which` enumerate
       The singular values to be sought (either largest or
       smallest).

setBV(self, BV V, BV U=None)

 
Associates basis vectors objects to the SVD solver.

Parameters
----------
V: BV
    The basis vectors context for right singular vectors.
U: BV
    The basis vectors context for left singular vectors.

setCrossEPS(self, EPS eps)

 
Associate an eigensolver object (`EPS`) to the singular value
solver.

Parameters
----------
eps: EPS
     The eigensolver object.

setCyclicEPS(self, EPS eps)

 
Associate an eigensolver object (`EPS`) to the singular value
solver.

Parameters
----------
eps: EPS
     The eigensolver object.

setCyclicExplicitMatrix(self, flag=True)

 
Indicate if the eigensolver operator ``H(A) = [ 0 A ; A^T 0
]`` must be computed explicitly.

Parameters
----------
flag: boolean
      True if ``H(A)`` is built explicitly.

setDimensions(self, nsv=None, ncv=None, mpd=None)

 
Sets the number of singular values to compute and the
dimension of the subspace.

Parameters
----------
nsv: int, optional
     Number of singular values to compute.
ncv: int, optional
     Maximum dimension of the subspace to be used by the
     solver.
mpd: int, optional
     Maximum dimension allowed for the projected problem.

Notes
-----
Use `DECIDE` for `ncv` and `mpd` to assign a reasonably good
value, which is dependent on the solution method.

The parameters `ncv` and `mpd` are intimately related, so that
the user is advised to set one of them at most. Normal usage
is the following:

 - In cases where `nsv` is small, the user sets `ncv`
   (a reasonable default is 2 * `nsv`).
 - In cases where `nsv` is large, the user sets `mpd`.

The value of `ncv` should always be between `nsv` and (`nsv` +
`mpd`), typically `ncv` = `nsv` + `mpd`. If `nsv` is not too
large, `mpd` = `nsv` is a reasonable choice, otherwise a
smaller value should be used.

setFromOptions(self)

 
Sets SVD options from the options database. This routine must
be called before `setUp()` if the user is to be allowed to set
the solver type.

Notes
-----
To see all options, run your program with the ``-help``
option.

Overrides: petsc4py.PETSc.Object.setFromOptions

setImplicitTranspose(self, mode)

 
Indicates how to handle the transpose of the matrix
associated with the singular value problem.

Parameters
----------
impl: boolean
      How to handle the transpose (implicitly or not).

Notes
-----
By default, the transpose of the matrix is explicitly built
(if the matrix has defined the MatTranspose operation).

If this flag is set to true, the solver does not build the
transpose, but handles it implicitly via MatMultTranspose().

setInitialSpaces(self, spaceright=None, spaceleft=None)

 
Sets the initial spaces from which the SVD solver starts to
iterate.

Parameters
----------
spaceright: sequence of Vec
   The right initial space.
spaceleft: sequence of Vec
   The left initial space.

setLanczosOneSide(self, flag=True)

 
Indicate if the variant of the Lanczos method to be used is
one-sided or two-sided.

Parameters
----------
flag: boolean
      True if the method is one-sided.

Notes
-----
By default, a two-sided variant is selected, which is
sometimes slightly more robust. However, the one-sided variant
is faster because it avoids the orthogonalization associated
to left singular vectors. It also saves the memory required
for storing such vectors.

setOperator(self, Mat A)

 
Sets the matrix associated with the singular value problem.

Parameters
----------
A: Mat
   The matrix associated with the singular value problem.

setOptionsPrefix(self, prefix)

 
Sets the prefix used for searching for all SVD options in the
database.

Parameters
----------
prefix: string
        The prefix string to prepend to all SVD option
        requests.

Notes
-----
A hyphen (-) must NOT be given at the beginning of the prefix
name.  The first character of all runtime options is
AUTOMATICALLY the hyphen.

For example, to distinguish between the runtime options for
two different SVD contexts, one could call::

    S1.setOptionsPrefix("svd1_")
    S2.setOptionsPrefix("svd2_")

Overrides: petsc4py.PETSc.Object.setOptionsPrefix

setTRLanczosOneSide(self, flag=True)

 
Indicate if the variant of the thick-restart Lanczos method to
be used is one-sided or two-sided.

Parameters
----------
flag: boolean
      True if the method is one-sided.

Notes
-----
By default, a two-sided variant is selected, which is
sometimes slightly more robust. However, the one-sided variant
is faster because it avoids the orthogonalization associated
to left singular vectors.

setTolerances(self, tol=None, max_it=None)

 
Sets the tolerance and maximum iteration count used by the
default SVD convergence tests.

Parameters
----------
tol: float, optional
     The convergence tolerance.
max_it: int, optional
     The maximum number of iterations

Notes
-----
Use `DECIDE` for `max_it` to assign a reasonably good value,
which is dependent on the solution method.

setType(self, svd_type)

 
Selects the particular solver to be used in the SVD object.

Parameters
----------
svd_type: `SVD.Type` enumerate
          The solver to be used.

Notes
-----
See `SVD.Type` for available methods. The default is CROSS.
Normally, it is best to use `setFromOptions()` and then set
the SVD type from the options database rather than by using
this routine.  Using the options database provides the user
with maximum flexibility in evaluating the different available
methods.

setUp(self)

 
Sets up all the internal data structures necessary for the
execution of the singular value solver.

Notes
-----
This function need not be called explicitly in most cases,
since `solve()` calls it. It can be useful when one wants to
measure the set-up time separately from the solve time.

setWhichSingularTriplets(self, which)

 
Specifies which singular triplets are to be sought.

Parameters
----------
which: `SVD.Which` enumerate
       The singular values to be sought (either largest or
       smallest).

view(self, Viewer viewer=None)

 
Prints the SVD data structure.

Parameters
----------
viewer: Viewer, optional
        Visualization context; if not provided, the standard
        output is used.

Overrides: petsc4py.PETSc.Object.view