published elec. at
www.maplesoft.com
A PROCEDURE FOR FINDING
THE kTH POWER OF A MATRIX
Branko Malešević
Faculty of Electrical Engineering,
University of Belgrade, Serbia
malesh@EUnet.yu
Ivana Jovović
PhD student, Faculty of Mathematics,
University of Belgrade, Serbia
ivana121@EUnet.yu
1 Introduction
This worksheet demonstrates the use of Maple in Linear Algebra.
We give a new procedure (PowerMatrix) in Maple for finding the kth power of n-by-n square
matrix A, in a symbolic form, for any positive integer k, k ≥ n. The algorithm is based on an
application of Cayley-Hamilton theorem. We used the fact that the entries of the matrix Ak satisfy
the same recurrence relation which is determined by the characteristic polynomial of the matrix A
(see [1]). The order of these recurrences is n− d, where d is the lowest degree of the characteristic
polynomial of the matrix A.
For non-singular matrices the procedure can be extended for k not only a positive integer.
2 Initialization
> restart:
with(LinearAlgebra):
2.1 Procedure Definition
2.1.1 PowerMatrix
Input data are a square matrix A and a parameter k. Elements of the matrix A can be numbers
and/or parameters. The parameter k can take numeric value or be a symbol. The output data is
the kth power of the matrix. The procedure PowerMatrix is as powerful as the procedure rsolve.
> PowerMatrix := proc(A::Matrix,k)
local i,j,m,r,q,n,d,f,P,F,C;
P := x->CharacteristicPolynomial(A,x);
n := degree(P(x),x);
d := ldegree(P(x),x);
http://arxiv.org/abs/0704.0755v2
http://www.maplesoft.com/
mailto:malesh@EUnet.yu
mailto:ivana121@EUnet.yu
F := (i,j)->rsolve(sum(coeff(P(x),x,m)*f(m+q),m=0..n)=0,seq(f(r)=(A^r)[i,j],
r=d+1..n),f);
C := q->Matrix(n,n,F);
if (type(k,integer)) then return(simplify(A^k)) elif (Determinant(A)=0 and
not type(k,numeric)) then printf("The %ath power of the matrix for %a>=%d:",
k,k,n) elif (Determinant(A)=0 and type(k,numeric)) then return(simplify(A^k)) fi;
return(simplify(subs(q=k,C(q))));
3 Examples
3.1 Example 1.
> A := Matrix([[4,-2,2],[-5,7,-5],[-6,6,-4]]);
4 −2 2
−5 7 −5
−6 6 −4
> PowerMatrix(A,k);
−2k + 2 · 3k 2(1+k) − 2 · 3k −2(1+k) + 2 · 3k
−5 · 3k + 5 · 2k 5 · 3k − 4 · 2k −5 · 3k + 5 · 2k
6 · 2k − 6 · 3k −6 · 2k + 6 · 3k −6 · 3k + 7 · 2k
> Determinant(A);
> B := A^(-1);
> PowerMatrix(B,k);
−2(−k) + 2 · 3(−k) 2(1−k) − 2 · 3(−k) −2(1−k) + 2 · 3(−k)
−5 · 3(−k) + 5 · 2(−k) 5 · 3(−k) − 4 · 2(−k) −5 · 3(−k) + 5 · 2(−k)
−6 · 3(−k) + 6 · 2(−k) −6 · 2(−k) + 6 · 3(−k) −6 · 3(−k) + 7 · 2(−k)
3.2 Example 2.
> A := Matrix([[1-p,p],[p,1-p]]);
1− p p
p 1− p
> PowerMatrix(A,k);
(1− 2 p)k
(1− 2 p)k
(1− 2 p)k
(1− 2 p)k
The example is from [4], page 272, exercise 19.
3.3 Example 3.
> A := Matrix([[a,b,c],[d,e,f],[g,h,i]]);
a b c
d e f
g h i
> PowerMatrix(A,k)[1,1];
R = RootOf( (gbf + hdc + iea − gce − hfa− idb) Z3 + (gc + hf + db − ie − ia − ea) Z2 + (i+ e + a) Z − 1 )
R2ie− R2hf − Re− Ri+ 1
(3 R2gbf + 3 R2hdc+ 3 R2iea
−3 R2gce− 3 R2hfa− 3 R2idb+ 2 Rgc+ 2 Rhf + 2 Rdb− 2 Rie− 2 Ria− 2 Rea+ i+ e+ a) R
# Warning!
In this example MatrixPower and MatrixFuction procedures cannot be done in real-time.
# MatrixPower(A,k)[1,1];
# MatrixFunction(A,v^k,v)[1,1];
3.4 Example 4.
> A := Matrix([[0,0,1,0,1],[1,0,0,0,1],[0,0,0,1,1],[0,1,0,0,1],[1,1,1,1,0]]);
0 0 1 0 1
1 0 0 0 1
0 0 0 1 1
0 1 0 0 1
1 1 1 1 0
> PowerMatrix(A,k)[1,5];
Replace ’:’ with ’;’ and see result!
> MatrixPower(A,k)[1,5]:
> assume(m::integer):simplify(MatrixPower(A,k)[1,5]):
The example is from [3], page 101.
3.5 Example 5. and Example 6.
Pay attention what happens for singular matrices.
3.5.1 Example 5.
> A := Matrix([[0,2,1,3],[0,0,-2,4],[0,0,0,5],[0,0,0,0]]);
0 2 1 3
0 0 −2 4
0 0 0 5
0 0 0 0
> PowerMatrix(A,2);
0 0 −4 13
0 0 0 −10
0 0 0 0
0 0 0 0
> PowerMatrix(A,3);
0 0 0 −20
0 0 0 0
0 0 0 0
0 0 0 0
> PowerMatrix(A,k);
The kth power of the matrix A for k ≥ 4:
0 0 0 0
0 0 0 0
0 0 0 0
0 0 0 0
> MatrixPower(A,k);
Error, (in LinearAlgebra:-LA_Main:-MatrixPower)
power k is not defined for this Matrix
> MatrixFunction(A,v^k,v);
Error, (in LinearAlgebra:-LA_Main:-MatrixFunction)
Matrix function vk is not defined for this Matrix
The example is from [2], page 151, exercise 23.
3.5.2 Example 6.
> A := Matrix([[1,1,1,0],[1,1,1,-1],[0,0,-1,1],[0,0,1,-1]]);
1 1 1 0
1 1 1 −1
0 0 −1 1
0 0 1 −1
> PowerMatrix(A,k);
> The kth power of the matrix for k ≥ 4:
2(−1+k) 2(−1+k)
(−1)(1+k) · 2k
5 · 2k
(−1)k · 2k
2(−1+k) 2(−1+k)
5 · 2k
5 · (−1)(1+k) · 2k
5 · (−1)k · 2k
0 0 (−1)k · 2(−1+k) (−1)(1+k) · 2(−1+k)
0 0 (−1)(1+k) · 2(−1+k) (−1)k · 2(−1+k)
> MatrixPower(A,k);
Error, (in LinearAlgebra:-LA_Main:-MatrixPower)
power k is not defined for this Matrix
> MatrixFunction(A,v^k,v);
Error, (in LinearAlgebra:-LA_Main:-MatrixFunction)
Matrix function vk is not defined for this Matrix
4 References
[1] BrankoMalešević: Some combinatorial aspects of the composition of a set of functions, NSJOM
2006 (36), 3-9, URLs: http://www.im.ns.ac.yu/NSJOM/Papers/36 1/NSJOM 36 1 003 009.pdf,
http://arxiv.org/abs/math.CO/0409287.
[2] JohnB. Johnston, G.BaleyPrice, Fred S.Van Vleck: Linear Equations and Matrices, Addi-
son-Wesley, 1966.
[3] Carl D.Meyer: Matrix Analysis and Applied Linear Algebra Book and Solutions Manual SIAM,
2001.
[4] Robert Messer: Linear Algebra Gateway to Mathematics, New York, Harper-Collins College
Publisher, 1993.
5 Conclusions
This procedure has an educational character. It is an interesting demonstration for finding the kth
power of a matrix in a symbolic form. Sometimes, it gives solutions in the better form than the
existing procedure MatrixPower (see example 4.). See also example 5. and example 6., where we
consider singular matrices. In these cases the procedure MatrixPower does not give a solution. The
procedure PowerMatrix calculates the kth power of any singular matrices. In some examples it is
possible to get a solution in the better form with using the procedure allvalues (see example 3.).
Legal Notice: The copyright for this application is owned by the authors. Neither Maplesoft nor the
author are responsible for any errors contained within and are not liable for any damages resulting
from the use of this material. This application is intended for non-commercial, non-profit use only.
Contact the author for permission if you wish to use this application in for-profit activities.
http://www.im.ns.ac.yu/NSJOM/Papers/36_1/NSJOM_36_1_003_009.pdf
http://arxiv.org/abs/math.CO/0409287
	Introduction 
	Initialization
	Procedure Definition
	PowerMatrix
	Examples
	Example 1.
	Example 2.
	Example 3.
	Example 4.
	Example 5. and Example 6.
	Example 5. 
	Example 6.
	References
	Conclusions

  We give a new procedure in Maple for finding the k-th power of a martix. The
algorithm is based on the article [1].

Introduction
This worksheet demonstrates the use of Maple in Linear Algebra.
We give a new procedure (PowerMatrix) in Maple for finding the kth power of n-by-n square
matrix A, in a symbolic form, for any positive integer k, k ≥ n. The algorithm is based on an
application of Cayley-Hamilton theorem. We used the fact that the entries of the matrix Ak satisfy
the same recurrence relation which is determined by the characteristic polynomial of the matrix A
(see [1]). The order of these recurrences is n− d, where d is the lowest degree of the characteristic
polynomial of the matrix A.
For non-singular matrices the procedure can be extended for k not only a positive integer.
2 Initialization
> restart:
with(LinearAlgebra):
2.1 Procedure Definition
2.1.1 PowerMatrix
Input data are a square matrix A and a parameter k. Elements of the matrix A can be numbers
and/or parameters. The parameter k can take numeric value or be a symbol. The output data is
the kth power of the matrix. The procedure PowerMatrix is as powerful as the procedure rsolve.
> PowerMatrix := proc(A::Matrix,k)
local i,j,m,r,q,n,d,f,P,F,C;
P := x->CharacteristicPolynomial(A,x);
n := degree(P(x),x);
d := ldegree(P(x),x);
http://arxiv.org/abs/0704.0755v2
http://www.maplesoft.com/
mailto:malesh@EUnet.yu
mailto:ivana121@EUnet.yu
F := (i,j)->rsolve(sum(coeff(P(x),x,m)*f(m+q),m=0..n)=0,seq(f(r)=(A^r)[i,j],
r=d+1..n),f);
C := q->Matrix(n,n,F);
if (type(k,integer)) then return(simplify(A^k)) elif (Determinant(A)=0 and
not type(k,numeric)) then printf("The %ath power of the matrix for %a>=%d:",
k,k,n) elif (Determinant(A)=0 and type(k,numeric)) then return(simplify(A^k)) fi;
return(simplify(subs(q=k,C(q))));
3 Examples
3.1 Example 1.
> A := Matrix([[4,-2,2],[-5,7,-5],[-6,6,-4]]);
4 −2 2
−5 7 −5
−6 6 −4
> PowerMatrix(A,k);
−2k + 2 · 3k 2(1+k) − 2 · 3k −2(1+k) + 2 · 3k
−5 · 3k + 5 · 2k 5 · 3k − 4 · 2k −5 · 3k + 5 · 2k
6 · 2k − 6 · 3k −6 · 2k + 6 · 3k −6 · 3k + 7 · 2k
> Determinant(A);
> B := A^(-1);
> PowerMatrix(B,k);
−2(−k) + 2 · 3(−k) 2(1−k) − 2 · 3(−k) −2(1−k) + 2 · 3(−k)
−5 · 3(−k) + 5 · 2(−k) 5 · 3(−k) − 4 · 2(−k) −5 · 3(−k) + 5 · 2(−k)
−6 · 3(−k) + 6 · 2(−k) −6 · 2(−k) + 6 · 3(−k) −6 · 3(−k) + 7 · 2(−k)
3.2 Example 2.
> A := Matrix([[1-p,p],[p,1-p]]);
1− p p
p 1− p
> PowerMatrix(A,k);
(1− 2 p)k
(1− 2 p)k
(1− 2 p)k
(1− 2 p)k
The example is from [4], page 272, exercise 19.
3.3 Example 3.
> A := Matrix([[a,b,c],[d,e,f],[g,h,i]]);
a b c
d e f
g h i
> PowerMatrix(A,k)[1,1];
R = RootOf( (gbf + hdc + iea − gce − hfa− idb) Z3 + (gc + hf + db − ie − ia − ea) Z2 + (i+ e + a) Z − 1 )
R2ie− R2hf − Re− Ri+ 1
(3 R2gbf + 3 R2hdc+ 3 R2iea
−3 R2gce− 3 R2hfa− 3 R2idb+ 2 Rgc+ 2 Rhf + 2 Rdb− 2 Rie− 2 Ria− 2 Rea+ i+ e+ a) R
# Warning!
In this example MatrixPower and MatrixFuction procedures cannot be done in real-time.
# MatrixPower(A,k)[1,1];
# MatrixFunction(A,v^k,v)[1,1];
3.4 Example 4.
> A := Matrix([[0,0,1,0,1],[1,0,0,0,1],[0,0,0,1,1],[0,1,0,0,1],[1,1,1,1,0]]);
0 0 1 0 1
1 0 0 0 1
0 0 0 1 1
0 1 0 0 1
1 1 1 1 0
> PowerMatrix(A,k)[1,5];
Replace ’:’ with ’;’ and see result!
> MatrixPower(A,k)[1,5]:
> assume(m::integer):simplify(MatrixPower(A,k)[1,5]):
The example is from [3], page 101.
3.5 Example 5. and Example 6.
Pay attention what happens for singular matrices.
3.5.1 Example 5.
> A := Matrix([[0,2,1,3],[0,0,-2,4],[0,0,0,5],[0,0,0,0]]);
0 2 1 3
0 0 −2 4
0 0 0 5
0 0 0 0
> PowerMatrix(A,2);
0 0 −4 13
0 0 0 −10
0 0 0 0
0 0 0 0
> PowerMatrix(A,3);
0 0 0 −20
0 0 0 0
0 0 0 0
0 0 0 0
> PowerMatrix(A,k);
The kth power of the matrix A for k ≥ 4:
0 0 0 0
0 0 0 0
0 0 0 0
0 0 0 0
> MatrixPower(A,k);
Error, (in LinearAlgebra:-LA_Main:-MatrixPower)
power k is not defined for this Matrix
> MatrixFunction(A,v^k,v);
Error, (in LinearAlgebra:-LA_Main:-MatrixFunction)
Matrix function vk is not defined for this Matrix
The example is from [2], page 151, exercise 23.
3.5.2 Example 6.
> A := Matrix([[1,1,1,0],[1,1,1,-1],[0,0,-1,1],[0,0,1,-1]]);
1 1 1 0
1 1 1 −1
0 0 −1 1
0 0 1 −1
> PowerMatrix(A,k);
> The kth power of the matrix for k ≥ 4:
2(−1+k) 2(−1+k)
(−1)(1+k) · 2k
5 · 2k
(−1)k · 2k
2(−1+k) 2(−1+k)
5 · 2k
5 · (−1)(1+k) · 2k
5 · (−1)k · 2k
0 0 (−1)k · 2(−1+k) (−1)(1+k) · 2(−1+k)
0 0 (−1)(1+k) · 2(−1+k) (−1)k · 2(−1+k)
> MatrixPower(A,k);
Error, (in LinearAlgebra:-LA_Main:-MatrixPower)
power k is not defined for this Matrix
> MatrixFunction(A,v^k,v);
Error, (in LinearAlgebra:-LA_Main:-MatrixFunction)
Matrix function vk is not defined for this Matrix
4 References
[1] BrankoMalešević: Some combinatorial aspects of the composition of a set of functions, NSJOM
2006 (36), 3-9, URLs: http://www.im.ns.ac.yu/NSJOM/Papers/36 1/NSJOM 36 1 003 009.pdf,
http://arxiv.org/abs/math.CO/0409287.
[2] JohnB. Johnston, G.BaleyPrice, Fred S.Van Vleck: Linear Equations and Matrices, Addi-
son-Wesley, 1966.
[3] Carl D.Meyer: Matrix Analysis and Applied Linear Algebra Book and Solutions Manual SIAM,
2001.
[4] Robert Messer: Linear Algebra Gateway to Mathematics, New York, Harper-Collins College
Publisher, 1993.
5 Conclusions
This procedure has an educational character. It is an interesting demonstration for finding the kth
power of a matrix in a symbolic form. Sometimes, it gives solutions in the better form than the
existing procedure MatrixPower (see example 4.). See also example 5. and example 6., where we
consider singular matrices. In these cases the procedure MatrixPower does not give a solution. The
procedure PowerMatrix calculates the kth power of any singular matrices. In some examples it is
possible to get a solution in the better form with using the procedure allvalues (see example 3.).
Legal Notice: The copyright for this application is owned by the authors. Neither Maplesoft nor the
author are responsible for any errors contained within and are not liable for any damages resulting
from the use of this material. This application is intended for non-commercial, non-profit use only.
Contact the author for permission if you wish to use this application in for-profit activities.
http://www.im.ns.ac.yu/NSJOM/Papers/36_1/NSJOM_36_1_003_009.pdf
http://arxiv.org/abs/math.CO/0409287
	Introduction 
	Initialization
	Procedure Definition
	PowerMatrix
	Examples
	Example 1.
	Example 2.
	Example 3.
	Example 4.
	Example 5. and Example 6.
	Example 5. 
	Example 6.
	References
	Conclusions
