Разработка стохастической модели в форме уравнения регрессии в пассивном эксперименте
Вариант №21
Выполнил: ст.гр. ТС 03-01
Валеев Я.Х.
Проверил: проф. Самойлов Н.А.
УФА 2007
Цель работы:
· рассмотреть принципы разработки стохастической модели на примере конкретной задачи;
· выполнить расчет коэффициентов уравнения регрессии;
· оценить адекватность полученного уравнения регрессии исходному эксперименту;
· закрепить методологию компьютерной обработки табличных материалов по алгоритму прямой задачи.
 |
Исходные данные:
Вариант №21
Зависимость константы скорости реакции K от температуры t приведена ниже. Рассмотреть возможность описания функции K=f(t) уравнением Аррениуса K=K0e-E/RT, где K0 – предэкспоненциальный множитель; E – энергия активации; T – абсолютная температура; R – универсальная газовая постоянная, и рассчитать коэффициенты уравнения Аррениуса.
| Температура, оС | |||||||||
| К, с-1 | 0,05 | 0,10 | 0,20 | 0,40 | 0,68 | 1,00 | 1,54 | 2,35 | 4,63 |
Для решения первой части (расчет константы скорости реакции, коэффициентов уравнения Аррениуса) поставленной задачи воспользуемся следующим алгоритмом:
1. Представим уравнение Аррениуса, прологарифмировав его, в виде: 
.
2. Рассмотрим возможность описания зависимости константы скорости реакции от температуры по методу наименьших квадратов (в виде линейной зависимости). Для этого построим график зависимости lnK = f(t).
3. При возможности описания указанной выше зависимости по методу наименьших квадратов произведем расчет lnK, а в дальнейшем и K.
4. Произведем расчет коэффициентов уравнения Аррениуса.
| Температура, оС | |||||||||
| ln(К) | -2,9957 | -2,3026 | -1,6094 | -0,9163 | -0,3857 | 0,0000 | 0,4318 | 0,8544 | 1,5326 |
График зависимости lnK = f(t).
Текст программы
После визуальной оценки графика зависимости сделан вывод, что, в первом приближении, возможно описание зависимости в форме уравнения регрессии вида y=a+bx. Для решения поставленной з 
адачи была написана программа в среде Turbo Pascal, представленная ниже.
program lab4;
uses crt;
const n=2;nn=9;tabl_stroki=13; tabl_stolb=12;
fisher: array [0..tabl_stroki,0..tabl_stolb] of real =
((0.0 ,1.0 ,2.0 ,3.0 ,4.0 ,5.0 ,6.0 ,7.0 ,8.0 ,9.0 ,10 ,12 ,16),
(1.0 ,161 ,200 ,216 ,225 ,230 ,234 ,237 ,239 ,241 ,242 ,244 ,246),
(2.0 ,18.5,19.0,19.1,19.2,19.3,19.3,19.3,19.3,19.3,19.3,19.4,19.4),
(3.0 ,10.1,9.55,9.28,9.12,9.01,8.94,8.88,8.84,8.81,8.78,8.74,8.69),
(4.0 ,7.71,6.94,6.59,6.39,6.26,6.16,6.09,6.04,6.00,5.96,5.91,5.84),
(5.0 ,6.61,5.79,5.41,5.19,5.05,4.95,4.88,4.82,4.78,4.74,4.68,4.60),
(6.0 ,5.99,5.14,4.76,4.53,4.39,4.28,4.21,4.15,4.10,4.06,4.00,3.92),
(7.0 ,5.59,4.74,4.35,4.12,3.97,3.87,3.79,3.73,3.68,3.63,3.57,3.49),
(8.0 ,5.32,4.46,4.07,3.84,3.69,3.58,3.50,3.44,3.39,3.34,3.28,3.20),
(9.0 ,5.12,4.26,3.86,3.63,3.48,3.37,3.29,3.23,3.18,3.13,3.07,2.98),
(10.0,4.96,4.10,3.71,3.48,3.33,3.22,3.14,3.07,3.02,2.97,2.91,2.82),
(12.0,4.75,3.88,3.49,3.26,3.11,3.00,2.92,2.85,2.80,2.76,2.69,2.60),
(14.0,4.60,3.74,3.34,3.11,2.96,2.85,2.77,2.70,2.65,2.60,2.53,2.44),
(16.0,4.49,3.63,3.24,3.01,2.85,2.74,2.66,2.59,2.54,2.49,2.42,2.33));
t:array [1..nn] of real =(100,150,200,250,300,350,400,450,500);
var i:integer;
sx,sx2,sy,sxy:real;
y,yy,w:array [1..nn] of real;
eps,sigma_y,sigma_ost,ysr,F:real;
stroka,stolb:integer;
a:array [1..n,1..n] of real;
b,xx:array [1..n] of real;
s,z,zz,ko,e:real;
j,k:integer;
begin
clrscr;
sx:=0;
sx2:=0;
sy:=0;
sxy:=0;
ko:=0;
e:=0;
z:=0;
zz:=0;
{Ввод данных}
for i:=1 to nn do begin
write ('t= ',t[i]:6:1,'(град) k = '); readln (w[i]);
y[i]:=ln(w[i]);
end;
{Метод наименьших квадратов}
for i:=1 to nn do begin
sx:=sx+t[i];
sx2:=sx2+t[i]*t[i];
sy:=sy+y[i];
sxy:=sxy+t[i]*y[i];
end;
writeln ('sx= ',sx:5:5);
writeln ('sx2= ',sx2:5:5);
writeln ('sy= ',sy:5:5);
writeln ('sxy= ',sxy:5:5);
readln; 
b[1]:=sy;
b[2]:=sxy;
a[1,1]:=nn; a[1,2]:=sx;
a[2,1]:=sx; a[2,2]:=sx2;
for k:=1 to n-1 do begin
for i:=(k+1) to n do begin
s:=-a[i,k]/a[k,k];
for j:=k to n do
a[i,j]:=a[i,j]+s*a[k,j];
b[i]:=b[i]+s*b[k];
end; end;
xx[n]:=b[n]/a[n,n];
for i:=(n-1) downto 1 do begin
s:=0;
for j:=(i+1) to n do
s:=s+a[i,j]*xx[j];
xx[i]:=(b[i]-s)/a[i,i];
end;
{Вывод результатов расчетов по методу наименьших квадратов}
writeln ('Зависимость имеет вид:');
writeln ('k=a+b*t, где');
writeln ('a= ', xx[1]:5:9);
writeln ('b= ', xx[2]:5:9);
readln;
writeln ('T (град) ','K (опыт) ','K (расч) ','Ошибка ');
for i:=1 to nn do begin
yy[i]:=xx[1]+xx[2]*t[i];
eps:=(yy[i]-y[i])*(yy[i]-y[i]);
z:=yy[1];
zz:=yy[9];
writeln (t[i]:7:2,w[i]:10:4,exp(yy[i]):15:4,eps:10:2);
end;
{Расчет коэффициентов уравнения Аррениуса}
E:=8.31*373.15*773.15*(zz-z)/(773.15-373.15);
Ko:=exp(zz)/(exp(-E/(8.31*773.15)));
writeln ('Предэкспоненциальный множитель Ko=',ko);
writeln ('Предэкспоненциальный множитель E=',e:8:5);
readln;
{Регрессионный анализ}
ysr:=sy/nn;
for i:=1 to nn do
begin
sigma_y:=sigma_y+sqr(y[i]-ysr);
sigma_ost:=sigma_ost+sqr(y[i]-yy[i]);
end;
sigma_y:=sigma_y/(nn-1);
sigma_ost:=sigma_ost/(nn-n);
writeln ('Дисперсия опытных данных равна ',sigma_y:10:2);
writeln ('Остаточная дисперсия равна ',sigma_ost:10:2);
F:=sigma_y/sigma_ost;
writeln ('Расчетный критерий Фишера равен F=',f:10:2);
stroka:=1;stolb:=1;
for i:=1 to tabl_stroki do
if (nn-n)=round (fisher[i,0]) then stroka:=i;
for i:=1 to tabl_stolb do
if (nn-1)=round (fisher[ 
0,i]) then stolb:=i;
writeln ('Табличный критерий Фишера равен F=',fisher[stroka,stolb]:5:2);
if fisher[stroka,stolb]<f then writeln ('Модель адекватна')
else writeln ('Модель неадекватна');
readln;
end.
Список идентификаторов
| Параметр | Расшифровка |
| w | Опытное значение константы скорости |
| y | Логарифм опытного значения константы скорости |
| yy | Логарифм расчетного значения константы скорости |
| sx | Сумма температур |
| sx2 | Сумма квадратов температур |
| sy | Сумма констант скоростей |
| sxy | Сумма величины ТК |
| eps | Величина квадрата отклонения |
| sigma_y | Дисперсия опытных данных |
| sigma_ost | Остаточная дисперсия |
| ysr | Вспомогательная величина |
| F | Расчетный критерий Фишера |
| stroka | Параметр для выбора табличного значения критерия Фишера |
| stolb | Параметр для выбора табличного значения критерия Фишера |
| xx[1] | Коэффициент а уравнения регрессии |
| xx[2] | Коэффициент b уравнения регрессии |
| ko | Предэкспоненциальный множитель |
| e | Энергия активации |
| a,b,s,j,k | Вспомогательные параметры |
| z | Величина логарифма расчетного значения константы скорости для температуры 100оС |
| zz | Величина логарифма расчетного значения константы скорости для температуры 500оС |
Бло 
к-схема
Результаты расчетов
В ходе расчетов получены следующие результаты.
Уравнение регрессии в виде зависимости lnK = f(t):
lnK = -3.8573 + 
0.0109 t.
| Температура, оС | |||||||||
| ln(К) | -2,9957 | -2,3026 | -1,6094 | -0,9163 | -0,3857 | 0,0000 | 0,4318 | 0,8544 | 1,5326 |
| ln(К) (расч.) | -2,7712 | -2,2281 | -1,6851 | -1,1420 | -0,5990 | -0,0559 | 0,4871 | 1,0301 | 1,5732 |
| Ошибка | 0,05 | 0,01 | 0,01 | 0,05 | 0,05 | 0,00 | 0,00 | 0,03 | 0,00 |
По полученным данным строится график зависимостей lnK = f(t).
График зависимости lnK = f(t). Пунктирной линией изображен график зависимости для опытных значений; сплошной – для расчетных.
По найденным значениям lnK рассчитываются значения K. По полученным данным строится график зависимости K = f(t).
| Температура, оС | |||||||||
| К, с-1 (опытн.) | 0,05 | 0,10 | 0,20 | 0,40 | 0,68 | 1,00 | 1,54 | 2,35 | 4,63 |
| К, с-1 (расч.) | 0,0626 | 0,1077 | 0,1854 | 0,3192 | 0,5494 | 0,9456 | 1,6276 | 2,8015 | 4,8220 |
| Ошибка | 0,05 | 0,01 | 0,01 | 0,05 | 0,05 | 0,00 | 0,00 | 0,03 | 0,00 |
График зависимости K = f(t). Пунктирной линией изображен график зависимости для опытных значений; сплошной – для расчетных.
По рассчитанным значениям рассчитываются величины энергии активации и предэкспоненциального множителя соответственно по формулам:
,
.
Приведенные выше формулы могут использоваться лишь в первом приближении для оценки величин значений, и не дают точного адекватного результата. Для более точных оценок не 
обходимо использовать графический метод (в работе не приводится).
Величина энергии активации: 
.
Предэкспоненциальный множитель: 
.
Качество разработанного уравнения регрессии, его адекватность, оценивалась по критерию Фишера.
Получены следующие результаты:
Дисперсия опытных данных 2.24
Остаточная дисперсия 0.03
Расчетный критерий Фишера 79.51
Табличный критерий Фишера 3.73
Так как, расчетный критерий Фишера на порядок превышает табличное значение, то можно говорить о высоком уровне адекватности уравнения регрессии.
Используя все полученные данные критерия Фишера и визуальную оценку графика зависимости K = f(t) для опытн 
ых и расчетных значений констант скоростей, можно сделать вывод, что уравнение найденное уравнение регрессии адекватно, программа работает правильно, и функция K = f(t) может описываться уравнением Аррениуса.
Выводы.
В ходе выполнения лабораторной работы были рассмотрены принципы разработки стохастической модели на примере расчета зависимости константы скорости реакции от температуры. Для получения уравнения регрессии использовалось 
прологарифмированное уравнение Аррениуса, представленное в виде:
.
После построения графика зависимости lnK = f(t) была выдвинута гипотеза, что зависимость может быть представлена в форме уравнения регрессии вида: y=a+bx. Адекватность модели оценивалась по величине критерия Фишера. Расчетное значение на порядок превышало табличное, следовательно, разработанная модель адекватна, программа работает правильно.
После расчета констант скорости реакции была построена зависимость K = f(t). При визуальной оценке графиков зависимостей для опытных и расчетных значений видно, что характер изменения практически совпадает для двух графиков.
Для рассчитанных величин были найдены значения энергии активации и предэкспоненциального множителя соответственно по формулам:
,
.
Значения составили:
Величина энергии активации .
Предэкспоненциальный множитель .