MAPLE is a commercial symbolic algebra system which allows computations in many fields of mathematics and engineering. A general introduction to MAPLE use and how MAPLE can be accessed at UMBC can be found in a separate document.
Materials for a linear algebra class taught completely with
MAPLE notebooks are available from
www.mapleapps.com.
Most of MAPLE's linear algebra functions are contained in
the linalg package. Thus, the first command you
need to execute when starting MAPLE to work on linear algebra
problems is:
> with(linalg);
This command loads the linear algebra package for use.
To define a matrix you use either the matrix function, declaring the number of rows and columns and giving the entries row by row:
> a:=matrix(2,2,[1,2,3,4]);
[1 2]
a := [ ]
[3 4]
> a:=matrix([[1,2],[3,4]]); # Or equivalently
[1 2]
a := [ ]
[3 4]
MAPLE does not normally understand how to do simple operations with
matrices. The evalm function allows these operations, though.
The one sticky point when working with evalm is that matrix multiplication
is given by the &* operator, not the usual *
operator. (Note: there must be spaces around the a times itself is:
> evalm(a &* a);
[ 7 10]
[ ]
[15 22]
You can pull out an individual element of a matrix
by specifying its column and row. Similarly, you can pull
out an individual row or column with the row and
col functions, respectively.
> a[2,1]; # the element in the second row and first column
3
> row(a,2); # the second row
[3, 4]
>> col(a,2); # the second column
[2, 4]
A common use of matrices is to express and solve a system of linear
equations. Consider the matrix equation a*x=b, where a is the 2x2 matrix
[1,2,3,4] and b is the vector [1,1]. This system can be solved with the
linsolve function. Thus, to find a
solution of this system use the command:
> b:=vector([1,1]);
b := [1, 1]
> linsolve(a,b);
[-1, 1]
We can confirm this result with the inverse operation, multiplying the matrix a against the vector we got, ans:
> evalm(a &* ");
[1, 1]
Another way to solve this linear system would be to compute the inverse
of the matrix a (if an inverse exists):
> inverse(a);
[-2 1 ]
[ ]
[3/2 -1/2]
> evalm(inverse(a) &* b);
[-1, 1]
In many cases a system of linear equations may not have a single
unique solution. In those cases linsolve will return
a parametrized description of the set of ALL solutions. Consider
the system A*x=b:
> A:=matrix([[1,2,3,4],[2,3,1,0]]);
[1 2 3 4]
A := [ ]
[2 3 1 0]
> b:=vector([1,2]);
b := [1, 2]
> linsolve(A,b);
[_t[1], _t[2], -2 _t[1] - 3 _t[2] + 2, 5/4 _t[1] + 7/4 _t[2] - 5/4]
Thus the complete set of solutions of the system A*x=b is all
vectors of the form:
[t1, t2, -2*t1 - 3*t2 + 2, 5/4*t1 + 7/4*t2 - 5/4]
where t1 and t2 range over all possible values.
MAPLE implements the LU decomposition with the command LUdecomp().
By default this command will perform an LU decomposition with partial
pivoting, producing three matrices, L, U and P such that if A is the
original matrix then L*U=P*A. The matrix U is printed out by the command
and, if the optional parameters are used, the matrices U and P are
saved as named matrices. (The permutation is the identity matrix if
the matrix entries are not decimal fractions.)
> a:=matrix([[1,2,1],[1,3,1],[2,1,1]]); # Matrix with no decimal entries
[1 2 1]
[ ]
a := [1 3 1]
[ ]
[2 1 1]
> LUdecomp(a, P='p', L='l', U='u' , U1='u1', R='r', rank='ran', det='d');
[1 2 1]
[ ]
[0 1 0]
[ ]
[0 0 -1]
> print(l);
[1 0 0]
[ ]
[1 1 0]
[ ]
[2 -3 1]
> print(u);
[1 2 1]
[ ]
[0 1 0]
[ ]
[0 0 -1]
> evalm(l &* u); # Check our work
[1 2 1]
[ ]
[1 3 1]
[ ]
[2 1 1]
When matrices have floating point entries (ie decimal fractions)
the LUdecomp() command will also perform permutations
or partial pivoting:
> b:=matrix([[1,2,1],[3,1,1],[2,3,1.0]]); # Matrix w/ floating point
[1 2 1 ]
[ ]
b := [3 1 1 ]
[ ]
[2 3 1.0]
> LUdecomp(b, P='p', L='l', U='u' , U1='u1', R='r', rank='ran', det='d');
[3. 1. 1. ]
[ ]
[0 2.333333333333333 .3333333333333334]
[ ]
[0 0 .4285714285714286]
> print(p);
[0 0 1.]
[ ]
[1. 0 0 ]
[ ]
[0 1. 0 ]
> print(l);
[ 1. 0 0 ]
[ ]
[.6666666666666666 1. 0 ]
[ ]
[.3333333333333333 .7142857142857143 1.]
> print(evalm(p &* l &* u)); # Confirm that b=p*l*u
[.9999999999 1.999999999 1.000000000]
[ ]
[ 3. 1. 1. ]
[ ]
[2.000000000 3.000000000 1.000000000]
The MAPLE commands eigenvalues() and eigenvects() will
compute the eigenvalues and eigenvectors of a matrix:
> A:=matrix([[0,0,0.33],[0.18,0,0],[0,0.71,0.94]]);
[ 0 0 .33]
[ ]
A := [.18 0 0 ]
[ ]
[ 0 .71 .94]
> eigenvalues(A);
-.02179636988 + .2059184805 I, -.02179636988 - .2059184805 I, .9835927398
> eigenvects(A); # List the eigenvalues, following each by its eigenvectors
[.9835927403, 1, {[.3822648106, .0699554439, 1.139372397]}], [
-.02179636983 + .2059184806 I, 1, {[-.6797453503 - .6620875268 I,
-.5101419741 + .6481856590 I, .4580364800 - .3804273484 I]}], [
-.02179636983 - .2059184806 I, 1, {[-.6797453503 + .6620875268 I,
-.5101419741 - .6481856590 I, .4580364800 + .3804273484 I]}]
You can also use the commands in MAPLE to compute the eigenvalues of a matrix. Here we use the power method to estimate the dominant (largest) eigenvalue of a matrix and its corresponding eigenvector. We choose a random start vector - in this case (1,1,1) and multiply the matrix A a large number of times (the more often the more accurate). The resulting vector is approximately the eigenvector. The eigenvalue is approximated by multiplying by A again and seeing the ratio of expansion of the vector (we compared first elements).
> v1:=evalm((A^30) &* vector([1,1,1]));
v1 := [.3478324912, .06365424009, 1.036743978]
> v2:=evalm(A &* v1);
v2 := [.3421255127, .06260984842, 1.019733850]
> v2[1]/v1[1]; # Estimate the eigenvalue
.9835927389