CC BY
235
К. С. Шефов, М. М. Степанова, А. Н. Макаров, А. Г. Стародубов
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ МОЛЕКУЛЯРНО-ДИНАМИЧЕСКИХ СИСТЕМ С ИСПОЛЬЗОВАНИЕМ ПАКЕТА LAMMPS*
Введение. В данной публикации представлены результаты начального этапа разработки методики численного моделирования нового физического эксперимента, направленного на исследование физических условий образования длинных нанотрубок.
В реальных условиях эксперимента используется аргоновая среда с примесью метана (1-5 %) при атмосферном давлении и комнатной температуре. В ней на скорости, равной средней скорости движения молекул в газе, пролетает мелкодисперсная металлическая (Ni) частица, размер которой около 500 A. На границе области разрежения (следа), возникающей при пролёте частицы, можно ожидать эффект образования новых микроструктур непосредственно за частицей. Диаметр частицы и средняя длина свободного пробега атомов газа являются величинами одного порядка. Это обстоятельство не позволяет рассматривать газ как сплошное вещество и использовать методы газодинамики. В данном случае нужно применять молекулярно-динамический подход для большой многочастичной системы, что требует ресурсоёмких расчётов на параллельном кластере.
Задача представленной работы состояла в проведении на кластере ряда численных экспериментов по моделированию процессов молекулярной динамики и исследовании некоторых эмпирических потенциалов межатомного взаимодействия.
Были проведены следующие расчёты с использованием эмпирических потенциалов Леннарда-Джонса, AIREBO и MEAM моделирование: 1) кристалла алмаза; 2) никеля в твёрдой фазе; 3) движения твёрдого никелевого шара в метане; 4) образования углеродных цепей из газовой фазы; 5) движения в газообразном аргоне твёрдой мелкодисперсной металлической частицы (в двумерии); 6) фуллеренов. Здесь представлены последние два из этих экспериментов.
Все вычисления производились на базе научного комплекса GRID технологий СПбГУ. Использовался параллельный вычислительный кластер из 10 восьмиядерных узлов под управлением ОС Scientific Linux с пакетом MPICH2, обеспечивающим поддержку MPI.
Специализированное ПО. Пакет LAMMPS [1] использовался для численного решения уравнений молекулярной динамики для многочастичной системы. Данный пакет работает только с классическими уравнениями Ньютона. Он включает в себя широкий набор потенциалов межатомного и межмолекулярного взаимодействия и позволяет моделировать химические реакции.
Программа написана на языке С+—Н с использованием стандарта MPI для параллельного программирования. Разбиение задачи на несколько параллельных производилось с помощью алгоритма пространственного разбиения [2]. Для вычисления использовался скоростной вариант алгоритма Верле [3].
Ввод исходных параметров модели осуществлялся через текстовый конфигурационный файл, содержащий последовательный набор команд и параметров. Результат
* По материалам доклада на юбилейном семинаре «Вычислительная физика» 29—30 октября 2009 г., С.-Петербург.
© К. С. Шефов, М. М. Степанова, А. Н. Макаров, А. Г. Стародубов, 2010
записывался также в текстовый файл в виде порядковых номеров, координат и прочей информации об атомах и молекулах (тип, скорость, энергия и др.).
ЬАММРБ является мощным инструментом для проведения сложных численных экспериментов для самых разнообразных моделей. Однако ЬАММРБ имеет недостаток, который выражается в невозможности обработки выходных данных без использования сторонних программ, а также в невозможности создания молекулярных систем средствами самого пакета. Для этих целей были использованы средства визуализации (VMD, GNUPlot) и программа для создания молекулярных систем (НурегСЪет).
Потенциалы взаимодействия. Потенциалом Леннарда-Джонса описывается взаимодействие неполярных атомов и молекул друг с другом. Такое взаимодействие характерно для атомов инертных газов, которые не участвуют в химических реакциях, например аргона, а также многих других подобных систем. Потенциальная энергия частицы, поведение которой описывается данным потенциалом, имеет следующий вид:
где r - расстояние между центрами частиц, а параметры е и о определяются из эксперимента и являются характеристиками вещества.
Потенциал AIREBO (adaptive intermolecular reactive empirical bond-order) [4] является эмпирическим потенциалом, который позволяет моделировать межмолекулярное взаимодействие и химические реакции в кристаллах, жидкостях и газообразных веществах, состоящих из углерода и водорода.
Данный потенциал базируется на эмпирическом потенциале Бреннера, описывающем непосредственное взаимодействие атомов через химические связи. Кроме того, AIREBO включает взаимодействие Леннарда-Джонса и взаимодействие, зависящее от двугранных углов, образованных атомами в молекуле. Это позволяет учитывать межатомное взаимодействие на расстояниях, превышающих длину химических связей, что улучшает качество моделирования молекул.
Потенциал MEAM (modified embedded-atom method) [5] использовался нами для моделирования поведения металла с кристаллической решёткой. С помощью потенциала MEAM можно моделировать не только металлы и сплавы с различными типами кристаллических решёток, но и ковалентные кристаллы (кремний, углерод). MEAM учитывает угловое взаимодействие (валентные углы между связями).
Ряд проведённых нами численных экспериментов с использованием перечисленных потенциалов подтверждает, что данные потенциалы позволяют достаточно правдоподобно моделировать поведение газов (потенциал Леннарда-Джонса, AIREBO) и твёрдых тел (AIREBO, MEAM), а также взаимодействие тех и других между собой (потенциал Леннарда-Джонса).
Обзор численных экспериментов.
Двумерное моделирование твёрдой мелкодисперсной металлической частицы, движущейся в газообразном аргоне. Задача состояла в моделировании движения абсолютно твёрдой никелевой шайбы, размер которой одного порядка со средней длиной свободного пробега атомов аргона при нормальных условиях. Моделирование производилось в двух измерениях.
Такое упрощение оправдано в связи с изотропностью трёхмерной модели по двум направлениям, что позволяет исключить одно из них. Движение частицы моделировалось при различных скоростях: звуковой, средней тепловой, двух- и пятикратной звуковой. Для моделирования взаимодействия никеля с аргоном и атомов аргона между собой
использовался потенциал Леннарда-Джонса. Параметр обрезания потенциала модели выбран 2,5 А из-за пренебрежимо малого вклада взаимодействия на расстояниях, больших 2,5 А, в общее взаимодействие. Скорость звука в аргоне при нормальных условиях - 325 м/с. Средняя тепловая скорость частиц - 432,6 м/с.
Опишем условия численного эксперимента. В двумерном прямоугольном ящике, наполненном аргоном (588 000 частиц), при атмосферном давлении и комнатной температуре (300 К), в направлении оси X движется шайба. Граничные условия в задаче - периодические по всем направлениям. Размер шайбы - 540 А, масса - 4,52 • 10_8 г/моль, средняя длина свободного пробега частиц аргона - 750 А. Иллюстрации модели при разных скоростях наночастицы - средней тепловой и пятикратной звуковой представлены на рис. 1 и 2 соответственно. На рисунках показаны также размеры системы.
Для расчёта использовались 32 ядра. Длина шага по времени составляла 1,5 фс. Временной шаг выбирается исходя из минимизации погрешности компьютерной математики и численного алгоритма. Общее число шагов - 3 • 106, реальное время движения моделируемой частицы - 4,5 нс. Время расчёта для одного значения скорости составило около 12 ч.
Как видно на приведённых иллюстрациях, чем быстрее движется наночастица, тем более выражен оставляемый ею след (тем больше разность плотностей аргона слева и справа от шайбы). Следует также отметить, что след почти полностью затухает на расстоянии порядка 10 средних длин свободного пробега (8000 А) от наночастицы, что согласуется с теорией. После перехода скорости шайбы через значение тепловой скорости атомов аргона картина существенно меняется. При скорости частицы до тепловой включительно видна относительно протяжённая область с увеличенной плотностью газа перед частицей, за частицей изменения плотности почти нет. При скоростях выше тепловой след приобретает форму овала, образованного газом относительно высокой плотности, с областью существенного разрежения внутри.
Моделирование фуллеренов с использованием потенциала AIREBO. Одним из наиболее распространенных методов получения фуллеренов (соединения большого числа углеродных атомов в структуры регулярной формы) на практике является метод возгонки и десублимации графита в электрической дуге, горящей между графитовыми электродами в потоке инертного газа (обычно гелия) [6]. В процессе образуется сажа, содержащая малый процент фуллеренов. Температура в дуге - около 4000 К. Давление гелия в камере - 106 кПа. Очень важный момент в процессе получения фул-леренов - это нагреть углерод, чтобы он сначала перешёл в газовую фазу.
Модель, описанная ниже, строилась с целью проверить, будут ли фуллерены С60 (рис. 3), у которых удалены несколько атомов, замыкаться до целого фуллерена в атмосфере газообразного углерода при температуре 4000 К. Температура фиксировалась
0 50000 А
Рис. 1. Частица в аргоне (тепловая скорость)
0 35000 А
Рис. 2. Частица в аргоне (пятикратная звуковая скорость)
термостатом Хувера [7]. Из каждого фуллерена было удалено по 9 атомов. Размеры системы:
320 х 130 х 130 А; полное число атомов - порядка 4000. Начальная конфигурация системы приведена на рис. 4: в центре - фуллерены, по бокам - газообразный углерод. Граничные условия - периодические по всем трём осям координат. Время вычисления составило 70 мин на 80 процессорах.
За время моделирования (250 пс) ни один фуллерен не замкнулся. Спустя четверть времени моделирования все фуллерены разрушились.
Однако к концу эксперимента число связей возросло, из газовой фазы образовались углеродные цепи и кластеры, состоящие из пяти- и шестиугольников и почти замкнутых шестиугольников (рис. 5). Образование соединений в форме многоугольников объясняется их относительно большой стабильностью (ими и образованы фуллерены). Данную модель, безусловно, следует рассматривать только как тестовую, для адекватного моделирования процесса формирования фуллеренов необходимо учитывать квантовые эффекты.
Рис. 3. Молекула фулерена Сво
Рис. 4. Моделирование фуллеренов: начальное положение
Рис. 5. Результат численного моделирования:
образуются связанные шестиугольники (показаны только химические связи)
Заключение. В результате исследований, во-первых, изучены возможности и принципы работы программного пакета ЬАММРЯ для численного моделирования молекулярно-динамических систем, а именно, сделаны успешные попытки применения данного пакета к решению ряда задач молекулярной динамики на параллельном вычислительном кластере, также проанализированы численные методы, используемые
ЬАММРБ для моделирования. Во-вторых, изучены и протестированы на физических моделях численные алгоритмы реализации некоторых эмпирических потенциалов взаимодействия, а именно потенциалов Леннарда-Джонса, АШЕВО и МЕАМ. Подобные физические модели являются основой для возможного создания в будущем численной модели эксперимента по получению новых наноструктур.
Литература
1. URL: http://lammps.sandia.gov/ (дата обращения 26.10.2009).
2. Plimpton S. // J. Comp. Phys. 1995. Vol. 117. P. 1-19.
3. Гулд Х, Тобочник Я. Компьютерное моделирование в физике. Ч. 1. М., 1990.
4. Stuart S. J. // J. Chem. Phys. 2000. Vol. 112. P. 6472-6486.
5. Baskes M. I. // J. Phys. (B). 1992. Vol. 46. P. 2727-2742.
6. Раков Э. Г. Нанотрубки и фуллерены. М., 2006.
7. Hoover W. G. // J. Phys. (A). 1985. Vol. 31. P. 1695-1697.
Статья поступила в редакцию 19 марта 2010 г.
