Сделай Сам Свою Работу на 5

НЕКОТОРЫЕ МЕТОДЫ ЧИСЛЕННОГО РЕШЕНИЯ СИСТЕМ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ.





Содержание работы и порядок выполнения

 

Задача, предлагаемая в курсовой работе – задача прямого моделирования – получение однозначного решения при условии, что математическая модель определена полностью, т. е. известны все параметры модели, начальные условия, задается:

1. Схема механизма химической реакции.

2. Константы скоростей отдельных стадий реакций.

3. Начальные концентрации компонентов.

4. Продолжительность реакции.

5. Метод численного решения дифференциальных уравнений кинетики.

Требуется:

1. Составить кинетическую модель данной химической реакции.

2. Выполнить на калькуляторе численное решение дифференциальных уравнений кинетики с целью получения кинетических зависимостей компонентов реакции в виде таблиц и графиков функций Сi (t) для пяти равноотстоящих значение t.

3. Выполнить программирование задачи для ЦВМ на алгоритмическом языке Паскаль.

4. Решить задачу на ЦВМ для 20 равноотстоящих значений t.

 

1.2. Содержание и оформления отчета.

Отчет по лабораторной работе должен включать формулировку задания, а также подробное описание порядка выполнения работы в соответствии с предыдущим разделом « Содержание и порядок выполнения ».



В отчете должны быть представлены:

система дифференциальных уравнений, представляющая кинетическую модель данной химической реакции;

описание метода численного решения системы дифференциальных уравнений;

результаты численного решения задачи на калькуляторе в виде таблиц, графиков;

Паскаль-программа решения задачи;

результаты решения задачи на ЦВМ в виде таблиц и графиков решений, построенных на миллиметровой бумаге.

Отчет оформляется лабораторной тетради.

 

ПОСТАНОВКА ЗАДАЧИ.

Маршрут химической реакции принято отражать с помощью стехиометрического уравнения, которое показывает, в каких соотношениях вещества вступают во взаимодействие.

Компоненты реакции обозначают буквами A,B,C,…

Например, А + В С – стехиометрическое уравнение.

Если стехиометрическое уравнение отражает механизм химической реакции, то оно служит основанием для составления кинетических уравнений.



Кинетические уравнения определяют связь между скоростью химической реакции концентрациями реагирующих веществ.

Кинетические уравнения записываются на основании закон действующих масс, в соответствии с которым, скорость химической реакции пропорциональна концентрациям взаимодействующих веществ и не зависит от концентрации продуктов.

Так, для приведенного выше стехиометрического уравнения имеем:

 

 

Здесь К – константа скорости химической реакции.

Решая систему дифференциальных уравнений с некоторыми заданными начальными условиями СА(О) = САО , СВ(О) = СВО , СС(О) = ССО , получим кривые зависимости СА(t) , СВ(t) , СС(t) – кинетические кривые, таким образом, можем рассчитать концентрации компонентов в определенный момент времени.

 

НЕКОТОРЫЕ МЕТОДЫ ЧИСЛЕННОГО РЕШЕНИЯ СИСТЕМ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ.

Математическое моделирование кинетики сложных химических реакций сводится к решению систем обыкновенных дифференциальных уравнений вида:

с начальными условиями yi(x0) = yio , где i = 1,…,К.

Как правило, система уравнений, описывающих кинетику изучаемой реакции, является нелинейной и поэтому не может быть решена аналитически. Возникает необходимость использования метод приближенного численного интегрирования.

Эти методы позволяют приближенно отыскать решения дифференциальных уравнений только на некотором конечном интервале [a; b].

Пусть имеем некоторое дифференциальное уравнение первого порядка y=f(x;y)

с начальными условиями y(x0)=y0. Будем искать решение этого уравнения на отрезке [x0;b].

Разобьем этот отрезок на n равных частей. Тогда получим систему равноотстоящих узлов



x0 , x1=x0+h, x2=x1+h, …, xn=b.

Здесь h=(b-x0)/n – шаг интегрирования.

Численные методы дают возможность найти в некотором числе точек x1, x2, …, xn приближения y1, y2, …, yn для значений точного решения y(x1), y(x2), …, y(x0) .

Наиболее простым методом решения дифференциальных и их систем является

 

3.1 МЕТОД ЭЙЛЕРА.

 

Пусть дано дифференциальное уравнение y'=f(x,y) (1), с начальными условиями y(x0)=y0 .

Пусть y=y(x) – искомое точное решение. Интегральная кривая проходит через точку (x0,y0).

Найдем приближенные значения функции в точках x1, x2, …, xn. Построим систему равноотстоящих точек x0, x1=x0+h, …, xn=b. Проведём прямые x=x0, x=x1, …, x=b.

y0
y1
y2

 

 

Рассмотрим отрезок [x0,x1]. На этом отрезке есть одна точка, которая принадлежит искомой кривой – это точка А (x0,y0). Заменим дугу искомой кривой y=y(x) на отрезке [x0,x1] касательной к ней, проведенной в точке (x0,y0). В качестве y1 возьмём ординату точки пересечения прямой x=x1 к касательной.

Очевидно y1=y0+h•tgα0. Но tgα0=y'(x0), т.е. y1=y0+h•y'(x0). Но из уравнения (1) следует, что

y'(x0)=f(x0,y0).

Итак, получаем y1=y0+h•f(x0,y0), x1=x0+h .

Предположим теперь, что точка (x1,y1) принадлежит искомой кривой. В этой точке опять проведём касательную к графику функции до пересечения с прямой x=x2 .

Тогда аналогично:

y2= y1+ h f(x1,y1);

x2= x1+ h .

Продолжая и так далее, получим систему значений y1,y2,…,yn, которые и будут приближенными значениями функции y=y(x) в точках x1,x2,…xn.

Итак, расчётные формулы метода Эйлера:

ym+1= ym+ h f(xm,ym);

xm+1= xm+ h

Для системы дифференциальных уравнений

y'I = fi(x,y1,y2,…yk)

yi(x0) = yi0 i = 1,…,K

расчетные формулы записываются аналогично

yim+1 = yim + hfi(x,y1m,y2m,…ykm);

xm+1 = xm+ h,

здесь i – номер уравнения в системе, n – номер шага.

Метод Эйлера является грубым методом, ошибка, которую мы допускаем, на каждом шаге, пропорциональна h, т.е. |y(xk) - yk|≈ h.

Чтобы повысить точность вычислений, используют некоторые усовершенствованные методы.

3.2. МЕТОД ЭЙЛЕРА-КОШИ

 

 

Пусть опять решаем уравнение y'= f(x,y),

y(x0)= y0.

Решение ищем на отрезке [x0,xn].

Пусть нам известны координаты некоторой точки, принадлежащей искомому решению (xm,ym). Найдём средний тангенс угла наклона касательной для двух точек: (xm,ym) и (xm+h,ym+h·ỹm).

Последняя точка, есть та самая, которую в методе Эйлера мы обозначаем (xm+1,ym+1), но здесь точка будет вспомогательной. ỹm+1

m+1   ym+1   ym

Итак, сначала по методу Эйлера находится точка А, лежащая на прямой L1, тангенс угла наклона которой tgα1= f(xm,ym). В этой точке снова вычисляется тангенс угла наклона касательной L2

tgα2= f(xm+h,ym+hỹm)

Затем через точку (xm,ym) проводим прямую L, тангенс угла наклона которой равен (tgα1+tgα2)/2. Точка, в которой L пересекается с прямой x=xm+1 и будет искомой (xm+1,ym+1). Таким образом, ym+1 есть искомое приближение значений функции на данном шаге интегрирования.

Расчетные формулы метода Эйлера-Коши следующие:

m+1= ym+hf(xm,ym);

ym+1= ym+(h/2)·[f(xm,ym)+f(xm+1,ỹm+1)];

xm+1= xm+h.

Аналогично, для системы дифференциальных уравнений:

yi'= fi(x,y1,y2,…yk), yi(x0)=yio, i=1,2,..k;

im+1= yim+ hfi(xm,y1m,y2m,…ykm);

yim+1= yim+ (h/2)·[fi(xm,y1m,y2m,…ykm)+fi(xm+1,ỹ1m+1,ỹ2m+1,ỹkm+1)];

xm+1= xm+h.

Здесь i – номер уравнения системы, m – номер шага.

 

ПРОГРАММА МЕТОДА ЭЙЛЕРА.

 

PROGRAM EILER;

USES CRT;

VAR

C, Z: ARRAY[1..4] OF REAL;

I: INTEGER;

T, TMAX: REAL;

BEGIN

CLRSCR;

WRITE (‘TMAX=’); READLN (TMAX);

T:=0;

WRITE (‘ВВЕДИТЕ С1, С2, С3, С4’);

READLN(C[1]);

READLN(C[2]);

READLN(C[3]);

READLN(C[4]);

WRITELN (‘T C1 C2 C3 C4’);

WHILE T<=TMAX DO

BEGIN

Z[1]:={ПРАВАЯ ЧАСТЬ ПЕРВОГО УРАВНЕНИЯ СИСТЕМЫ};

Z[2]:={ПРАВАЯ ЧАСТЬ ВТОРОГО УРАВНЕНИЯ СИСТЕМЫ};

Z[3]:={ПРАВАЯ ЧАСТЬ ТРЕТЬЕГО УРАВНЕНИЯ СИСТЕМЫ};

Z[4]:={ПРАВАЯ ЧАСТЬ ЧЕТВЁРТОГО УРАВНЕНИЯ СИСТЕМЫ};

T:=T+0.5;

WRITE (T:5:1);

FOR I:=1 TO 4 DO

WRITE (Z[1]:10:5);

WRITELN;

FOR I:=1 TO 4 DO

C[I]:=Z[I];

END;

READKEY;

END.

 

Program Eiler;

Uses Crt;

Type Vec=Array[1..6] of Real;

 

Var y, k, f:Vec;

i, KolKom, KolRec, KolForPrint,count:Byte;

x, xk, h, s:Real;

Fout:text;

NameOut:string[12];

Procedure f_dir(y:vec; var f:vec);

var r:vec;

begin

r[1]:=k[1]*y[1];

r[2]:=k[2]*y[2];

r[3]:=k[3]*y[2];

r[4]:=k[4]*y[3];

r[5]:=k[5]*y[2];

f[1]:=-r[1]+r[5] ;

f[2]:=r[1]-r[2]-r[3]+r[4]-r[5];

f[3]:=r[2]-r[4];

f[4]:=r[3]

end;

 

BEGIN

ClrScr;

x:=0; xk:=15; h:=0.1;count:=0;s:=0;

 

write('введите количество компонентов ');readln(KolKom);

For i:=1 to KolKom do

begin

write(введите концентрации компонентов');readln(y[i])

end;

write('введите количество химических реакций ');readln(KolRec);

For i:=1 to KolRec do

begin

write(введите константы химической реакции');readln(k[i])

end;

write('введите число точек для вывода значения на экран');

readln(KolForPrint);

write('введите имя файлов для результатов');

readln(NameOut);

ClrScr;

assign(fout,NameOut);

rewrite(fout);

write('время');write(fout,'время');

For i:=1 to KolKom do begin

write(' y[',i,']');

write(fout,' y[',i,']'); end;

writeln(' сумма '); writeln(fout,' сумма ');

write(x:5:1);write(fout,x:5:1);

for i:=1 to KolKom do

begin s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5) end;

writeln(s:8:2);writeln(fout,s:8:2);

 

repeat

f_dir(y,f);

For i:=1 to KolKom do y[i]:=y[i]+h*f[i];

x:=x+h;

count:=count+1;

if count=round(xk/h/KolForPrint) then

begin

write(x:5:1);write(fout,x:5:1);s:=0;

for i:=1 to KolKom do

begin

s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5);

end;

writeln(s:8:2);writeln(fout,s:8:2);count:=0

end;

until x>xk;

close(fout);

readkey

end.

ПРОГРАММА МЕТОДА ЭЙЛЕР-КОШИ

Program Eiler_Koshi;

Uses Crt;

Type MAS=Array[1..6] of Real;

Var

y, y0,k, f0,f:MAS;

i, KolKom, KolRec, KolForPrint, count:Byte;

x

xk,

h, s:Real;

Fout:text;

NameOut:string[12];

 

Procedure f_dir(y:vec; var f:vec);

var r:vec;

begin

r[1]:=k[1]*y[1];

r[2]:=k[2]*y[2];

r[3]:=k[3]*y[2];

r[4]:=k[4]*y[3];

r[5]:=k[5]*y[2];

f[1]:=-r[1]+r[5] ;

f[2]:=r[1]-r[2]-r[3]+r[4]-r[5];

f[3]:=r[2]-r[4];

f[4]:=r[3]

end;

 

BEGIN

ClrScr;

x:=0; xk:=15; h:=0.1;count:=0;s:=0;

 

write('введите количество компонентов ');readln(KolKom);

For i:=1 to KolKom do

begin

write(‘введите концентрацию компонента');readln(y[i])

end;

write('введите количество химических реакций ');readln(KolRec);

For i:=1 to KolRec do

begin

write('введите константы химической реакции');readln(k[i])

end;

write('введите число точек для вывода значений на экран ');

readln(KolForPrint);

write('введите имя файла для результатов');

readln(NameOut);

ClrScr;

assign(fout,NameOut);

rewrite(fout);

 

write(' время');write(fout,' время');

For i:=1 to KolKom do begin

write(' y[',i,']');

write(fout,' y[',i,']'); end;

writeln(' сумма '); writeln(fout,' сумма ');

write(x:5:1);write(fout,x:5:1);

for i:=1 to KolKom do

begin s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5) end;

writeln(s:8:2);writeln(fout,s:8:2);

 

repeat

f_dir(y,f);

For i:=1 to KolKom do

begin

y0[i]:=y[i];

f0[i]:=f[i];

y[i]:=y0[i]+h*f[i];

end;

f_dir(y,f);

For i:=1 to KolKom do y[i]:=y0[i]+0.5*h*(f0[i]+f[i]);

x:=x+h;

count:=count+1;

if count=round(xk/h/KolForPrint) then

begin

write(x:5:1);write(fout,x:5:1);s:=0;

for i:=1 to KolKom do

begin

s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5);

end;

writeln(s:8:2);writeln(fout,s:8:2);count:=0

end;

until x>xk;

close(fout);

readkey

end.

 

3.3. МЕТОД РУНГЕ-КУТТА 2-ГО ПОРЯДКА

 

Пусть имеем дифференциальное уравнение y'=f(x,y) с начальными условиями у(x0)=y0.

Ищем решение на отрезке [x0,xn].

Пусть имеем точку (xm,ym), принадлежащему искомому решению. Для того, чтобы найти следующую точку проведем касательную к кривой в точке (xm,ym) до пересечения с прямой x=xm+0.5, где xm+0.5=xm+h/2. Тогда получим координату (по формуле Эйлера)

ym+0.5=ym+h/2·f(xm,ym).

Теперь найдём тангенс угла наклона касательной в т.В (xm+0.5,ym+0.5), (прямая L). Через точку А проведём прямую L'||L. Ординату точки пересечения прямых L' и x=xm+1 возьмём в качестве ym+1.

Таким образом

ym+0.5= ym+h/2·f(xm,ym), xm+0.5=xm+h/2

ym+1= ym+h·f(xm+0.5,ym+0.5), xm+1=xm+h

Для системы дифференциальных уравнений

yi'= fi(x,y1,y2,…yk),

yi(x0)= yi0, i= 1,2,…k

расчетные формулы имеют вид:

yim+0.5= yim+h/2·fi(xm,y1m,y2m,…ykm), xm+0.5= xm+h/2;

yim+1= yim+h·fi(xm+0.5,y1m+0.5,y2m+0.5,…ykm+0.5);

xm+1= xm+h.

 

ПРОГРАММА РУНГЕ-КУТТА 2-ГО ПОРЯДКА

Program Runge_Kutte_Method;

Uses Crt;

Type Vec=Array[1..6] of Real;

Var y, y0, k, f:Vec;

i, KolKom, KolRec, KolForPrint, count:Byte;

x, xk, h, s:Real;

Fout:text;

NameOut:string[12];

 

Procedure f_dir(y:vec; var f:vec);

var r:vec;

begin

r[1]:=k[1]*y[1];

r[2]:=k[2]*y[2];

r[3]:=k[3]*y[2];

r[4]:=k[4]*y[3];

r[5]:=k[5]*y[2];

f[1]:=-r[1]+r[5] ;

f[2]:=r[1]-r[2]-r[3]+r[4]-r[5];

f[3]:=r[2]-r[4];

f[4]:=r[3]

end;

 

BEGIN

ClrScr;

x:=0; xk:=15; h:=0.1;count:=0;s:=0;

 

write('введите количество компонентов ');readln(KolKom);

For i:=1 to KolKom do

begin

write(i,’введите концентрации компонентов ');readln(y[i])

end;

write('введите количество химических реакций ');readln(KolRec);

For i:=1 to KolRec do

begin

write(i,'введте константы химических реакции');readln(k[i])

end;

write('введите число точек для вывода значения на экран ');

readln(KolForPrint);

write('введите имя файлов для результатов');

readln(NameOut);

ClrScr;

assign(fout,NameOut);

rewrite(fout);

 

write(' время ');write(fout,' время ');

For i:=1 to KolKom do begin

write(' y[',i,']');

write(fout,' y[',i,']'); end;

writeln(' сумма '); writeln(fout,' сумма ');

write(x:5:1);write(fout,x:5:1);

for i:=1 to KolKom do

begin s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5) end;

writeln(s:8:2);writeln(fout,s:8:2);

 

repeat

f_dir(y,f);

For i:=1 to KolKom do

begin

y0[i]:=y[i];

y[i]:=y0[i]+0.5*h*f[i];

end;

f_dir(y,f);

For i:=1 to KolKom do y[i]:=y0[i]+h*f[i];

x:=x+h;

count:=count+1;

if count=round(xk/h/KolForPrint) then

begin

write(x:5:1);write(fout,x:5:1);s:=0;

for i:=1 to KolKom do

begin

s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5);

end;

writeln(s:8:2);writeln(fout,s:8:2);count:=0

end;

until x>xk;

close(fout);

 

readkey

end.

 

3.4. МЕТОД РУНГЕ-КУТТА 4-ГО ПОРЯДКА

Этот метод один из самых распространённых методов интегрирования дифференциальных уравнений.

Для одиночного дифференциального уравнения y'= f(x,y), y(x0)=y0 расчётные формулы имеют следующий вид:

ym+1= ym+h/6(k1+2k2+2k3+k4),

где k1= f(xm,ym);

k2= f(xm+ ,ym+ );

k3= f(xm+ ,ym+ );

k4= f(xm+h,ym+hk3);

xm+1= xm+h.

Для системы дифференциальных уравнений

yi'= fi(x,y1,y2,…yk), yi(x0)=yi0, i=1,2,..k,

расчетные формулы запишутся следующим образом:

yim+1= yim+ (k1i+2k2i+2k3i+k4i), i=1,2,…k,

где k1i= fi(xm,y1m,y2m,…ykm);

k2i= fi(xm+h/2,y1m+ k11,y2m+ k12,…ykm+ k1k);

k3i= fi(xm+ ,y1m+ k21,y2m+ k22,…ykm+ k2k);

k4i= fi(xm+h,y1m+hk31,y2m+hk32,…ykm+hk3k).

PROGRAM RYNGE_KYTT;

USES CRT;

VAR

X,X0,Y0,XMAX:REAL;

YT,DY,K1,K2,K3,K4:REA;

H:REAL;

I,N:LONGINT;

 

FUNCTION F(X,Y:REAL):REAL;

BEGIN

F:=(-2*Y/X+EXP(-X*X)/X;

END;

 

BEGIN

WRITE(‘H=’); READLN(H);

WRITE(‘X0=’); READLN(X0);

WRITE(‘XMAX=’); READLN(XMAX);

WRITE(‘Y0=’); READLN(Y0);

N:=ROUND((XMAX-X0)/H);

YT:=-EXP(-X0*X0)/(2*X0*X0)+15/(X0*X0);

WRITELN(‘X0=’,X0:8:4,’ Y0=YT= ’,YT:8:4);

FOR I:=1 TO N DO

BEGIN

K1:=F(X0,Y0);

K2:=F(X0+H/2,Y0+H*K1/2);

K3:=F(X0+H/2,Y0+H*K2/2);

K4:=F(X0+H,Y0+H*K3);

Y0:=Y0+H/6*(K1+2*K2+2*K3+K4);

X0:=X0+H;

YT=-EXP(-X0*X0)/(2*X0*X0)+15/(X0*X0);

DY:=ABS(YT-Y0);

WRITELN(‘X=’,X0:8:4,’YT=’,YT:8:4,’YR=’,YO:8:4,’DY=’,DY:10:8);

END;

READKEY;

END.

 

Program Runge_Kutte_4_Method;

Uses Crt;

Type =Array[1..6] of Real;

Var y, y0, k, ka, f:Vec;

i, KolKom, KolRec, KolForPrint, count:Byte;

x, xk, h, s:Real;

Fout:text;

NameOut:string[12];

 

 

Procedure f_dir(y:vec; var f:vec);

var r:vec;

begin

r[1]:=k[1]*y[1];

r[2]:=k[2]*y[2];

r[3]:=k[3]*y[2];

r[4]:=k[4]*y[3];

r[5]:=k[5]*y[2];

f[1]:=-r[1]+r[5] ;

f[2]:=r[1]-r[2]-r[3]+r[4]-r[5];

f[3]:=r[2]-r[4];

f[4]:=r[3]

end;

 

BEGIN

ClrScr;

x:=0; xk:=15; h:=0.1;count:=0;s:=0;

 

write('введите количество компонентов ');readln(KolKom);

For i:=1 to KolKom do

begin

write(i,'введите концентрации компонентов');readln(y[i]);

y0[i]:=y[i]

end;

write('введите количество химических реакций ');readln(KolRec);

For i:=1 to KolRec do

begin

write(i,'введите константы химической реакции ');readln(k[i])

end;

write(' введите число точек для вывода значения на экран ');

readln(KolForPrint);

write('введите имя файла для результатов ');

readln(NameOut);

ClrScr;

assign(fout,NameOut);

rewrite(fout);

 

write('время');write(fout,' время ');

For i:=1 to KolKom do begin

write(' y[',i,']');

write(fout,' y[',i,']'); end;

writeln(' сумма '); writeln(fout,' сумма');

write(x:5:1);write(fout,x:5:1);

for i:=1 to KolKom do

begin s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5) end;

writeln(s:8:2);writeln(fout,s:8:2);

 

repeat

f_dir(y,f);

For i:=1 to KolKom do

begin

ka[i]:=h*f[i];

y[i]:=y0[i] + ka[i]/2;

end;

f_dir(y,f);

For i:=1 to KolKom do

begin

ka[i]:=ka[i] + 2*h*f[i];

y[i]:=y0[i] + h*f[i]/2;

end;

f_dir(y,f);

For i:=1 to KolKom do

begin

ka[i]:=ka[i] + 2*h*f[i]; y[i]:=y0[i] + h*f[i];

end;

f_dir(y,f);

For i:=1 to KolKom do

begin

y[i]:=y0[i] + (ka[i] + h*f[i])/6;

y0[i]:=y[i]

end;

x:=x+h;

count:=count+1;

if count=round(xk/h/KolForPrint) then

begin

write(x:5:1);write(fout,x:5:1);s:=0;

for i:=1 to KolKom do

begin

s:=s+y[i];write(y[i]:8:5);write(fout,y[i]:8:5);

end;

writeln(s:8:2);writeln(fout,s:8:2);count:=0

end;

until x>xk;

close(fout);

readkey

end.

 

4. МОДЕЛИРОВАНИЕ КИНЕТИКИ ХИМИЧЕСКОЙ РЕАКЦИИ

Пример.

Задано:

1. Схема механизма химической реакции

 

 

К1 К2

А В С

К5 К3 К4

 

Д

 

2. Константы скоростей отдельных стадий реакции

К1 = 0,76 с-1

К2 = 0,9 с-1

К3 = 0,5 с-1

К4 = 0,45 с-1

К5 = 0,6 с-1

3. Начальные концентрации компонентов

СА(0) = 0,7 моль/л

СВ(0) = 0

СС(0) = 0

СД(0) = 0

4. Продолжительность реакции 15 секунд

5. Метод численного решения – Эйлера.

Выполнение работы:

1. Система дифференциальных уравнений, представляющая кинетическую модель данной химической реакции:

2. Расчетные формулы метода Эйлера:

CAm+1= CAm+h(-K1CAm+K5CBm);

CBm+1= CBm+h(-K2CBm+K1CAm+K4CCm-K5CBm-K3CBm);

CCm+1= CCm+h(K2CBm-K4CCm);

CDm+1= CDm+h·K3CBm;

tm+1= tm+h.

3. Результаты численного решения системы дифференциальных уравнений на калькуляторе. h = 0.2; 5 шагов по времени.

Пример вычислений:

 

СA1 = 0.7+0.2(-0.76*0.7+0.6*0)=0.5936;

CB1 = 0+0.2(-0.2*0+0.76*0.7+0.45*0-0.6*0-0.5*0)=0.01064

CC1 = 0+0.2(0.9*0-0.45*0)=0;

CD1=0+0.2(0.5*0)=0;

t=0+0.2=0.2

Аналогично вычисляем СА2, СВ2, СС2, СД2,…СА5, СВ5, СС5, СД5.

 

 

Таблица результатов расчёта

t CA CB CC CD
0,7
0,2 0,5936 0,1064
0,4 0,5161 0,1241 0,0192 0,0106
0,6 0,4562 0,1726 0,0452 0,026
0,8 0,4076 0,177 0,0722 0,0433
0,3669 0,1747 0,0976 0,061

 

Графики кинетических кривых

4. ПАСКАЛЬ-программа решения системы дифференциальных уравнений

PROGRAM EILER;

USES CRT;

VAR

C, Z: ARRAY[1..4] OF REAL;

I: INTEGER;

T, TMAX: REAL;

BEGIN

CLRSCR;

WRITE (‘TMAX=’); READLN (TMAX);

T:=0;

WRITE (‘ВВЕДИТЕ С1, С2, С3, С4’);

READLN(C[1]);

READLN(C[2]);

READLN(C[3]);

READLN(C[4]);

WRITELN (‘T C1 C2 C3 C4’);

WHILE T<=TMAX DO

BEGIN

Z[1]:=C[1]+0.2*(-0.76*C[1]+0.6C[2]);

Z[2]:=C[2]+0.2*(-0.9*C[2]+0.76*C[1]+0.45*C[3]-0.6*C[2]-0.5*C[2]);

Z[3]:=C[3]+0.2*(0.9*C[2]-0.45*C[3]);

Z[4]:=C[4]+0.2*0.5*C[2];

T:=T+0.2;

WRITE (T:5:1);

FOR I:=1 TO 4 DO

WRITE (Z[1]:10:5);

WRITELN;

FOR I:=1 TO 4 DO

C[I]:=Z[I];

END;

READKEY;

END.

 

 








Не нашли, что искали? Воспользуйтесь поиском по сайту:



©2015 - 2024 stydopedia.ru Все материалы защищены законодательством РФ.