Многопиковые колебания в трехмерной кинетической модели реакции окисления оксида углерода на наночастицах палладия Текст научной статьи по специальности «Математика»

УДК 517.928.4:517.929.5

МНОГОПИКОВЫЕ КОЛЕБАНИЯ В ТРЕХМЕРНОЙ КИНЕТИЧЕСКОЙ МОДЕЛИ РЕАКЦИИ ОКИСЛЕНИЯ ОКСИДА УГЛЕРОДА НА НАНОЧАСТИЦАХ ПАЛЛАДИЯ*)

Е, А. Лашина, Н, А. Чумакова, Г, А. Чумаков

Настоящее исследование является продолжением работы [1], в которой изучается двумерная модель гетерогенной каталитической реакции окисления СО и найдены параметры модели, при которых в одпопараметрическом семействе двумерных динамических подсистем существует гистерезис стационарных состояний и на одной из ветвей узлов и фокусов происходит потеря устойчивости стационарного состояния.

Здесь мы предприняли попытку дальнейшего развития сценария возникновения многопиковых и хаотических автоколебаний, предложенного в [2], и рассмотрели в качестве примера трехмерную динамическую систему, моделирующую кинетические автоколебания в реакции окисления СО. Основной является идея блочного построения сложной динамики системы с двумя быстрыми и одной медленной переменными на основе анализа фазовых портретов однопараметрического семейства подсистем меньшей размерности. Целью нашей работы является получение достаточных условий существования сложных автоколебания в динамической системе трех нелинейных автономных обыкновенных дифференциальных уравнений.

Работа выполнена при поддержке Сибирского отделения РАН (Междисциплинарный проект 80).

@ 2013 Лашина Е. А., Чумакова Н. А., Чумаков Г. А.

1. Бифуркации периодических решений

Рассмотрим однопараметрическое семейство динамических систем с параметром £0 ^ = (0,1):

х = кх( 1 - х - у) - к— (у)х - к3(у)ху£0 = Р(х,у,£0), ^ У = к2(у)(1 - х - у)2£0 - к3(у)ху£0 = £оЯ(х,у),

где (х,у) £ О = {(х, у) £ Е : х ^ 0, у ^ 0, х + у ^ 1} х и у — безразмерные концентрации СО и кислорода, адсорбированных на поверхности катализатора (см. [1]). Будем предполагать, что

где постоянные Дк, и к^ положительны, параметры 6 £ (0,0.5) и ус £ (6,1 - 6) определяют ширину и расположение интервала значений у, та котором в(у) монотонно убывает от +1 до -1. В численных примерах будем использовать следующие значения параметров [3], если не оговорено другое:

Заметим, что значения параметров к\ и пока не заданы.

к

ция Андронова — Хопфа при варьировании параметра [4-9].

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

О

ствует стационарное состояние (х^, у^) такое, что матрица Якоби, вычисленная в этой точке, имеет собственные значения вида ±гш, ш / 0:

к,-(у) = Дк,- • *(у) + к-, з = -,2,3,

Дк,- = (к, - Щ)/2, к? = (кГ + к,)/2,

,алг

к—1 = 49.82, кг_1= 0Л9, к12 = 105, Щ = 2 • 104, к13 = 1.18 ^О5, кг3 = 11.8, ус = 0.37, 6 = 0Л.

(3)

Р(хь,уь, £о) = 0, Я{хъ,,у^ = 0,

(4)

р' (xh ,Vh ,b) + bQ'y( xh ,Vh)= о, (5)

{PLQ'y - РУQ'x) U=*h,y=yh >o. (6)

Уравнение главной изоклины P(x, y, £о) = 0 можно записать в виде

x = k1{\ — y)/[h + k— (y) + k3(y)y£0],

а второй главной изоклины Q(x, у) = 0 — как x = (1 — y)f(y), где f(y) определяется формулой

f, ч v/4(l-y) + A(y)- уЩу) к3(у)

= -Ттгт' А(у> = ТТ\У• 7)

v/4(l-y) + A(y) + vM(y) Mf)

yh

+ + = /Ы' (8)

xh

xh - yh f yh . След матрицы Якоби, вычисленной в (xh, yh)> можно записать в

виде

р' (xh,yh,£o) +&Q'( xh ,yh)

= -кг - к-\(уь) - ь (кг{ун)ун + д{ун) + ~~j~Xhyh

где

g(yh) = 2k2(yh)(l -xh- yh) - dk2jyh^ (i-Xh- yhf + k3(yh)xh. (10)

dy

Отсюда видно, что для yh вне интервала (yc — ö, yc + 5) след матрицы Якоби отрицателен. Следовательно, условие (5) эквивалентно равенству

h + k-i{yh) + k3(yh)yh£о + д(ун)Ь = - dkz\Vh^ XhVh£,о, (11)

dy

причем yh G (yc — ö, у^ + ö).

Таким образом, из равенств (8) и (11) следует, что при выполнении условий (3), параметры к\ = к° и £о = Со > ПРИ которых существует стационарное состояние (хо, уо) с чисто мнимыми корнями, удовлетворяют системе уравнений

¿к = _к-!(ун)_

(<?Ы + ^¡Г1^) (ДЫ - 1) - МЫг/й'

(12)

/•1 ./'(.'// ) ^.'/(/// ) ~ )

где уо € (ус — д, Ус + д) и значения хо, /(уо) и д(уо) определены формулами (9), (7) и (10) соответственно. Кроме того, необходимо, чтобы выполнялось неравенство (6) и, более того, к° >0и0<£° < 1-

Рис. 1. Бифуркационные кривые, определяющие параметры (£о>&1)> при которых система (1) имеет негру бое стационарное состояние типа седло-узел (кривая 1). неустойчивый сложный фокус (кривая 2) или полуустойчивый предельный цикл (кривая 3).

Численно показано, что, используя формулы (12) и (6), значение к\ = к° можно выбрать так, что при некотором £о = Со система (1) име-

ет негрубое стационарное состояние (хн, ун), являющееся сложным фокусом [8]. Например, на рис. 1 показана кривая 2 па плоскости (£0, к) такая, что для каждой пары значений к\ и на этой кривой в системе существует неустойчивый сложный фокус. Устойчивость фокуса определяется с помощью первой ляпуновской величины I (см. [6,9]). В рассматриваемом случае I > 0, и при малом изменении параметров в системе (1) происходит бифуркация Андронова Хопфа и зарождается грубый неустойчивый предельный цикл (рис. 2).

0.6 0.5 0.4

0.2 0.1 0

0 0.2 0.4 0.6 0.8

Рис. 2. Зависимость y-коордпнаты стационарных состояний от параметра CD и KB\{K} — устойчивые стационарные состояния, BC\({B} и {C}) — седлa, AK — неустойчивые фокусы и узлы. Кривые 1 (и 2) — максимальные и минимальные значения переменной y на устойчивых (и неустойчивых) периодических решениях системы (1) при = 4000.

2. Максимальные семейства периодических решений

На рис. 3 слева кривая AKBCD является максимальным семейством стационарных состояний одноиараметрического семейства систем (1) с параметром £0- Проекция этой кривой па плоскость (£o,y)

BC

ционарным состояниям системы, имеющим тип седло-узел. Точкой К обозначено негрубое стационарное состояние (жо, уо) пр и Со = 0.370108, где жк = 0.301507 и ук = 0.408646. Кривые СВ\{С}ш КВ\{К} состоят из устойчивых стационарных состояний, кривая АК — из неустойчивых стационарных состояний, имеющих тип узел или фокус. Стационарные состояния, лежащие на кривой ВС\{В}, являются седловыми. Для определения типа и устойчивости стационарного состояния рассмотрены собственные значения матрицы Якоби системы, линеаризованной в его окрестности.

При Со = 0.370108 в системе (1) происходит бифуркация Андронова — Хопфа и зарождается неустойчивое периодическое решение периода Т « 3.8 • 10-4.

Продолжая грубый предельный цикл по параметру Со па максимально возможный интервал Со < Со < Со> получим максимальное семейство грубых предельных циклов одпопараметрического семейства динамических систем (1). При Со = Со и & = Со в системе (1) происходит бифуркация периодического решения, так что Со = Со или Со = Со- Построение максимальных семейств периодических решений будем проводить с помощью метода продолжения по параметру, основанного на определении неподвижной точки функции последования (см. [8,10]).

Численный эксперимент показал, что (1) имеет максимальное семейство грубых неустойчивых предельных циклов 7й (Со) для Со € (0.348538,0.370108). Кроме того, существует максимальное семейство грубых устойчивых предельных циклов 7® (Со) для Со € (0.348538,1) С.

слияния устойчивого и неустойчивого предельного циклов. На рис. 2

у

С

Таким образом, при рассматриваемых значениях параметров изме-

С

С.

бально устойчивое стационарное состояние — узел (рис. 4 слева). При С.

рых — устойчивые узел и фокус, и одно — седло. На рис. 5 слева

Рис. 3. Слева. Максимальные семейства ЛБСО стационарных состояний и йа устойчивых предельных циклов. Справа. Семейства устойчивых (сплошная линия) и неустойчивых (пунктирная линия) предельных циклов одиопараметрического семейства системы (1).

показан фазовый портрет при £0 = 0.365. В этом случае существуют три стационарных состояния: устойчивые узел и фокус и неустойчивое седло. Устойчивый фокус окружают два предельных цикла: кривые 1 и 2 это устойчивая и неустойчивая периодические траектории. Наконец, при £о = 0.94 система (1) имеет одно стационарное состояние (неустойчивый фокус) и устойчивый предельный цикл.

Рис. 4. Фазовые портреты системы (1) для = 4000 и £о = 0.001 (слева) и 5о = 0.2 (справа).

Рис. 5. Фазовые портреты системы (1) для = 4000 и £о = 0.365 (слева) и £о = 0.94 (справа).

Таким образом, на плоскости (£о, можно выделить область параметров и к, при которых система имеет два грубых предельных цикла, устойчивый и неустойчивый. На параметрическом портрете (см. рис. 1) эта область ограничена кривыми 2 и 3.

3. Динамика трехмерной системы

Рассмотрим систему уравнений

х = к1{\ - х - у) - к— (у)х - к3{у)ху£,

у = к2(у)(1 - х - у)2£ - кг{у)ху£, V = к±(\ - V) - к_4V, £ = ф(1 - £) - о.уС),

содержащую малый параметр е > 0. Как показано в [1], для любых начальных данных (х(0), у(0), £(0)) = (щ, уо^о,£о) существует момент времени Ь = т такой, что для всех Ь > т траекторию системы (1) можно приблизить решением трехмерной динамической системы

х = кг(1 - х - у) - к— {у)х - кг(у)ху£, у = к2(у)(1 - х - у)2£ - к3(у)ху£, (14)

£ = ф^(1 - О - ау0,

где ^Ь) « vst и Vst = кА/(кА + к_4).

Многопиковые автоколебания. Предположим, что значения параметров системы (1) выбраны так, что соответствующее однопара-метрическое семейство двумерных динамических систем относительно переменных x и y с параметром 0 < £ < 1 имеет гистерезис ABCD стационарных состояний и семейства (ys(£)} и {7"(£)} устойчивых и неустойчивых периодических решений. Отметим, что в фазовом пространстве G х множество ABCD = AK U KB U BC U CD является гладкой кривой и существуют две поверхности Ss = {7я (£)} и Su = {7й (£) }, образованные периодическими решениями. Рассмотрим множество

г = {(x, y, О е G х vat(l — £) - ay£ = 0}

точек фазового пространства G х в которых £ = 0. Можно показать, что существуют такие значения параметра а, что если (x,y,£) е CD, то £ > 0, и если (x,y,£) е AB или (x,y,£) е Ss, то £ < 0. Предполо-

ABCD

O = (x*,y*,£*), O е BC. Тогда система (14) имеет единственное стационарное состояние (x*, y* ,£*), неустойчивое при малых значениях е. (При рассматриваемом наборе значений параметров имеем е = 0.415.) Кроме того, предположим, что поверхность Г не пересекает поверхно-стн Ss = {7s(£)}и Su = {7u(£)}•

Тогда система (14) имеет многопиковый колебательный режим. Действительно, пусть начальная точка (щ, yo,£o) лежит в окрестности KB CD

яний системы (1). Согласно классическим результатам А. П. Тихонова [11] о решениях систем ОДУ с малым параметром при малых значе-

е > x x y y £ £

ние (x, y, £) = (x(t), y(t), £(t)) системы (14) на некотором интервале времени [0, T] лежит в окрестности кривой KB (или CD). В силу выбора а компонента £(t) решения уменьшается (или увеличивается) с увеличением времени, и при t = T точка (x(T), y(T), £(T)) приходит в окрест-BC

точка B называется точкой срыва. Существует момент времени t = T\, Ti > T, такой, что (x(Ti), y(Ti), £(Ti)) попадает в окрестность кривой CD и |£(Ti) — £(T)| = O(e). Отметим, что C также является точкой

срыва, и если (х(Т),у(Т), £(Т)) лежит в окрестности этой точки, то существует момент времени Ь = Т > Т такой, что {х{Т2),у{Т2),£{Т2)) находится в окрестности поверхности и |£(Т2) - £{Т) 1 = О(е).

Предположим, что точка (щ, уо,£о) лежит в окрестности поверхности 5я. Тогда при малых значениях е и начальньж условиях х(0) = хд, у(0) = уо, £(0) = & решение {х,у,£) = {х{г), у{г), £{г)) системы (14) на некотором интервале времени [0, Тя] проходит в окрестности поверхности 5При Ь = Тя точка (х(Тя),у(Тя),£(Тя)) попадает в окрестность полуустойчивого предельного цикла. При дальнейшем увеличении Ь существует Ь = Т{ > Тя такое, что, пройдя в окрестности дуги КБ, точка (х(Т{) ,у{Т(), £(ТЯ)) попадает в окрестность кривой СБ.

е

стема (14) генерирует многопиковые автоколебания ^е), причем колебания соответствующего периодического решения являются релаксационными. Колебания компонент х(Ь) и у(€) имеют многопиковый характер: чередуются колебания, имеющие большую амплитуду, и колебания с малой амплитудой.

Исходная система (13) является малым возмущением системы (14). Поскольку многопиковые колебания системы (14) и стационарное состояние vs = к4/(к4 + к_4) уравнения V = к4(1 - V) - к_4V являются

е

же имеет многопиковые автоколебания такие, что их проекция в пространство О х Ц близка к многопиковым колебаниям ^е) системы (14), а компонента ^Ь) соответствующего колебательного режима лежит в окрестности точки ^ = vs.

к к_

бания в полной системе (1) из [1] и vst = 0.107143 (рис. 6).

Основной результат этой работы состоит в том, что с использованием принципа блочного построение сложной динамики (см. [2]) найдены параметры модели, при которых в трехмерной динамической системе с быстрыми и медленными переменными существует класс решений, позволяющий моделировать многопиковые автоколебания, полученные в химическом эксперименте [3,13].

Рис. 6. Слева. Многопиковые автоколебания, предельный цикл. Справа. Зависимости компонент периодического решения системы (1) в [1] от времени. Значения параметров определены в (3), a = 5.37, е = 0.415, fc4 = 1500, fc-4 = 5000.

ЛИТЕРАТУРА

1. Лашина, Е. А., Чумакова Н. А., Чумаков Г: А. О двумерной кинетической модели реакции окисления оксида углерода на наночастицах палладия // Мат. заметки ЯГУ. 2013. Т. 20, вып. 2. С. 271-285.

2. Чумаков Г: А., Слинько М. Г. Кинетическая турбулентность (хаос) скорости реакции окисления водорода на металлических катализаторах // Докл. АН СССР. 1982. Т. 266, № 5. С. 1194-1197.

3. Lasbina Е. A., Slavinskava Е. M., Cbumakova N. A., Stonkus О. A., Gulvaev R. V., Stadnicbenko А. L, Cbumakov G. A., Boronin A. L, Demidenko G. V. Self-oscillations in СО oxidation on РсЮ/АЬОз catalyst // Chem. Eng. Sci. 2012. V. 83. P. 149-158.

4. Андронов A. A., Леонтовнч E. A., Гордон И. И., Майер А. Г. Качественная теория динамических систем второго порядка. М.: Наука, 1966.

5. Kuznetsov Yu. A. Eléments of applied bifurcation theory. New York: Springer-Verl., 2004.

6. Хэссард В., Казарннов H., Вэн И. Теория и приложения бифуркации рождения цикла. М.: Мир, 1985.

7. ¡Пильняков Л. П., Шильников А. Л., Тураев Д. В., Чуа Л. Методы качественной теории в нелинейной динамике. Ч. 2. М.; Ижевск: НИЦ Регулярная и хаотическая динамика, Ин-т. компьютерных исследований, 2009.

8. Андронов А. А., Леонтович Е. А., Гордон И. И., Майер А. Г. Теория бифуркаций динамических систем на плоскости. М.: Наука, 1967.

9. Баутин H. Н., Леонтович Е. А. Методы и приемы качественного исследования динамических систем на плоскости. М.: Наука, 1976.

10. Чумаков Г: А., Лашина Е. А., Чумакова, Н. А. Максимальные семейства периодических решений кинетической модели гетерогенной каталитической реакции // Вестн. НГУ. Сер. Математика, Механика, Информатика. 2005. Т. 5, вып. 4. С. 3-20.

11. Тихонов А. Н. О зависимости решений дифференциальных уравнений от малого параметра. Мат. сб. 1948. Т. 22, № 2. С. 193-204.

12. Мищенко Е. Ф., Колесов Ю. С'., Колесов А. Ю., Розов H. X. Периодические решения и бифуркационные процессы в сингулярно-возмущенных системах. М.: Физматлит, 1995.

13. Slavinskava Е. M., Stonkus О. A., Gulvaev R. V., Ivanova A. S., Zaikovskïi V. L, Kuznetsov P. A., Boronin A. I Structural and chemical states of palladium in Pd/AhOe catalysts under self-sustained oscillations in reaction of CO oxidation // Appl. Catalysis A: General. 2011. V. 401. P. 83-97.

г. Новосибирск

15 октября 2013 г.