УДК 517.9
ЧИСЛЕННОЕ ИССЛЕДОВАНИЕ ПРОЦЕССОВ С ПЕРЕМЕННЫМ РЕАКЦИОННЫМ ОБЪЕМОМ В УСЛОВИЯХ НЕОПРЕДЕЛЕННОСТИ КИНЕТИЧЕСКИХ ДАННЫХ
© В. А. Вайтиев*, С. А. Мустафина
Башкирский государственного университета, Стерлитамакский филиал Россия, Республика Башкортостан, 453103 г. Стерлитамак, пр. Ленина, 47а. Тел./факс: +7 (347) 343 10 39.
E-mail: [email protected]
На основе методов интервального анализа получены границы решения прямой задачи химической кинетики, обусловленного интервальным представлением кинетических параметров.
Ключевые слова: олигомеризация а-метилстирола, кинетическая модель, прямая задача, кинетические параметры, интервальный анализ, оптимальное двустороннее решение.
Введение
При исследовании математических моделей с помощью численных методов возникает проблема решения задач, в которых часть параметров моделей задана неточно. К таким задачам относится и прямая задача химической кинетики, точность решения которой зависит от точности кинетических параметров. Как правило, кинетические параметры, к которым относятся константы скорости и энергии активации, определяются экспериментально либо при решении обратной задачи. Это означает, что кинетические параметры содержат погрешность, величина которой может достигать более 20%. Поэтому при решении прямой задачи целесообразно задавать кинетические параметры не числами, а интервалами. Отсюда возникает задача построения гарантированных оценок приближенных решений, а также разработки алгоритмов и комплекса программ, позволяющих проводить анализ чувствительности полученного решения к вариации кинетических параметров.
В работе построен метод и алгоритм учета влияния неопределенности в кинетических параметрах на решение прямой задачи. На примере реакции олигомеризации а-метилстирола проведен вычислительный эксперимент, позволяющий оценить степень чувствительности решения прямой задачи к вариации констант скоростей.
Построение кинетической модели реакции олигомеризации а-метилстирола
Продукты реакции олигомеризации а-метилстирола являются ценным нефтехимическим сырьем. Насыщенный циклический димер - 1,1,3-триметил-3-фенилиндан представляет интерес в качестве смазочного или изоляционного материала, стойкого к радиолизу теплоносителя, реактивного топлива. Фенилиндан применяют в качестве пластификатора оргстекла и сцинтилляционных счетчиков ядерных излучений. Линейные ненасыщенные димеры - 4-метил-2,4-дифенилпентен-1 и 4-метил-2,4-дифенилпентен-2 используются как растворители для лаков, диэлектрические жидкости, основа смазочных масел, материалы для получения фрикционных жидкостей, пластификаторы каучуков.
Для построения модели введем следующие обозначения: А1 - а-метилстирол; А2 - 4-метил-2,4 дифенилпентен-1; А3 - 4-метил-2,4 дифенилпентен-2; А4 - 1,1,3-триметил-3-фенилиндан; А5 - триммеры.
Совокупность химических превращений, описывающих реакцию олигомеризации а-метилстирола с учетом введенных обозначений представляется следующей схемой стадий [1]:
2А1 А2,
2А1 А3,
2А1 ~^А2<
А2 — А4,
А3 —— А4,
А1 +А2 -—А5>
А1 +А3 —А5>
А1 +А4 —— А5,
А2 0А3.
Согласно закону действующих масс кинетические уравнения, соответствующие схеме химических превращений, можно выразить уравнениями:
61 = к1 х1 — кю%2,
62 = к2 Х12 — кцХз,
6 = кз Х1 ,
6)4 = к4 Х2,
6 = к5 Хз,
6 = кб Х2 Х1,
6 = ку Х1 Хз, б)^ = к$ Х1Х4,
69 = кд Х2 — кпХз,
где б(1,Х,к) - скорость г-й стадии (моль/(л-ч)), г'=1,2,...,9; 0,1] - время реакции (ч); х=(х1,х2,...,х5)т
- вектор концентраций компонентов (мольные доли); к=(к1,к2,...,к12)Т - вектор параметров - кинетических констант скоростей у-й реакции, причем к^=к]0'вхр(-Е/(ЯТ)), ]=1,2,...,12, Еу - энергии активации (Дж/моль), Я - универсальная газовая постоянная (Дж/(моль-К)), Т=Т(г) - температура (К).
Матрица стехиометрических коэффициентов реакции имеет вид, представленный в табл. 1.
Экспериментальные данные получены в тер-мостатируемом лабораторном реакторе с мешалкой периодического действия. Такой реактор достаточно корректно может быть описан в приближении идеального смешения.
* автор, ответственный за переписку
Таблица 1
Матрица стехиометрических коэффициентов реакции олигомеризации а-метилстирола
СЯ1 | Я?. 1 с | Я4 1 с | Сйб | Я? 1 Ян 1 (Яд
X1 -2 -2 -2 0 0 -1 -1 -1 0
Х2 1 0 0 -1 0 -1 0 0 -1
Хз 0 1 0 0 -1 0 -1 0 1
х4 0 0 1 1 1 0 0 -1 0
Х5 0 0 0 0 0 1 1 1 0
£ -1 -1 -1 0 0 -1 -1 -1 0
При разработке математического описания данного реактора учитывается изменение числа молей (или реакционного объема) в ходе протекания химических реакций. Это отражается стехиометрическими коэффициентами последней строки табл. 1. Тогда уравнения материального баланса процесса имеют вид [2]:
= р , где р = £ у 6 , (1)
йг г 1 11 1
с начальными условиями при г=0: хг=хг0. Здесь N - относительное изменение числа молей реакционной среды, Уу - стехиометрические коэффициенты веществ, участвующих в реакции. Систему уравнений (1) замыкает условие нормировки по
5
компонентам жидкой
фазы Т хі = 1 [3].
і=1
Продифференцировав уравнение (1), получим
(2)
йг йг 1 йг
Просуммировав обе части уравнения (2), по лучим
N ^ і) +—х. = Е. или — = Е’
“г “г 1=1 ‘ ІҐ ‘ “г п
где
Е=Тс Ек}.
і=і ‘=1
Таким образом, система уравнений (2) примет
вид:
Т = = М-X. * )■ = 5'
“г N
N
“г
(3)
= ¥п = /6(г, х, *)>
с начальными условиями: хг(0)=хг , 1=1,...,5, N=1, где хг - концентрация г-го компонента (мольные доли).
Полученная система дифференциальных уравнений представляет собой кинетическую модель сложной (многостадийной) реакции олигомеризации а-метилстирола, учитывающую изменение числа молей в ходе ее проведения. Правые части уравнения (2) с учетом матрицы стехиометрических коэффициентов следующие:
¥1=-2 ю1-2 6-2 63- 66- бу- 68, р2=б1-б4-бб-бд, р3= 66-6-6+бд, р4=б3+б4+б5-б8, ¥5=6+67+6 FN=-6-66-66-б-6Ь-6.
(4)
Интервальный метод анализа чувствительности решения прямой задачи к погрешности кинетических параметров
Основная идея метода заключается в анализе частных производных правых частей дифференциальных уравнений, входящих в запись кинетической модели (з) исследуемой реакции, по константам скоростей ку и по искомым величинам хг. Для реализации этого метода используется аппарат интервального анализа.
Введем обозначение для относительного изменения числа молей реакционной среды, участвующего в записи кинетической модели реакции олигомеризации а-метилстирола: ^х^г=хб. Тогда система дифференциальных уравнений (з) с начальными условиями запишется в виде: х/=/(г,х,к), 1=1,...,6,
х(0)= х0, 1=1,...,б,
пП *-> ^
где хек - вектор значений концентраций веществ, к - вектор параметров (вектор кинетических констант), ге[0,1] - время реакции, х(0)= х0еЯ -вектор начальных значений. Далее будем считать, что х=х(г,к,х0).
Будем искать решение системы (4) в виде х = (х1,...,х6)Т е X , где
х1 = [хг] = |х; е Я|хг < х1 < Хг}, х, Хг - нижние и
верхние границы компонентов вектора неизвестных соответственно [4]. Аналогичное представление характерно для вектора констант скоростей реакции к = (к1,...,к12)Т е к.
В работе [5] представлены утверждения, позволяющие построить алгоритм интервального решения задачи Коши с заданными начальными условиями. Воспользуемся данными формулировками для решения прямой задачи химической кинетики.
Оценим х<г) - верхнюю границу х(г) по г-й
_(‘') ^ тт
координате, х > х. • Для оценки х рассмотрим систему однородных дифференциальных уравне-
х' = /(г,х,~), ~ є X, ~ є К . ~(0) = ~0 є Х0.
Здесь
\ку, если Xк (г) < 0, ~у = и у, если хк (г) > о, [к у, если о е хк (г),
(5)
9 5
и
x0j, если X 0 (t) < 0, x 0 j, если X 0. (t) > 0, x если 0 e X0 (t),
(7)
гДе Xk (t) и x 0 (t) - интервальные расширения dx /Эк,, и 3xi /Эх0. соответственно.
(g)
Интервальные функции хк (г), X 0. (г) можно
определить, одновременно решая систему (4) и системы однородных дифференциальных уравнений
хк = £ (г, х, к) хЦ + ^ (г, х, к )>
I=1 дх1 дХк
хк (0) = 0, г = 1.....6, у = 1,...,12,
х0 = Ц ^ х к) х1'
1=1 дх„ (9)
х" (0) = 8у, I,у = I,...,6,
где - символ Кронекера. При выполнении условий
0(0,x0,k)^ 0g
f
dxr
(0, x0 ,k),
(10)
Эх.
ку --к
можно говорить о существовании такого момента времени г0>0, при котором 0 £ Xк. (г),
0 £ X 0 (г). До момента времени гд границы полу-
ченного множества решений являются оптимальными. При 0 е Xк (г), 0 е X 0 (г) система (5) содержит в себе интервальные параметры, в этом случае решение системы (5) можно получить, к примеру, двусторонним методом, идея которого основывается на понятиях изотонности и антитон-ности функции по переменным и на использовании двусторонних неравенств [5-6].
Введем в рассмотрение вспомогательные вещественные функции 0(г,у1,у2,...уп,11,12,...,1п), определяющие границы двустороннего решения, 1=1,2, такие что при у < у и £ < % выполнены следующие неравенства
0‘ (г, у, £) < 0‘ (г, у у], £),
l = 1,2, i = 1,..., n,
(11)
где у[г‘] = (у19..., у._1, уг-+1уп ). Кроме того, выполнены неравенства
Gг1 (г, х, х) < fi (г, х, к) < Gг (г, х, х).
Любое решение х системы (4) с начальным условием х0 < х(0) < Х0 удовлетворяет оценкам Х < х(г) < Х1 в том случае, если вектор-функции х,хе Яп удовлетворяют соотношениям
х_<G1(г,х,х)> х >G2(г,х,х), х(0) < Х0, х(0) > Х0.
Заметим, что при этом
G1(t, x, x) < inf f (t, x, к),
G2(t, x, x) > sup f (t, x, к), и в качестве G можно
взять границы монотонного по включению интервального расширения функции f - правых частей дифференциальных уравнений в записи кинетической модели реакции.
Пусть функция f(x,y) изотонная по x и анти-тонная по y, тогда интервальное расширение f(x,y) имеет вид
f(X, У) = f (x, y), f (X, y) = f (x, y).
Систему (11) можно перезаписать таким образом
x_t < f.(t,X[Xi],k), x i > (t,X[Xi],k), (2
—l - - (12) x(0) < x0, x(0) > x0*
Решение построенной системы дает в общем случае более широкое, чем оптимальное, двустороннее решение.
При выполнении следующих условий
f ф 0’ i ф j > i, j = 1,•••,n>
dx,
(13)
м§п^ = СОп$г, у: (к у) * °> Vk е к,
йку (14)
Ух е X,
систему (12) можно представить в виде
Х = / (г, х, к'), хг = / (г, х, к2),
х(0) = Х0, х(0) = Х0, где (к.) - ширина интервальной величины
(15)
к ек; Х и Х являются некоторыми частными решениями исходной системы и, следовательно, они оптимальны.
Результаты вычислительного эксперимента и их обсуждение Вычислительный эксперимент реакции олигомеризации а-метилстирола проводился при начальных концентрациях веществ х1(0)=1, х(0)=0,
1=2,...,5, Х}(0)=1 и времени реакции ге [0,5] ч.
При решении прямой задачи интервальным методом были использованы численные значения кинетических констант (табл. 2), определенные при решении обратной кинетической задачи на базе математического описания процесса олигомеризации а-метилстирола и экспериментальных данных, полученных в термостатируемом лабораторном реакторе с мешалкой периодического действия при температуре 80 °С (353.15 К) и 10%-ной концентрации катализатора «Цеокар-10».
Для системы уравнений, описывающей кинетическую модель реакции олигомеризации а-метилстирола, как и для большинства уравнений химической кинетики, выполняются условия (1з). Условия (14), как правило, наоборот не выполняются, что подтверждается в том числе и в нашем случае. Следовательно, решая систему (з), приведенную к виду (12), неравенства в которой заменяются част-
0
ным случаем равенств, обыкновенным двусторонним методом, не удается получить оптимальное решение [7]. Двустороннее решение получается шире оптимального. Значительного уменьшения диапазона решения при использовании двустороннего метода можно добиться за счет введения дифференциального уравнения, учитывающего дополнительную информацию о протекании химической реакции. Однако выполнение условия (10) говорит о существовании такого момента времени г0> 0, до которого возможно построение оптимальных границ множества решений системы обыкновенных дифференциальных уравнений, соответствующих кинетической модели исследуемой реакции.
На рис. 1, 2 представлено графическое решение прямой кинетической задачи для реакции олигомеризации а-метилстирола. Интервальное решение получено двумя методами. В качестве кинетических констант подобраны интервалы, содержащие 5%-ную погрешность от экспериментальных значений.
Предложенный метод интервального анализа чувствительности при 5%-ной погрешности от экспериментальных значений позволяет построить границы интервального решения, наиболее близкие к границам оптимального решения (в сравнении с двусторонним методом). При этом можно утверждать, что до момента времени г0=0.6 ч. полученное решение само является оптимальным.
Таблица 2
Численные значения кинетических констант при температуре 80 °С и 10% концентрации катализатора
| к | к2 | к3 | к4 | к5 | к6 | 1.222 0.095 0.212 0.066 0.010 0.368
к7 к8 к9 к10 к11 к12
0.024 0.861 0.105 10-5 0.000 0.114
время, ч.
Рис. 1. Изменение концентрации продукта 4-метил-2,4-дифенилпентен-1: 1 - точное решение, 2 - двустороннее решение, 3 - решение, полученное методом интервального анализа чувствительности.
время, ч.
Рис. 2. Изменение концентрации продукта 4-метил-2,4-дифенилпентен-2: 1 - точное решение, 2 - двустороннее решение, 3 - решение, полученное методом интервального анализа чувствительности.
Анализ результатов, полученных при изменении погрешности кинетических параметров в диапазоне от 5% до 10%, позволяет сделать вывод, что погрешность экспериментальных данных не изменяет динамику кривых, при этом выход основных продуктов чувствителен к погрешности не более чем на 15%-30%. Кроме того, при увеличении погрешности в кинетических параметрах уменьшается временной интервал, на котором возможно построение оптимальных границ множества решений прямой задачи предложенным методом.
ЛИТЕРАТУРА
1. Байтимерова А. И., Мустафина С. А., Спивак С. И. Поиск оптимального управления в каскаде реакторов для процессов с переменным реакционным объемом // Системы управления и информационные технологии. 2008. №2(32). С.38-42.
2. Мустафина С. А., Байтимерова А. И., Степашина Е. В. О свойствах решения задач моделирования каталитических процессов с переменным реакционным объемом // Труды Средневолжского математического общества. 2010. Т. 12. №3. С. 145-150.
3. Мустафина С. А., Валиева Ю. А., Давлетшин Р. С. Оптимальные технологические решения для каталитических процессов и реакторов // Кинетика и катализ. 2005. №5. С.
749-757.
4. Шокин Ю. И. Интервальный анализ. Новосибирск: Наука, 1981. 113 с.
5. Добронец Б. С. Интервальная математика. Красноярск: КГУ, 2004. 216 с.
6. Курпель Н. С., Шувар Б. А. Двусторонние неравенства и их приложения. Киев: Наук.думка, 1980. 268 с.
7. Вайтиев В. А., Мустафина С. А. Чувствительность решения прямой задачи химической кинетики к вариации кинетических констант // Математические методы в технике и технологиях - ММТТ-25. Сб. трудов XXV Международной научной конференции в 10 т. Т. 7. Волгоград. Волгоградский гос. техн. Ун-т. 2012. С. 8-9.
Поступила в редакцию 31.10.2012 г.