Advanced matrix stuff

Ed Bueler's Matrices in Matlab: Advanced

Here are some more advanced operations in Matlab. Start with the matrix

» A=[2 -3 2; -4 2 -6; 2 2 4]

A =

     2    -3     2
    -4     2    -6
     2     2     4

"LU" decomposition is a single command, but since Matlab does the full algorithm with pivoting, you get a surprise if you write:

» [L, U]=lu(A)

L =

  -0.50000000000000  -0.66666666666667   1.00000000000000
   1.00000000000000                  0                  0
  -0.50000000000000   1.00000000000000                  0


U =

  -4.00000000000000   2.00000000000000  -6.00000000000000
                  0   3.00000000000000   1.00000000000000
                  0                  0  -0.33333333333333

That is, L isn't lower triangular! A row permutation of it is lower triangular, and we will recover that. First, though, it is easy to check that we have a decomposition of A:

» L*U

ans =

     2    -3     2
    -4     2    -6
     2     2     4

The right way is easy: ask the "lu" command for the permutation matrix P:

» [L,U,P]=lu(A)

L =

   1.00000000000000                  0                  0
  -0.50000000000000   1.00000000000000                  0
  -0.50000000000000  -0.66666666666667   1.00000000000000

U =

  -4.00000000000000   2.00000000000000  -6.00000000000000
                  0   3.00000000000000   1.00000000000000
                  0                  0  -0.33333333333333

P =

     0     1     0
     0     0     1
     1     0     0

The matrix P says that R1 -> R2, R2 -> R3, and R3 -> R1. Check that we have a decomposition: PA = LU.

» P*A

ans =

    -4     2    -6
     2     2     4
     2    -3     2

» L*U

ans =

    -4     2    -6
     2     2     4
     2    -3     2

Matlab does eigen stuff. Here is a symmetric matrix:

» C=[1.59 1.69 2.13; 1.69 1.31 1.72; 2.13 1.72 1.85]

C =

   1.59000000000000   1.69000000000000   2.13000000000000
   1.69000000000000   1.31000000000000   1.72000000000000
   2.13000000000000   1.72000000000000   1.85000000000000

Easy check for symmetry:

» C'-C

ans =

     0     0     0
     0     0     0
     0     0     0

Eigenvalues and vectors can be approximated by built in routines:

» eig(C)

ans =

  -0.13654626528583
  -0.42131129760583
   5.30785756289166

» [V,D]=eig(C)

V =

   0.23382978369234  -0.77264233596480   0.59020966862060
  -0.84501101513945   0.13876277780700   0.51643128854502
   0.48091581025309   0.61949068771866   0.62044441433418

D =

  -0.13654626528583                  0                  0
                  0  -0.42131129760583                  0
                  0                  0   5.30785756289166

It should be clear that the first command just gets the eigenvalues, and the second gets both eigenvalues and eigenvectors.