Письма в ЖЭТФ, том 111, вып. 2, с. 101 - 106

Моделирование взаимодействия графена с поверхностью меди

с помощью модифицированного потенциала Морзе

С. В. Колесников¹⁾, А. В. Сидоренков, А. М. Салецкий

Московский государственный университет им. М. В. Ломоносова, 119899 Москва, Россия

Поступила в редакцию 19 ноября 2019 г.

После переработки 5 декабря 2019 г.

Принята к публикации 5 декабря 2019 г.

Предложен новый потенциал взаимодействия углерод-медь, позволяющий моделировать муаровую

структуру графена на поверхности меди. Показано, что получающаяся в результате моделирования

муаровая структура качественно согласуется с изображениями, полученными с помощью сканирующе-

го туннельного микроскопа. Толщина муаровой структуры и энергия связи графена с поверхностью в

пределах погрешности совпадают с известными на сегодняшний день экспериментальными данными.

Предложенный потенциал может быть использован и для моделирования диффузии атомов меди по

поверхности графена. Исследована диффузия атома и димера меди в широком интервале температур.

Обнаружено, что при моделировании диффузии необходимо учитывать вклад колебательной свободной

энергии атомов меди.

DOI: 10.31857/S0370274X20020095

Графен известен своими уникальными свойства-

кластеров необходимо моделирование системы боль-

ми [1-5], благодаря которым он является одним из

шого числа атомов, что на сегодняшний день воз-

наиболее перспективных материалов для техниче-

можно только в классическом приближении [10, 15,

ского применения. Одним из наиболее перспектив-

16]. При этом точность моделирования существен-

ных методов получения графена высокого качества

но зависит от выбора потенциала взаимодействия

является метод химического осаждения паров на ме-

углерод-медь.

таллических подложках [6, 7], в частности, на по-

В данной работе мы предлагаем использовать мо-

верхности меди [8]. Графен на поверхности Cu(111)

дифицированный потенциал Морзе для моделирова-

обладает муаровой структурой [9-11], появление ко-

ния взаимодействия графена с атомами меди. Ниже

торой связано с небольшим различием в межатом-

мы покажем, что использование данного потенциа-

ных расстояниях Cu-Cu и C-C [6]. Период и толщина

ла при МД моделировании позволяет адекватно вос-

муаровой структуры графена зависит от угла пово-

произвести как муаровую структуру графена, так и

рота Θ, который определяется как угол между на-

диффузию атомов меди на поверхности графена.

правлением [110] в кристалле меди и направлением

Моделирование графена на поверхности Cu(111)

“зигзаг” в графене. Экспериментально обнаружены

проводилось методом классической МД с использо-

муаровые структуры с углами поворота 0^◦, 7^◦ [9] и

ванием цепочки термостатов Нозе-Гувера [17, 18].

10.4^◦ [10]. При напылении на графен атомов металла

Для моделирования взаимодействия атомов уг-

в углублениях муаровой структуры могут образовы-

лерода между собой был использован потенциал

ваться нанокластеры, выстроенные в упорядоченную

Терсоффа-Бреннера

[19]. В этом приближении

сверхрешетку [12].

энергия связи атомов углерода равна

Для изучения взаимодействия графена с поверх-

∑∑

E_C-C =

(V_R(r_ij ) -

B_ij · V_A(r_ij)),

(1)

ностью металлов широко используются методы ком-

i j>i

пьютерного моделирования, такие как метод молеку-

лярной динамики (МД). Для исследования свойств

√

A^C-C

небольших участков графена можно использовать

V_R(r_ij) =

e^-

2Sβ(rij -r0-C)f_c(r_ij ),

(2)

S-1

первопринципную МД, основанную на теории функ-

√

ционала плотности (ТФП) [13, 14]. Однако для иссле-

A^C-CS

V_A(r_ij) =

e^-

2/Sβ(rij -r0-C)f_c(r_ij ),

(3)

дования муаровой структуры графена и роста нано-

S-1

где r_ij

- расстояние между атомами i и j, а A^C-C,

¹⁾e-mail: kolesnikov@physics.msu.ru

rC-C0, S и β - подгоночные параметры. Для ограни-

Письма в ЖЭТФ том 111 вып. 1 - 2

2020

101

102

С. В. Колесников, А. В. Сидоренков, А. М. Салецкий

чения области взаимодействия атомов используется

и 2704 атома углерода), а для угла Θ = 10.4◦ -

функция обрезания

93.9Å × 162.7Å (8 слоев по 2702 атома меди и 5844

⎧

атома углерода). Атомы двух нижних слоев меди бы-

⎪1

⎪

rij < R⁽¹⁾,

ли зафиксированы. На движение остальных атомов

⎪

⎨

{

были наложены периодические граничные условия

^[π(r_ij - R⁽¹⁾)^]}

f_c(r_ij)=

в плоскости поверхности. Для моделирования диф-

1 + cos

R⁽¹⁾<r_ij<R⁽²⁾,

⎪

(R⁽²⁾ - R⁽¹⁾)

⎪

фузии атома и димера меди по поверхности графе-

⎪

на была использована вычислительная ячейка разме-

^⎩0

rij > R⁽²⁾,

ром 34.4Å× 34.1Å (448 атомов углерода). На движе-

(4)

ние атомов были наложены периодические гранич-

где R⁽¹⁾ и R⁽²⁾ - радиусы обрезания. Многочастич-

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

ный характер взаимодействия между атомами угле-

Использование центрально-симметричных потен-

рода учитывается с помощью функции

B_ij = (B_ij +

циалов для моделирования взаимодействия углерод-

+ B_ji)/2, где

медь приводит к неправильной муаровой структу-

[

]_-δ

∑

ре графена [15] на поверхности Cu(111). Муаровая

B_ij =

1 + G(θ_ijk)f_c(r_ik)

(5)

структура оказывается инвертированной по отноше-

нию к той, что наблюдается в эксперименте [9, 11],

где θ_ijk - угол между связями i - j и i - k. Функция

т.е. вместо “холмов” на поверхности графена образу-

G(θ_ijk ) задается формулой

ются “впадины”. Это связано с тем, что атому угле-

(

)

рода оказывается энергетически выгодно иметь как

c²⁰

G(θ_ijk) = a₀

(6)

можно больше ближайших соседей меди. Для то-

d²⁰

+ [1 + cos θ_ijk]²

го, чтобы при моделировании получилась правиль-

ная муаровая структура графена, необходимо, что-

Параметры потенциала Терсоффа-Бреннера взяты

из статьи [19]: A^C-C = 6.325 эВ, S = 1.29, β = 1.5Å-1,

бы атому углерода было энергетически выгодно на-

rC-C0 = 1.315Å, R⁽¹⁾ = 1.7Å, R⁽²⁾ = 2Å, δ = 0.805,

ходиться в положении над атомом меди [6]. Это озна-

a₀ = 0.0113, c₀ = 19.0, d₀ = 2.5.

чает, что взаимодействие графена с поверхностью

В качестве потенциала взаимодействия атомов

Cu(111) осуществляется в основном за счет элек-

меди между собой использован потенциал [20], по-

трона, находящегося на p_z орбитали атома углеро-

лученный в приближении сильной связи

да и не участвующего в образовании связей углерод-

углерод. Поскольку волновая функция электрона на

E_Cu-Cu =

(7)

p_z орбитали зависит от угла φ = (

,n) между радиус-

вектором электрона r и нормалью к поверхности гра-

⎡

(

))

∑

фена n, как cosφ, то притягивающая часть потенци-

rij

⎣_ACu-Cu

exp

-p

-1

f_c(r_ij) -

ала взаимодействия углерод-медь должна тоже за-

r^Cu-Cu

j>i

висеть от cosφ.

⎤

⎛

⎞

1/2

(

))

Если центрально-симметричный потенциал име-

∑

rij

⎥

ет вид_ij(V_A(r_ij)+ V_R(r_ij)), где слагаемые V_A и V_R

-ξ^⎝

exp

-2q

−1

f_c(r_ij)^⎠

⎦,

r^Cu-Cu

отвечают за притяжение и отталкивание атомов, то

мы можем модифицировать его следующим образом:

где A^Cu-Cu, ξ, p, q, rCu-Cu0 - параметры потенциа-

∑∑(

)

ла. Функция обрезания f_c(r_ij ) была взята в форме

E_C-Cu =

V_A(r_ij)| cos(φ_ij)|^δ + V_R(r_ij)

(8)

(4) с радиусами обрезания, равными R⁽¹⁾

= 6.5Å

i j>i

и R⁽²⁾ = 7.5Å. Параметры потенциала Cu-Cu взя-

ты из работы [20]: A^Cu-Cu = 0.0854 эВ, ξ = 1.224 эВ,

где r_ij - расстояние между атомами i и j, φ_ij - угол

p = 10.94, q = 2.28, rCu-Cu0 = 2.556Å. Использование

между нормалью к поверхности графена и радиус-

потенциалов (7) хорошо зарекомендовало себя при

вектором, соединяющем атомы i и j. Параметр δ дол-

моделировании поверхностных свойств меди [21-23].

жен быть подобран таким образом, чтобы графен

Для моделирования муаровой структуры графе-

имел правильную муаровую структуру. При этом

на на поверхности Cu(111) были использованы вы-

толщина муаровой структуры Δz при различных уг-

числительные ячейки следующих размеров. Вычис-

лах поворота Θ и энергия связи графена с поверхно-

лительная ячейка для угла Θ

= 0^◦ имела раз-

стью меди E должны совпадать с экспериментально

мер 63.9Å × 110.7Å (8 слоев по 1250 атомов меди

измеренными величинами (см. табл. 1).

Письма в ЖЭТФ том 111 вып. 1 - 2

2020

Моделирование взаимодействия графена с поверхностью меди. . .

103

Рис. 1. (Цветной онлайн) Муаровые структуры графена на поверхности Cu(111), полученные с использованием мо-

дифицированного потенциала Морзе: (а), (c) - Θ = 0^◦, (b), (d) - Θ = 10.4^◦. Цвета на рисунках (a), (b) обозначают

расстояния в ангстремах между графеном и поверхностью меди в соответствии с цветовой шкалой слева от рисунков.

На рисунках (c), (d) показаны профили графена вдоль линий, показанных на рисунках (a), (b)

Таблица 1. Сравнение вычисленных значений энергии связи

качестве потенциала взаимодействия углерод-медь

графена с поверхностью меди E, толщины муаровой структу-

модифицированный потенциал Морзе

ры графена Δz, энергии связи графена с атомом меди E₁ и

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

∑∑

(

сти графена E_D с экспериментальными данными и результа-

E_C-Cu =

A^C-Cu e-2α(rij-r0-Cu) -

тами ТФП расчетов. В скобках указан угол поворота муаровой

i j>i

)

структуры Θ

- 2e-α(rij-r0-Cu) |cos(φ_ij

)|^δ

f_c(r_ij),

(9)

Вычисленное

Данные из

Величина

значение

литературы

где AC-Cu = 0.120 эВ, rC-Cu0

= 2.054Å, α = 1.05,

110 (0^◦)

E, мэВ

110 ± 11 [24]

δ = 6.75, а функция обрезания f_c(r_ij) имеет такие

110 (10.4^◦)

же параметры, как в случае взаимодействия медь-

0.351 (0^◦)

0.35 ± 0.10 [11]

Δz,Å

медь. Из таблицы 1 и рис. 1 видно, что при данных

0.149 (10.4^◦)

0.15 ± 0.05 [10]

значениях параметров потенциала величины энергии

E₁, мэВ

182.1

227.3 [12]

связи графена с поверхностью меди E и толщины му-

E_D, мэВ

1.81

3.7 [12]

аровой структуры графена Δz совпадают в пределах

погрешности с экспериментально измеренными вели-

чинами, а вид муаровой структуры и профили гра-

Попытка модифицировать таким образом потен-

фена вдоль линий, соединяющих “холмы” муаровой

циал Леннарда-Джонса привела к неудовлетвори-

структуры, качественно согласуются с СТМ изобра-

тельному результату: “холмы” муаровой структуры

жениями графена [9-11].

оказываются слишком узкими по сравнению с тем,

Итак, модифицированный потенциал Морзе (9)

что наблюдается в эксперименте [9, 11]. Более адек-

позволяет правильно моделировать муаровую струк-

ватным для описания ковалентной связи является

туру графена на поверхности Cu(111). Теперь мож-

потенциал Морзе. Поэтому далее мы используем в

но попробовать расширить область его применения

Письма в ЖЭТФ том 111 вып. 1 - 2

2020

104

С. В. Колесников, А. В. Сидоренков, А. М. Салецкий

и проверить, можем ли мы его использовать для мо-

делирования диффузии атомов меди по поверхности

графена. Ниже рассматривается диффузия одиноч-

ного атома и димера меди по поверхности свободно-

го графена. Данное исследование необходимо для де-

монстрации возможности использования потенциала

(9) при решении более сложных задач, связанных с

диффузией атомов.

Из соображений симметрии очевидно, что равно-

весные положения атома меди на поверхности графе-

на могут находиться в точках одного из трех типов:

над атомом углерода (T), над серединой ковалент-

ной связи C-C (B) или над центром шестиугольни-

ка, образованного атомами углерода (H)²⁾. Энергия

связи атома меди с графеном максимальна в поло-

Рис. 2. (Цветной онлайн) Зависимость натурального

[K]

жении (T) и равна E = 182.1 мэВ. В положениях

логарифма величины α = 4 · 108D

от обрат-

[m²/s] T

(B) и (H) энергия связи оказывается ниже на 1.81

ной температуры T-1 для диффузии атома меди по

и 13.44 мэВ соответственно. Таким образом, диффу-

поверхности свободного графена. На вставках изобра-

зионный барьер для прыжка атома меди по поверх-

жены траектории атома меди при температурах 7 и

ности графена равен E_D = 1.81 мэВ. Эти результаты

50 K

хорошо согласуются с результатами ТФП расчетов

ма меди при температурах 7 и 50 K. Видно, что точ-

(см. табл. 1).

ки, отвечающие прыжковому механизму диффузии

Характер диффузии атома меди по поверхно-

атома лежат вдоль прямой линии, проведенной ме-

сти графена существенно зависит от температу-

тодом наименьших квадратов. Из угла наклона этой

ры. Согласно приведенным выше величинам энер-

прямой следует, что ΔF = 1.28 ± 0.04 мэВ, откуда

гии связи мы ожидаем, что при температуре ниже

Δf_vib = -0.53 мэВ. Мы видим, что разность колеба-

1.81 мэВ/k_B ≈ 21 K диффузия будет иметь прыжко-

тельных свободных энергий Δf_vib оказывается одно-

вый характер. В этом случае коэффициент диффу-

го порядка величины с диффузионным барьером E_D,

зии атома может быть вычислен по формуле [25]

поэтому ей нельзя пренебрегать при моделировании

(

)

диффузии атомов меди по поверхности графена³⁾.

k_BT nl

ΔF

exp

(10)

Аналогичные вычисления были выполнены для

k_BT

диффузии димера меди по поверхности графена. Ре-

где n = 3 - количество направлений диффузии ато-

зультаты представлены на рис. 3. Движение диме-

ма, l = 2.46Å - расстояние между устойчивыми по-

ра складывается из перемещения его центра масс

ложениями атома, ΔF = E_D + Δf_vib - разница сво-

и вращения атомов вокруг центра масс. При низ-

бодных энергий атома в равновесном положении и в

ких температурах (3-7 K) преобладает вращательное

седловой точке, f_vib - колебательная свободная энер-

движение димера. При повышении температуры до

гия атома.

15-100 K происходит диффузия почти свободно вра-

Коэффициент диффузии атома может быть вы-

щающегося димера. Поэтому на зависимости вели-

числен методом МД согласно формуле Эйнштейна

чины ln α от обратной температуры 1/T можно вы-

делить два прямолинейных участка, отвечающих по-

Δx(t)² + Δy(t)

ступательному и вращательному движению димера.

(11)

Из линейной аппроксимации этих участков зависи-

мости находим величины ΔF для поступательного

где Δx(t) и Δy(t) - смещения атома относительно

4.32 ± 0.13 мэВ и вращательного 0.33 ± 0.04 мэВ дви-

начального положения в плоскости графена. Усред-

жения димера.

нение проводится по различным траекториям ато-

ма меди. Зависимость величины ln α ≡ ln(4D/T ) от

³⁾Отметим, что грубая оценка по формуле Аррениуса без

обратной температуры 1/T представлена на рис. 2.

учета вклада колебательной свободной энергии приводит

На вставках изображены примеры траекторий ато-

к неправильной температурной зависимости коэффициента

диффузии. Кроме того, как видно из рис. 2, существенные от-

²⁾Обозначения T, B, H происходят от английских слов “top”,

клонения от закона Аррениуса начинаются при температуре

“ bond” и “hollow”.

выше 25 К.

Письма в ЖЭТФ том 111 вып. 1 - 2

2020

Моделирование взаимодействия графена с поверхностью меди. . .

105

При выполнении работы были использованы вы-

числительные ресурсы Научно-исследовательского

вычислительного центра Московского государ-

ственного университета М.В.Ломоносова (НИВЦ

МГУ) [26].

K. S. Novoselov, A.K. Geim, S. V. Morozov, D. Jiang,

Y. Zhang, S. V. Dubonos, I. V. Grigorieva, and

A. A. Firsov, Science 306, 666 (2004).

K. S. Novoselov, A.K. Geim, S. V. Morozov, D. Jiang,

M. I. Katsnelson, I. V. Grigorieva, S. V. Dubonos, and

A. A. Firsov, Nature 438, 197 (2005).

A. H. Castro Neto, F. Guinea, N. M. R. Peres,

K. S. Novoselov, and A.K. Geim, Rev. Mod. Phys. 81,

109 (2009).

Рис. 3. (Цветной онлайн) Зависимость натурального

M. A. H. Vozmediano, M. I. Katsnelson, and F. Guinea,

[K]

логарифма величины α = 4 · 108D

от обратной

Phys. Rep. 496, 109 (2010).

[m²/s] T

температуры T-1 для диффузии димера меди по по-

E. Nakhmedov, E. Nadimi, S. Vedaei, O. Alekperov,

верхности свободного графена. На вставках изображе-

F. Tatardar, A. I. Najafov, I. I. Abbasov, and

ны траектории атомов меди при температурах 3, 7 и

A. M. Saletsky, Phys. Rev. B 99, 125125 (2019).

25 K

J. Wintterlin and M.-L. Bocquet, Surf. Sci. 603, 1841

(2009).

H. Tetlow, J. P. de Boer, I. J. Ford, D. D. Vvedensky,

J. Coraux, and L. Kantorovich, Phys. Rep. 542, 195

Итак, мы показали, что использование модифи-

(2014).

цированного потенциала Морзе (9) приводит к пра-

И. В. Антонова, УФН 183, 10 (2013).

вильному моделированию муаровой структуры гра-

L. Gao, J. R. Guest, and N. P. Guisinger, Nano Lett. 10,

фена на поверхности меди. Муаровая структура ка-

3512 (2010).

чественно согласуется с изображениями, получен-

10.

P. Süle, M. Szendrö, C. Hwang, and L. Tapasztó, Carbon

ными экспериментально с помощью сканирующего

77, 1082 (2014).

туннельного микроскопа [9-11], а толщина муаровой

11.

L. Zhao, K. T. Rim, H. Zhou, R. He, T. F. Heinz,

структуры и энергия связи графена с поверхностью

A. Pinczuk, G. W. Flynn, and A. N. Pasupathy, Solid

меди в пределах погрешности совпадают с экспери-

State Comm. 151, 509 (2011).

ментальными данными [10, 11, 24]. Кроме того, мо-

12.

X. Liu, Y. Han, J. W. Evans, A. K. Engstfeld,

R. J. Behm, M. C. Tringides, M. Hupalo, H.-Q. Lin,

дифицированный потенциал Морзе дает адекватные

L. Huang, K.-M. Ho, D. Appy, P. A. Thiel, and C.-

значения энергии связи и диффузионного барьера

Z. Wang, Prog. Surf. Sci. 90, 397 (2015).

для атома меди на поверхности графена. При этом

13.

I. Hamada and M. Otani, Phys. Rev. B 82, 153412

оказывается, что энергия связи атома меди с гра-

(2010).

феном максимальна, если он находится в положе-

14.

T. Olsen, J. Yan, J. J. Mortensen, and K. S. Thygesen,

нии над атомом углерода (T). Согласно результатам

Phys. Rev. Lett. 107, 156401 (2011).

ТФП вычислений [12] положение (T) энергетически

15.

X. Shi, Q. Yin, and Y. Wei, Carbon 50, 3055 (2012).

выгодно также для атомов свинца и золота. Поэто-

16.

A. V. Sidorenkov, S.V. Kolesnikov, and A. M. Saletsky,

му мы ожидаем, что потенциал (9) с другими пара-

Eur. Phys. J. B 89, 220 (2016).

метрами может быть использован для моделирова-

17.

W. G. Hoover, Phys. Rev. A 31, 3 (1985).

ния диффузии этих атомов по поверхности графена.

18.

G. J. Martyna, M. E. Tuckerman, D. J. Tobias, and

Другим важным результатом нашей работы яв-

M. L. Klein, Mol. Phys. 87, 5 (1996).

ляется оценка разностей свободных энергий ΔF для

19.

D. W. Brenner, Phys. Rev. B 42, 9458 (1990).

диффузии атома и димера меди по поверхности гра-

20.

F. Cleri and V. Rosato, Phys. Rev. B 48, 22 (1993).

фена. Показано, что в случае диффузии одиночного

21.

T. Siahaan, O. Kurnosikov, H. J. M. Swagten,

атома вклад от разности колебательных свободных

B. Koopmans, S. V. Kolesnikov, A. M. Saletsky,

энергий Δf_vib оказывается одного порядка с диффу-

and A. L. Klavsyuk, Phys. Rev. B 94, 195435 (2016).

зионным барьером E_D. Поэтому при моделировании

22.

S. A. Dokukin, S. V. Kolesnikov, A.M. Saletsky, and

диффузии величиной Δf_vib нельзя пренебрегать.

A. L. Klavsyuk, J. Alloys Compd. 763, 719 (2018).

Письма в ЖЭТФ том 111 вып. 1 - 2

2020