Matrices

  1. Creating matrices
  2. Basic functionality
  3. Solving sets of linear equations
  4. Matrix decomposition
  5. Echelon forms and Null Space
  6. Linear algebra
  7. Special functionality
  8. Random matrices

Creating matrices

Matrices may be created as follows:


B=(2, 3, 4 ; 4, 6, 6)
2  3  4
4  6  6

If the rest of values in a row are zero, you may use a ; and begin entering the next row. Zeros will be added:


C=(1 ; 0, 1 ; 0, 0, 3)
1  0  0
0  1  0
0  0  3

The following functions also create matrices:

Matrices can also be created by creating any entry:


D[2,3]=-2
-2
D
0  0  0
0  0 -2

You may use functions joinright() and joinbottom() to augment matrices. A typical session follows.


A = ( 1, 0 ; 3, 2 )
1  0
3  2
B = ( 5 ; 7 )
5
7
C = joinright( A, B )
1  0  5
3  2  7
D = ( 9 8 2 )
 9  8  2
E = joinbottom( C, D )
1  0  5
3  2  7
9  8  2

Use function submatrix() to build a smaller matrix from an existing one. Functions row() and column() are special cases of function submatrix().

Also note that column matrices and row matrices are treated like arrays too.

Basic functionality

Use the + and - operators provided the matrices have the same dimensions.

Operator * (matrix multiplication) works as long as the inner dimensions of the matrices are the same.

The following functions can be used with matrix arguments:

Solving sets of linear equations

Use function solve() to solve a set of linear equations of the form:

A*x = B

where both A and B are matrices with equal number of rows.
There are three cases to consider. Assuming matrix A has r rows and c columns, we show each case with an example.

Case r = c

Given A and B:


A=(1,2,3;5,6,8;0,-2,3)
  1  2  3
  5  6  8 
  0 -2  3 
B=(3;5;7)
   3 
   5 
   7

We solve this systems as follows:


X = solve(A,B)
   -14/13 
   -19/26 
    24/13 

Calcugator uses an LU Decomposition to solve the system of equations. Matrices A and B can be real or complex.

Case r<c

Let's create a matrix with 3 rows and 5 columns with random integers between 0 and 10:


A=randommatrix(3,5,10)
6  5  2  8  8 
2  1  7  1  7 
6  2  8  8  5 
B=(7;5;6)

The systems is solved using the same function as above:


X = solve(A,B)
0.1559
0.1811
0.1572
0.1387
0.4667

To solve the system Calcugator uses the pseudo inverse of matrix A computed with a Singular Value Decomposition. Both matrices A and B must be real.

Case r>c

Given A and B:


A=(3,1 ; 6,6 ; 3,3 ; 0, 4)
3  1
6  6
3  3
0  4
B=(10 ; 9 ; 0 ; 6)
10
 9
 0
 6

Again, use function solve:


X = solve(A,B)
0.9425
0.6896 

In this case, Calcugator finds the least squares solution. This solution minimizes the norm of A*X-B. Both matrices A and B must be real.

Matrix decomposition

Three types of matrix decompositions are supported.

QR Decomposition

A matrix A (square or rectangular matrix with more rows than columns) that has full rank, can be decomposed into the product:

A = Q*R

where Q is an orthonormal matrix and R is upper triangular. Use function qrgetQ or function qrgetR to compute those matrices.

Example:


A=(3,1 ; 6,6 ; 3,3 ; 0, 4)
3  1
6  6
3  3
0  4
Q=qrgetQ(A)
0.4082  -0.3790 
0.8164   0.1516
0.4082   0.0758 
0.0      0.9097 
R=qrgetR(A)
7.3484  6.5319 
0.0     4.3969
Q*R
3.0000  1.0000 
6.0000  6.0000
3.0000  3.0000 
0.0     4.0000

For properties of matrices Q and R, consult a book on linear Algebra.

Singular Value Decomposition

A square matrix A (or rectangular with more rows than columns), can be decomposed into the product:

A = U*S*VT

where U and V are orthonormal and S is diagonal. Use functions svdgetU(), svdgetS() and svdgetV() to compute those matrices.

Example:


A=(3,1 ; 6,6 ; 3,3 ; 0, 4)
3  1
6  6
3  3
0  4
U=qrgetU(A)
0.2685  -0.4880 
0.8228  -0.1124
0.4114  -0.0562
0.2856   0.8636 
S=qrgetS(A)
10.3037  0.0 
 0.0     3.1358
V=qrgetV(A)
0.6771  -0.7358 
0.7358   0.6771
U*S*transpose(V)
3.0000      1.0000 
5.9999      5.9999
2.9999      2.9999 
1.1102E-15  4.0000
normf(ans-A)
2.4825E-15

The norm returned by function normf() returns the Frobenius norm of the matrix.

Eigenvalue Decomposition

A square matrix A can be decomposed into the following:

A = S*L*S-1

where S is a matrix such that its columns are the eigenvectors and L is a diagonal matrix containing the eigenvalues. Use functions eigenmatrix() and eigenvalues() to compute those matrices.

Example:


A=(2,4,6 ; 4,5 ; 0, -3, 5)
2   4  6  
4   5  0  
0  -3  5 
S=eigenmatrix(A)
-0.8436   0.0885 + 0.8650 i    0.0885 - 0.8650 i 
 0.4918   0.9228 + 0.5488 i    0.9228 - 0.5488 i 
 0.2150  -0.9174 + 0.3833 i   -0.9174 - 0.3833 i 
L=eigenvalues(A)
-1.8613  0                 0 
 0       6.9306+2.6011 i   0
 0       0                 6.9306-2.6011 i
S*L*inverse(S)
 2.0000+1.1102E-16 i   4.0000   5.9999+1.102E-15 i
 3.9999                4.9999  -2.1094E-15+4.4408E-16 i 
-3.0531E-16           -3.0000   5.0000+4.4408E-16 i

Notice all imaginary values are small; let's take the real part of the matrix:


Re(ans)
 2.0000        4.0000   5.9999
 3.9999        4.9999  -2.1094E-15
 -3.0531E-16  -3.0000   5.0000

LU Decomposition

Any matrix A can be decomposed as the multiplication of a lower triangular matrix L times and upper triangular matrix U. If row permutations were required, then a permutation matrix P needs to premultiply the product L*U in order to retrieve A. To obtain matrices L, U and P use functions lugetL(), lugetU() and lugetP() respectively.


A = ( 0, 2, 3, 0 ; 1, 1, 0, 5 ; 1, 2, 3, 1 )
0  2  3  0
1  1  0  5
1  2  3  1
L = lugetL( A )
1   0   0
0   1   0
1  1/2  1
U = lugetU( A )
1  1   0    5
0  2   3    0
0  0  3/2  -4
P = lugetP( A )
0  1  0
1  0  0
0  0  1
P*L*U
0  2  3  0
1  1  0  5
1  2  3  1

Echelon forms and Null Space

To find the echelon form of a matrix, use function echelon().


A = ( 0, 2, 3, 0 ; 1, 1, 0, 5 ; 1, 2, 3, 1 )
0  2  3  0
1  1  0  5
1  2  3  1
echelon( A )
1  1   0    5
0  2   3    0
0  0  3/2  -4

To find the reduced echelon form of a matrix, use function echelonr()


echelonr( A )
1  0  0   1
0  1  0   4
0  0  1  -8/3

Function nullspace returns a matrix that spans the null space of a matrix.


A = ( 1, 2, 2, 3 ; 8, -8, -2, 6 ; -2, 4, 2, 0 )
 1   2   2   3
 8  -8  -2   6
-2   4   2   0
N = nullspace( A )
-1/2  -3/2
-3/4  -3/4
 1     0
 0     1

Function solven() solves exactly for a matrix equation A*X=B for the case that the number of columns of matrix A is bigger than the number of rows. Function solven() returns a solution such that adding a linear combination of columns of the null space, the equality still hods.


b = ( 1 ; -1 ; 1 )
 1
-1
 1
X = solven( A, b )
1/4
3/8
 0
 0
X2 = X + 2*column(N,1) - 7*column(N,2)
39/4
33/8
 2
 -7
A*X2
 1
-1
 1

Linear algebra

The following functions are implemented:

Special functionality

When testing for singularity while finding an inverse, Calcugator uses a default tolerance equal to 1.e-40.

To change this value do:

Where tol is the new tolerance to be used in the check for singularity.

Example:


set("matrixSingularityTolerance", 1e-20)
Property "matrixSingularityTolerance" was set to 1.0E-20.

sets the tolerance to 1e-20.

Other functions to use with matrices are:

Random matrices

Matrices with random values can be created using the randommatrix() function.

Example:

A random matrix with real values between 0. and 5.:


a=randommatrix(4,4,5.)
4.8761  0.5769  0.6039  4.5528   
2.2761  4.1217  1.3126  0.0923   
0.6550  4.8586  1.3032  4.7012   
3.2179  0.6328  1.2388  2.8118   

A random matrix with integer values between 3 and 7:


a=randommatrix(4,4,3,7)
7  6  3  5
6  7  5  6
4  4  4  7
3  7  3  5