 Hands-On Exercises

# Exercise 4: Singular Value Decomposition

In this exercise we turn our attention to the Singular Value Decomposition (SVD). Remember that in real symmetric (or complex Hermitian) matrices, singular values coincide with eigenvalues, but in general this is not the case. The SVD is defined for any matrix, even rectangular ones. Singular values are always non-negative real values.

This example works also by reading a matrix from a file. In particular, the matrix to be used is related to a 2D reaction-diffusion model. More details about this problem can be found at Matrix Market.

## Compiling

Copy the file ex14.c [plain text] to your directory and add these lines to the makefile

```ex14: ex14.o
\${RM} ex14.o
```

Note: In the above text, the blank space in the 2nd and 3rd lines represents a tab.

Build the executable with the command

```\$ make ex14
```

## Running the Program

In order to run the program, type the following

```\$ ./ex14 -file \$SLEPC_DIR/share/slepc/datafiles/matrices/rdb200.petsc
```

You should get an output similar to this

```Singular value problem stored in file.

Reading REAL matrix from a binary file...
Number of iterations of the method: 3
Solution method: cross

Number of requested singular values: 1
Stopping condition: tol=1e-08, maxit=100
SVD solve converged (1 singular triplet) due to CONVERGED_TOL; iterations 3
---------------------- --------------------
sigma            relative error
---------------------- --------------------
35.007519           1.40809e-10
---------------------- --------------------
```

## Source Code Details

The way in which the SVD object works is very similar to that of EPS. However, some important differences exist. Examine the source code of the example program and pay attention to the differences with respect to EPS. After loading the matrix, the problem is solved by the following sequence of function calls:

SVDCreate`(MPI_Comm comm,SVD *svd);`
SVDSetOperator`(SVD svd,Mat A);`
SVDSetFromOptions`(SVD svd);`
SVDSolve`(SVD svd);`
SVDGetConverged`(SVD svd, int *nconv);`
SVDGetSingularTriplet`(SVD svd,int i,PetscReal *sigma,Vec u,Vec v);`
SVDDestroy`(SVD svd)`;

First, the singular value solver (SVD) context is created and the matrix associated with the problem is specified. Then various options are set for customized solution. After that, the program solves the problem, retrieves the solution, and finally destroys the SVD context.

Note that the singular value, `sigma`, is defined as a `PetscReal`, and that the singular vectors are simple `Vec`'s.

## SVD Options

Most of the options available in the EPS object have their equivalent in SVD. A full list of command-line options can be obtained by running the example with the option `-help`.

To show information about the SVD solver, add the `-svd_view` option:

```\$ ./ex14 -file \$SLEPC_DIR/share/slepc/datafiles/matrices/rdb200.petsc -svd_view
```

Note: All the command-line options related to the SVD object have the `-svd_` prefix.

Your output will include something like this

```SVD Object: 1 MPI processes
type: cross
implicit matrix
EPS Object: (svd_cross_) 1 MPI processes
type: krylovschur
50% of basis vectors kept after restart
using the locking variant
problem type: symmetric eigenvalue problem
selected portion of the spectrum: largest real parts
number of eigenvalues (nev): 1
number of column vectors (ncv): 16
maximum dimension of projected problem (mpd): 16
maximum number of iterations: 100
tolerance: 1e-09
convergence test: relative to the eigenvalue
BV Object: (svd_cross_) 1 MPI processes
type: svec
17 columns of global length 200
vector orthogonalization method: classical Gram-Schmidt
orthogonalization refinement: if needed (eta: 0.7071)
block orthogonalization method: GS
doing matmult as a single matrix-matrix product
DS Object: (svd_cross_) 1 MPI processes
type: hep
solving the problem with: Implicit QR method (_steqr)
ST Object: (svd_cross_) 1 MPI processes
type: shift
shift: 0.
number of matrices: 1
using a shell matrix
transpose mode: explicit
selected portion of the spectrum: largest
number of singular values (nsv): 1
number of column vectors (ncv): 16
maximum dimension of projected problem (mpd): 16
maximum number of iterations: 100
tolerance: 1e-08
convergence test: relative to the singular value
```

The output shows all the options that are susceptible of being changed, either from the command line or from the source code of the program: the method, the portion of the spectrum (largest or smallest singular values), the number of singular values (`nsv`), etc.

Try to change some of the values, for instance:

```\$ ./ex14 -file \$SLEPC_DIR/share/slepc/datafiles/matrices/rdb200.petsc
-svd_nsv 10 -svd_ncv 40 -svd_smallest
```

The "transpose mode" refers to whether the transpose of matrix A is being built explicitly or not (see SVDSetImplicitTranspose for an explanation).

Note that in the sample output above, the SVD object contains an EPS object. This only happens in some SVD solver types, as detailed below.

## Changing the Singular Value Solver

SLEPc provides several solvers for computing the SVD, which can be selected in the source code with the function SVDSetType, or at run time:

```\$ ./ex14 -file \$SLEPC_DIR/share/slepc/datafiles/matrices/rdb200.petsc
-svd_type trlanczos
```

The following table shows the list of SVD solvers available in SLEPc.

Solver Command-line Name Parameter
Cross Product cross SVDCROSS
Cyclic Matrix cyclic SVDCYCLIC
Lanczos with Explicit Restart lanczos SVDLANCZOS
Lanczos with Thick Restart trlanczos SVDTRLANCZOS
Lapack lapack SVDLAPACK

Note: The Lapack solver is not really a full-featured singular value solver but simply an interface to some LAPACK routines. These routines operate in dense mode with only one processor and therefore are suitable only for moderate size problems. This solver should be used only for debugging purposes.

Note: The default solver is `cross`.

The first two solvers, `cross` and `cyclic`, are not real methods implemented in the SVD module, but are two convenient ways of solving the SVD problem by making use of the eigensolvers available in the EPS module. In those two cases, the SVD object manages an EPS object internally, whose parameters can be set as desired (typically only the method). For example:

```\$ ./ex14 -file \$SLEPC_DIR/share/slepc/datafiles/matrices/rdb200.petsc
-svd_type cyclic -svd_cyclic_eps_type lanczos
```
[Previous exercise] [Index] [Next exercise]