next | previous | forward | backward | up | top | index | toc | Macaulay2 website
Macaulay2Doc :: SVD

SVD -- singular value decomposition of a matrix

Synopsis

Description

If Sigma is the diagonal m by n matrix whose (i,i) entry is the i-th element of S, then M = U Sigma Vt. This is the singular value decomposition of M. The entries of S are (up to roundoff error) the eigenvalues of the Hermitian matrix M * (conjugate transpose M)

M may also be a MutableMatrix in which case the returned values U and Vt are also mutable matrices.

If M is over CC, then U and Vt are unitary matrices over CC. If M is over RR, U and Vt are orthogonal over RR.

i1 : printingPrecision=2;
i2 : M = map(RR^3, RR^5, (i,j) -> (i+1)^j * 1.0)

o2 = | 1 1 1 1  1  |
     | 1 2 4 8  16 |
     | 1 3 9 27 81 |

                3          5
o2 : Matrix RR    <--- RR
              53         53
i3 : (S,U,V) = SVD(M)

o3 = ({88 }, | -.016 -.4  -.92  |, | -.014 -.038 -.11 -.32 -.94  |)
      {3.9}  | -.21  -.89 .4    |  | -.28  -.41  -.57 -.59 .29   |
      {.8 }  | -.98  .2   -.068 |  | -.74  -.41  .066 .51  -.15  |
                                   | .084  .33   -.81 .47  -.081 |
                                   | -.61  .74   .077 -.27 .062  |

o3 : Sequence
i4 : M' = (transpose U) * M * (transpose V)

o4 = | 88      1.5e-14 2.1e-14  -3.1e-15 1e-14    |
     | 1e-14   3.9     -5.3e-15 4e-16    2.2e-15  |
     | 7.5e-15 2.2e-14 .8       -1.6e-16 -1.4e-15 |

                3          5
o4 : Matrix RR    <--- RR
              53         53
We can clean the small entries from the result above with clean.
i5 : e = 1e-10;
i6 : clean_e M'

o6 = | 88 0   0  0 0 |
     | 0  3.9 0  0 0 |
     | 0  0   .8 0 0 |

                3          5
o6 : Matrix RR    <--- RR
              53         53
i7 : clean_e norm (1 - U * transpose U)

o7 = 0

o7 : RR (of precision 53)
Alternatively, if the only issue is display of the matrix, we may set the printing accuracy.
i8 : printingAccuracy = 2

o8 = 2
i9 : M'

o9 = | 88 0   0  -0 0  |
     | 0  3.9 -0 0  0  |
     | 0  0   .8 -0 -0 |

                3          5
o9 : Matrix RR    <--- RR
              53         53
Now let's try the divide and conquer algorithm and compare answers.
i10 : (S', U', V') = SVD(M, DivideConquer => true)

o10 = ({88 }, | -.02 -.4  -.92 |, | -.01 -.04 -.11 -.32 -.94 |)
       {3.9}  | -.21 -.89 .4   |  | -.28 -.41 -.57 -.59 .29  |
       {.8 }  | -.98 .2   -.07 |  | -.74 -.41 .07  .51  -.15 |
                                  | .08  .33  -.81 .47  -.08 |
                                  | -.61 .74  .08  -.27 .06  |

o10 : Sequence
i11 : norm \ ({S', U', V'}-{S, U, V})

o11 = {0, 0, 0}

o11 : List
The SVD routine calls on the SVD algorithms in the lapack library.

See also

Ways to use SVD :

For the programmer

The object SVD is a method function with options.