function laba1

  % task 1

  P = [0.173, 0.287, 0.540; 0, 0.587, 0.413; 0, 0.159, 0.841];
  n_min = 100;
  delta = [];
  n = 2;
  printf(' P(1) = \n'), disp(P)

  while n <= n_min

    P(:, :, n) = P(:, :, n-1) * P(:, :, 1); % 3-merny massiv, kotory soderjit vse matricy stepeney ot 1 do n_min
    delta(n) = max(max(abs( P(:, :, n) - P(:, :, n-1) ))); % nahodim maksimalnoe otklonenie sredi vseh veroyatnostey mejdu nastoyacshey matricey i predyduschey

    if delta(n) < 0.00001
      n_min = n;
    endif

    n += 1;
  endwhile

  for n = 2:n_min

%    printf('\n\n P(%i) = \n', n);
%    disp(P(:, :, n));
%    printf(' \ndelta%i = %i\n', n, delta(n));

  endfor



  % task 2

  % q - stolbets
  % qP = q => qP - q = 0 => q(P - E) = 0, gde E - edinichnaya matrica, prichem v nashem kode eto eye(n)
  % q = (r1, r2, r3)
  % eye(n) - edinichnaya matrica n x n

  A = P(:, :, 1)' - eye(3); % vychitaem edinichnuyu matricu
  A = [A, zeros(3, 1)]; % pripisyvaem 4-y stolbets s nulevymi svobodnymi chlenami
  A(3,:) = ones(1, 4); % t.k. vectory lineyno zavisimye, to vycherkivaem poslednyuyu stroku
  % pri etom u nas est' uslovie normirovki: r1 + r2 + r3 = 1 -- stavim ego na mesto vycherknutoy stroki

  A = gauss(A); % reshaem SLAU

  printf('\n\n');
  q = A(:, 4)' % q - naydennoe stacionarnoe raspredelenie
  qP = q*P(:, :, 1) % Prover 04ka:

  if q == qP
    disp(' => q = qP')
  else
    disp(' => q =/= qP')
  endif



  % task 3

  % p(0) = ( p1(0), p2(0), p3(0) )
  % p(n) = p(0)*P^n
  % prichem v nashem kode P^n - eto "P(:, :, n)"

  states = eye(3); % potomu chto sostoyaniya (1, 0, 0), (0, 1, 0), (0, 0, 1) - eto stolbtsy edinichnoy matricy
  for p0 = states
    printf('\n\n');

%    p0 = p0'
%    delta0 = max(abs( p0 - q) )
%
%    ProbDistr(p0);

  endfor



  % task 4

  % orientiruyus' po kartinke na doske

  Seq = [];
  for i = 1:3

    Seq = [Seq; gen(i)];

  endfor
  %Seq(:, 1:100)



  size(Seq)
  % task 5
  all_R = R = [];
  %for state = 1:3
    for i = 1:3
    sum(i == Seq(state,:)
     % R = [R, sum(i == Seq(state,:)];
    endfor
    %all_R = [all_R; R]
  %endfor
R








  %_____________________________________________________________________________


  function ProbDistr(p0) % task 3

    n = 1;
    m_min = 100;
    while n <= m_min

      p(n, :) = p0 * P(:, :, n);
      delta_(n) = max(abs( p(n, :) - q ));

      if delta_(n) < 0.00001
        m_min = n;
      endif

      n += 1;
    endwhile

    for n = 1:m_min

      printf('\n\n  p(%i) = \n', n);
      disp(p(n, :));
      printf(' \n  delta%i = %i\n', n, delta_(n));

    endfor


   endfunction



   function I = gen(i) % task 4
      I = [];
      E = 0;

      for n = 1:1000
        E = rand(1);

        if E < P(i, 1, 1)
           i = 1;
         elseif E < P(i, 1, 1) + P(i, 2, 1)
           i = 2;
         else
           i = 3;
         endif

         I(n) = i;
      endfor

   endfunction


end






% for task 4
%        switch i
%
%          case 1
%            if E < P(1, 1, 1)
%              i = 1;
%            elseif E < P(1, 1, 1) + P(1, 2, 1)
%              i = 2;
%            else
%              i = 3;
%            endif
%            I(n) = i;
%
%
%          case 2
%            if E < P(2, 1, 1)
%              i = 1;
%            elseif E < P(2, 1, 1) + P(2, 2, 1)
%              i = 2;
%            else
%              i = 3;
%            endif
%            I(n) = i;
%
%
%          case 3
%            if E < P(3, 1, 1)
%              i = 1;
%            elseif E < P(3, 1, 1) + P(3, 2, 1)
%              i = 2;
%            else
%              i = 3;
%            endif
%            I(n) = i;