CC BY
40
УДК 519.876.5
МЕТОД ПОИСКА ИНТЕРВАЛОВ НЕОПРЕДЕЛЕННОСТИ КИНЕТИЧЕСКИХ КОНСТАНТ ХИМИЧЕСКОЙ РЕАКЦИИ
© Э. Р. Ахматсафина*, И. М. Губайдуллин, С. И. Спивак
Башкирский государственный университет Россия, Республика Башкортостан, 450074 г. Уфа, ул. Заки Валиди, 32.
Тел./факс: +7 (34 7) 273 6 7 78.
Институт нефтехимии и катализа РАН Россия, Республика Башкортостан, 450075 г. Уфа, пр. Октября, 141.
Тел./факс: +7 (347) 231 35 44.
E-mail: elzaahmat@mail. ru
В работе решается задача поиска интервалов неопределенности для кинетических констант химической реакции циклоалюминирования олефинов с использованием параллельных вычислительных систем. Создан, описан и обоснован метод поиска интервалов. Приведены результаты расчетов.
Ключевые слова: интервалы неопределенности, кинетические константы, циклоалюмини-рование олефинов, параллельные вычисления.
Введение
Традиционно алюминийорганические соединения широко применяются для осуществления таких реакций как термическое гидроалюминиро-вание олефинов, диенов и ацетиленов, восстановление карбонильных соединений, сложных эфиров и нитрилов. В ряду этих исследований к числу наиболее крупных достижений отечественной и мировой науки последних лет в области химии алюми-нийорганических соединений следует отнести открытие реакции каталитического циклоалюмини-рования непредельных соединений (реакция Дже-милева) [1] с помощью простейших триалкил- и алкилгалогеналанов с участием Ть и Zr-содержащих комплексных катализаторов.
Изучение механизмов сложных реакций металлокомплексного катализа, к которым относятся реакции циклоалюминирования, представляет огромный практический интерес. При этом возникают как физико-химические, так и математические проблемы. Физико-химические проблемы сводятся к тому, что очень трудно и подчас практически не представляется возможным получить данные о промежуточных соединениях реакции с использованием существующего оборудования, из чего следует математическая неоднозначность решения обратных задач определения кинетических констант. Математические проблемы предполагают значительный объем вычислений при поиске этой многомерной области для констант. Для повышения эффективности решения данной задачи целесообразно использовать высокопроизводительные параллельные вычислительные системы.
Целью исследования является создание метода поиска интервалов неопределенности кинетических констант с использованием многопроцессорных технологий, поиск интервалов неопределенности констант скоростей для реакции циклоалюминиро-вания олефинов, физико-химическая интерпретация полученных результатов.
Среди интервальных методов оценивания параметров физико-химических моделей серьезное развитие получил метод [2], сочетающий принцип выравнивания по П. Л. Чебышеву с методологией Л. В. Канторовича и успешно примененный в задачах обработки данных кинетических измерений.
Первые работы по приложению нестатистических интервальных методов в аналитической химии опубликовали Г. И. Рузайкин [3] и В. М. Белов [4]. Среди работ, имеющих отношение к нестатистическому оцениванию параметров эмпирических зависимостей, выделим работы А. П. Вощинина [5].
Классическая постановка задач математической обработки эксперимента состоит в поиске таких значений параметров математического описания, при которых минимизируется какой-либо критерий соответствия расчета измерениям [6]. Определение точки, соответствующей глобальному минимуму критерия - достаточно сложная вычислительная задача. Возникает вопрос: чем константы, соответствующие минимуму критерия, предпочтительнее констант, соответствующих некоторой окрестности критерия, если при этом сохраняется достаточная точность описания эксперимента, сопоставимая с его погрешностью? Ответ на этот вопрос об определении области в пространстве параметров, каждая точка внутри которой достаточно хорошо описывает измерения, был поставлен Л. В. Канторовичем [7] как решение системы неравенств:
I л
х..
(1)
1=1 і=1
Величины є, в системе неравенств (1) есть характеристики предельно допустимой погрешности эксперимента, информация о величине которой, как правило, присутствует. И тогда выполнение условий (1) означает, что модель описывает измерения в пределах, обусловленных величиной максимально допустимой погрешности измерений, что совершенно естественно.
* автор, ответственный за переписку
Рассмотрение многомерной области для констант скоростей начнем с одномерного и двумерного случаев. Для одномерного случая область представляет собой интервал неопределенности. Определим интервал неопределенности как отрезок, вариация константы внутри которого сохраняет совместность системы неравенств (1). Для определения интервала по К-й константе необходимо решить две задачи математического программирования: найти тт К и тах К при выполнении ограничений (1).
В начальный момент времени нам известна начальная точка, принадлежащая искомому интервалу - точка, полученная в решении обратной задачи. Необходимо найти интервал по каждой из констант, беря начальную точку за исходную в поиске. Предлагается следующий алгоритм. Пусть ищется интервал неопределенности по К константе. Чтобы найти границы искомого интервала, нужно зафиксировать все константы, кроме К . Для этой константы выберем произвольный отрезок, (например, [0; 2К]), содержащий исходную точку, покроем его равномерной сеткой с определенным малым шагом (например, К/1000). В узлах сетки получим различные значения К константы, в каждом из узлов проверим выполнение системы неравенств (1), Если неравенства верны, то данный узел принадлежит искомому интервалу. Определяя наименьшее и наибольшее значения принадлежащих интервалу узлов, получаем значения границ интервала неопределенности по К константе. Аналогичным образом поступаем для поиска границ двумерной области [8].
Оценим объем вычислений. Для определенности будем рассматривать двумерные области. Количество вычислений, которое необходимо произвести для нахождения всех областей неопределенности для всех наборов кинетических констант, можно вычислить по следующей комбинаторной формуле: р = Ст ■ д , где п - общее количество
констант, т - количество констант, по набору которых получаем срез области, д - количество узлов сетки, которой покрывается п-мерный куб. Например, количество констант п=10, сетка двумерная 1000 на 1000 узлов, тогда потребуется р = С120 ■ 1000000 = 45000000 вычислений [8]. Нами было подсчитано, что в среднем на одну прямую задачу (одно вычисление) затрачивается примерно 0.01 с. Тогда для поиска двумерной области придется затратить приблизительно 450000 секунд. Это составляет 125 часов или 5.2 дня. При увеличении количества узлов сетки время расчетов, естественно, увеличивается. Кроме того, существует реакции с большим количеством констант, также нельзя исключать возможность увеличения размерности области. Все это приведет к тому, что время расчетов существенно увеличится. В связи с этим для ускорения расчетов было решено использовать многопроцессорные вычисления.
Обсуждение результатов Объектом исследования является разработанная на основании экспериментально обоснованной схемы (рис. 1) процесса кинетическая модель реакции циклоалюминирования олефинов триэтилалю-минием в присутствии Ср^гС12 [9].
Рис. 1. Схема протекания реакции циклоалюминирования олефинов.
Для этой реакции схема химических превращений имеет вид (2), соответствующая ей система уравнений - (3)
лг + Ъ А2 5 Ъ Аз + Ъ А4,
А1 + Аз А4 + А12,
Аз + А4 А1 + Ъ Аб + Ъ А13,
А3 ~>А5 + A13, (2)
А1 + А5 А13 + А 8>
А5 + А9 ~^А1,
А1 + А10 А3 + А11
А1 + Ъ Аб А4 + Ъ А7,
А7 ->Аз + А5, где А1 = А1(С2Н5)з, А2 = Ср 2^тС 12,
А3 = СР2^г(С2Н5)С1 ■A^(C2H5)з, А4 = ClAl(C2H5)2,
А5 = Ср21гСН2СН2А!(С1)(С2Н5)2,
Аб = (С1)Ср21гСН2СН21гСр2(С1) -2[СШ(СН5)21,
А7 = (С1)Ср21тСН:СН21гСр2(С1) ^[МСН)],
А8 = Ср21г(С1)СН2СН[А!(С2Н5)2]2, Ад = СН£НК, А10 = Ср21г(С1)СН2СНКСН2СН2А1(С2Н5)2, А11 = (СНМЧСН^зРНК, А12 = Ср21г(С2Н5)2 А1(СН5), А13 = С2Нб, Ср= С5Н5, Я — олефин.
ах1
= -Ш1 -Ш2 + Ш3 -Ш5 -Ш7 -W8;
dt
dX
2 1
- = - 2
dt
dX 3 1
------ = - W1 - W2 - W3 - W4 + W7 - W9;
dt 2
dX 4 1
------- = - Wj + W2 - W3 + W8;
dt
dX
2
5 1
-------- = - W3 + W4 + W ■
dt 2 3 4 5
dX-6
dt
dX
= W4 - W5 - W6 + W9;
7 1 1
W3 - W8; dt 2 3 2 8 dX8 = 1W8 - W9; dt 2 8 9
(3)
dX 9 dt
dX 10 dt
dX jj dt
dX 12 dt
dX 13 dt
= W,
- = -W,
= W,
=W
где W, = klXlX1 - kwxfX4^, W2 = k2X 2X 3, W3 = k3 X3 X4,
W4 = k4X3, W5 = k5X!X5, W = k6X6X9, W7 = k7X1X10,
W8 = kgX^X^2, W9 = k9X 7 - скорости стадий;
X1,...,X13 - концентрации веществ A1,..., A13.
Были найдены следующие приведенные значения кинетических констант реакции для соотношения Cp2ZrC12:A1Et3=1:10 (табл. 1)
В нашем случае измеряемыми компонентами являются Х2, Х3, Х5, Х7, Х8, Х12, поэтому система неравенств (1) приобретает вид:
1 N I 1 N. .
VIIХг--Х'1 ^МIIХ7- -Х-^
М -=1
N
j=1
N
L7j
1 N . 1 N
—у \X;i - Xр. < e3, — Ях8э
N 3j 3j| 35 N 8
1 N 1 N
—У Xэ. - x5р <e5, — У x1
Njg^l 5j 5jl 5’ N ' 1
э - X р < 8
12j 12j — 12
(4)
j=1
где 8t задаются из условий анализа экспериментальных данных, но не более 15% от экспериментальных значений.
С использованием созданного программного комплекса найдены области неопределенности кинетических констант реакции циклоалюминирова-ния олефинов (табл. 2).
Можно заметить, что начальная точка, вокруг которой находился интервал, не является геометрическим центром этого интервала. На рис. 2 приведены наиболее типичные фазовые плоскости по некоторым константам при разных значениях s.
Программная реализация
Для решения поставленной задачи использовалась разработанная в Институте нефтехимии и катализа РАН информационно-аналитическая система обратных задач химической кинетики (ИАС ОЗХК) [10].
Вся информация, относящая к эксперименту, сведена в базу данных, в которой можно выделить: входную информацию (начальные (исходные) данные): количество веществ, стадий, констант, точность численного решения прямой кинетической задачи методом Кутта-Мерсона [11-12], начальные данные концентраций химических веществ, участвующих в реакции; выходную информацию: рассчитанные константы скоростей стадий, энергия активации, левые и правые границы интервалов для пар констант. Для хранения экспериментальных данных кинетических экспериментов в лабораториях ИНК РАН используется формат DAT.
Прямая кинетическая задача решается методом Кутта-Мерсона. Алгоритм метода поиска интервалов неопределенности кинетических констант был описан выше. Распараллеливание вычислений осуществлялось средствами библиотеки MPI (Message Passing Interface) [13-14], с использованием языка программирования Fortran. Для проведения расчетов использовался суперкомпьютер Уфимского государственного авиационно-технического
университета, имеющий классическую для подобных систем массивно-параллельную архитектуру.
Таблица 1
Найденные значения кинетических констант
і 1 2 3 4 5 6 7 8 9 10
кі 165.556 1.121 26.478 4.061 4.843 90.3435 373.5228 26.792 1.986 0.179
Рис. 2. Фазовые плоскости по некоторым константам при разных е.
Таблица 2
Интервалы неопределенности кинетических констант к
і Значение Область при є = 5% Область при є = 10%
1 165.555 [164.451; 169.97] [154.518; 176.592]
2 1.121 [1.0909; 1.2029] [1.068; 1.225]
3 26.478 [25.595; 33.715] [25.595; 36.363]
4 4.061 [3.953; 4.359] [3.845; 4.576]
5 4.843 [4.617; 5.134] [4.553; 5.263]
6 90.34 [85.43; 100.06] [82.3; 101.4]
7 373.523 [370.56; 381.22] [365.29; 382.44]
8 26.792 [17.326; 29.65] [10.002; 70.554]
9 1.986 [1.655; 2.105] [1.43; 2.25]
10 0.179 [0.166; 0.185] [0.162; 0.185]
Параллельная обработка данных основана на следующем принципе: один процессор системы выделяется для сбора данных со всех остальных процессоров и сохранения их в базе данных, также этот процессор (и соответствующий ему процесс), назовем его нулевым, производит считывание всей необходимой входной информации из базы данных и первоначальную рассылку данных остальным процессам для обработки. В результате алгоритм представляет собой следующую последовательность шагов:
- нулевой процесс считывает из базы все необходимые исходные данные и распределяет их по остальным рабочим процессам. Распределение данных между процессорами производится при помощи функций МР1_8епё/МР1_Кесу, МР1_Бсаз^ МР1_8сайег/МР1_8сайегу библиотеки МР1;
- каждый рабочий процесс рассчитывает область неопределенности для своей пары кон-
стант (или одной константы в случае поиска интервала) по описанной выше схеме. В результате получается массив граничных значений интервала или области неопределенности;
- каждый рабочий процесс по мере окончания счета производит передачу найденного массива нулевому процессу при помощи функций МР1_8епё/МР1_Яесу, МР1_ОаЛег/МР1_ОаЛегу библиотеки МР1;
- нулевой процесс, в свою очередь, производит запись полученных данных в базу.
Заключение Определены кинетические параметры реакции, наилучшим образом описывающие экспериментальные данные, полученные при различных начальных концентрациях компонентов и температурах протекания реакции. Создан метод и описан алгоритм распараллеливания в обратных задачах химической кинетики. Разработан комплекс прикладных программ для расчета интервалов и областей неопределенности кинетических параметров с использованием многопроцессорных технологий. Проведен вычислительный эксперимент для анализа данных реакции циклоалюминирования олефи-нов, построены фазовые плоскости для констант исследуемой реакции.
ЛИТЕРАТУРА
1. Ибрагимов А. Г. // Химия в интересах устойчивого развития. 2008. Т. 16. С. 715-719.
2. Спивак С. И. Оценка значимости влияния измерения на кинетическую модель химической реакции // Математические проблемы химии. Ч. 2. Новосибирск: ВЦ СОАН СССР, 1973. С. 3-9.
3. Рузайкин Г. И. Метод построения по экспериментальным данным сглаженной кривой с контролем за выбором для нее погрешности // Математические методы и ЭВМ в аналитической химии. М.: Наука, 1986. С. 54.
4. Белов В. М. // Изв. вузов. Химия и хим. Технология. 1993. Т. 36. №5. С. 46-50.
5. Вощинин А. П. // Заводская лаборатория. 1990. Т. 56. №7.
С. 76-81.
6. Царева З. М., Орлова Е. А. Теоретические основы химтех-нологии. Киев: Высшая школа, 1986. 271 с.
7. Канторович Л. В. // Сиб. мат. журн. 1962. Т. 3. № 5. С. 701-709.
8. Ахматсафина Э. Р., Губайдуллин И. М., Спивак С. И. Ме-
тод поиска интервалов неопределенности кинетических констант химической реакции с использованием параллельных вычислительных систем. // Параллельные вычислительные технологии (ПаВТ’2010): Труды международ-
ной научной конференции (Уфа, 29 марта - 2 апреля 2010 г.) Челябинск: изд-во ЮУрГУ, 2010. С. 387-394.
9. Балаев А. В., Парфенова Л. В., Русаков С. В., Губайдуллин И. М., Спивак С. И., Халилов Л. М., Джемилев У. М. // Докл. АН. 2001. Т. 381. С. 364-367.
10. Губайдуллин И. М., Спивак С. И. // Системы управления и информационные технологии. 2008. №1.1(31). С. 150.
11. Самарский А. А., Гулин А. В. Численные методы: учеб. пособие для вузов. М.: Наука, 1989. 432 с.
12. Арушан О. Б., Залеткин С. В. Численное решение обыкновенных дифференциальных уравнений на Фортране. М.: изд-во МГУ, 1990. 336 с.
13. Воеводин В. В., Воеводин Вл. В. Параллельные вычисления. СПб: БХВ-Петербург, 2002. 608 с.
14. Гергель В. П. Теория и практика параллельных вычислений. М: Интуит, 2007. 423 с.
Поступила в редакцию 01.03.2010 г.
