Этот метод один из самых распространенных методов интегрирования дифференциальных уравнений.
Для одиночного дифференциального уравнения расчетные формулы имеют следующий вид:
где
Для системы дифференциальных уравнений
Для системы дифференциальных уравнений
расчетные формулы запишутся следующим образам:
Моделироваие кинетики химической реакция
Пример.
Задано:
I. Схема механизма химической реакции
2. Константы скоростей отдельных стадий реакции
3. Начальные концентрации компонентов
4. Продолжительность реакции 15 сек.
5. Метод численного решения - Эйлера.
Выполнение работы
I. Система, дифференциальных уравнения, представлявшая кинетическую модель данной химической реакции:
2. Расчетные формулы метода Эйлера:
3. Результаты численного решения системы дифференциальных уравнений на калькуляторе. h= 0,2; 5 шагов по времени.
Таблица результатов расчета
№ | СА | СВ | СС | СД | |
Q | 0,7 | ||||
I | 0,2 | 0,5936 | 0,1064 | ||
0,4 | 0,5161 | 0,1541 | 0,0192 | 0,0106 | |
0.6 | 0,4562 | 0,1726 | 0,0452 | 0,0260 | |
0,8 | 0,4076 | 0,1770 | 0,0722 | 0,0433 | |
I | 0,3669 | Q.I747 | 0,0976 | 0,0610 |
5. ПРОГРАММА МЕТОДОВ:
1. Эйлера;
2. Эйлера-Коши;
3. Рунге-Кутта 2-го порядка;
4. Рунге-Кутта 4-го порядка.
1)
Program EULER;
const
eps = 1e-5;
x0 = 3;
y0 = 3;
type
TFunc = function(x,y:extened): extended;
function func(x,y: extended): extended; far;
begin
func:= y*cos(x) - 2*sin(2*x);
end;
function Euler(f: TFunc; x0_,x_end,y0_: extended; n: word): extended;
{ где x0_ и y0_ - начальное условие
x_end - точка, в которой необходимо вычислить результат
n - количество шагов для вычисления результата }
var
i: word; { счетчик цикла }
x,h: extended; { текущая точка и длина шага }
res: extended; { переменная для накопления конечного результата функции}
begin
h:= (x_end - x0_)/n; { Находим длину шага }
res:= y0_; { устанавливаем начальные значения}
x:=x0_;
for i:=1 to n do
begin { вычисляем результат по методу Эйлера }
res:=res+h*f(x,res);
x:=x+h; { переходим к следующей точке }
end;
Euler:=res; { присваиваем конечный результат функции }
end;
begin
writeln('Численное решение дифференциальных уравнений:');
writeln(#10,' y'' = y*cos(x) - 2*sin(2*x); y(0)=3; x_end=',(5*x0+3.5):5:5);
writeln(#10,'Метод Эйлера:');
writeln('n=5 000, result: ',Euler(func,x0,5*x0+3.5,y0,5000):5:5);
writeln('n=10 000, result: ',Euler(func,x0,5*x0+3.5,y0,10000):5:5);
writeln('n=25 000, result: ',Euler(func,x0,5*x0+3.5,y0,25000):5:5);
write(#10,'Press Enter to continue.');
readln;
end.
2)
Program EULER-KOSHI;
const
eps = 1e-5;
x0 = 3;
y0 = 3;
type
TFunc = function(x,y:extened): extended;
function func(x,y: extended): extended; far;
begin
func:= y*cos(x) - 2*sin(2*x);
end;
function Euler2(f: TFunc; x0_,x_end,y0_: extended; n: word): extended;
{ где x0_ и y0_ - начальное условие
x_end - точка, в которой необходимо вычислить результат
n - количество шагов для вычисления результата }
var
i: word; { счетчик цикла }
x,h: extended; { текущая точка и длина шага }
res: extended; { переменная для накопления конечного результата функции}
begin
h:= (x_end - x0_)/n; { Находим длину шага }
res:= y0_; { устанавливаем начальные значения}
x:=x0_;
for i:=1 to n do
begin { вычисляем результат по исправленному методу Эйлера }
res:=res+h*(f(x,res)+f(x+h,res+h*f(x,res)))/2;
x:=x+h; { переходим к следующей точке }
end;
Euler2:=res; { присваиваем конечный результат функции }
end;
begin
writeln(#10,'Эйлера-Коши:');
writeln('n=50, result: ',Euler2(func,x0,5*x0+3.5,y0,50):5:5);
writeln('n=100, result: ',Euler2(func,x0,5*x0+3.5,y0,100):5:5);
writeln('n=250, result: ',Euler2(func,x0,5*x0+3.5,y0,250):5:5);
write(#10,'Press Enter to continue.');
readln;
end.
3)
Program RK2;
const
eps = 1e-5;
x0 = 3;
y0 = 3;
type
TFunc = function(x,y:extened): extended;
function func(x,y: extended): extended; far;
begin
func:= y*cos(x) - 2*sin(2*x);
end;
function Euler3(f: TFunc; x0_,x_end,y0_: extended; n: word): extended;
{ где x0_ и y0_ - начальное условие
x_end - точка, в которой необходимо вычислить результат
n - количество шагов для вычисления результата }
var
i: word; { счетчик цикла }
x,h: extended; { текущая точка и длина шага }
res: extended; { переменная для накопления конечного результата функции}
begin
h:= (x_end - x0_)/n; { Находим длину шага }
res:= y0_; { устанавливаем начальные значения}
x:=x0_;
for i:=1 to n do
begin { вычисляем результат по модифицированному методу Эйлера }
res:=res+h*f(x+h/2,res+(h/2)*f(x,res));
x:=x+h; { переходим к следующей точке }
end;
Euler3:=res; { присваиваем конечный результат функции }
end;
begin
writeln(#10,'Метод Рунге-Кутта 2:');
writeln('n=50, result: ',Euler3(func,x0,5*x0+3.5,y0,50):5:5);
writeln('n=100, result: ',Euler3(func,x0,5*x0+3.5,y0,100):5:5);
writeln('n=250, result: ',Euler3(func,x0,5*x0+3.5,y0,250):5:5);
write(#10,'Press Enter to continue.');
readln;
end.
4)
Program RK4;
const
eps = 1e-5;
x0 = 3;
y0 = 3;
type
TFunc = function(x,y:extened): extended;
function func(x,y: extended): extended; far;
begin
func:= y*cos(x) - 2*sin(2*x);
end;
function RungeKutt(f: TFunc; x0_,x_end,y0_: extended; n: word): extended;
{ где x0_ и y0_ - начальное условие
x_end - точка, в которой необходимо вычислить результат
n - количество шагов для вычисления результата }
var
i: word; { счетчик цикла }
x,h: extended; { текущая точка и длина шага }
res: extended; { переменная для накопления конечного результата функции }
k1,k2,k3,k4: extended; { вспомогательные переменные вычисления результата }
begin
h:= (x_end - x0_)/n; { Находим длину шага }
res:= y0_; { устанавливаем начальные значения}
x:=x0_;
for i:=1 to n do
begin { вычисляем результат по методу Рунге-Кутта 4го порядка }
k1:=f(x,res);
k2:=f(x+h/2,res+h*k1/2);
k3:=f(x+h/2,res+h*k2/2);
k4:=f(x+h,res+h*k3);
res:=res+h*(k1+2*k2+2*k3+k4)/6;
x:=x+h; { переходим к следующей точке }
end;
RungeKutt:=res; { присваиваем конечный результат функции }
end;
begin
writeln(#10,'Метод Рунге-Кутта 4:');
writeln('n=50, result: ',RungeKutt(func,x0,5*x0+3.5,y0,50):5:5);
writeln('n=100, result: ',RungeKutt(func,x0,5*x0+3.5,y0,100):5:5);
writeln('n=250, result: ',RungeKutt(func,x0,5*x0+3.5,y0,250):5:5);
write(#10,'Press Enter to continue.');
readln;
end.
ЗАДАНИЯ
Схема химической реакции | Продолжительность реакции | Константы скоростей отдельных стадий реакции | Начальное значение концентрации реагентов | Метод численного решения | |||||||||
К1 | К2 | К3 | К4 | К5 | К6 | СА | СВ | СС | СД | Диф. Уравнений кинетики | |||
1. | 0,3 | 0,2 | 0,1 | - | - | - | 0,7 | Метод Эйлера | |||||
2. | 0,7 | 0,65 | 0,5 | 0,3 | - | - | 0,5 | Эйлера-Коши | |||||
3. | 0.65 | 0.6 | 0.55 | - | - | - | 0.3 | Рунге-Кутта 2-го порядка | |||||
4. | 0,2 | 0,3 | 0,4 | - | - | - | 1,0 | Эйлера | |||||
5. | 0,7 | 0,6 | 0,5 | - | - | - | 1,0 | - | Эйлера-Коши | ||||
6. | 0,37 | 0,3 | - | - | - | - | 0,96 | - | Рунге-Кутта 2-го порядка | ||||
7. | 0,5 | 0,36 | 0,5 | 0,3 | 0,25 | 0,1 | 1,1 | - | Эйлера | ||||
8. | 0,95 | 0,8 | 0,63 | 0,45 | - | - | 0,3 | Эйлера-Коши | |||||
9. | 0,85 | 0,9 | 0,85 | 0,6 | - | - | 0,2 | 0,1 | Рунге-Кутта 2-го порядка | ||||
10. | 0,09 | 0,08 | 0,1 | 0,06 | - | - | 0,85 | Эйлера |
11. | 0,2 | 0,16 | 0,1 | - | - | - | 1,1 | - | Эйлера | ||||
12. | 0,9 | 0,95 | 0,66 | - | - | - | 0,8 | - | Рунге-Кутта 2-го порядка | ||||
13. | 0,26 | 0,09 | 0,22 | 0,15 | 0,18 | - | 1,3 | - | Эйлера | ||||
14. | 0,5 | 0,5 | 0,4 | 0,3 | - | - | 0,45 | Эйлера-Коши | |||||
15. | 0,09 | 0,09 | 0,08 | 0,07 | - | - | 0,85 | Рунге-Кутта 2-го порядка | |||||
16. | 0,2 | 0,3 | 0,4 | 0,1 | - | - | 0,7 | - | Эйлера | ||||
17. | 0,95 | 0,8 | 0,46 | 0,39 | - | - | 0,87 | - | Эйлера-Коши | ||||
18. | 0,24 | 0,2 | 0,15 | - | - | - | 1,3 | - | Рунге-Кутта 2-го порядка | ||||
19. | 0,33 | 0,45 | 0,4 | 0,3 | 0,21 | - | 0,4 | Эйлера | |||||
20. | 0,76 | 0,9 | 0,5 | 0,45 | 0,6 | - | 0,7 | Эйлера-Коши | |||||
21. | 0,95 | 0,9 | 0,85 | 0,6 | - | 0,2 | Рунге-Кутта 2-го порядка | ||||||
22. | 0,82 | 0,3 | 0,15 | 0,1 | - | - | 0,65 | - | Эйлера | ||||
23. | 0,6 | 0,48 | - | - | - | - | 1,35 | - | Эйлера-Коши | ||||
24. | 0,76 | 0,76 | 0,94 | 0,35 | 0,2 | 0,2 | 0,85 | Рунге-Кутта 2-го порядка | |||||
25. | 0,09 | 0,1 | 0,08 | 0,06 | 0,05 | - | 0,3 | Эйлера | |||||
26. | 0,45 | 0,4 | 0,3 | 0,2 | 0,1 | - | 1,3 | Эйлера-Коши | |||||
27. | 0,5 | 0,5 | 0,4 | 0,3 | 0,2 | - | 0,8 | Рунге-Кутта 2-го порядка | |||||
28. | 0,8 | 0,8 | 0,7 | 0,5 | 0,4 | - | 0,85 | Эйлера | |||||
29. | 0,55 | 0,5 | 0,64 | 0,3 | 0,27 | - | 1,5 | Эйлера-Коши | |||||
30. | 0,93 | 0,93 | 0,55 | 0,4 | - | - | 1,0 | Рунге-Кутта 2-го порядка | |||||
31. | 0,93 | 0,8 | 0,76 | 0,4 | 0,5 | 0,5 | 0,8 | Эйлера | |||||
32. | 0,25 | 0,2 | 0,3 | 0,1 | 0,09 | - | 1,2 | - | Эйлера-Коши | ||||
33. | 0,85 | 0,7 | 0,53 | 0,35 | - | - | 1,35 | Рунге-Кутта 2-го порядка | |||||
34. | 0,1 | 0,25 | 0,3 | 0,2 | 0,18 | - | 0,7 | Эйлера | |||||
35. | 0,9 | 0,8 | 0,7 | 0,6 | 0,5 | 0,4 | 1,0 | Эйлера-Коши | |||||
36. | 0,09 | 0,08 | 0,07 | 0,06 | 0,05 | - | 1,25 | 0- | Рунге-Кутта 2-го порядка | ||||
37. | 0,62 | 0,35 | 0,5 | 0,25 | - | - | 0,9 | Эйлера | |||||
38. | 0,9 | 0,7 | 0,7 | 0,6 | 0,5 | 0,4 | 0,85 | Эйлера-Коши | |||||
39. | 0,57 | 0,42 | 0,3 | 0,25 | - | - | 1,6 | - | Рунге-Кутта 2-го порядка | ||||
40. | 0,09 | 0,06 | 0,05 | - | - | - | 1,2 | - | Эйлера |
Литература
1. Л. И. Турчак. Основы численных методов: Учебное пособие.-М.: Наука, 1987.
2. Методические указания. Выполнение лабораторных работ по вычислительной технике, 1979, №96.
3. Электронные вычислительные машины: Учебное пособие для вузов, Кн.3. Алгоритмизация и основы программирования.-М.: Высшая школа, 1987.
Редактор Н.В. Калачаева