УДК 544
Ю. В. Ерандаева, Е. С. Воробьев, Ф. И. Воробьева
РАСЧЕТ СКОРОСТИ ХИМИЧЕСКИХ РЕАКЦИЙ В MS EXCEL
Ключевые слова: Скорость реакции, кинетика, Рунге-Кутт, программа.
В работе представлена программа для решения обратной задачи кинетики химической реакции с поиском скорости реакции на основании экспериментальных кинетических кривых и надстройки Excel «Поискрешения».
Keywords: Сomposition-property diagrams, program, coefficients, equation, diagrams.
The paper presents a program to generate a plan-matrix composition-property with subsequent calculation of the coefficients of the model and the charting of natural and normalized coordinates.
Нахождение скорости реакции на основании экспериментальных данных является актуальной задачей, с которой сталкиваются многие исследователи [1]. Использование графических методов давно уже устарели. Применение упрощенных формул скорости и приемов линеаризации данных тоже не всегда дают хорошие результаты. Широкое развитие вычислительной техники позволяет решать обратные задачи химической кинетики и вычислять скорости химических реакций, подбирая различные схемы их протекания.
Использование MS Excel для решения данной задачи требует создания дополнительной подпрограммы функции, которая реализует решение системы дифференциальных уравнений методом Рунге-Кутта. Анализ возможных реализаций данной функции показал, что наилучшие результаты может дать модифицированный метод Рунге-Кутта-Мерсона, который позволяет изменять шаг при решении задачи с обеспечением заданной точности решения [2]. Для получения массива расчетных данных, который может иметь произвольный набор в зависимости от эксперимента временных отсечек, используем функцию массив. Функция имеет четыре формальных входных параметра:
- Вектор времен из набора экспериментальных данных.
- Вектор начальных концентраций компонентов реакций.
- Вектор подбираемых констант скоростей реакций.
- Имя подпрограммы для вычисления значений производных по кинетической схеме реакции.
Для использования представленной функции необходимо в среде VBA записать подпрограмму для вычисления значений производных, которая имеет следующие входные параметры:
- Текущие концентрации компонентов.
- Набор констант, которые подбираются во время решения.
- Массив для хранения выходных результатов.
Пример данной подпрограммы для двух последовательных реакций показан ниже Public Function Sys3(C, K, s)
S(1) = (-K(1) * C(1))
S(2) = (-K(2) * C(2) Л 2 + K(1) * C(1)) s(3) = K(2) * C(2) Л 2 Sys3 = S End Function
Таких подпрограмм в таблице может быть несколько, что даст возможность их быстро выбирать, меняя их имена в формальном параметре функции массиве. Допустимо выполнение нескольких расчетов с одними и теми же данными, но по различным схемам протекания реакции.
Пример реализации решения показан на рис.1.
А Б С D Е F G Н
1 2 Таблица расчетных и экспериментальных данных Начальные концентрации
Время С2экс С1 С2 СЗ С1 = 1
3 0 0 1 0 0 С2= 0
4 5 0 5 0 216436 0.769839 0 013726 сз= 0
£ 10 0 69 0 046844 0 902696 0 05046 Константы реакций
Ё 15 0.73 0 010139 0 898441 0 09142 к1= 0 3
7 20 0,72 0 002194 0 867305 0 1305 к2= 0 01
8 30 0.62 0 000103 0 799876 0 200021 Функция Svs3 I
9 40 0.52 4 82Е-06 0 740666 0.259329
10
11 |Минимум 0 249172
/ ► ♦ < ► 1
/ ♦ <
♦ С2экс С2
О 10 20 30 40
Время
Рис. 1 - Пример реализации поиска скорости реакции
В ячейки А3:В9 вносятся экспериментальные данные (Время и концентрация ключевого компонента), в ячейки G2:H4 начальные концентрации компонентов реакций, в ячейки G6:H7 начальные приближения констант скоростей реакций, в ячейку H8 имя подпрограммы, где записаны правые части системы кинетических уравнений, которую используем в расчете. Для подбора минимизации констант используем стандартную надстройку MS Excel «Поиск решения» с критерием оптимизации в виде суммы квадратов разностей экспериментальных и расчетных значений концентрации ключевого компонента. Для этого в ячейку D11 вносим встроенную функцию вычисления суммы квадратов разностей =СУММКВРАЗН(В3:В9; D3:D9). Остается записать обращение к функции-массиву для реализации решения. Выделяем весь диапазон расчетных концентраций реакций С3:Е9. Активной должна быть ячейка С3, в нее вставляем нашу функцию и выбираем в мастере построения функции нужные входные параметры. В результате получаем следующую запись - =R_KA(A4:A9;C3:E3;H6:H7;H8), где R_KA - имя функции и далее в скобках набор входных параметров. Ввод функции завершаем комбинацией клавиш [Ctrl+Alt+Enter] чем сообщаем системе, что это функция массив, иначе в ответ будет выведено только одно число в ячейке С3.
Для визуализации результатов расчета строим график зависимости концентраций от времени. Теперь можно искать константы с использованием надстройки «Поиск решения». Устанавливаем курсор в ячейку D11 и вызываем надстройку по команде Сервис - Поиск решения. В открывшемся окне задаем основные настройки для поиска решения:
- Установить целевую ячейку - D11, где записана функция суммы квадратов разностей.
- Равной - выбираем «минимальному значению».
- Изменяя - указываем на константы скорости реакции G6:H7.
Переходим в настройку «Параметры» и устанавливаем следующие параметры:
- Разности - квадратичная.
- Оценки - центральные.
- Метод поиска - Ньютона, если возникают проблемы с нахождением решения, можно воспользоваться методом сопряженных градиентов.
- Автоматическое масштабирование - ставим галочку.
- Неотрицательные значения - ставим галочку, если все константы должны быть положительными.
Завершаем настройку кнопкой «Ок», и нажимаем конку «Выполнить». Компьютер начинает поиск констант. Когда поиск завершается и выводится окно с сообщением, нажимаем «Ок». В результате в таблице выводятся найденные константы, вычисляются концентрации по всем компонентам и строится график (рис.2).
Рис. 2 - График экспериментальных и расчетных значений концентрации после подбора констант скоростей реакций
Для построения графиков с более хорошей детализацией данных можно также воспользоваться функцией массивом Я_КЛ (рис.3).
Рис. 3 - Построение графика с равномерным шагом по времени
Готовим таблицу, в столбец А определяем интервал от 0 до нужного значения времени и заполняем его с помощью команды Правка - Заполнить - Прогрессия. В ячейки Б30:Б44 вносим функцию Я_КЛ, как было описано выше и завершаем ввод комбинацией клавиш [С1х1-Л11;-Еп1ег]. На основании полученных данных строим график.
Литература
1. Абрамов, А.Г. Кинетические параметры реакции каталитической дегидратации фенилэтанолов / Д. Т. Мухамадиев, В. А. Васильев, Э. А. Каралин // Вестник Казан. технол. ун-та. - 2010. - Т.11. -С.30-35.
2. Коробов, В.И. Химическая кинетика: введение с МаШса^Мар1е/МС8 / В. Ф. Очков. - М.: Горячая линия - Телеком, 2009. - 384 с.
© Ю. В. Ерандаева - магистр КГТУ, you.1.40nok@mai1.ru; Е. С. Воробьев - канд. техн. наук, доц. каф. общей химической технологии КГТУ, Vorobiev@kstu.ru; Ф. И. Воробьева - канд. хим. наук, доц. той же кафедры.