Code covered by the BSD License  

Highlights from
polyvalm2 - A faster matrix polynomial evaluator

5.0

5.0 | 1 rating Rate this file 0 Downloads (last 30 days) File Size: 4.47 KB File ID: #25780
image thumbnail

polyvalm2 - A faster matrix polynomial evaluator

by James Tursa

 

09 Nov 2009 (Updated 11 Nov 2009)

polyvalm2 evaluates a polynomial with a matrix argument faster than the MATLAB function polyvalm.

| Watch this File

File Information
Description

polyvalm2 evaluates a polynomial with a square matrix argument faster than the MATLAB built-in functions polyvalm or mpower.

Y = polyvalm2(P,X), when P is a vector of length N+1 whose elements are the coefficients of a polynomial, is the value of the polynomial evaluated with matrix argument X. X must be a square matrix.

       Y = P(1)*X^N + P(2)*X^(N-1) + ... + P(N)*X + P(N+1)*I

Class support for inputs P, X:
      float: double, single

The polyvalm2 speed improvements come from the following:

1) The MATLAB built-in function polyvalm uses Horner's method. polyvalm2 uses a binary decomposition of the matrix powers to do the calculation more efficiently, reducing the total number of matrix multiplies used to calculate the answer.

2) polyvalm calculates the product of a scalar times a matrix as the product of diag(scalar*ones(etc))*matrix ... i.e. it does a matrix multiply. polyvalm2 will calculate this more efficiently as the simple product scalar*matrix.

3) polyvalm does all of the calculations shown above, even if the coefficient P(i) is zero. polyvalm2 does not do calculations for P(i) coefficients that are zero.

4) polyvalm converts sparse matrix inputs into full matrices to do the calculations, whereas polyvalm2 keeps the intermediate calculations and the answer sparse.

The trade-off is that polyvalm2 uses more memory for intermediate variables than polyvalm, so for very large matrices polyvalm2 can run out of memory. In these cases polyvalm2 will abandon the efficient calculation method and just call the built-in polyvalm. For sparse matrix inputs, however, polyvalm2 will typically be more memory efficient than the MATLAB polyvalm function.

Caution: Since polyvalm2 uses different calculations to form the matrix powers, the end result may not match polyvalm exactly. Also, for the case where only the leading coefficient is non-zero, polyvalm2 may not match mpower exactly. But the answer will be just as accurate. And if there are inf's or NaN's involved, then the end result will not, in general, match polyvalm or mpower results. This should not be a great drawback to using polyvalm2, however, since even the MATLAB built-in functions polyvalm and mpower will not match each other in these cases. (By reordering the calculations, the NaN's propagate differently)

MATLAB release MATLAB 7.4 (R2007a)
Tags for This File  
Everyone's Tags
Tags I've Applied
Add New Tags Please login to tag files.
Comments and Ratings (3)
10 Nov 2009 Bruno Luong

I wish James could do a treatment for SPARSE matrix, because initialize an intermediate matrix with EYE is sometime not an option with large (sparse matrix).

11 Nov 2009 James Tursa

Bruno: Done. Thanks for the suggestion. polyvalm2 now keeps all the intermediate calculations and the answer as sparse for sparse inputs.

12 Nov 2009 Bruno Luong

Neatly coded. How often do you see a Matlab stock matrix function get beaten in speed?

Please login to add a comment or rating.
Updates
11 Nov 2009

Modified the code so that if an input matrix is sparse, the intermediate calculations and final answer will also be sparse. Note that this is different behavior from the MATLAB intrinsic polyvalm.

Tag Activity for this File
Tag Applied By Date/Time
polyvalm James Tursa 09 Nov 2009 10:33:48
mpower James Tursa 09 Nov 2009 10:33:48
polynomial James Tursa 09 Nov 2009 10:33:48
matrix James Tursa 09 Nov 2009 10:33:48

Contact us at files@mathworks.com