Практическая часть. Расчет кинетики реакции

Расчет кинетики реакции

Исходные данные:

A B

C E F G

D

G-?

Начальные концентрации: C_A0= 200г/л; C_B0=0г/л; C_C0= 130г/л; C_D0= 0г/л; C_E0= 0г/л; C_F0= 10г/л; C_G0= 0г/л;

Константы скоростей: k₁=0.003 c^-1; k₂=0.09 c^-1; k₃=0.2 c^-1; k₄=0.07 c^-1; k₅=0.055 c^-1; k₆=0.02 c^-1; k₇=0.035 c^-1;

переделать

На рисунке 3 приведена блок-схема решения задачи в общем виде, где концентрации компонентов Сi, представленных как А, В, С, D, E, F, G для предыдущих (n-x) расчетов и как А1, В1,С1, D1, E1, F1 и G1 для последующих (n+1-x) расчетов. При достижении максимальной величины концентрации G, равной СG mах, решение задачи прекращается.

Найдем основной компонент G методом дифференциальных уравнений:

Код программы

CLS

ca0 = 300: cb0 = 0: cc0 = 130: cd0 = 0: ce0 = 50: cf0 = 10: cg0 = 0

k1 =.003: k2 =.09: k3 =.02: k4 =.07: k5 =.055: k6 =.02: k7 =.035

INPUT "dt="; dt

t = 0

a = ca0: b = cb0: c = cc0: d = cd0: e = ce0: f = cf0: g = cg0

30 a1 = a + (-k1 * a - k2 * a - k3 * a) * dt

b1 = b + k1 * a * dt

c1 = c + (k2 * a - k4 * c) * dt

d1 = d + k3 * a * dt

e1 = e + (k4 * c - k5 * e + k6 * f) * dt

f1 = f + (k5 * e - k6 * f - k7 * f) * dt

g1 = g + k7 * f * dt

IF g1 > g THEN 10 ELSE 20

10 a = a1: b = b1: c = c1: d = d1: e = e1: f = f1: g = g1: t = t + dt

GOTO 30

20 PRINT "g=", g, "t=", t

END

РЕЗУЛЬТАТЫ ВЫЧИСЛЕНИЯ ПРОГРАММЫ

dt=? 1

g= 428.9373 t= 638

Поскольку точность решения задачи по методу Эйлера зависит от величины шага интегрирования (чем меньше шаг, тем точнее решение, но больше продолжительность расчета), рассмотрим принципы автоматического обеспечения правильности расчета, задав погрешность определения концентрации СG.

В число исходных данных вводим произвольно малое значение Gmax, например, Gmax=0. После определения лучшей концентрации компонента G в ходе расчета кинетики с шагом Δτ и сравнения G с Gmax расчет полностью повторяют с уменьшенным в 2-10 раз шагом интегрирования приемлемого совпадения G и Gmax. Как видно из рисунка 2 получение верного решения достигается дополнением блок-схемы принципиального решения (рисунок 2) всего двумя блоками (8 и 9).

Рисунок 2 – Блок-схема расчета кинетики реакции с автоматическим обеспечением правильности расчета

Код программы

CLS

ca0 = 300: cb0 = 0: cc0 = 130: cd0 = 0: ce0 = 50: cf0 = 10: cg0 = 0

k1 =.003: k2 =.09: k3 =.02: k4 =.07: k5 =.055: k6 =.02: k7 =.035

INPUT "dt="; dt

60 t = 0

a = ca0: b = cb0: c = cc0: d = cd0: e = ce0: f = cf0: g = cg0

30 a1 = a + (-k1 * a - k2 * a - k3 * a) * dt

b1 = b + k1 * a * dt

c1 = c + (k2 * a - k4 * c) * dt

d1 = d + k3 * a * dt

e1 = e + (k4 * c - k5 * e + k6 * f) * dt

f1 = f + (k5 * e - k6 * f - k7 * f) * dt

g1 = g + k7 * f * dt

IF g1 > g THEN 10 ELSE 20

10 a = a1: b = b1: c = c1: d = d1: e = e1: f = f1: g = g1: t = t + dt

GOTO 30

20 IF ABS(g - gmax) <.01 * g THEN 40 ELSE 50

50 gmax = g: dt = dt / 2: PRINT "dt="; dt

GOTO 60

40 PRINT "g=", g, "t=", t

END

РЕЗУЛЬТАТЫ ВЫЧИСЛЕНИЯ ПРОГРАММЫ

dt=? 1

dt=.5

g= 428.9369 t= 609.5

Для обеспечения вывода информации только в реперной точке, когда фактическое время реакции совпадает с позицией реперной точки, можно воспользоваться двумя приемами.

1Способ:

Если продолжительность реакции примерно известна, то приняв Δτrep, можно сделать следующую вставку между блоками 4 и 5; во вставке j – время по оси реперных точек:

Код программы

CLS

ca0 = 300: cb0 = 0: cc0 = 130: cd0 = 0: ce0 = 50: cf0 = 10: cg0 = 0

k1 =.003: k2 =.09: k3 =.02: k4 =.07: k5 =.055: k6 =.02: k7 =.035

INPUT "dt="; dt

1 t = 0

a = ca0: b = cb0: c = cc0: d = cd0: e = ce0: f = cf0: g = cg0

30 a1 = a + (-k1 * a - k2 * a - k3 * a) * dt

b1 = b + k1 * a * dt

c1 = c + (k2 * a - k4 * c) * dt

d1 = d + k3 * a * dt

e1 = e + (k4 * c - k5 * e + k6 * f) * dt

f1 = f + (k5 * e - k6 * f - k7 * f) * dt

g1 = g + k7 * f * dt

FOR j = 0 TO 600 STEP 60

IF j = t + dt THEN 5 ELSE 6

5 PRINT a1; b1; c1; d1; e1; f1; g; j

6 NEXT j

IF g1 > g THEN 10 ELSE 20

10 a = a1: b = b1: c = c1: d = d1: e = e1: f = f1: g = g1: t = t + dt

GOTO 30

20 IF ABS(g - gmax) <.01 * g THEN 40 ELSE 50

50 gmax = g: dt = dt / 2: PRINT "dt="; dt: GOTO 1

40 PRINT "g=", g, "t=", t

END

РЕЗУЛЬТАТЫ ВЫЧИСЛЕНИЯ ПРОГРАММЫ

428.9316 540

1.703303E-29 7.964601 9.31933E-17 53.09734 5.97196E-04 9.903376E-04

428.9361 600

dt=.5

.2793506 7.957183 9.956594 53.04789 93.76899 119.8573 203.0227 60

2.601226E-04 7.964593.1460681 53.09728 24.36211 39.61463 364.1142 120

2.42218E-07 7.964597 2.03864E-03 53.09731 6.493336 10.75564 411.4967

2.255458E-10 7.964597 2.836082E-05 53.09731 1.738846 2.883373 424.2647

2.100213E-13 7.964597 3.944601E-07 53.09731.4657882.7724198 427.6861

1.955653E-16 7.964597 5.486319E-09 53.09731.1247737.2069138 428.6027

1.821043E-19 7.964597 7.630598E-11 53.09731 3.342399E-02 5.542741E-02

428.8481 420

1.695699E-22 7.964597 1.061295E-12 53.09731 8.95351E-03 1.484772E-02

428.9138 480

1.578982E-25 7.964597 1.476093E-14 53.09731 2.398437E-03 3.977357E-03

428.9314 540

1.470299E-28 7.964597 2.05301E-16 53.09731 6.424856E-04 1.065442E-03

428.9363 600

g= 428.9369 t= 609.5

Способ.

Если τр неизвестна, то вставка между блоками 4 и 5 примет вид:

Код программы

CLS

ca0 = 300: cb0 = 0: cc0 = 130: cd0 = 0: ce0 = 50: cf0 = 10: cg0 = 0

k1 =.003: k2 =.09: k3 =.02: k4 =.07: k5 =.055: k6 =.02: k7 =.035: dtr = 60

INPUT "dt="; dt

j = dt

60 t = 0: j = dtr

a = ca0: b = cb0: c = cc0: d = cd0: e = ce0: f = cf0: g = cg0

30 a1 = a + (-k1 * a - k2 * a - k3 * a) * dt

b1 = b + k1 * a * dt

c1 = c + (k2 * a - k4 * c) * dt

d1 = d + k3 * a * dt

e1 = e + (k4 * c - k5 * e + k6 * f) * dt

f1 = f + (k5 * e - k6 * f - k7 * f) * dt

g1 = g + k7 * f * dt

IF j = t + dt THEN 150 ELSE 160

150 PRINT a1; b1; c1; d1; e1; f1; g1; j

j = j + dtr

160 GOTO 170

170 IF g1 > g THEN 10 ELSE 20

10 a = a1: b = b1: c = c1: d = d1: e = e1: f = f1: g = g1: t = t + dt

GOTO 30

20 IF ABS(g - gmax) <.01 * g THEN 40 ELSE 50

50 gmax = g: dt = dt / 2: PRINT "dt="; dt

GOTO 60

40 PRINT "g=", g, "t=", t

END

РЕЗУЛЬТАТЫ ВЫЧИСЛЕНИЯ ПРОГРАММЫ

2.269211E-26 7.964601 7.251169E-15 53.09734 2.245724E-03 3.724113E-03

428.9317 540

1.703303E-29 7.964601 9.31933E-17 53.09734 5.97196E-04 9.903376E-04

428.9361 600

dt=.5

.2793506 7.957183 9.956594 53.04789 93.76899 119.8573 205.1326 60

2.601226E-04 7.964593.1460681 53.09728 24.36211 39.61463 364.8149 120

2.42218E-07 7.964597 2.03864E-03 53.09731 6.493336 10.75564 411.687 180

2.255458E-10 7.964597 2.836082E-05 53.09731 1.738846 2.883373 424.3157

2.100213E-13 7.964597 3.944601E-07 53.09731.4657882.7724198 427.6998

1.955653E-16 7.964597 5.486319E-09 53.09731.1247737.2069138 428.6063

1.821043E-19 7.964597 7.630598E-11 53.09731 3.342399E-02 5.542741E-02

428.849 420

1.695699E-22 7.964597 1.061295E-12 53.09731 8.95351E-03 1.484772E-02

428.9141 480

1.578982E-25 7.964597 1.476093E-14 53.09731 2.398437E-03 3.977357E-03

428.9315 540

1.470299E-28 7.964597 2.05301E-16 53.09731 6.424856E-04 1.065442E-03

428.9363 600

g= 428.9369 t= 609.5

Блок-схема расчета кинетики реакции

A B

C E F G

D

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

Код программы

CLS

ca0 = 300: cb0 = 0: cc0 = 130: cd0 = 0: ce0 = 50: cf0 = 10: cg0 = 0

k1 =.003: k2 =.09: k3 =.02: k4 =.07: k5 =.055: k6 =.02: k7 =.035

n = 200

dtr = 39

DIM a0(n), b0(n), c0(n), d0(n), e0(n), f0(n), g0(n), a(n), b(n), c(n), d(n), e(n), f(n), g(n)

INPUT "dt="; dt

FOR i = 1 TO n

a0(i) = 0: b0(i) = 0: c0(i) = 0: d0(i) = 0: e0(i) = 0: f0(i) = 0: g0(i) = 0

NEXT i

150 t = 0: i = 0: k = 0: j = dtr

a = ca0: b = cb0: c = cc0: d = cd0: e = ce0: f = cf0: g = cg0

a(i) = a: b(i) = b: c(i) = c: d(i) = d: e(i) = e: f(i) = f: g(i) = g

50 a1 = a + (-k1 * a - k2 * a - k3 * a) * dt

b1 = b + k1 * a * dt

c1 = c + (k2 * a - k4 * c) * dt

d1 = d + k3 * a * dt

e1 = e + (k4 * c - k5 * e + k6 * f) * dt

f1 = f + (k5 * e - k6 * f - k7 * f) * dt

g1 = g + k7 * f * dt

IF j = t THEN GOTO 10 ELSE 20

10 i = i + 1: k = k + 1

a(i) = a: b(i) = b: c(i) = c: d(i) = d: e(i) = e: f(i) = f: g(i) = g

j = j + dtr

20 IF g1 > g THEN GOTO 30 ELSE GOTO 40

30 a = a1: b = b1: c = c1: d = d1: e = e1: f = f1: g = g1

t = t + dt

GOTO 50

40 tOPT = t

gmax = g

FOR i = 1 TO k

IF ABS(a0(i) - a(i)) <.01 * a(i) THEN GOTO 70 ELSE GOTO 130

70 IF ABS(b0(i) - b(i)) <.01 * b(i) THEN GOTO 80 ELSE GOTO 130

80 IF ABS(c0(i) - c(i)) <.01 * c(i) THEN GOTO 90 ELSE GOTO 130

90 IF ABS(d0(i) - d(i)) <.01 * d(i) THEN GOTO 100 ELSE GOTO 130

100 IF ABS(e0(i) - e(i)) <.01 * e(i) THEN GOTO 110 ELSE GOTO 130

110 IF ABS(f0(i) - f(i)) <.01 * f(i) THEN GOTO 120 ELSE GOTO 130

120 IF ABS(g0(i) - g(i)) <.01 * g(i) THEN GOTO 60 ELSE GOTO 130

60 NEXT i

GOTO 140

130 dt = dt / 2: PRINT "dt="; dt

FOR i = 1 TO n

a0(i) = a(i): b0(i) = b(i): c0(i) = c(i): d0(i) = d(i): e0(i) = e(i): f0(i) = f(i): g0(i) = g(i)

NEXT i

GOTO 150

140 j = dtr

FOR i = 1 TO k

PRINT a(i); b(i); c(i); d(i); e(i); f(i); g(i); j

j = j + dtr

NEXT i

PRINT "tOPT="; tOPT; "gmax="; g

END

РЕЗУЛЬТАТЫ ВЫЧИСЛЕНИЯ ПРОГРАММЫ

dt=.125

dt=.0625

dt=.03125

dt=.015625

dt=.0078125

dt= 3.90625E-03

3.653958 7.867608 41.76391 52.45067 142.5556 132.2051 109.5032 39

4.450447E-02 7.963397 3.128239 53.08921 63.13806 91.61901 271.0163 78

5.420607E-04 7.964242.208884 53.0944 26.42174 42.72292 359.5831 117

6.602178E-06 7.964242 1.367833E-02 53.0944 11.21387 18.51671 399.1925

8.041361E-08 7.964242 8.924915E-04 53.0944 4.780559 7.922224 416.2341

9.79429E-10 7.964242 5.819509E-05 53.0944 2.039757 3.38219 423.5142 234

1.192925E-11 7.964242 3.794126E-06 53.0944.870446 1.443447 426.6267

1.452971E-13 7.964242 2.473593E-07 53.0944.3714626.6159993 427.9482

1.769688E-15 7.964242 1.612655E-08 53.0944.1585222.2628793 428.5064

2.155456E-17 7.964242 1.051368E-09 53.0944 6.764965E-02.1121842

428.8111 390

tOPT= 390.2383 gmax= 428.8129

Поиск по сайту: