## `numeric::singularvalues` -- numerical singular values of a matrix

#### Introduction

`numeric::singularvalues``(A)` returns numerical singular values of the matrix `A`.

#### Call(s)

`numeric::singularvalues(A)`

#### Parameters

`A` | - | a numerical matrix of domain type `DOM_ARRAY` or of category `Cat::Matrix`.

#### Returns

an ordered list of real floating point values

#### Side Effects

The function is sensitive to the environment variable `DIGITS`, which determines the numerical working precision.

#### Related Functions

`linalg::eigenvalues`, `linalg::eigenvectors`, `numeric::eigenvalues`, `numeric::eigenvectors`, `numeric::singularvectors`, `numeric::spectralradius`

#### Details

- The singular values of an m x n matrix A are the p=min(m,n) real nonnegative square roots of the eigenvalues of A^H*A (for p=n) or of A*A^H (for p=m). The Hermitean transpose A^H is the complex conjugate of the transpose of A.
- `numeric::singularvalues` returns a list of real singular values [d[1],..,d[p]] sorted by `numeric::sort`, i.e., 0.0 <= d[1] <= .. <= d[p].
- All entries of `A` must be numerical. Numerical expressions such as `exp(PI), sqrt(2)` etc. are accepted and converted to floats. Non-numerical symbolic entries lead to an error.
- `Cat::Matrix` objects, i.e., matrices `A` of a matrix domain such as `Dom::Matrix``(..)` or `Dom::SquareMatrix``(..)`, are internally converted to arrays over expressions via `A::dom::expr(A)`.

- Singular values are approximated with an absolute precision of 10^(-DIGITS)*r, where r is the spectral radius of `A` (i.e., r is the absolute value of the largest eigenvalue). Consequently, large singular values should be computed correctly to `DIGITS` decimal places. The numerical approximations of the small singular values are less accurate.

- Singular values may also be computed via

      map(numeric::eigenvalues(A*A^H), sqrt);

  or

      map(numeric::eigenvalues(A^H*A), sqrt);

  respectively. The use of `numeric::singularvalues` avoids the costs of the matrix multiplication. Further, the eigenvalue routine requires about twice as many `DIGITS` to compute small singular values with the same precision as `numeric::singularvalues`. Cf. example 2.

#### Example 1

The singular values of A and A^H coincide:

    >> A := array(1..3, 1..2, [[1, 2*I], [2, 3],[3, PI]]):

    >> numeric::singularvalues(A)

                            [1.503668692, 5.882906158]

The Hermitean transpose B=A^H:

    >> B := array(1..2, 1..3,[[1, 2, 3], [-2*I, 3, PI]]):

    >> numeric::singularvalues(B)

                            [1.503668692, 5.882906158]

    >> delete A, B:

#### Example 2

We use `numeric::eigenvalues` to compute singular values:

    >> M := Dom::Matrix():

    >> A := M([[1, 2*I], [PI, 312689/49766*I], [2, 4*I]]):

The Hermitean transpose B=A^H can be computed by the methods `conjugate` and `transpose` of the matrix domain:

    >> B := M::conjugate(M::transpose(A)):

Note that A^H*A is positive semi-definite and cannot have negative eigenvalues. However, computing small eigenvalues is numerically ill-conditioned and a small negative value occurs due to round-off:

    >> numeric::eigenvalues(B*A)

                          [-8.67361738e-19, 74.34802201]

Consequently, an illegal imaginary singular value is computed:

    >> map(%, sqrt)

                      [0.0000000009313225746 I, 8.622529908]

We have to increase `DIGITS` in order to compute this value more accurately:

    >> DIGITS := 20: map(numeric::eigenvalues(B*A), sqrt)

             [0.000000000015115433585415141592, 8.6225299075259371493]

With `numeric::singularvalues` the standard precision suffices:

    >> DIGITS := 10: numeric::singularvalues(A)

                          [1.511542232e-11, 8.622529908]

    >> delete M, A, B:

#### Background

- The function implements standard numerical algorithms from the Handbook of Automatic Computation by Wilkinson and Reinsch.

#### Changes

- Internal conversion of `Cat::Matrix` objects now uses the method `"expr"` of the matrix domain.
- The singular values are now sorted by `numeric::sort` from small values to large values.

Do you have questions or comments?
Copyright © SciFace Software GmbH & Co. KG 2000