#
#
# Written by:
# -- 
# John L. Weatherwax                2009-04-21
# 
# email: wax@alum.mit.edu
# 
# Please send comments and especially bug reports to the
# above email address.
# 
#-----

# fix the random.seed, so I'll always get the same answer:
#set.seed(10131985)

M = matrix(runif(4),nrow=2,ncol=2)
M

# sum the columns:
sc = apply(M,1,sum)

# divided each row by the corresponding number
DM = matrix(sc,nrow=2,ncol=2) # a "divide" matrix
Ms = M / DM

apply(Ms,1,sum) # verify rows sum to 1 ... they do

MP = Ms # start with power of 1 (doing nothing)
print(apply(MP,1,sum)) # display the row sum 
max_power = 10000 # now iterate
for( ii in 1:max_power ){
  MP = MP %*% Ms
  #print(apply(MP,1,sum)) # display the row sum 
}
print(apply(MP,1,sum)) # display the row sum