Page 22 - 4861
P. 22

a(N+j)=lamda+N*mu+j*nu;
                   end
                   for i=1:n
                      b(N+i)=N*mu+i*nu;
                   end
                   a(n+N)=N*mu+n*nu;
                   %Обчислення інфінітезимальної матриці
                   A(1,1)=-lamda;A(1,2)=b(1);
                   for i=3:N+n+1 A(1,i)=0; end
                   m=0;
                   for i=2:N+n
                           if i>2
                               m=m+1;
                             for ii=1:m
                                 A(i,ii)=0;
                             end
                               end
                   A(i,m+1)=lamda;A(i,m+2)=-a(i-1);A(i,m+3)=b(i);
                   if i<N+n
                          for q=m+4:N+n+1
                           A(i,q)=0;
                       end
                   end
                   end
                   for j=1:N+n-2 A(N+n+1,j)=0; end
                   A(N+n+1,N+n)=lamda;A(N+n+1,N+n+1)=-a(N+n);
                   %Виклик солвера для файл-функції, початкового і кінцевого моментів часу
             і
                   %вектора початкових умов
                   options=odeset('RelTol',E)
                   [T,X]=ode45(@rksol,[0 tk],X0,options,A);
                   %Вивід графіків роз'вязку моделі
                   subplot(2,2,1)
                   plot(T,X(:,1),'r',T,X(:,2),'c',T,X(:,3),'g');
                   %Пояснення до графіка
                   grid on
                   legend('P0','P1','P2')
                   subplot(2,2,2)
                   plot(T,X(:,4),'r',T,X(:,5),'c',T,X(:,6),'g');
                   %Пояснення до графіка
                   grid on
                   legend('P3','P4','P5')
                   subplot(2,2,3)
                   plot(T,X(:,7),'r',T,X(:,8),'c',T,X(:,9),'g');
                   %Пояснення до графіка
                   grid on
                   legend('P6','P7','P8')
                   subplot(2,2,4)
                   plot(T,X(:,10),'r',T,X(:,11),'c');
                   %Пояснення до графіка
                   grid on
                   legend('P9','P10')
                   %ОБЧИСЛЕННЯ ЕРГОДИЧНИХ ЙМОВІРНОСТЕЙ СТАНІВ СИСТЕМИ МО
                   alfa=lamda/mu;
                   hama=nu/mu;
                   %Обчислення N!
                   F=factorial(N);
                   %Обчислення Р0


                                                               21
   17   18   19   20   21   22   23   24   25   26   27