CHAPTER 2: MATRIX COMPUTATIONS AND LINEAR ALGEBRA

 

Lecture 2.1: Two-dimensional arrays (matrices) and matrix operations

 

In linear algebra, a matrix is a rectangular array of elements (numbers). If a matrix A = [ ai,j ] has m rows and n columns, it has m-by-n elements ai,j  where the first index i  runs for rows of A and the second index j runs for columns of A. MATLAB matrices are direct and convenient numerical realization of matrices in linear algebra.

 

   A = [ 0.1 0.2 0.3 0.4 ; -0.1,-0.2,-0.3,-0.4 ; 0 1 0 1 ]

% Notice that rows in MATLAB matrices are separated by the sign ";"

% Notice that elements at each row can be separated by spaces and commas  

 

   A =

      0.1000    0.2000    0.3000    0.4000

   -0.1000   -0.2000   -0.3000   -0.4000

            0    1.0000         0    1.0000  

 

 A1 = [ 0.1 0.2 0.3 0.4 ; 0 1 0 1 0] ;

  % Number of elements in each row must be identical  

 

??? Error using ==> vertcat

All rows in the bracketed expression must have the same

number of columns.  

 

The colon sign ":" represents all elements of a row or a column of a MATLAB matrix.

 

   A(:,1), A(:,4) % columns of A are column-vectors     

 

   ans =

        0.1000

       -0.1000

             0

ans =

       0.4000

     -0.4000

    1.0000  

 

   A(1,:), A(3,:) % rows of A are row-vectors  

 

  ans =

        0.1000    0.2000    0.3000    0.4000

ans =

          0     1     0     1  

 

  A(1,2), A(3,3) % individual elements of A are scalar variables  

 

   ans =

        0.2000

 ans =

         0   

 

 

 

 

Matrix operations:

 

·         all vector operations, such as pointwise additions, subtractions, multiplications, divisions, and powers of matrices of the same size

 

  A = [ 1 2 3 4 ; -1,-2,-3,-4 ];

B = [ 4 3 2 1 ; -4,-3,-2,-1 ];

C1 = A + B, C2 = A.*B, C3 = A./B, C4 = A.^B

% Pointwise matrix operations are preferable compared to double loops in i,j

% because of (i) clear mathematical notations, (ii) computational efficiency

  

 

C1 =

     5     5     5     5

    -5    -5    -5    -5

C2 =

     4     6     6     4

     4     6     6     4

C3 =

    0.2500    0.6667    1.5000    4.0000

    0.2500    0.6667    1.5000    4.0000

C4 =

    1.0000    8.0000    9.0000    4.0000

    1.0000   -0.1250    0.1111   -0.2500  

 

·         size: outputs numbers of rows and columns

 

  [nRow,nCol] = size(A)  

 

nRow =     2

nCol =     4  

 

·         transpose: A', outputs a transposed matrix AT

 

   B = A', size(B)  

 

B =  1    -1

     2    -2

     3    -3

     4    -4

ans =

     4     2  

 

·         elementary matrices:

 

   zeros(2,6) % a null 2-by-6 matrix

 zeros(3) % a null 3-by-3 matrix  

 

ans =    0     0     0     0     0     0

          0     0     0     0     0     0

ans =    0     0     0

     0     0     0

          0     0     0  

 

 

  ones(2,6) % a 2-by-6 matrix of ones

  ones(3) % a 3-by-3 matrix of ones  

 

ans =

     1     1     1     1     1     1

     1     1     1     1     1     1

ans =

     1     1     1

     1     1     1

     1     1     1  

 

 

eye(2,6), eye(3)  

 

ans =

     1     0     0     0     0     0

     0     1     0     0     0     0

ans =

     1     0     0

     0     1     0

     0     0     1  

 

 

  rand(2,6) % random numbers are uniformly distributed in the interval [0,1]

  randn(3) % random numbers are normally distributed with mean zero and variance one  

 

ans =

    0.9501    0.6068    0.8913    0.4565    0.8214    0.6154

    0.2311    0.4860    0.7621    0.0185    0.4447    0.7919

ans =

   -0.4326    0.2877    1.1892

   -1.6656   -1.1465   -0.0376

    0.1253    1.1909    0.3273  

 

·         special matrices:

 

 

C1 = hilb(6) % Hilbert matrices are examples of ill-conditioned matrices  

 

C1 =   1.0000    0.5000    0.3333    0.2500    0.2000

       0.5000    0.3333    0.2500    0.2000    0.1667

  0.3333    0.2500    0.2000    0.1667    0.1429

       0.2500    0.2000    0.1667    0.1429    0.1250

       0.2000    0.1667    0.1429    0.1250    0.1111  

 

 

 

 

C2 = pascal(6)

  % Pascal matrices produce a table of binomial coefficients

  % Coefficients at the k-th anti-diagonal of Pascal matrices produce

     % coefficients of the binomial expansion: (x + y)^(k-1)  

 

C2 =      1     1     1     1     1     1

      1     2     3     4     5     6

       1     3     6    10    15    21

       1     4    10    20    35    56

       1     5    15    35    70   126

       1     6    21    56   126   252  

 

 

   C3 = vander([1:6])  

 

C3 =       1           1           1           1           1           1

          32          16           8           4           2           1

         243          81          27           9           3           1

        1024         256          64          16           4           1

        3125         625         125          25           5           1

        7776        1296         216          36           6           1  

 

r = [ 1 2 3 4 5 6 ]; s = [ -1,-2,-3,-4, -5];

C4 = toeplitz(r,s) % r and s are the column and row vectors of C4

C5 = toeplitz(r)   % the symmetric Toeplitz matrix

rInv = r(length(r):-1:2); rInv = [r(length(r)),rInv];

C6 = toeplitz(r,rInv)' % the circulant matrix   

 

Warning: Column wins diagonal conflict.

> In c:\progra~1\matlab6\toolbox\matlab\elmat\toeplitz.m at line 18

C4 =

     1    -2    -3    -4    -5

     2     1    -2    -3    -4

     3     2     1    -2    -3

     4     3     2     1    -2

     5     4     3     2     1

     6     5     4     3     2

C5 =

     1     2     3     4     5     6

     2     1     2     3     4     5

     3     2     1     2     3     4

     4     3     2     1     2     3

     5     4     3     2     1     2

     6     5     4     3     2     1

Warning: Column wins diagonal conflict.

> In c:\progra~1\matlab6\toolbox\matlab\elmat\toeplitz.m at line 18

C6 =

     1     2     3     4     5     6

     6     1     2     3     4     5

     5     6     1     2     3     4

     4     5     6     1     2     3

     3     4     5     6     1     2

     2     3     4     5     6     1  

·         matrix-vector multiplications

 

A = [ 1 2 3 1 ; 0 1 4 2 ; 3 0 2 3 ], x = [ 1; 2; 1; 1; ]  

 

A =  1     2     3     1

     0     1     4     2

     3     0     2     3

x =   1

       2

       1

       1  

 

  % the number of columns of A must be equal to the number of rows of x

  y1 = A*x  % MATLAB matrix multiplication operator

 

% a row-oriented loop:

[n,m] = size(A);

for k = 1:n,     y2(k) = A(k,:)*x; end

y2

 

% a column-oriented loop

y3 = zeros(n,1);

for j = 1 : m,   y3 = y3 + A(:,j)*x(j); end

y3  

 

y1 =  9

       8

       8

y2 =     9     8     8

y3 =   9

       8

       8  

 

·         matrix-matrix mutiplications:

 

If A is a n-by-m matrix, B is a m-by-q matrix, then the matrix multiplication makes sense and C = A*B is a n-by-q matrix with elements: ci,j = ai,k bk,j . The other product B*A is defined if n = q. For a square matrices of the same size, both products are defined but they are not equal in a general case: A*B B*A.

 

A = [ -2 1 0 ; 1,-2,1 ; 0,1,-2];

B = [1 0 1 ; 0 1 0; 1 0 1];

C1 = A*B, C2 = B*A  

 

C1 =

    -2     1    -2

     2    -2     2

    -2     1    -2

 

C2 =

    -2     2    -2

     1    -2     1

    -2     2    -2  

 

 

n = 3;

for k = 1 : n

  for j = 1 : n

     C3(k,j) = A(k,:)*B(:,k); % the dot product algorithm

end

end

C3

 

C4 = zeros(n,n);

for j = 1 : n

  for p = 1 : n

     C4(:,j) = C4(:,j) + A(:,p)*B(k,j); % the column-oriented algorithm

  end

end

C4

 

C5 = zeros(n,n);

for p = 1 : n

  C5 = C5 + A(:,k)*B(k,:); % the outer-product oriented algorithm

end

C5

 

C3 =

    -2    -2    -2

    -2    -2    -2

    -2    -2    -2

 

C4 =

    -1     0    -1

     0     0     0

    -1     0    -1

 

C5 =

     0     0     0

     3     0     3

    -6     0    -6