CC BY
31
ФИЗИКА КОНДЕНСИРОВАННОГО СОСТОЯНИЯ ВЕЩЕСТВА
Суперкомпьютерное моделирование процесса напыления тонких пленок диоксида кремния с использованием программы LAMMPS
А. А. Горох, Ф. В. Григорьевa, Е.В. Каткова, А. В. Сулимов, С. А. Шарапова
Московский государственный университет имени М. В. Ломоносова, Научно-исследовательский вычислительный центр. Россия, 119991, Москва, Ленинские горы, д. 1, стр. 4.
E-mail: afedor.grigoriev@gmail.com
Статья поступила 16.07.2015, подписана в печать 02.11.2015.
Методом классической молекулярной динамики проведено суперкомпьютерное моделирование процесса высокоэнергетического напыления (ion beam sputtering) тонких пленок диоксида кремния. Обсуждаются особенности используемого в рамках программы LAMMPS метода напыления. Проведен структурный анализ полученной пленки (плотность в зависимости от толщины, радиальная функция распределения, среднее расстояние связи Si—O, концентрация точечных дефектов различных видов).
Ключевые слова: суперкомпьютерное моделирование, молекулярная динамика, тонкие пленки, диоксид кремния, напыление, LAMMPS.
УДК: 539.231. PACS: 81.15.Aa.
Введение
Напыление тонких пленок осаждением на подложку высокоэнергетических атомов кремния (10-100 эВ) в настоящее время является одним из перспективных методов, позволяющих получить плотную и однородную пленку [1], что важно с точки зрения ее оптических свойств. Параметры процесса напыления (энергия и угловое распределение осаждаемых атомов, температура подложки и др.) во многом определяют структурные свойства пленок, важные в практическом использовании. Из-за малой, около одного микрона, толщины пленок и неупорядоченной структуры экспериментальное исследование такого влияния затруднено. Вследствие этого использование математического моделирования для изучения влияния параметров процесса напыления на свойства получаемых пленок представляется целесообразным.
Основными проблемами при моделировании процесса напыления являются малая скорость роста пленки — около 0.3 нм/с [2] и большие размеры атомистических кластеров, необходимых для исследования значимых неоднородностей структуры с характерным размером в десятки нм. Вследствие этого плотность потока осаждаемых атомов при моделировании существенно выше, чем в эксперименте. В рамках схемы, разработанной в [3, 4] для пленок диоксида кремния, процесс напыления организован как последовательность коротких (длительность 6-10 пс) молекулярно-динамических запусков. В начале каждого запуска осаждаемые атомы кремния и кислорода в количественном соотношении 1 : 2 располагаются в верхней части области моделирования и имеют начальные скорости, направленные вертикально вниз. Как следует из приведенных в [3, 4] оценок, временного интервала в 6-10 пс достаточно
для достижения атомами пленки, образования связей кремний-кислород и температурной релаксации. Количество атомов, напыляемых за один запуск, достигало 330.
С точки зрения адекватности схемы моделирования важно проверить, зависят ли свойства пленки от способа построения потока осаждаемых атомов. В настоящий работе, в отличие от [3, 4], рассматривается метод непрерывного введения атомов кремния и кислорода в область моделирования, реализованный средствами молекулярно-динамической программы LAMMPS [5]. Ранее эта программа уже применялась для моделирования процесса осаждения пленок на подложку [6], в том числе оксидных пленок [7]. Нами проведено сравнение структурных свойств пленки, полученной в [4] и в настоящей работе, сделаны выводы о преимуществах и недостатках двух методов построения потока осаждаемых атомов.
1. Метод моделирования
Рассмотрим процедуру напыления пленки с использованием молекулярно-динамической программы LAMMPS. Версия LAMMPS v. 18 Sep. 2014 была скомпилирована на суперкомпьютере «Ломоносов» суперкомпьютерного центра МГУ им. М. В. Ломоносова [8] с дополнительными библиотеками USER-CUDA (позволяет использовать видеокарты в расчетах) и USER-MISC (описывает команду для вставки новых атомов в область симуляции fix_de-posit). В качестве компилятора был использован mpicxx (OpenMPI).
Область симуляции представляет собой параллелепипед размерами 28 х 23 х 30 нм (первые два размера по горизонтали, последний — по вертикали) с периодическими граничными условиями по всем
направлениям. Размеры этой области не меняются во время вычислений. Используется ансамбль ЫУТ (постоянные число частиц, объем, температура Т = 300 К). Для поддержания в области моделирования постоянной температуры применялся термостат Берендсена [9], в рамках которого через определенное число шагов скорости атомов подложки и напыленной пленки масштабируются так, чтобы средняя кинетическая энергия в системе соответствовала заданной температуре. Соблюдение эргодичности в ходе моделирования обеспечивается использованием стандартных молекулярно-динами-ческих методик (термостат Берендсена не нарушает эргодичность для систем с большим (десятки тысяч) числом атомов при моделировании относительно коротких по времени (до нескольких наносекунд) процессов), а также небольшой величиной плотности потока осаждаемых атомов, что позволяет системе релаксировать после каждого взаимодействия с высокоэнергетическим атомом кремния.
На нижней границе области симуляции размещается подложка толщиной 2.5 нм, подготовленная в соответствии с процедурой, описанной в [3]. Общее число атомов кремния и кислорода в подложке составляет 9 • 104. Горизонтальные размеры подложки совпадают с горизонтальными размерами области моделирования. Для описания потенциальной энергии межатомного взаимодействия используется силовое поле ОЕБ1Ь [4].
При моделировании напыления атомы кремния и кислорода поочередно (после одного кремния два
кислорода) инжектируются в область моделирования командой fix_deposit. Атомы вставляются на высоте 4 нм от поверхности пленки с энергиями 10 эВ (кремний) и 0.05 эВ (кислород). Начальные координаты атомов выбираются случайно, с равномерным распределением, независимо друг от друга. Инжек-ция атомов кремния происходит каждые 40 фс, кислорода — 20 фс (рис. 1).
В результате взаимодействия с высокоэнергетическими атомами кремния нижние атомы подложки с вертикальной координатой, близкой к нулю, могут пересечь границу области моделирования. Тогда, в соответствии с периодическими граничными условиями, они появятся в верхней части области. Чтобы избежать этого эффекта, в плоскости z = 0 была установлена стенка с потенциалом Леннард-Джонса (команда fix wall/lj93 zlo EDGE 1.0 1.0 2.5). Сверху от этой стенки и размещалась подложка.
Чтобы избежать накопления заряженного облака атомов кислорода над пленкой, на высоте 5 нм от пленки устанавливается еще одна стенка (командой fix oneway). Тогда те атомы кислорода, которые отразились при первом столкновении с пленкой, снова попадают на нее через время порядка 10 пс. Во время напыления высота пленки постепенно растет, поэтому стенка равномерно сдвигается вверх со скоростью 2.3 нм/нс, что на 50% превышает скорость роста пленки: 1.5 нм/с. Данная величина была рассчитана для режима, когда частота вставки новых атомов равна 75 пс_1, а плотность пленки 2.2 г/см3.
Рис. 1. Атомистическая
структура подложки и напыленной пленки. Показан инжектируемый атом и направление его начальной скорости
ФИЗИКА КОНДЕНСИРОВАННОГО СОСТОЯНИЯ ВЕЩЕСТВА
93
2. Результаты и обсуждение
Вид радиальной функции распределения для подложки и пленки показан на рис. 2, общее число атомов в пленке 142 000. Зависимости близки друг к другу, только высота первого пика, соответствующего длине связи 51—0 0.164 нм, для подложки несколько выше, при этом ширина пика практически не изменилась. Положение первого пика соответствует экспериментальному [10].
Рис. 2. Радиальная функция распределения Т(г) напыленной пленки (сплошная линия) и подложки (точечная линия) в зависимости от расстояния от центра между атомами
Профили плотности подложки, напыленной пленки и суммарной плотности показаны на рис. 3. Плотность напыленной пленки несколько ниже 2.5 г/см3, размер переходных областей пленка-вакуум около 1.5 нм, пленка-подложка — около 2 нм, что соответствует результатам, полученным с использованием методики пошагового напыления [4]. Ширина переходной области пленка-подложка соответствует средней глубине проникновения атома кремния
Рис. 3. Зависимость плотности р от толщины Н пленки: подложка (кривая 1), пленка (2), суммарная плотность (3)
с начальной кинетической энергией 10 эВ в подложку [3]. Как показано в [3], кинетическая энергия атома кремния после столкновения с подложкой снижается до величины, соответствующей поддерживаемой в системе температуре, в течение 1 пс (характерное время термализации). Плотность потока осаждаемых атомов выбиралась так, чтобы средний промежуток времени между попаданием атома кремния в точки подложки, расположенные на расстоянии менее 1 нм, в несколько раз превышал время термализации. Это условие обеспечивает температурную релаксацию подложки и напыленных слоев пленки в области их взаимодействия с высокоэнергетическим атомом кремния.
Основные структурные характеристики подложки и пленки приведены в таблице. Плотность пленки несколько выше, чем плотность подложки, что соответствует полученным ранее результатам. Средняя
Структурные характеристики подложки и пленки
Объект Подложка Пленка
Плотность, г/см3 2.16 2.44
Среднее расстояние, нм Si-0 0.166 0.166
0-0 0.270 0.270
Процентное содержание основного состояния кремния $¡4 и дефектов (нижний индекс указывает координацию), % Si3 0.43 0.79
Si4 99.45 98.79
Si5 0.12 0.42
Процентное содержание основного состояния кислорода Оо и дефектов (нижний индекс указывает координацию), % Oi 2.51 6.43
o2 96.16 90.46
Оз 1.32 2.95
Средний угол, град O-Si-O 109.3 109.3
Si-O-Si 144.8 142.0
длина связи Б1—0 одинакова для подложки и пленки, ее отличие от экспериментального значения 0.164 нм незначительно. Расстояние между ближайшими атомами кислорода (второй пик РФР) также совпадает для подложки и пленки.
Процентное содержание трехкоординированных атомов кремния Б13 в пленке почти вдвое выше, чем в подложке, а процентное содержание пятико-ординированных атомов 5 выше почти в четыре раза. Однако главным видом точечных дефектов и в подложке и в пленке является терминальный атом кислорода 01. В целом число атомов кислорода в пленке с координацией, отличной от двойной, достигает 10% и почти втрое превышает аналогичное значение для подложки. Значения валентных углов для подложки и пленки почти совпадают.
Отметим, что концентрация точечных дефектов каждого из перечисленных видов в пленке, полученной пошаговым напылением [4], лежит в пределах 1 %. Таким образом, пленка, напыленная изложенным в настоящей работе методом, содержит существенно больше точечных дефектов и поэтому является менее равновесной. Возможно, такое различие обусловлено случайным распределением начальных координат атомов кислорода и кремния при равномерном напылении, в то время как в методе пошагового напыления атомы организованы в стехиометрические группы Э102 (не обязательно молекулы), начальные координаты которых также случайны с некоторыми ограничениями (подробнее описано в [4]). В первом случае напыляемые атомы Б1 и О формируют химические связи Б1-0 не друг с другом, а с атомами подложки и ранее напыленными слоями пленки в области соударения, что повышает вероятность формирования дефектов в сравнении с вариантом, когда связи Б1-0 образуются напыляемыми атомами, локализованными достаточно близко друг от друга. Поскольку в рамках силового поля ОЕБ1Ь атомы Б1 и О заряжены (заряд кремния положителен и согласно стехиометрии вдвое больше по величине заряда кислорода) локализация атомов в различных областях подложки и пленки (непрерывное напыление со случайными координатами) увеличивает флуктуацию зарядовой плотности, что способствует образованию дефектов.
Таким образом, за исключением концентрации точечных дефектов, метод непрерывного напыления воспроизводит основные структурные характеристики пленки (профили плотности, РФР), полученные развитым ранее методом пошагового напыления. Поскольку экспериментальные значения концентра-
ция дефектов все же существенно ниже тех, что приведены в таблице, то в настоящий момент предпочтительнее выглядит метод пошагового напыления, поскольку дает менее дефектную структуру. Возможно, увеличение промежутка времени между последовательным введением в область моделирования осаждаемых атомов в рамках метода непрерывного напыления приведет к уменьшению концентрации дефектов вследствие более длительной релаксации.
Заключение
В работе с использованием молекулярно-дина-мической программы LAMMPS реализован метод атомистического моделирования процесса высокоэнергетического напыления тонких пленок диоксида кремния. В отличие от развитого ранее подхода [4] с пошаговым введением групп осаждаемых атомов в область моделирования, предложенный метод использует непрерывное введение осаждаемых атомов. Получено, что структурные характеристики пленки (плотность, радиальная функция распределения) близки для обоих подходов, однако концентрация точечных дефектов при непрерывном напылении существенно выше.
Работа выполнена при финансовой поддержке Российского научного фонда (грант 14-11-00409).
Список литературы
1. Piegari A., Flory F. Optical Thin Films and Coatings. Cambridge, 2013.
2. Kaiser N., Pulker H.K. Optical Interference Coatings.
B.; Heidelberg, 2003.
3. Григорьев Ф.В., Сулимов В.Б., Кондакова О.А. и др. // Вестн. Моск. ун-та. Физ. Астрон. 2013. 68, № 3.
C. 80. (Grigoriev F.V., Sulimov V.B., Kondakova O.A. et al. // Moscow University Phys. Bull. 2013. 68, N 3. P. 259.)
4. Grigoriev F.V., Sulimov A.V., Kochikov I.V. et al. // Int. J. of High Perf. Comp. Appl. 2015. 29, N 2. P. 184.
5. http://lammps.sandia.gov/.
6. Xie L., Brault P., Bauchire J.-M. et al. // J. Phys. D: Appl. Phys. 2014. 47. 224004.
7. Bahramian A. // Surf. Interface Anal. 2013. 45, N 11-12. P. 1727.
8. http://parallel.ru/cluster.
9. Berendse H.J.C., Postma J.P.M., van Gunsteren W.F. et al. // J. Chem. Phys. 1984, 81, N 8. P. 3684.
10. Johnson P.A.V., Wright A.C., Sinclair R.N. // J. Non-Cryst. Solids. 1983. 58, N 1. P. 109.
High-performance modeling of the deposition of a silicon dioxide thin Him using the LAMMPS program
A. A. Gorokh, F. V. Grigorieva, E.V. Katkova, A. V. Sulimov, S.A. Sharapova
Research Computing Center, Lomonosov Moscow State University, Moscow 119991, Russia. E-mail: a fedor.grigoriev@gmail.com.
Supercomputer modeling of the process of high-energy deposition (ion-beam sputtering) of thin films of silicon dioxide using the molecular dynamics (MD) approach was carried out. The deposition method based on the facilities of the LAMMPS MD program was compared with another method that is known from the literature. An analysis of the structure of the deposited film (density, radial distribution function, concentration of defects) was carried out.
Keywords: high-performance simulation, molecular dynamics, thin-film growth, deposition process, silicon dioxide, LAMMPS. PACS: 81.15.Aa. Received 16 July 2015.
English version: Moscow University Physics Bulletin. 2016. 71, No. 1. Pp. 114-117.
Сведения об авторах
1. Горох Артур Анатольевич — специалист; тел.: (495) 939-2346, e-mail: arturgorokh@yahoo.com.
2. Григорьев Федор Васильевич — канд. хим. наук, вед. науч. сотрудник; тел.: (495) 939-32-53, e-mail: fedor.grigoriev@gmail.com.
3. Каткова Екатерина Владимировна — канд. физ.-мат. наук, науч. сотрудник; тел.: (495) 939-32-53, e-mail: katkova@dimonta.com.
4. Сулимов Алексей Владимирович — системный программист; тел.: (495) 939-32-53, e-mail: as@dimonta.com.
5. Шаропова Светлана Анатольевна — мл. науч. сотрудник; тел.: (495) 939-2346, e-mail: svet.sharapova@gmail.com.
