% Nettoyage
clf;
clear;

% Obtention et test de la matrice de transition.
P=rmarmat(3) % Création et affichage d'une matrice de transition 3x3.
sum(P,2)     % On vérifie que c'est bien une matrice de transition.
P>0          % Test de condition suffisante d'irréductibilité.
abs(eig(P))  % Spectre de P 
             % Si RI, alors Apériodique ssi une seule vp de module 1.

% Tracé d'un histogramme pour un temps grand fixé (loi de X_n).
% Cela revient juste à rdist puisque la loi cherchée est P^n(:,X_0).
n=500;       % temps à considérer (i.e. X_500)
m=10000;     % Nombre de réalisations.
xzero=1;     % Point initial.
Q=P^n;
for (i=1:m),
  realis(i)=rdist([1,2,3],Q(:,xzero));
end
figure(1);
hist(realis,m/100);

% Mesure invariante mu et tracé de l'erreur relative 
% de son approximation. mu est bien approximée par 
% P^n(:,X_0) quand n est grand et X_0 quelconque.
Q=P^1000;           
mu=Q(:,1);
max((mu*P-mu)./mu) % Maximum des erreurs relatives pour m=1000.
R=P^100;           
Q=R;
for (i=100:100:10000), 
 Q=Q*R;
 mu=Q(:,1);
 erreur(i)=max((mu*P-mu)./mu);
end
figure(2);
plot(erreur,'b-') % Tracé de l'erreur relative par pas de 100.

% Tracé d'une trajectoire particulière (une seule réal. de X(1),...,X(m)).
m=500;       % temps maximum à considérer 
x(1)=1;      % Point initial.
for (i=1:m-1),
  x(i+1)=rdist([1,2,3],P(:,x(i)));
end
figure(3);
plot(x,'b-');