CC BY
8
УДК 541.128.1
ХИМИЧЕСКИЕ НЕУСТОЙЧИВОСТИ ПРИ ОКИСЛЕНИИ 1,4-НАФТОДИОЛА В ГОМОГЕННОЙ СРЕДЕ. II. ТЕРМОДИНАМИЧЕСКИЙ АНАЛИЗ И МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
У.Г. Магомедбеков, Х.М. Гасанова, У.Г. Гасангаджиева, К.А. Абдулхамидов, С.Ш. Муцалова
(Дагетанский государственный университет, кафедра общей и неорганической химии; e-mail: ukhgmag@mail. ru)
Рассмотрен возможный механизм окисления 1,4-нафтодиола в колебательном режиме и составлена математическая модель, описывающая кинетику. На основе термодинамической функции Ляпунова показано, что причинами возникновения химических неустойчивостей является реализация автокаталитических и обратных динамических связей. Проведен качественный анализ и численное решение системы дифференциальных уравнений, моделирующих кинетику реакций. Определены характер стационарного состояния и возможность бифуркации. Показано, что математическая модель удовлетворительно описывает протекающие в системе процессы.
1. Термодинамический анализ процесса
окисления 1,4-нафтодиола в колебательном режиме
Одной из важных задач при изучении колебательных химических реакций является установление возможности и причин возникновения критических явлений на основе принципов неравновесной нелинейной термодинамики. Для определения термодинамической устойчивости системы в сильно неравновесных условиях применяют подход, связанный с оценкой производной второй вариации энтропии, которая для стационарного состояния является функцией Ляпунова [1]. На основании этого, если
(д/дг) рГ5 > о, то стационарное состояние вдали от равновесия устойчиво, а в
2
случае, когда
(д/дг) рГ 5 < 0 — оно неустойчиво и возможно возникновение диссипативных структур (р — средняя плотность).
Для оценки устойчивости системы использовано уравнение для производной от второй вариации энтропии в виде [1]:
I (д/дг)р^ЭйУ = |§(IVТ)8х^сР, + I —ЪТ^Щ +
Fs
Fs
J 5w5(AJ7) dV = J8( i/7) bxkvndFs + J 5w5( AJT) dV ■
Fs
KjFs( 5T) 2/f + pC Tv( 5T) 2/T,
(1)
где АЦ1 — сродство химической реакции, — скорость реакции, V — рабочий объем реактора, р — поверх-
ность реактора, ^ — тепловой поток через поверхность, !к — химический потенциал к-го компонента, V — объемный расход раствора, Т — температура,СТ — теплоемкость, Кт— коэффициент теплопередачи.
При рассмотрении процессов, протекающих в непроточном реакторе, член 15( 1к/Т)
ха-
Fs
рактеризующий обмен массой с окружающей средой, можно считать равным нулю (§Хк = 0), и уравнение для производной от второй вариации энтропии преобразуется к виду:
| 8и>8(А^/Т)dV + Кр( §Т21Т- + рСТч(§ТгТ (2)
V
В соответствии с соотношением (2) для оценки величины ( д/дг) р82SdV необходимо провести вычисление произведения 8^ 8 ( А к/Т) dV для определенных кинетических схем элементарных стадий протекания сложной реакции.
Для обоснования механизма процесса окисления 1,4-нафтодила (ОН2) в колебательном режиме в соответствии с литературными данными [2—4] были учтены следующие факты: 1) процесс проводили при рН ~ 7,77—8,14; в этих условиях 1,4-нафтодиол находится преимущественно в форме ОН - [2]; 2) окислителем выступает молекулярный кислород, связанный обратимо со смешаннолигандными комплексами кобальта (II) с БМО и Ру [3]; в) на образование частиц типа О. Н02" и НО- указано в [4].
+
С учетом этих соображений схема протекающих процессов представляется в виде:
1. <<н- + Со2ко24+ ^ 2Сок2+ + д- + но2;
2. <- + Н02' ^ < + И02",
3. < + <Н- ^ 2<- + Н+,
4. дн- + но2' ^ д- + н2о2,
5. д + но2" ^ д- + но2.
6. но2" + ho2" ^ h2o2 + о2,
7. 2Cok2+ + O2 ^ Co2kO24+,
(3)
где QH ", Q - и Q - различные формы превращений формы 1,4-нафтодиола, Co2kO24+ и Co 2 - оксигени-рованные и неоксигенированные формы комплексов кобальта (II). При введении обозначений [QH "]исх = a; [Co2O+] = b; [Q= х; [HO2] = y; [Q] = z; [QH ~] = a - х - z = h производная от термодинамической функции Ляпунова (избыточного производства энтропии) для приведенной кинетической схемы в явном виде имеет вид (процедура вычисления bw и b( Aw/T) показана в [5]):
bwß(AJT) dV = V{[( Ьх/х) 2Rw2+ (by/y)2 R( w2 + + w4 + 4w6) + (bz/z)2 R( w3 + w5) + (bh/h)2 R( w1 + + w3 + w5) + (Ъх/х) (by/y) R( 2w2- w4) + (Ъх/х) X (bz/z) R(- w2 - 2w3 - w5) + (Ъх/х) (bh/h) R( -w1 -
- 2w3 - Wj) + (by/y) (bz/z) R( -w2 - w5) + (by/y) X (bh/h) R( -w1+ 2w4) + (bz/z) (bh/h) 2Rw3] +
+ (bT/T2) (Ьх/х) [ -w1E1 - w2( Q2 - E2) - 2w3E3 -
- W4E4 - W5E5) ] + (bT/T2) (by/y) [-W1E1 - w( Q2 -
- E2) - w( Q - Ea) - W5E5 - w6( Q6 + E,) ] +
+ (bT/T2) (bz/z) [-w2E2 - w3( Q3 - E3) - w5( Q5 - E5) ] + + (bT/T2) (bh/h) [ - w1( Q1 - E1) - w3( Q3 - E3) -
W
4( Q4 - E4) ] + [(bT2/T)(-1/RT2) ■L_wtEiQi +
KTF/V + pvCjJV]}.
(4)
Анализ соотношения (4) показывает, что слагаемые, отвечающие за протекание прямых реакций,
2 2 равны: (Ъх/х) Rw2+ (Ъу/у) Я(w2 + w4 + 4w() +
(Ъх/х) 2R( w3 + w5) + (ЪН/Н) 2R( w1+ w3 + w5), и они будут положительными, что способствует устойчивости системы. В таких системах могут возникать тепловые неустойчивости в том случае, когда член (-1/RT2) £wiEiQi, характеризующий тепловой эффект химической реакции, становится больше члена связанного с теплоотводом, т.е. при довольно высоких тепловых эффектах реакции и недостаточном отводе тепла.
Кинетические схемы показывают, что автокаталитическое образование одного из компонентов - Q -связано с протеканием стадий 2 и 3. В этом случае в выражении (4) появляются отрицательные слагаемые:
[( Ьх/х) (bz/z) R( w2 + 2w3) + 2 (Ьх/х) (bh/h) Rw3 + + (by/y) (bz/z) RW2].
Если вклад этих членов будет преобладать, то процесс может становиться неустойчивым и возможно возникновение критических явлений.
Из кинетических схем (3) следует, что обратная связь осуществляются за счет образования Q - в результате взаимодействия QH2 с HO2' (стадии 1 и 4). В этом случае в выражении для производной второй вариации энтропии появляются слагаемые
-R[ (Ьх/х) (by/y) w4 + (Ьх/х) (bh/h) (w1 + w4) + + (by/y) (bh/h) W1],
которые имеют отрицательный знак. При преобладании вклада данных членов также вероятен переход системы в нестабильное состояние.
При рассмотрении слагаемых в производной термодинамической функции Ляпунова, отвечающих за возникновение термокинетических колебаний и связанных с автокатализом, получаем, что они имеют вид:
-{[(bT/T2) (Ьх/х) [w2( Q2- E2) + 2w3E3] +
+ (bT/T2) (by/y) [w2( Q2-E2) ] + (bT/T2) (bz/z) X
X[w2E2 + w3( Q3 - E3) ] + (bT/T2) (bh/h)
X W( Q3-E3) ]}.
Если повышение температуры будет способствовать увеличению скорости стадий 2 и 3, то в системе возможно накопление анион-радикала Q' (bT > 0 ^ Ьх > 0), и тогда рассматриваемые члены могут становиться более отрицательными, в результате чего будет усиливаться дестабилизация системы. Наличие обратных связей также может привести систему к нестабильности. С реализацией обратных связей связаны члены:
-{(bT/T2) (Ьх/х) [w1E1 + w4E4)] + (bT/T2) (by/y) [w1E1 + + w4( Q4-E4) ] + (bT/T2) (bh/h) [w1( Q1 - E1) + + w4( Q4 - E4) ]}.
В том случае, когда вклады этих слагаемых становятся преобладающими, могут проявляться критические явления.
Таким образом, проведенный термодинамический анализ на основе оценки производной второй вариации, как функции Ляпунова, позволяет заключить, что возникновение химических осцилляций в исследуемой системе является следствием наличия автокаталитических стадий, а также обратных динамических и термокинетических связей. Наряду с этим необходимо подчеркнуть, что предложенная кинетическая схема может быть рассмотрена как одна из альтернативных возникновению колебательного режима в исследуемой гомогенной системе.
2. Составление и анализ математической модели
Математические модели представляют собой системы обыкновенных дифференциальных уравнений нелинейного типа, полученных на основе рассмотрения механизмов протекания химических процессов. В настоящей работе в соответствии с кинетическими схемами окисления 1,4-нафтодиола в присутствии ок-сигенированных комплексов кобальта (II), после некоторых упрощений (учет размерности химической системы, постоянства pH и растворимости кислорода в воде, а также использование принципа квазистацинар-ности по отношению к некоторым стадиям механизма), математическая модель представлена в виде системы из трех нелинейных дифференциальных уравнений:
dCXldt' = к1( CA - CX - CZ) CB - k2CXCY + 2 k3( CA -cz) ^Z + k'4( Cz) CY + ksCz;
dCyldt ~k^( Ca ^^z) ^B ^^Y k,4( ^A.
- Cx - CZ) CY + k^Cz - k6^y ;
dCz/dt kC^y k'3 (^A ^^z) ^Z k5 C^z, (5)
edxldt = (|a -p-1x b - xy + 2(|£pa —x - |£2pz)z + + (црф о-ф x-$ipq>z)y + z,
e'dyldt = (|a -p-1x - |l £z) b - xy - (|рфо - ф x -
- |^z)y + z - oy2,
dz/dt = xy - 2(|£pa-£x -|£2pz) z - z, (6)
где ffl = kjk^, П = k2lk3; e = k5lk2; 5 = k5; | = k1lk2^> p = k3lk1; £ = k2lk5; ф = k4ik2; o = k2 k6 lk2; e'= k3 k5 lk2.
Для качественного анализа систем нелинейных уравнений необходимо определить стационарное состояние. Для этого приравняем левые части уравнений к нулю и решаем систему алгебраических уравнений (1):
(|a - p-1x £z) b - xy + 2(|£pa-£x - |£2pz) z + (|pqa -
- ф x-$pqz)y + z = 0,
(|ia - p-1x - b - xy - (дофа-ф x-|£pфz)y + z -
- oy2 = 0,
xy - 2(|£pa - £x pz) z - z = 0.
(1)
Можно показать, что координаты неподвижной (особой) точки соответствуют:
х = ацр; у = 0; I = 0.
Матрица Раусса-Гурвица, при помощи которой можно проанализировать устойчивость данного стационарного состояния, запишется в виде:
-p-1b-Ä
-a|p +1
-p-1b -a|p-A +1
a|p
-1-A,
= 0.
(8)
где [д"-], Су= [Н02-], с2 = [д], сА = сЛ(ИСХ), св = скат („с^ СК(исХ) = [0Щ + [01 + [д], I' - время.
В соответствии с существующей процедурой для установления характерных особенностей динамики протекающих процессов был проведен качественный анализ рассматриваемой математической модели и ее численное решение.
2.1. Качественный анализ математической модели
Для решения поставленной задачи вначале система была приведена к следующему безразмерному виду [6]:
Отсюда получаем, что коэффициенты характеристического многочлена
X3 + А1Х + А2Х + А3 = 0 (9)
равны:
А1 = (р-1Ь + а|р + 1); А2 = (р-1Ь + оЬ|Др); А3 = 0. (10)
Тогда корни характеристического уравнения (10) соответствуют: Х1 = 0, Х2 Ф 0 и Х3 Ф 0.
В том случае, когда один из корней уравнения обращается в нуль, имеет место смена качественных структур в окрестности стационарного состояния при переходе через границу области Раусса-Гурвица и реализуется двойная особая точка [7]. Анализ выражения А* = -А2А2 + 4А33 при помощи проце-
Рис. 1. Зависимость изменения концентрации частиц X,, У1 иот времени (а). Фазовые портреты систем в координатах X- У, У1 -X;Х1 -(б) и X;- У1 (в) (а = 3,010~3, Ъ = 1Д10-4)
дуры, рассмотренной в [7], приводит к тому, что ционарное состояние может стать неустойчивым,
^ > \4Аз3
при А* < 0 (или \А2А2\ > \4А3'\) рассматриваемое ста- причем особая точка в этом случае относится к
Рис. 2. Зависимость изменения концентрации частиц X,, У{ и от времени (а). Фазовые портреты систем в координатах X - У, У1 - г; Х1 - ^ (6) и X; - У1 (в) (а = 1,010~2, Ъ = 1Д10-4)
Рис. 3. Зависимость изменения концентрации частиц X; У1 и от времени (а). Фазовые портреты систем в координатах
X- - у, у1 - г; х¡- X (б) и х - у1 - X («) (а=з,010-3, ъ = 1,0-ю_3)
типу "пространственное седло - фокус"; из данного стационарного состояния возможна бифуркация в предельный цикл, т.е. возможна реализация колебательного режима.
2.2. Численное интегрирование математической модели
Основной задачей численного решения математической модели является нахождение области значений параметров, при которых возможно возникновение химических неустойчивостей [8]. Нами при выполнении работы для реализации решений систем обыкновенных дифференциальных уравнений была использована программа Mathcad с подпрограммой Rkadapt, которая внутренне адаптирована к шагу интегрирования [9]. При этом основное внимание уделялось определению концентрационных пределов реагента и катализатора возникновения осцилляций. На рис. 1 приведены результаты расчета в виде траекторий динамической системы (а) и ее фазовых портретов (б и в) для случая, когда а = 3,010 3 и Ь = 1,0-Ю"4 (число точек принято равным 600, а величина шага соответствует 1). Результаты численного эксперимента пока-
зывают, что в указанных условиях реализуются незатухающие колебания, а фазовые портреты имеют вид предельного цикла (условия приведены в подписях под рисунками; Х1 = [О-]; У1= [НО^]; %= [0]).
Данные, представленные на рис. 2 и 3, иллюстрируют характер эволюции динамики процессов, протекающих в исследуемой системе, в зависимости от изменения значений управляющих параметров - концентраций реагента (а) и катализатора (Ь). Эти результаты свидетельствуют о том, что при значениях параметров a и Ь, отличных от вышеуказанных, возможно возникновение затухающих колебаний.
Таким образом, результаты проведенного численного эксперимента позволяют заключить, что в целом математическая модель удовлетворительно описывает процессы окисления 1,4-нафтодиола в колебательном режиме. Вместе с этим следует отметить, что для более точного количественного описания наблюдаемых эффектов, видимо, потребуется усложнить математическую модель, а также учесть численные значения констант скоростей определенных стадий протекающих процессов.
Работа выполнена при финансовой поддержке РФФИ (проект № 06-03-96621).
СПИСОК ЛИТЕРАТУРЫ
1. Кольцова Э.М., Третьяков Ю.Д., Гордеев Л. С., ВертегелА.А. Нелинейная динамика и термодинамика необратимых процессов в химии и химической технологии. М., 2001.
2. МецлерД. Биохимия. Химические реакции в живой клет-
ке. М., 1980. Т. 2. С. 383.
3. Братушко Ю.И. Координационные соединения 3ё пере-
ходных металлов с молекулярным кислородом. Киев,
1987.
4. Сычев А.Я., Исак В. Г. // Усп. хим. 1995. 64. № 12. С. 1183.
5. Магомедбеков У.Г. // Вестн. Моск. ун-та. Сер. 2. Химия.
2001.42. С. 75.
6. ХолодниокМ., Клич А., КубинекМ., Марек М. Методы анализа нелинейных математических моделей. М., 1991. 1. Баутин Н.Н. Поведение динамических систем вблизи границ области устойчивости. М., 1984.
8. МалинецкийГ.Г. Хаос. Структуры. Вычислительный экспе-
римент: Введение в нелинейную динамику. М., 2000.
9. Pojman J.A. Studying Nonlinear Chemical Dynamics with
Numerical Experiments. Department of Chemistry & Biochemistry. University of Southern Mississippi, 1991. P. 339.
Поступила в редакцию 26.10.06
CHEMICAL INSTABILITIES DURING 1,4-NAPHTODIOL OXIDATION IN HOMOGENEOUS MEDIUM. II. THERMODYNAMIC ANALYSIS AND MATHEMATICAL MODELING
U.G. Magomedbekov, Kh.M. Gasanova, U.G. Gasangadzhieva, K.A. Abdulkhamidov, S.Sh. Mutsalova
(Daghestan State University, Department of general and inorganic chemistry)
A possible mechanism of the processes of 1,4-naphtodiol oxidation in oscillatory regime was proposed and a kinetic model was developed. On the base of analysis of Lyapunov thermodynamic function had been shown that reasons of the origin of chemical instabilities were the realization of autocatalytic stages and reversible dynamics bonds. The system of differential equations describing the kinetics of the reaction were analyzed qualitatively and solved numerically. The character of the stationary state and the possibility of bifurcations were examined. The proposed mathematical model satisfactorily describes the processes occurring in the system.
