Questao 4 - Lista 8

markov.py

import numpy as np

# matriz de transicoes
P = np.array([[0.0, 1.0/3, 1.0/3, 1.0/3], [0.0, 0.0, 0.0, 1.0], [1.0/2, 0.0, 0.0, 1.0/2], [0.0, 0.0, 1.0, 0.0]])

def sistema_linear():
  # declarando matriz original
  A_original = np.array([[1.0,-(1.0/2),0.0,0.0], [-1.0/3,1.0,0.0,0.0], [-1.0/3,0.0,1.0,-1.0], [-1.0/3, -1.0, -1.0/2, 1.0]])

  # substituindo a primeira linha por [1, 1, 1, 1] para podermos resolver o sistema
  A = A_original.copy()
  np.put(A, [0, 1, 2, 3], [1, 1, 1, 1])

  b = np.array([1, 0, 0, 0]).T

  print "O vetor PI e:"
  print np.linalg.solve(A,b)

def calculaPk(k):
  return np.linalg.matrix_power(P, k)

def teste_k():
  erro = 0.0001

  i = 2
  
  P_k = calculaPk(i)
  P_k1 = calculaPk(i+1)
  
  while True:
    R = P_k1 - P_k
    diff_maxima = abs(R.max())
    diff_minima = abs(R.min())

    if diff_minima < erro and diff_maxima < erro:
      break
    
    i = i + 1
    P_k = P_k1
    P_k1 = calculaPk(i+1)
  
  print "Para uma diferenca de %f entre uma matriz e outra, K = %d" % (erro, i)
  print "A matriz de transicoes no caso estacionario eh: "
  print P_k

  return i

print "### Solucao da questao 2"
sistema_linear()
print

print "### Solucao da questao 3"
k = teste_k()
print

print "### Solucao da questao 4"
P_0 = np.array([ 0.3, 0.2, 0.4, 0.1 ]) # chutando pi_0
P_K = calculaPk(k)
print np.dot(P_0, P_K)
print