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