MATLAB Answers

# Why does matrix power give wrong answer?

3 views (last 30 days)
Joseph on 26 Feb 2014
Commented: per isakson on 1 Mar 2014
I am trying to do a MCMC simulation and need to calculate a small frequency event over a long time. To do this I need to raise
A=[1-1.01E-11,1.01E-11,0;1-3.024E-11,0,3.024E-11;0,0,1];
to a large power (10^10). I am particularly interested in the (1,3) element of the result (the probability of transferring between states 1 and 3 after 10 billion iterations). After 10^9 the (1,3) element is 3.05423999691491e-13. However after 10^10 it is -3.05423999991124e-22. Which is impossible because I am multiplying matrices with strictly non-negative values. As a check I tried (A^(10^9))^10, and this gave a (1,3) entry of 3.05423999966373e-12. Doing the calculation (A^(10^10) using repeated squaring and multiplication) in Excel gave 3.05424E-12.
I am using MATLAB 2013b on windows 7 64-bit.
Why does A^(10^10) fail but (A^(10^9))^10 work?
##### 1 CommentShowHide None
John D'Errico on 26 Feb 2014
Time to learn about floating point arithmetic. When your numbers are only held accurately to roughly 16 digits, what do you expect?

Sign in to comment.

### Answers (1)

per isakson on 27 Feb 2014
Edited: per isakson on 27 Feb 2014
##### 4 CommentsShowHide 3 older comments
per isakson on 1 Mar 2014

Sign in to comment.

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!