Параллельная реализация асинхронного клеточного автомата, моделирующего реакцию окисления Co на палладии Текст научной статьи по специальности «Химические технологии»

2011 Дискретные модели реальных процессов №1(11)

ДИСКРЕТНЫЕ МОДЕЛИ РЕАЛЬНЫХ ПРОЦЕССОВ

УДК 004.942, 54-44

ПАРАЛЛЕЛЬНАЯ РЕАЛИЗАЦИЯ АСИНХРОННОГО КЛЕТОЧНОГО АВТОМАТА, МОДЕЛИРУЮЩЕГО РЕАКЦИЮ ОКИСЛЕНИЯ CO НА ПАЛЛАДИИ

В. П. Маркова, А. Е. Шарифулина

Институт вычислительной математики и математической геофизики СО РАН,

г. Новосибирск, Россия

E-mail: [email protected], [email protected]

Для моделирования каталитической реакции окисления CO на поверхности металлов платиновой группы используется асинхронный клеточный автомат (КА) с вероятностными правилами переходов. При КА-моделировании кинетики поверхностных реакций КА требуется использовать массивы больших размеров и проводить вычисления в течение длительного времени. Следовательно, моделирование таких процессов в реальном времени может проводиться только с помощью распараллеливания задач на суперкомпьютере. Приводится параллельная реализация блочно-синхронного КА, аппроксимирующего кинетический КА.

Ключевые слова: каталитическая реакция окисления, кинетический клеточный автомат, синхронный режим, асинхронный режим, блочно-синхронный режим, параллельная реализация блочно-синхронного клеточного автомата.

Введение

Каталитическое окисление CO на платиновых металлах является классической модельной реакцией гетерогенного катализа, которая помимо фундаментального интереса имеет важное прикладное значение в связи с экологической проблемой очистки выбросов отходящих газов от примесей окиси углерода.

Экспериментальные исследования процессов на поверхности катализатора требуют значительных материальных затрат, поэтому особое значение приобретает компьютерное моделирование. Все процессы, протекающие на поверхности катализатора — адсорбция, десорбция, диффузия частиц по поверхности и взаимодействие частиц, — происходят асинхронно. Традиционные методы моделирования, основанные на решении дифференциальных и алгебраических уравнений, позволяют описать интегральные зависимости элементарных физико-химических процессов, но не учитывают возможности изменения атомарной структуры поверхности под воздействием реакционной среды [1].

Наиболее эффективным для описания пространственно-временной динамики процессов на поверхности катализатора, структура которой может изменяться в ходе реакции, является асинхронный КА [2]. Асинхронный клеточный автомат с вероятностными правилами переходов, использующийся для моделирования кинетических процессов на поверхности катализатора, называется кинетическим КА [2, 3], который известен ещё как кинетический метод Монте-Карло. Для исследования каталитических

процессов необходимо использовать большие клеточные массивы и проводить вычисления в течение длительного времени, поэтому при КА-моделировании таких задач требуется использовать эффективные алгоритмы распараллеливания. Но распараллеливание асинхронных КА в отличие от синхронных сопряжено с определёнными трудностями, так как межпроцессорный обмен должен выполняться всякий раз, когда изменяется состояние хотя бы одной граничной клетки. Аппроксимация исходного асинхронного КА блочно-синхронным позволяет достичь более высокой эффективности распараллеливания [3].

Целью работы является параллельная реализация и исследование эволюции кинетического клеточного автомата, моделирующего каталитическую реакцию окисления монооксида углерода (СО) на поверхности Рёцо, и сравнение асинхронного и блочносинхронного режимов работы КА. В п. 1 представлено описание каталитической реакции окисления СО, п. 2 посвящён формальному описанию КА-модели каталитической реакции. В п. 3 описан блочно-синхронный КА, аппроксимирующий кинетический КА, и приведено сравнение результатов моделирования реакции окисления с помощью кинетического и блочно-синхронного КА. В п. 4 рассмотрена параллельная реализация блочно-синхронного КА. В п. 5 приведён гистерезис скорости реакции, полученный в результате КА-моделирования.

1. Реакция окисления СО на поверхности палладия

Каталитическое окисление СО на платиновых металлах сопровождается колебаниями скорости образования СО2 и концентраций адсорбированных веществ. Эти колебания обусловлены сравнительно медленным процессом образования и расходования приповерхностного кислорода, который приводит к изменению каталитических и адсорбционных свойств поверхности. Реакция окисления СО на поверхности Рёцо описывается следующими элементарными стадиями [4]:

1) О2(§ав) + ** —? 20^ — адсорбция и диссоциация кислорода;

2) СОёав + * СО^ — адсорбция и десорбция СО;

3) СОаЙ8+ Оай8 — СОг^ав) + ** — реакция между СОадв и Оаав;

4) Оаав ОвиЬ — образование приповерхностного кислорода;

5) СОаав+ ОвиЬ СОг^ав) + ** — реакция между СОаав и Овиъ;

6) СОёав + ОвиЬ I—У [СОаав*ОвиЬ]—образование комплекса СОаав и ОвиЬ;

7) [СОаав*Овиъ] ^ СОг^ав) + * — реакция в комплексе [СОаав*Овиъ];

8) СОаав + * *+ СОаав — диффузия СОаав по поверхности;

9) СОаав + ОвиЬ *+ [СОаав*Овиъ] — диффузия СОаав с образованием СОаав;

10) [СОаёв*ОвиЬ] + ОвиЬ ^ ОвиЬ+ [СОаёв*ОвиЬ] диффузия комплекса.

Символ «*» обозначает свободный активный центр поверхности катализатора; символ «**»—два соседних свободных активных центра; кг и к-г — константы скорости прямых и обратных элементарных стадий реакции.

Все элементарные стадии, кроме третьей, реализуются с заданной вероятностью рг = кг/ £ к,, где кг — константа скорости данной стадии, г Е {1, 2, -2, 4, 5, 6, -6, 7}. з

Вероятность реализации третьей стадии равна единице, так как молекулы СОаав и Оаав, оказавшиеся в соседних клетках, немедленно вступают в реакцию. Константа скорости стадий 8-10, описывающих диффузию, вычисляется следующим образом:

к^ = М^ц х кг, где М^ц — параметр интенсивности диффузии. Параметр М^ц соответствует скорости перемещения реагентов по поверхности катализатора.

В ходе реакции окисления на поверхности катализатора в режиме автоколебаний происходит смена адсорбционных покрытий Оа^ ^ СОа^. В начальный момент времени из газовой фазы на чистую поверхность катализатора (активное состояние) адсорбируются монооксид углерода СО и кислород О2 (рис. 1).

□ СО.*; □ Оаас; ЕИ ОзиЬ; □ свободный активный центр

Рис. 1. Смена адсорбционных покрытий в ходе реакции

Поскольку парциальное давление кислорода в газовой фазе превышает парциальное давление СО, на поверхности накапливается Оа^. В результате взаимодействия СОа^ с О^в скорость образования СО2 увеличивается. Одновременно Оа^ частично преобразуется в приповерхностный кислород (ОвиЬ). Когда концентрация ОвиЬ достигает критического значения, адсорбция кислорода блокируется и на поверхности накапливается СО^в. Это соответствует неактивному состоянию поверхности, скорость образования СО2 достигает минимального значения. Молекулы СОа^ диффундируют по поверхности катализатора к участкам, занятым ОвиЬ, и вступают с ними в реакцию. В результате реакции освобождаются активные центры для адсорбции кислорода, и поверхность вновь становится активной. Колебательный цикл повторяется.

Все элементарные стадии реакции окисления (адсорбция, десорбция, взаимодействие адсорбированных реагентов и диффузия) на поверхности катализатора выполняются асинхронно, т. е. активные центры на поверхности катализатора выбираются случайно, стадии для выбранного активного центра выбираются с заданной вероятностью. Такая модель наиболее эффективно может быть описана кинетическим клеточным автоматом.

2. Кинетический клеточный автомат для моделирования реакции окисления

2.1. О п р е д е л е н и е к и н е т и ч е с к о г о к л е т о ч н о г о а в т о м а т а

Кинетический клеточный автомат определяется множеством Ма = {А, М, Т, 0), где А — алфавит состояний клеток; М — множество имён клеток; Т — множество шаблонов; 0 — множество подстановок; а — асинхронный режим функционирования. Алфавитом состояний клеток является множество А = {*, СО^, О^, ОвиЬ, [СО^ * ОвиЬ]}, элементы которого обозначают центры на поверхности, атомы и молекулы. Множество имён М = {(г,]) : г = 0,1,..., I; ] =0,1,... ,3} представляет собой поверхность катализатора. Каждой клетке (а, (г,])) соответствует конечный автомат с именем (г]) и состоянием а. Множество клеток образует клеточный массив П(£) = {(а, (г,]))} С С А х М, где а — это состояние клетки с именем (г,]) в момент времени £.

На множестве имён определены именующие функции ^ : М — М. Множество именующих функций для клетки с именем (г,]) Е М задает шаблон соседства:

Т (г,]) = {(г,]),<Р1(г,]),...,(Рд(г,] )}-

Клетка (г,]) является центральной клеткой шаблона. Наиболее распространенным шаблоном соседства является «крест» Т(г,]) = {(г,]), (г — 1,]), (г,] + 1), (г + 1,]), (г,] — 1)}. Множество клеток с именами из шаблона

Я (г,] ) = {(^0 , (г,] )), (V1, ^1(г,] )), ..., (Уд , ^д (г,] ))}

называется локальной конфигурацией. Клетка (г>0, (г,])) — центральная клетка конфигурации Я (г,]). Две конфигурации Я (г,]) и Я '(г,]), Я '(г,]) = {(и0, (г,])), (и1, ^1(г,])),... , (ид, <£>д(г,]))}, с одной и той же центральной клеткой образуют подстановку

0(г,]): Я(г,]) А Я'(г,]). (1)

Ниже будем использовать подстановки, для которых Т(г,]) = Т'(г,]). Подстановка (1) применима к случайно выбранной клетке (г,]), если Я (г,]) С П. Если это условие не выполняется, попытка применения подстановки считается неудачной. Применение 0(г,]) к клетке с именем (г,]) состоит в следующем. Сначала вычисляются значения ик = fk (ь0,ь1,... ,уд), к = 0,1,...,д. Затем состояния клеток (ьк ,рк (г,])) Е Б'(г,]) заменяются с некоторой вероятностью р на состояния ик. Применение подстановки 0(г,]) к одной клетке выполняется за определенный отрезок дискретного времени т, называемый тактом. Применение 0(г,]) ко всем клеткам массива П(£) приводит к изменению его глобального состояния П(£) А П(£ + 1) и называется итерацией. Она состоит из |М| х М^ц тактов. Последовательность П(0), П(1),... , П(£),..., где П(0) — состояние клеточного массива в начальный момент времени, называется эволюцией.

В КА-модели реакции окисления предполагается, что взаимодействие СО^в и О^в, находящихся в соседних клетках, происходит мгновенно. Поэтому при наличии соответствующей молекулы реакция должна происходить сразу же после адсорбции молекул СОадв и Оадв и после диффузии СОадв.

Например, рассмотрим адсорбцию кислорода и последующее взаимодействие Оаав с СОаав. При применении подстановки

01(г,]) : {(* (г,])), (*,^1(г,]))} ^ {(Оаёв, (г,])), (0adв,^1(г,]))}, 1 = 1, 2, 3, 4,

для клетки с именем (г, ] ) по шаблону «крест» выбирается одна из четырёх соседних клеток ^1(г,]), сразу же после адсорбции к клеткам (г,]) и ^1(г,]) применяется подстановка

03(г,]) : { (СОадв, (г,])), (Оадв ,^1(г,]))} А {(*, (г,])), (*,^1(г,]))}, 1 = 1 2, 3, 4,

которая также по шаблону «крест» выбирает одну из четырёх соседних клеток. Следовательно, для применения подстановок 01(г,]), 03(г,]) могут понадобиться состояния тринадцати клеток (рис. 2). Обозначим полученный шаблон моделирования за В(г,]).

2.2. Результаты КА-моделирования реакции окисления

Для исследования реакции окисления СО на палладии проводились вычислительные эксперименты для КА размером 400 х 400 клеток и при следующих значениях констант скорости элементарных стадий: к1 = 1, к2 = 1, к-2 = 0,2, к4 = 0,03, к5 = 0,01,

Рис. 2. Шаблон моделирования B(i, j)

k6 = 1, k-6 = 0,5, k7 = 0,02, = 50. В исходном состоянии все клетки массива

активные, граничные условия периодические.

Динамика реакции представлена колебаниями концентраций реагентов (Oads, COads и Osub) и скорости реакции. Концентрация реагентов определяется долей клеток в массиве, находящихся в соответствующем состоянии, скорость реакции — количеством сформированных молекул CO2 в одной клетке за одну итерацию.

Рассмотрим динамику реакции в течение одного периода (рис. 3). С момента времени tl (он соответствует минимальному значению Osub) начинается смена покрытий Oads о COads (рис. 3, a), в результате которой возрастает скорость реакции. При достижении максимального значения скорости реакции (t2, рис. 3, с) происходит перераспределение концентрации Oads: Oads ^ Osub. Положение максимума на кривой Osub определяет начало снижения скорости реакции и накопление адсорбированного монооксида углерода, который медленно взаимодействует с Osub. В момент времени t3 концентрация COads достигает максимального значения и сохраняет его в течение времени t3 -t4. Уменьшение концентрации Osub до критического создает условия для последующей адсорбции кислорода с повторением автоколебательного цикла.

Полученный характер динамики реакции окисления соответствует экспериментальным данным и результатам вычислительных экспериментов, представленных в работе [1].

3. Моделирование реакции окисления

с помощью блочно-синхронного КА

3.1. Трудности распараллеливания асинхронных КА

Для изучения пространственно-временной динамики реакции окисления CO на поверхности Pdll0 требуется проводить моделирование с использованием больших клеточных массивов (1012) в течение нескольких сот тысяч итераций (106). Распараллеливание таких задач позволяет получить решение задачи за значительно меньшее время и существенно снизить необходимый для решения объём ресурсов. Однако распаралле-

к

5

со

а

I-

I

0)

X

о

О

к

со

I-

X

0)

X

о

Рис. 3. Колебания концентрации реагентов и скорости реакции: а — С0а^й и 0а^й;

Ь — 08иь; с — скорость реакции

ливание КА с асинхронным режимом функционирования сопряжено с определёнными трудностями [3].

Клетки выбираются случайным образом, поэтому в любой момент времени может быть выбрана граничная клетка, а её состояние сразу же должно быть обновлено в клеточном массиве соседнего процессора. Следовательно, при асинхронном режиме работы межпроцессорный обмен должен выполняться сразу же после изменения состояния хотя бы одной граничной клетки. Кроме того, при распараллеливании КА должно выполняться условие равноправия всех клеток, т. е. вероятность выбора всех клеток должна быть одинаковой для всех процессоров, между которыми распределён клеточный массив.

Решением проблемы является аппроксимация исходного кинетического КА блочносинхронным КА.

3.2. Алгоритм построения блочно-синхронного КА

Алгоритм преобразования асинхронного КА в блочно-синхронный КА заключается в следующем [2, 3, 5].

1) На множестве имен М определяется шаблон блока

Тв(г,і) = {(^ І ),Рl(i, І ),Р2 (і,з),...,Рг (І,І)}.

Здесь именующие функции р1(І,І ),р2(І,І),..., Рг (І, І) перечисляют г соседей центральной клетки с именем (і,і), включенным в подстановки из множества 0. Обозначим блок Б1рк(^), где Рк(І,І) —имя центральной клетки, через Вк. Блок Вк характеризуется следующими свойствами:

— Т(і,і) С Твк, где Т(і,і) —основной шаблон подстановки из множества 0;

— на множестве имён М блок В к определяет множество разбиений П = {П1, П2,... ,

Пг}, каждое из которых состоит из С = -—- блоков Пк = {В^, В\,... , В^ }. Так как Пк Є П — разбиения, для них выполняются следующие соотношения:

а

и Вд = М; Вд П В£ для всех д = к, д,к Є {1, 2,..., С}.

9=1

Второе соотношение является условием корректности. Оно требует, чтобы состояние клетки не изменялось двумя подстановками одновременно.

2) Каждая итерация разбивается на г синхронных шагов. На каждом к-м шаге, к = 1, 2,... ,г, подстановки применяются синхронно к к-й центральной клетке во всех блоках к-го подмножества.

В КА-модели реакции окисления шаблон блока Тв(^) состоит из 13 клеток и равен шаблону моделирования В(і,і), представленному выше на рис. 2. На рис. 4 показано подмножество Пб. Оно включает блоки, для которых клетка с номером 6 является центральной.

11 1 2 3 12 4 5 6 7 8

4 5 6 7 8 0 9 10 11 1

° I9 10 11 1 2 3 12 4 5

2 3 12 4 5 6 7 8 0 9

6 7 8 0 9 10 11 1 2 3

10 11 1 2 3 12 4 5 6 7

12 4 5 6 7 8 0 9 10 11

8 0 9 10 11 1 2 3 12 4

Рис. 4. Подмножество Пб

3.3. Сравнение блочно-синхронного и асинхронного

р е ж и м о в р а б о т ы К А

В экспериментах кинетический КА и блочно-синхронный КА демонстрируют одинаковый характер поведения реакции окислении. На рис. 5 показаны колебания концентрации 08иь для обоих режимов работы клеточного автомата размером 400 х 400 (—да = 50). При асинхронном и блочно-синхронном режимах работы используются одни и те же генераторы случайных чисел, поэтому последовательности применения под-

становок совпадают. Получены следующие среднеквадратичные разности между концентрациями основных реагентов (Oads, Osub, COads, CO2) для асинхронного и блочносинхронного режимов работы клеточного автомата: dCO2 = 0,0085, dCOads = 0,0278, dOads = 0,0154, dOsub = 0,0086, d[COads*Osub] = 0,0052, что показывает, что аппроксимация исходного кинетического КА блочно-синхронным корректна.

Сг

к

s

и

со

I-

0

О

Число итераций

Рис. 5. Концентрация Osub при асинхронном (АКА) и блочно-синхронном (БСКА) режимах работы

Блочно-синхронный режим работы КА показывает меньшую временную сложность по сравнению с асинхронным. Это объясняется тем, что за итерацию случайным образом выбирается 13 клеток при блочно-синхронном режиме работы и |М | клеток при асинхронном. Для нашего эксперимента время вычислений для асинхронного КА в 1,57 раз больше, чем для блочно-синхронного КА.

4. Параллельная реализация блочно-синхронного клеточного автомата

Распараллеливание блочно-синхронного клеточного автомата основано на методе декомпозиции области и заключается в следующем.

1) КА-массив размером |М | разрезается на равные непересекающиеся части (домены Дот). Домены распределяются между п процессорами суперкомпьютера. В памяти каждого процессора хранится массив размера (К + 4) х (Ь + 4), где |К| х |Ь| = |Дот| =

|М |

=-------размер домена. На рис. 6 представлены два домена в соседних процессорах

п

для подмножества П6. В результате декомпозиции клеточного массива на прямоугольные домены граничные клетки блоков хранятся в памяти разных процессоров. И тогда часть подстановок не может быть применена к граничным (рис. 6, (/+1)-й процессор) и приграничным (рис. 6, /-й процессор) клеткам. Для того чтобы применить подстановки к граничным клеткам с номером 6, необходимо передать из /-го в (/ + 1)-й процессор клетки с номерами 4, 1, 5, 9 и из (/ + 1)-го в /-й процессор — клетку с номером 8.

2) Итерация эволюции блочно-синхронного КА состоит из 13 синхронных шагов. На каждом шаге случайным образом выбирается одно из подмножеств Пк, к = 1, 2,...,г, для всех процессоров, затем во всех доменах к клеткам с номером к синхронно применяются подстановки. После применения подстановок каждый процессор пересылает новые значения граничных и приграничных клеток соседним процессорам. В результа-

процессор I процессор ( I + 1)

Рис. 6. Структура доменов в соседних процессорах

те на каждом шаге у доменов есть значения всех клеток, необходимых для выполнения следующего шага. Пакет пересылаемых данных составляет 8К клеток.

Распараллеливание блочно-синхронного КА выполнено на суперкомпьютере МВС-100К (МСЦ РАН), каждый вычислительный узел которого содержит два четырёх-ядерных микропроцессора 1П;е1 Хеоп 5365. Скорость обмена данными между узлами — 1400 Мбайт/с, время пересылки одного байта Ьь = 0,007мкс. Латентность составляет = 3,2 мкс. Время обработки одной клетки равно т = 0,058 мкс. Исходя из известного условия эффективного распараллеливания т • |Дот| > к • (Ьы + V • Ьь), где к =13 — количество обменов за итерацию, размер домена должен удовлетворять условию

|£от| > 1200 х 1200. (2)

В таблице для КА размеров 9000 х 9000 и 12000 х 12000 приведены значения эф-

Т

фективности распараллеливания Q(n) = —-----, где Т —время эволюции КА на одном

ТП • п

процессоре; Тп — время эволюции КА на п процессорах суперкомпьютера. Эксперименты показали, что эффективность параллельной реализации блочно-синхронного КА достигает 90 % при условии, что размер домена не меньше указанного в формуле (2).

Время Тп и эффективность Q(n) распараллеливания

I х 3 Параметры п

1 4 16 32 64 128

9000 х 9000 Тп, с 568,87 17,18 4,37 4,36 1,00 0,81

Q{v) 1 0,97 0,95 0,95 0,93 0,89

12000 х 12000 Тп, с 1005,54 29,61 7,53 3,74 1,76 0,87

1 0,98 0,94 0,93 0,93 0,92

5. Гистерезис покрытия поверхности ОЭаай в каталитической реакции окисления

В каталитической реакции окисления СО могут возникать различные критические явления: множественность стационарных состояний, автоколебания, хаос, гистерезис.

При этом на поверхности катализатора наблюдаются различные пространственно-временные структуры: спирали, кольца, турбулентности. Кинетический КА позволяет моделировать и изучать динамику реакции окисления.

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

В КА-модели реакции окисления гистерезис возникает при изменении константы скорости кислорода к1 от 0,7 до 1 и постоянных значениях остальных коэффициентов к: к2 = 1, к-2 = 0,2, к4 = 0,03, к5 = 0,01, к6 = 1, к-6 = 0,5, к7 = 0,02. При моделировании использовался клеточный массив размером 1000 х 1000 клеток, параметр диффузии первые 3000 итераций задавался равным = 100, затем уменьшался до Мда = 20. Начиная с 3000-й итерации, при = 20 константа скорости кислорода к1 сначала уменьшалась от 1 до 0,7 (режим 1), а затем увеличивалась от 0,7 до 1 (режим 2). При пошаговом уменьшении к1 в режиме 1, а затем увеличении в режиме 2 реакция окисления демонстрирует различный характер колебаний концентрации СОаав (рис. 7). Период и амплитуда колебаний в режимах 1 и 2 существенно отличаются друг от друга. В режиме 2 колебания носят более регулярный характер, при этом наблюдается увеличение амплитуды и периода.

Режим 1 Режим 2

Число итераций

Рис. 7. Гистерезис концентрации СОаав при пошаговом уменьшении и увеличении константы скорости кислорода

Полученные результаты совпадают с результатами численных экспериментов, представленных в работе [1].

Заключение

Представлены результаты исследования кинетического клеточного автомата, моделирующего реакцию окисления СО на платиновых металлах. Показано, что асинхронный и блочно-синхронный режимы работы КА имеют одинаковый характер поведения

реакции, причём временная сложность эволюции блочно-синхронного КА меньше временной сложности эволюции асинхронного КА. Построена параллельная реализация блочно-синхронного КА с высокой эффективностью распараллеливания. С помощью КА-модели получен гистерезис скорости реакции и покрытий поверхности адсорбированными реагентами. Результаты моделирования согласуются с работами [І, 4].

ЛИТЕРАТУРА

1. Elokhin V. I., Latkin E. I., Matveev A. V., and Gorodetskii V. V. Application of statistical lattice models to the analysis of oscillatory and autowave processes in the reaction of carbon monoxide oxidation over platinum and palladium surfaces // Kinet. Catalys. 2003. V. 44. No. Б. P. 173-180.

2. Bandman O. L. Synchronous versus asynchronous cellular automata for simulating nanosystems kinetics // Bulletin November Computer Center. Novosibirsk: NCC Publisher, 2006. V.25. P. 1-12.

3. Bandman O. L. Parallel simulation of asynchronous cellular automata evolution // Proc. ACRI-2006. LNCS. 2007. V.4173. P.41-48.

4. Elokhin V. I., Latkin E. I., Matveev A.V., and Gorodetskii V. V. Manifestation of the adsorbed CO Diffusion anisotropy caused by the structure properties of the Pd(110) — (1x2) surface on the oscillatory behavior during co oxidation reaction — Monte-Carlo model // Chemist. Sustainab. Developm. 2003. No. 11. P. 173-180.

Б. NedeaS.V., LukkienJ.J., Jansen A. P. J., and Hilbers P. A. J. Methods for parallel simulations of surface reactions // arXiv:physics/0209017. 2002. V. 1. No. 4.