S3 |
λ13=0,4 |
|
|
λ31=3,333 |
λ14=0,2 |
λ32=1,333 |
S1 |
S4 |
|
λ21=0,333 |
|
|
|
λ41=4 |
S2 |
λ12=0,1 |
|
Рис. 2
2) Матрица интенсивностей в нашем случае имеет вид
( 12 13 14) |
12 |
13 |
14 |
|
|
|
|
21 |
( 21) |
0 |
0 |
|
|
|
|
|
||||
|
31 |
32 |
( 31 32) |
0 |
|
|
|
41 |
0 |
0 |
|
|
|
|
( 41) |
|
||||
0,7 |
0,1 |
0,4 |
0,2 |
|
|
|
|
|
|
|
|
0,333 |
0,333 |
0 |
0 |
. |
|
3,333 |
1,333 |
4,666 |
0 |
|
|
|
4 |
0 |
0 |
|
|
|
4 |
||||
3) Составим систему дифференциальных уравнений ЧепменаКолмогорова
d P(t) P(t) : dt
24
d |
P (t) ( |
|
|
|
|
)P (t) |
|
|
P (t) |
|
P (t) |
||||||||||||||||
|
|
|
|
21 |
31 |
||||||||||||||||||||||
dt |
1 |
12 |
|
13 |
14 |
1 |
|
2 |
|
3 |
|||||||||||||||||
|
|
|
|
|
|
d |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
||||
|
|
|
|
|
|
dt P2(t) 12P1(t) 21P2(t) 32P3(t), |
|||||||||||||||||||||
|
|
|
|
|
|
|
|||||||||||||||||||||
|
|
|
|
|
|
|
|
d |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
||
|
|
|
|
|
|
|
|
|
P (t) P (t) ( |
|
)P (t), |
||||||||||||||||
|
|
|
|
|
|
|
|
dt |
|
||||||||||||||||||
|
|
|
|
|
|
|
3 |
|
|
|
13 |
1 |
|
|
31 |
|
32 |
3 |
|
||||||||
|
|
|
|
|
|
|
|
|
|
|
|
d |
P (t) |
P (t) ( |
|
)P (t), |
|
||||||||||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
||||||||||||||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
||||||||||||
|
|
|
|
|
|
|
|
|
|
|
|
dt |
4 |
|
14 |
1 |
|
|
|
41 |
4 |
|
|
||||
|
|
|
|
|
|
|
|
|
|
|
|
P (t) P (t) P (t) P (t) 1. |
|
||||||||||||||
|
|
|
|
|
|
|
|
|
|
|
1 |
|
|
2 |
|
3 |
|
|
|
4 |
|
|
|
||||
|
|
d |
P (t) 0,7P (t) 0,333P (t) 3,333P (t) 4P |
||||||||||||||||||||||||
|
|
|
|
||||||||||||||||||||||||
|
|
dt |
1 |
|
|
|
|
|
|
|
|
|
|
|
1 |
|
|
2 |
|
|
|
|
3 |
|
4 |
||
|
|
|
|
d |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|||||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
||||||
|
|
|
|
dt P2 (t) 0,1P1(t) 0,333P2 (t) 1,333P3(t), |
|||||||||||||||||||||||
|
|
|
|
|
|||||||||||||||||||||||
|
|
|
|
|
|
|
|
|
|
|
d |
|
|
|
|
|
|
|
|
|
|
|
|
|
|
||
|
|
|
|
|
|
|
|
|
|
|
P (t) 0,4P (t) 4,666P (t), |
|
|
||||||||||||||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|||||||||||||
|
|
|
|
|
|
|
|
|
|
dt |
3 |
|
|
1 |
|
|
|
|
3 |
|
|
||||||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
d |
|
P (t) 0,2P (t) 4P (t), |
|
|
|||||||||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|||||||||||
|
|
|
|
|
|
|
|
|
|
|
|
dt |
|
|
|
||||||||||||
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
4 |
|
1 |
|
|
|
4 |
|
|
|
|||
|
|
|
|
|
|
|
|
|
|
|
P (t) P (t) P (t) P (t) 1. |
|
|
||||||||||||||
|
|
|
|
|
|
|
|
|
|
|
1 |
|
|
|
|
2 |
|
3 |
|
|
4 |
|
|
|
|
||
41P4 (t),
.
(t),
.
4) Предполагается, что в начальный момент времени вероятности состояний системы равны
P1 (0) 1, P2 (0) P3 (0) P4 (0) 0.
Для аналитического решения задачи Коши в пакете Maple можно ввести набор следующих строк:
sis:=diff(p1(t),t)+0.7*p1(x)-0.333*p2(t)-3.333*p3(t)-4*p4(t), diff(p2(t),t)-0.1*p1(t)+0.333*p2(t)-1.333*p3(t), diff(p3(t),t)-0.4*p1(t)+4.666*p3(t), diff(p4(t),t)-0.2*p1(t)+4*p4(t);
>fens:={p1(t),p2(t),p3(t),p4(t)}:
>F:=dsolve({sis, p1(t)=1, p2(t)=0, p3(t)=0, p4(t)=0}, fens );
25
Опыт показывает, что полученные решения больших систем дифференциальных уравнений слишком громоздки и не представляют практического интереса.
5) Для решения задачи Коши численно с помощью метода Эйлера воспользуемся программой, написанной на Pascal.
Program lab_eiler; uses Crt;
var
k,n,i: integer; a,b,h,p10,p20,p30,p40 : real; p1,p2,p3,p4 : array[0..100] of real;
function f1(p1,p2,p3,p4:real):real; begin
f1:= 0.7*p1+0.333*p2+3.333*p3+4*p4; end;
function f2(p1,p2,p3,p4:real):real; begin
f2:= 0.1*p1 0.333*p2+1.333*p3; end;
function f3(p1,p2,p3,p4:real):real; begin
f3:= 0.4*p1 4.666*p3; end;
function f4(p1,p2,p3,p4:real):real; begin
f4:= 0.2*p1 4*p4; end;
begin ClrScr;
write(‘vvedite a:’); readln(a); write(‘vvedite b:’); readln(b);
26
write(‘vvedite n:’); readln(n); write(‘vvedite p10:’); readln(p10); write(‘vvedite p20:’); readln(p20); write(‘vvedite p30:’); readln(p30); write(‘vvedite p40:’); readln(p40);
p1[0]:=p10; p2[0]:=p20; p3[0]:=p30; p4[0]:=p40; h:=(b-a)/n;
k:=0;
for i:=1 to n do begin k:=k+1;
p1[i]:=p1[i 1]+h*f1(p1[i 1],p2[i 1],p3[i 1],p4[i 1]); p2[i]:=p2[i 1]+h*f2(p1[i 1],p2[i 1],p3[i 1],p4[i 1]); p3[i]:=p3[i 1]+h*f3(p1[i 1],p2[i 1],p3[i 1],p4[i 1]); p4[i]:=p4[i 1]+h*f4(p1[i 1],p2[i 1],p3[i 1],p4[i 1]);
if (k=4) then begin
writeln(‘p1[‘ ,i:2, ‘]=’, p1[i]:7:3, ‘p2[‘ ,i:2, ‘]=’, p2[i]:7:3, ‘p3[‘ ,i:2, ‘]=’, p3[i]:7:3, ‘p4[‘ ,i:2, ‘]=’, p4[i]:7:3); k:=0;
end;
end;
readln;
end.
Предполагается, а=0, b=8, n=64, p10=1, p20= p30= p40=0.
Результаты вычислений записываются в таблицу и оформляются в
виде графиков зависимостей P1(t), |
P2 (t) , P3 (t), |
P4 (t) . |
||||
|
|
|
|
|
|
Таблица 4 |
t |
P1(t) |
P2 (t) |
|
P3 (t) |
|
P4 (t) |
0 |
1 |
0 |
|
0 |
|
0 |
0.5 |
0,816 |
0,072 |
|
0,071 |
|
0,041 |
1 |
0,753 |
0,142 |
|
0,067 |
|
0,039 |
1.5 |
0,707 |
0,195 |
|
0,062 |
|
0,036 |
2 |
0,672 |
0,235 |
|
0,059 |
|
0,034 |
|
|
27 |
|
|
|
|
Продолжение таблицы 4
2.5 |
0,646 |
0,265 |
0,056 |
0,033 |
3 |
0,625 |
0,288 |
0,054 |
0,032 |
3.5 |
0,610 |
0,306 |
0,053 |
0,031 |
4 |
0,599 |
0,319 |
0,052 |
0,030 |
4.5 |
0,590 |
0,330 |
0,051 |
0,030 |
5 |
0,583 |
0,337 |
0,050 |
0,029 |
5.5 |
0,578 |
0,343 |
0, 050 |
0,029 |
6 |
0,574 |
0,348 |
0,049 |
0,029 |
6.5 |
0,571 |
0,351 |
0, 049 |
0,029 |
7 |
0,569 |
0,354 |
0, 049 |
0,029 |
7.5 |
0,567 |
0,355 |
0, 049 |
0,028 |
8 |
0,566 |
0,357 |
0, 049 |
0,028 |
Для графического оформления табличных данных воспользуемся командой plot.
> p1:=[[0,1],[0.5,0.816],[1,0.753],[1.5,0.707],[2,0.672],[2.5,0.646], [3,0.625],[3.5,0.610],[4,0.599],[4.5,0.590],[5,0.583],[5.5,0.578],[6,0.574 ],[6.5,0.571],[7,0.569],[7.5,0.567],[8,0.566]]: >p2:=[[0,0],[0.5,0.072],[1,0.142],[1.5,0.195],[2,0.235],[2.5,0.265], [3,0.288],[3.5,0.306],[4,0.319],[4.5,0.330],[5,0.337],[5.5,0.343],[6,0.348 ],[6.5,0.351],[7,0.354],[7.5,0.355],[8,0.357]]:
>p3:=[[0,0],[0.5,0.071],[1,0.067],[1.5,0.062],[2,0.059],[2.5,0.056],
[3,0.054],[3.5,0.053],[4,0.052],[4.5,0.051],[5,0.050],[5.5,0.050],[6,0.049
],[6.5,0.049],[7,0.049],[7.5,0.049],[8,0.049]]:
>p4:=[[0,0],[0.5,0.041],[1,0.039],[1.5,0.036],[2,0.034],[2.5,0.033],
[3,0.032],[3.5,0.031],[4,0.030],[4.5,0.030],[5,0.029],[5.5,0.029],[6,0.029
],[6.5,0.029],[7,0.029],[7.5,0.028],[8,0.028]]:
>plot({p1,p2,p3,p4});
28