MATRICES AND MATRIX OPERATIONS: Unit 21

Dr. Wlodzislaw Kostecki

The Papua New Guinea University of Technology (PNGUT)

Department of Electrical and Communication Engineering

Lae, Morobe Province

Papua New Guinea

Copyright © 2000 by Wlodzislaw Kostecki

All rights reserved

-------------------------------------------------------------------

(21) Eigenvalues and eigenvectors of matrices

OBJECTIVES :

• To provide a concise analysis of solution of a homogeneous system of algebraic equations.

• To introduce the concepts of the characteristic matrix , characteristic determinant , characteristic polynomial , and characteristic equation of a square matrix.

• To introduce the concepts of eigenvalues or characteristic roots of a square matrix and of eigenvectors or characteristic vectors corresponding to matrix eigenvalues.

• To provide a step-by-step example of computing eigenvalues and eigenvectors of a matrix with numerical elements.

• To introduce the function eigenvects and analyse the structure returned by the function.

• To provide selection procedures for extracting in a unique way eigenvalues and eigenvectors from the structure returned by the eigenvects function.

• To introduce the concept of the Frobenius norm of a vector and show how to use it for computing the norm of eigenvectors.

• To show how to perform normalization of eigenvectors using the normalize function.

• To specify and illustrate properties of eigenvalues and eigenvectors.

• To give the Cayley-Hamilton theorem and provide a numerical example for verification of the theorem.

• To suggest an application of the Cayley-Hamilton theorem to computing functions of matrices.

> restart : interface(warnlevel=0) : with(linalg, charmat, charpoly, det, diag, eigenvals, eigenvects, inverse, norm, normalize, trace, transpose, vector) :

A. Definitions and properties

It was stated in Unit (19) that a system of linear algebraic inhomogeneous equations may be written in the matrix form

[ A ] [ X ] = [ K ] (1)

where [ A ] is an ( × ) coefficient matrix, [ X ] is a ( × ) solution column matrix (vector), and [ K ] is a ( × ) inhomogeneous column matrix (vector).

It follows from matrix equation (1) that [ K ] is proportional to [ X ], or

[ K ] = [ X ] (2)

where is some scalar multiplier.

Therefore, equation (1) may be rewritten in the form

[ A ] [ X ] = [ X ] (3)

which is the matrix form of the system of homogeneous equations

( [ A ] – [ U ] ) [ X ] = [ 0 ] (4)

where [ U ] is the unit matrix and [ 0 ] is the zero column matrix (vector).

There is a theorem, which states that the homogeneous system of equations

[ A ] [ X ] = [ 0 ] (5)

always has the trivial solution [ X ] = [ 0 ] and it has a non-trivial solution only when the determinant of matrix [ A ] is . If [ X ] is a non-trivial solution, so also is [ X ], where is an arbitrary scalar.

According to this theorem, equation (4) can only have a non-trivial solution when

Det( [ A ] – [ U ] ) = [ 0 ] (6)

When expanded, this determinant gives rise to an algebraic polynomial equation of degree in of the form

+ ... + (7)

The matrix [ U ] – [ A ] is called the characteristic matrix of [ A ].

The determinant Det( [ A ] – [ U ] ) is called the characteristic determinant associated with [ A ], and equation (7) is called the characteristic equation of [ A ]. The left-hand side of equation (7) is called the characteristic polynomial of [ A ]. The characteristic equation of [ A ] has roots, , , ..., , each of which is called either an eigenvalue , a characteristic root , or, in some texts, a latent root of [ A ]. The roots may be real , complex , or multiples of each other.

Once the eigenvalues have been determined, any of them can be substituted into equation (3) to find a corresponding solution vector [ X ]. Thus, setting will yield a solution vector , which, because of the aforementioned theorem, will only be determined to within an arbitrary scalar multiplier.

This vector is called either an eigenvector , a characteristic vector or, a latent vector of [ A ] corresponding to .

Eigenvalues of a matrix may be , whereas an eigenvector may not be the zero vector.

For example, consider a ( × ) matrix [ A ] given as

> A := matrix(2, 2, [1, 2, 3, 0]) : A = matrix(A) ;

and find its characteristic equation, eigenvalues, and eigenvectors.

Step 1 . Obtain the matrix [ A ] – [ U ]:

The appropriately sized unit matrix corresponding to this case is

> U := diag(1, 1) : U = matrix(U) ;

Therefore,

> l := lambda : `A - lU` := evalm(A - l*U) : A - l*U = matrix(`A - lU`) ;

* * *

N.B. The matrix [ A ] – [ U ] may also be obtained using the charmat function, viz.

> `A - lU` := evalm(-charmat(A, l)) : A - l*U = matrix(`A - lU`) ;

* * *

Step 2 . Find the determinant Det( [ A ] – [ U ] ) , or the characteristic polynomial of [ A ]:

> `det(A - lU)` := det(`A - lU`) : Det(A - l*U) = `det(A - lU)` ;

* * *

N.B. The characteristic polynomial may be obtained in Maple directly , using the charpoly function. For this particular case,

> `charpoly(A)` := charpoly(A, l) : char_poly(A) = `charpoly(A)` ;

* * *

Step 3 . Form the relevant characteristic equation of [ A ]:

> char_eq(A) := `det(A - lU)` = 0 : char_eq(A) ;

or, using direct computation of the characteristic polynomial,

> char_eq(A) := `charpoly(A)` = 0 : char_eq(A) ;

Step 4 . Solve the characteristic equation of [ A ]:

> solution(char_eq(A)) := solve({char_eq(A)}, {l}) : solution(char_eq(A)) ;

Denoting the roots by and and extracting either of them from the solution sequence yields

> l[`1_n`] := lambda[1] : l[`2_n`] := lambda[2] : l[`1_v`] := subs(solution(char_eq(A))[1], lambda) :

> l[`2_v`] := subs(solution(char_eq(A))[2], lambda) : l[`1_n`] = l[`1_v`] ; l[`2_n`] = l[`2_v`] ;

The roots and are the eigenvalues or characteristic roots of [ A ].

* * *

N.B. The eigenvalues may be obtained in Maple directly , using the eigenvals function. For this particular case,

> charroots(A) := eigenvals(A) : char_roots(A) = charroots(A) ;

Denoting the characteristic roots by and and extracting each from the above sequence yields

> l[`1_v`] := charroots(A)[1] : l[`2_v`] := charroots(A)[2] : root1(A) := l[`1_v`] : root2(A) := l[`2_v`] :

> l[`1_n`] = l[`1_v`] ; l[`2_n`] = l[`2_v`] ;

Alternatively, the inert Eigenvals function may be used, but it must be evaluated with evalf . For this case,

> charroots(A) := evalf(Eigenvals(A)) : char_roots(A) = charroots(A) ;

The returned output is an array structure, from which the characteristic roots may be extracted using subscript notation, viz.

> charroots(A)[1] ; charroots(A)[2] ;

* * *

Step 5 . Determine the eigenvector of [ A ] corresponding to :

The zero column matrix (vector) appropriately sized for this case is

> `0` := matrix(2,1, [0, 0]) : `0` = matrix(`0`) ;

(a) The eigenvector is

> x[1][l1_n] := x[1][l[`1_n`]] : x[2][l1_n] := x[2][l[`1_n`]] : X[l1_n] := X[lambda[1]] :

> X[l1_sv] := matrix(2, 1, [x[1][l1_n], x[2][l1_n]]) : X[l1_n] = matrix(X[l1_sv]) ;

(b) The matrix [ A ] – [ U ] is

> `f(A - l1U)` := subs(lambda=l[`1_n`], matrix(`A - lU`)) : A - l[`1_n`]*U = matrix(`f(A - l1U)`) ;

or, upon substitution of the numerical value of ,

> `v(A - l1U)` := subs(l[`1_n`]=l[`1_v`], matrix(`f(A - l1U)`)) : A - l[`1_n`]*U = matrix(`v(A - l1U)`) ;

(c) The matrix equation [ A ] – [ U ] = [ 0 ] is

> matrix(`v(A - l1U)`) * matrix(X[l1_sv]) = matrix(`0`) ;

(d) The corresponding system of homogeneous equations is

> s_l1 := evalm(matrix(`v(A - l1U)`) &* matrix(X[l1_sv])) :

> Eq_1_l1 := s_l1[1,1] = 0 : Eq_1_l1 ; Eq_2_l1 := s_l1[2,1] = 0 : Eq_2_l1 ;

(e) The solution of the system is obtained as follows:

Assign an arbitrary non- value to , e.g.

> x[1][l1_v] := 1 : x[1][l1_n] = x[1][l1_v] ;

Substitute this value in, say, the second equation and solve it for the other unknown, i.e.

> x[2][l1_v] := solve(subs(x[1][l1_n]=1, Eq_2_l1)) : x[2][l1_n] = x[2][l1_v] ;

(f) Substitute the roots and in the symbolic expression for the eigenvector :

> X[l1_nv] := subs(x[1][l1_n]=x[1][l1_v], x[2][l1_n]=x[2][l1_v], matrix(X[l1_sv])) :

> X[l1_n] = matrix(X[l1_nv]) ;

(g) Check the resultant of the matrix operation [ A ] – [ U ]:

> A - l[`1_n`] * U = matrix(`v(A - l1U)`) * matrix(X[l1_nv]) ;

or

> A - l[`1_n`] * U = evalm(matrix(`v(A - l1U)`) &* matrix(X[l1_nv])) ;

* * *

N.B. For any arbitrary scalar , the product vector is also the eigenvector of [ A ], i.e.

> A - l[`1_n`] * U = mu[1] * X[l1_n] ; A - l[`1_n`] * U = mu[1] * matrix(X[l1_sv]) ;

since the resultant of the matrix equation

> A - l[`1_n`] * U = matrix(`v(A - l1U)`) * matrix(evalm(mu[1] * X[l1_nv])) ;

is

> A - l[`1_n`] * U = evalm(matrix(`v(A - l1U)`) &* matrix(evalm(mu[1]*X[l1_nv]))) ;

Notice, therefore, that the eigenvector of a matrix is not unique.

* * *

Step 6 . Determine the eigenvector of [ A ] corresponding to :

(a) The eigenvector is

> x[1][l2_n] := x[1][l[`2_n`]] : x[2][l2_n] := x[2][l[`2_n`]] : X[l2_n] := X[lambda[2]] :

> X[l2_sv] := matrix(2, 1, [x[1][l2_n], x[2][l2_n]]) : X[l2_n] = matrix(X[l2_sv]) ;

(b) The matrix [ A ] – [ U ] is

> `f(A - l2U)` := subs(lambda=l[`2_n`], matrix(`A - lU`)) : A - l[`2_n`]*U = matrix(`f(A - l2U)`) ;

or, upon substitution of the numerical value of ,

> `v(A - l2U)` := subs(l[`2_n`]=l[`2_v`], matrix(`f(A - l2U)`)) : A - l[`2_n`]*U = matrix(`v(A - l2U)`) ;

(c) The matrix equation [ A ] – [ U ] = [ 0 ] is

> matrix(`v(A - l2U)`) * matrix(X[l2_sv]) = matrix(`0`) ;

(d) The corresponding system of (identical) homogeneous equations is

> s_l2 := evalm(matrix(`v(A - l2U)`) &* matrix(X[l2_sv])) :

> Eq_1_l2 := s_l2[1,1] = 0 : Eq_1_l2 ; Eq_2_l2 := s_l2[2,1] = 0 : Eq_2_l2 ;

(e) The solution of the system is obtained as follows:

Assign an arbitrary non- value to , e.g.

> x[1][l2_v] := 1 : x[1][l2_n] = x[1][l2_v] ;

Substitute this value in, say, the second equation and solve it for the other unknown, i.e.

> x[2][l2_v] := solve(subs(x[1][l2_n]=1, Eq_2_l2)) : x[2][l2_n] = x[2][l2_v] ;

(f) Substitute the roots and in the symbolic expression for the eigenvector :

> X[l2_nv] := subs(x[1][l2_n]=x[1][l2_v], x[2][l2_n]=x[2][l2_v], matrix(X[l2_sv])) :

> X[l2_n] = matrix(X[l2_nv]) ;

(g) Check the resultant of the matrix operation [ A ] – [ U ]:

> A - l[`2_n`] * U = matrix(`v(A - l2U)`) * matrix(X[l2_nv]) ;

or

> A - l[`2_n`] * U = evalm(matrix(`v(A - l2U)`) &* matrix(X[l2_nv])) ;

* * *

N.B. For any arbitrary scalar , the product vector is also the eigenvector of [ A ], i.e.

> A - l[`2_n`] * U = mu[2] * X[l2_n] ; A - l[`2_n`] * U = mu[2] * matrix(X[l2_sv]) ;

since the resultant of the matrix equation

> A - l[`2_n`] * U = matrix(`v(A - l2U)`) * matrix(evalm(mu[2] * X[l2_nv])) ;

is

> A - l[`2_n`] * U = evalm(matrix(`v(A - l2U)`) &* matrix(evalm(mu[2]*X[l2_nv]))) ;

* * *

N.B. The eigenvectors may be obtained in Maple directly , using the eigenvects function. This function computes the eigenvalues and eigenvectors of a matrix. The result returned is a sequence of lists . For the matrix [ A ] analysed above, the eigenvects function yields

> `roots&vectors(A)` := eigenvects(A) : roots_and_vectors(A) = `roots&vectors(A)` ;

However, the eigenvects function returns an unordered sequence of its elements, i.e. lists may appear in the above or reverse order. This implies that the first eigenvalue in the sequence returned by the eigenvects function may not be the same as the first eigenvalue in the sequence returned by the eigenvals function, which returns eigenvalues in a unique, reproducible order.

On the other hand, it is essential that every eigenvalue be used together with the eigenvector corresponding to it. Therefore, it is necessary to develop a procedure that will extract pairs of lists in a unique way. This is ensured by the selection procedure that follows.

For the case under consideration, two different lists of eigenvalues and eigenvectors are extracted from the above solution in the desired order, viz.

> e[1] := `roots&vectors(A)`[1][1] : e[2] := `roots&vectors(A)`[2][1] :

• list containing the first eigenvalue:

> if e[1] = root1(A) then l := root1(A) : list1 := `roots&vectors(A)`[1] else l := e[2] : list1 := `roots&vectors(A)`[2] fi : list[1](roots_and_vectors(A)) = list1 ;

• list containing the second eigenvalue:

> if e[2] = root2(A) then l := root2(A) : list2 := `roots&vectors(A)`[2] else l := e[1] : list2 := `roots&vectors(A)`[1] fi : list[2](roots_and_vectors(A)) = list2 ;

In either list, the first number is an eigenvalue and the second value is its multiplicity , indicating how many times this eigenvalue appears in the solution to the characteristic equation of the matrix. This may be verified by inspection of the characteristic polynomial in its factored form. For this case,

> `charpoly(A)` := factor(charpoly(A, lambda)) : char_poly(A) = `charpoly(A)` ;

It is evident from the above that multiplicity of either eigenvalue is , or either eigenvalue is a simple root of the characteristic equation of [ A ].

The last object in the solution list is a set containing an eigenvector (or, in the case of matrices of a higher order, several eigenvectors separated with the comma) corresponding to the eigenvalue of a given list.

Extracting the eigenvalue and its corresponding eigenvector from the first list gives

> l1 := list1[1] : charvector[1] := list1[3][1] : v1 := charvector[1] :

> lambda[1] = l1 ; char_vector[lambda[1]](A) = eval(charvector[1]) ;

or, as a column matrix (vector),

> X[l1] := convert(eval(charvector[1]), matrix) : X[lambda[1]] = matrix(X[l1]) ;

Extracting the eigenvalue and its corresponding eigenvector from the second list gives

> l2 := list2[1] : charvector[2] := list2[3][1] : v2 := charvector[2] :

> lambda[2] = l2 ; char_vector[lambda[2]](A) = eval(charvector[2]) ;

or, as a column matrix (vector),

> X[l2] := convert(eval(charvector[2]), matrix) : X[lambda[2]] = matrix(X[l2]) ;

Both eigenvectors are equal to the eigenvectors obtained earlier "manually" and either of them is associated with a proper eigenvalue.

* * *

Since and are also eigenvectors of [ A ], this fact is often used to normalize eigenvectors of a matrix by setting the scalars to be equal to the the reciprocal of the Frobenius norm of the pertinent vector.

The Frobenius norm of a vector [ X ] having elements is defined to be the square root of the sum of the squares of the magnitudes of each element of the vector. In other words, it is the length (or magnitude) of the vector, which is given by

For the case under analysis, the Frobenius norm or length of eigenvectors and may be computed as follows.

• for eigenvector corresponding to the first eigenvalue:

> `abs(v1)` := sqrt(v1[1]^2 + v1[2]^2) : abs(char_vector[lambda[1]](A)) = `abs(v1)` ;

or, directly by using the function norm together with the norm name frobenius , viz.

> `norm(v1)` := norm(v1, frobenius) : abs(char_vector[lambda[1]](A)) = `norm(v1)` ;

• for eigenvector corresponding to the second eigenvalue:

> `abs(v2)` := sqrt(v2[1]^2 + v2[2]^2) : abs(char_vector[lambda[2]](A)) = `abs(v2)` ;

or, directly,

> `norm(v2)` := norm(v2, frobenius) : abs(char_vector[lambda[2]](A)) = `norm(v2)` ;

Consequently, the scalars assume the following values:

> mu[1] := 1/`abs(v1)` : mu[2] := 1/`abs(v2)` : 'mu[1]' = mu[1] ; 'mu[2]' = mu[2] ;

The normalized characteristic vectors of [ A ] are

> `v1[N]` := evalm(mu[1] * v1) : `v2[N]` := evalm(mu[2] * v2) : mu[1] := 'mu[1]' : mu[2] := 'mu[2]' :

> char_vector[lambda[1]][N](A) = op(`v1[N]`) ; char_vector[lambda[2]][N](A) = op(`v2[N]`) ;

Alternatively, normalization of vectors may be performed in Maple directly , using the normalize function. In this particular case, it gives

> `v1[N]` := normalize(v1) : `v2[N]` := normalize(v2) :

> char_vector[lambda[1]][N](A) = op(`v1[N]`) ; char_vector[lambda[2]][N](A) = op(`v2[N]`) ;

Note that the magnitude of normalized vectors is . Verification of this fact for the analysed case gives

> `abs(v1[N])` := sqrt(`v1[N]`[1]^2 + `v1[N]`[2]^2) : `abs(v2[N])` := sqrt(`v2[N]`[1]^2 + `v2[N]`[2]^2) :

> abs(char_vector[lambda[1]][N](A)) = `abs(v1[N])` ; abs(char_vector[lambda[2]][N](A)) = `abs(v2[N])` ;

* * *

N.B. If one of the eigenvalues of a matrix is , its determinant is and the matrix is singular.

For example, consider a ( × ) matrix [ A ] given as

> A := matrix(3, 3, [1, 2, 5, 3, 1, 5, -5, 0, -5]) : A = matrix(A) ;

(a) The eigenvalues of [ A ] are

> charroots(A) := eigenvals(A) : char_roots(A) = charroots(A) ;

(b) The determinant of [ A ] is

> `det(A)` := det(A) : Det(A) = `det(A)` ;

* * *

N.B. If is an eigenvalue of a matrix [ A ], then an eigenvalue of [ A ] is .

For example, consider a ( × ) matrix [ A ] given as

> A := matrix(3, 3, [-2, 2, 0, 2, 1, 0, -2, -1, -1]) : A = matrix(A) ;

(a) The eigenvalues of [ A ] are

> charroots(A) := eigenvals(A) : char_roots(A) = charroots(A) ;

(b) The eigenvalues of [ A ] are

> charroots(kA) := eigenvals(evalm(k*A)) : char_roots(k*A) = charroots(kA) ;

* * *

N.B. If , , ..., are all eigenvectors of a matrix [ A ] corresponding to the same eigenvalue , then any non- linear combination of these vectors is also an eigenvector of [ A ] corresponding to . This implies that the vector

+ ... +

is a solution vector to the equation

([ A ] – [ U ]) = [ 0 ]

with the selected eigenvalue .

Coefficients , , ..., are arbitrary scalars.

For example, consider a ( × ) matrix [ A ] given as

> A := matrix(3, 3, [-2, 2, 3, 2, 1, 6, 3, 6, 6]) : A = matrix(A) ;

(a) The eigenvalues of [ A ] are

> charroots(A) := eigenvals(A) : char_roots(A) = charroots(A) ;

from which two distinct roots are extracted, viz.

> root1(A) := charroots(A)[1] : root2(A) := charroots(A)[2] :

> char_root[1](A) = root1(A) ; char_root[2](A) = root2(A) ;

The eigenvalue is a double root (root of multiplicity two), as clearly seen from inspection of the characteristic polynomial of [ A ] in its factored form

> `charpoly(A)` := factor(charpoly(A, lambda)) : char_poly(A) = `charpoly(A)` ;

(b) The sequence of lists containing eigenvalues and eigenvectors of [ A ] is

> `roots&vectors(A)` := eigenvects(A) : roots_and_vectors(A) = `roots&vectors(A)` ;

It can be seen from the above sequence that only one eigenvalue ( ) is associated with two eigenvectors. This eigenvalue together with its multiplicity and both corresponding eigenvectors are contained in one of the two lists of the sequence.

Since only and its corresponding eigenvectors are of interest in this case, it is necessary to develop a selection procedure, which will extract these items from the sequence irrespective of the actual order of the lists in it. This is performed under (c).

(c) The eigenvalue involved and the sequence containing eigenvectors corresponding to it are

> e[1] := `roots&vectors(A)`[1][1] : e[2] := `roots&vectors(A)`[2][1] :

> if e[1] = root1(A) then l := root2(A) : vectors(A) := op(`roots&vectors(A)`[2][3]) else l := e[1] : vectors(A) := op(`roots&vectors(A)`[1][3]) fi :

> l := l : lambda = l ; char_vectors[lambda](A) = vectors(A) ;

(d) Assigning the names and to the eigenvectors of the above sequence, extracting them, and converting into column matrices give

> X[1][l] := convert(vectors(A)[1], matrix) : X[2][l] := convert(vectors(A)[2], matrix) :

> X[1][lambda] = matrix(X[1][l]) ; X[2][lambda] = matrix(X[2][l]) ;

(e) Let a linear combination of the two eigenvectors be of the form

> X[lambda] = mu[1] * X[1][lambda] + mu[2] * X[2][lambda] ;

Substituting both eigenvectors with their respective numerical elements yields

> X[lambda] = mu[1] * matrix(X[1][l]) + mu[2] * matrix(X[2][l]) ;

Evaluation of gives

> X[l] := evalm(mu[1]*X[1][l] + mu[2]*X[2][l]) : X[lambda] = matrix(X[l]) ;

(f) Let the unit matrix [ U ] appropriately sized for this case be

> U := diag(1,1,1) : U = matrix(U) ;

(g) Substituting the respective matrices and eigenvalue into ( [ A ] – [ U ] ) gives

> '(A - lambda*U)' * X[lambda] = (matrix(A) - l * matrix(U)) * matrix(X[l]) ;

Evaluation of ( [ A ] – [ U ] ) yields the zero column matrix (vector)

> (matrix(A) - l * matrix(U)) * matrix(X[l]) = evalm((A - l * U) &* X[l]) ;

irrespective of numerical values of the coefficients and .

* * *

N.B. Raising a matrix to a positive or negative integer power does not change its eigenvectors, but raises its eigenvalues to the same power.

For example, consider a ( × ) matrix [ A ] given as

> A := matrix(2, 2, [5, 8, 1, 7]) : A = matrix(A) ;

Let the positive and negative exponents be

> p := 3 : n := -2 : 'p' = p ; 'n' = n ;

(a) The eigenvalues of [ A ] are

> charroots(A) := eigenvals(A) : char_roots(A) = charroots(A) ;

whence

> root1(A) := charroots(A)[1] : root2(A) := charroots(A)[2] :

> char_root[1](A) = root1(A) ; char_root[2](A) = root2(A) ;

(b) The sequence of lists containing eigenvalues and eigenvectors of [ A ] is

> `roots&vectors(A)` := eigenvects(A) : roots_and_vectors(A) = `roots&vectors(A)` ;

(c) The individual eigenvalues and corresponding eigenvectors are

> e[1] := `roots&vectors(A)`[1][1] : e[2] := `roots&vectors(A)`[2][1] :

• pair consisting of the first eigenvalue and the corresponding eigenvector:

> if e[1] = root1(A) then l := root1(A) : `vector(A)` := op(`roots&vectors(A)`[1][3]) else l := e[2] : `vector(A)` := op(`roots&vectors(A)`[2][3]) fi :

> l1 := l : lambda[1] = l1 ; char_vector[lambda[1]](A) = op(`vector(A)`) ;

or, as a column matrix (vector),

> X[l1] := convert(`vector(A)`, matrix) : X[lambda[1]] = matrix(X[l1]) ;

• pair consisting of the second eigenvalue and the corresponding eigenvector:

> if e[2] = root2(A) then l := root2(A) : `vector(A)` := op(`roots&vectors(A)`[2][3]) else l := e[1] : `vector(A)` := op(`roots&vectors(A)`[1][3]) fi :

> l2:=l : lambda[2]=l2 ; char_vector[lambda[2]](A) = op(`vector(A)`) ; l:='l' : e[1]:='e[1]' : e[2]:='e[2]':

or, as a column matrix (vector),

> X[l2] := convert(`vector(A)`, matrix) : X[lambda[2]] = matrix(X[l2]) ;

(d) Raising the matrix [ A ] to the power yields the following matrix:

> `A^3` := evalm(A^3) : A^3 = matrix(`A^3`) ;

(e) The eigenvalues of [ A ] ^ are

> charroots(A^3) := eigenvals(`A^3`) : char_roots(A^3) = charroots(A^3) ;

whence

> root1(A^3) := charroots(A^3)[1] : root2(A^3) := charroots(A^3)[2] :

> char_root[1](A^3) = root1(A^3) ; char_root[2](A^3) = root2(A^3) ;

For a comparison, the eigenvalues of [ A ] raised to the power are

> [char_root[1](A)]^p = root1(A)^p ; [char_root[2](A)]^p = root2(A)^p ;

It can be seen from this comparison that either eigenvalue of [ A ] ^ is the cube of the corresponding eigenvalue of [ A ].

(f) The sequence of lists containing eigenvalues and eigenvectors of [ A ] ^ is

> `roots&vectors(A^3)` := eigenvects(`A^3`) : roots_and_vectors(A^3) = `roots&vectors(A^3)` ;

It follows from the inspection of the above that the eigenvectors of [ A ] ^ are the same as those of matrix [ A ].

(g) Raising the matrix [ A ] to the power yields the following matrix:

> `A^(-2)` := evalm(A^(-2)) : A^` -2` = matrix(`A^(-2)`) ;

(h) The eigenvalues of [ A ] ^( ) are

> charroots(A^(-2)) := eigenvals(`A^(-2)`) : char_roots(A^` -2`) = charroots(A^(-2)) ;

whence

> root1(A^(-2)) := charroots(A^(-2))[1] : root2(A^(-2)) := charroots(A^(-2))[2] :

> char_root[1](A^` -2`) = root1(A^(-2)) ; char_root[2](A^` -2`) = root2(A^(-2)) ;

For a comparison, the eigenvalues of [ A ] raised to the power are

> '[char_root[1](A)]'^` -2` = root1(A)^n ; '[char_root[2](A)]'^` -2` = root2(A)^n ;

It can be seen from this comparison that either eigenvalue of [ A ] ^( ) is the square of the reciprocal of the corresponding eigenvalue of [ A ].

(i) The sequence of lists containing eigenvalues and eigenvectors of [ A ] ^( ) is

> `roots&vectors(A^(-2))` := eigenvects(`A^(-2)`) : roots_and_vectors(A^` -2`) = `roots&vectors(A^(-2))` ;

It follows from the inspection of the above that the eigenvectors of [ A ] ^( ) are the same as those of matrix [ A ].

* * *

N.B. If is an eigenvalue of a matrix [ A ] with corresponding eigenvector , then for any scalar , is an eigenvector of [ A ] – [ U ] corresponding to the eigenvalue .

For example, consider a ( × ) matrix [ A ] given as

> A := matrix(2, 2, [-1, -5, 2, 6]) : A = matrix(A) ;

(a) The eigenvalues of [ A ] are

> charroots(A) := eigenvals(A) : char_roots(A) = charroots(A) ;

whence

> root1(A) := charroots(A)[1] : root2(A) := charroots(A)[2] :

> char_root[1](A) = root1(A) ; char_root[2](A) = root2(A) ;

(b) The sequence of lists containing eigenvalues and eigenvectors of [ A ] is

> `roots&vectors(A)` := eigenvects(A) : roots_and_vectors(A) = `roots&vectors(A)` ;

(c) Choosing and extracting, for instance, the smaller eigenvalue of [ A ] give

> e[1](A) := `roots&vectors(A)`[1][1] : e[2](A) := `roots&vectors(A)`[2][1] :

> l(A) := min(e[1](A), e[2](A)) : lambda[A] = l(A) ;

(d) The corresponding eigenvector of [ A ] is

> if e[1](A) = l(A) then `vector(A)` := op(`roots&vectors(A)`[1][3]) else `vector(A)` := op(`roots&vectors(A)`[2][3]) fi : char_vector[lambda](A) = op(`vector(A)`) ;

or, as a column matrix (vector),

> X[l(A)] := convert(`vector(A)`, matrix) : X[lambda[A]] = matrix(X[l(A)]) ;

(e) Let the unit matrix [ U ] appropriately sized for this case be

> U := diag(1,1) : U = matrix(U) ;

(f) Let a ( × ) matrix [ B ] be defined as

> B := A - mu*U : 'B' = B ;

Substituting [ A ] and [ U ] with their numerical values and evaluating [ B ] yield

> B := evalm(B) : B = matrix(B) ;

(g) The eigenvalues of [ B ] are

> charroots(B) := eigenvals(B) : char_roots(B) = charroots(B) ;

whence

> root1(B) := charroots(B)[1] : root2(B) := charroots(B)[2] :

> char_root[1](B) = root1(B) ; char_root[2](B) = root2(B) ;