potential

 potential.m computes the expected number of visits to j starting at i.  It infers locations where R(i,j) should be infinity by checking whether P^5000 has entries significantly larger than 0.  If so, R(i,j) is set to infinity.

0001 % potential.m computes the expected number of visits to j starting at i.  It infers locations where R(i,j) should be infinity by checking whether P^5000 has entries significantly larger than 0.  If so, R(i,j) is set to infinity.
0002 
0003 % define the transition matrix P before running potential.m
0004 
0005 R = eye(size(P));                 % an identity matrix
0006 A = P;                            % keep track of powers of P
0007 
0008 for n=1:5000,
0009   R = R + A;                      % cumulative sum of powers of P
0010   A = A * P;                      % calculate the next power of P
0011 end
0012 
0013 [s,t] = size(A);
0014 
0015 for i = 1:s,
0016   for j = 1:t,
0017     if A(i,j) > 0.000001,
0018       R(i,j) = Inf;
0019     end
0020   end
0021 end
0022 
0023 showmatrix(R);
0024 
0025 F = [];
0026 
0027 for i = 1:s,
0028   for j = 1:t,
0029     if j == i,
0030       F(i,j) = 1 - 1/R(i,i);
0031     else
0032       F(i,j) = R(i,j) / R(j,j);
0033     end
0034   end
0035 end
0036 
0037 showmatrix(F);