УДК 539.3, 539.412
Численное моделирование деформирования и разрушения графена при одноосном растяжении методом молекулярной динамики
С.П. Киселев, Е.В. Жиров
Институт теоретической и прикладной механики им. С.А. Христиановича СО РАН, Новосибирск, 630090, Россия
В работе представлены результаты численного моделирования методом молекулярной динамики деформирования и разрушения графена при одноосном растяжении. Получены зависимости модуля Юнга, критической силы и деформации разрушения от скорости деформации, температуры и угла между направлением растяжения и решеткой графена. Исследовано влияние дефектов на процесс разрушения графена.
Ключевые слова: графен, молекулярная динамика, дефекты, деформация, разрушение
Molecular dynamics simulation of deformation and fracture of graphene
under uniaxial tension
S.P. Kiselev and E.V Zhirov
Khristianovich Institute of Theoretical and Applied Mechanics SB RAS, Novosibirsk, 630090, Russia
The paper reports on molecular dynamics simulation of deformation and fracture of graphene under uniaxial tension. Dependences of Young’s modulus, critical force and fracture strain on the strain rate, temperature and angle between the tension direction and the graphene lattice are derived. The effect of defects on fracture of graphene is studied.
Keywords: graphene, molecular dynamics, defects, strain, fracture
1. Введение
Графен — монослой атомов углерода с гексагональной двумерной структурой. Его можно представить как одну плоскость графита, отделенную от объемного кристалла. Графен обнаружен в 2004 г., в настоящее время его свойства интенсивно изучаются как экспериментально, так и теоретически [1-5]. Эксперименты показали, что графен имеет уникальные электронные и механические свойства. В частности, обладает большой механической жесткостью (модуль Юнга Е — порядка 1 ТПа) и хорошей теплопроводностью (коэффициент теплопроводности — порядка 5-103 Вт/(м-К ). Высокая подвижность носителей заряда делает его перспективным материалом для использования в самых различных приложениях, в частности, как будущую основу наноэлектроники и возможную замену кремния в интегральных микросхемах.
В работе [ 1] приведены результаты эксперимента по прокалыванию круглых графеновых мембран, подве-
шенных на кремниевой подложке, с помощью атомного силового микроскопа. Предположив толщину графено-вого листа равной 0.335 нм, авторы получили значение для модуля Юнга 1.0 ± 0.1 ТПа. Кроме того, в работе найдено критическое напряжение разрушения 130 ± ±10 ГПа, что для данной толщины листа соответствует силе 42 ± 4 Н/м. Отметим, что эксперименты для графита дают близкое значение для модуля Юнга, равное 1.02 ± ± 0.03 ТПа.
В [2, 3] представлены результаты численного моделирования процессов деформации и разрушения листа графена на основе квантовомеханического расчета электронной плотности. Получены следующие характеристики графена: модуль Юнга — 1.0 ТПа, сила разрушения — 43 Н/м. Растяжение производилось последовательными смещениями краевых атомов с промежуточной релаксацией системы при нулевой температуре до момента разрушения. Методом функционала электронной плотности в работе [4] были получены следующие
© Киселев С.П., Жиров Е.В., 2012
характеристики графена: модуль Юнга — 1.05 ТПа, коэффициент Пуассона — 0.186, критические напряжения разрушения в направлении «кресло» — 110 ГПа (36.85 Н/м) и «зигзаг» — 121 ГПа (40.54 Н/м).
В работе [5] моделировался процесс деформации и разрушения графена методом молекулярной динамики с использованием пакета LAMMPS [6]. Взаимодействие между атомами углерода рассчитывалось с помощью потенциала АШЕВО [7]. Получены следующие характеристики: модуль Юнга — 1.01 ± 0.03 ТПа, коэффициент Пуассона—0.21 ± 0.01, критические напряжения разрушения в направлении «кресло» — 90 ГПа (30 Н/м), в направлении «зигзаг» — 107 ГПа (36 Н/м). Здесь растяжение также производилось последовательными смещениями краевых атомов с промежуточной релаксацией системы. Скорость деформации графена — 0.001/пс = = 109 с-1, температура — 300 К.
В упомянутых выше работах [2-5] особое внимание уделялось изучению влияния границы листа графе-на на его прочностные свойства. Были рассмотрены свободные границы графена типа «кресло» и «зигзаг», а также случай «безграничного» листа, который моделировался периодическими граничными условиями. Показано, что прочность листа графена со свободной границей «зигзаг» выше, чем у графена с границей «кресло» и приближается к прочности «безграничного» листа.
В данной работе изучается влияние скорости деформации, температуры и дефектов, расположенных на границе и поверхности листа графена, на механические свойства графена.
2. Методика расчета
2.1. Метод молекулярной динамики
Моделирование процессов деформации и разрушения графена проводилось с помощью метода молекулярной динамики. Все расчеты выполнены в пакете моделирования LAMMPS [6]. Метод молекулярной динамики заключается в численном интегрировании уравнений движения
дН . дН
Р] = дх} ’ Х = др у ’ где Н— полный гамильтониан системы; р/ — импульс; х — координата.
Для поддержания постоянной температуры графена применялись термостат Берендсена [8] и прямое масштабирование скоростей. Алгоритм Берендсена основан на введении знакопеременного трения. Отклонение температуры Тот ее равновесного значения Т0 корректируется согласно уравнению dT^)/dt = (Т0 - Т(^)/т. Известно, что термостат Берендсона плохо работает для относительно небольших систем (до сотни атомов) и на длинных траекториях из-за неравномерного распределения энергии по степеням свободы. Но в данном
(1)
случае эта проблема не возникала, так как система насчитывала порядка тысячи атомов, а траектории частиц локализованы вблизи узлов решетки. Начальная температура задавалась посредством задания тепловой скорости атомов согласно распределению Максвелла.
2.2. Выбор потенциала
Гамильтониан Н, входящий в (1), представляет собой сумму кинетической и потенциальной энергии взаимодействия атомов (потенциал). В данной работе были опробованы следующие потенциалы: МЕАМ [9], ТегеоН" [10] и АЖЕВО [7].
МЕАМ [9] является наиболее универсальным потенциалом. Имеется параметризация для значительной части таблицы Менделеева. Однако тестовые расчеты показали, что потенциал МЕАМ плохо описывает графен. Гексагональная структура в данном случае нестабильна даже при низких температурах.
Потенциал Тегеой- [10] был разработан для описания близкодействующих ковалентных связей атомов кремния и углерода. Однако тестовые расчеты показали, что он неправдоподобно описывает поведение графена при растяжении. Численный расчет, выполненный с помощью ТегеоГГ-потенциала приводит к поперечному расширению графена при его растяжении в продольном направлении (отрицательный коэффициент Пуассона).
Потенциал А1ИЕВО [7] создан для моделирования сложных процессов (включая химические реакции) с углеводородами. В рамках данного потенциала имеется также возможность моделировать взаимодействие углерода с водородом. Потенциал состоит из трех частей:
Е = ^11
ЕТ0+е" + II Е1
кг/1
(2)
где Е?-ЕВ0 = к/ + КУгА — потенциал REBO, является
у у и и
потенциалом типа ТегеоГС, описывает близкодействующие ковалентные связи в атомах углерода, где —
тгА ~ ^
отталкивающий член; V/ — притягивающий член; Ъу — коэффициент пропорциональности, зависящий от углов. Отталкивающий член в (2) выглядит следующим образом: V/ = ^ ^) [1 + Ц /гу ] 4е'у , где Ц/, 4у и аг/ зависят от типов взаимодействующих атомов i и j. Обрезающая функция (г) обеспечивает плав-
ное выключение части ИЕВО при достижении максимального радиуса действия и имеет вид [7]:
1, гу < гс
°.5(1 + с°5(п(Гу - гсуА), гс < Гу < гс + А, (3)
0, гу > гс +А,
где гс — радиус обрезания; А — ширина «размазывания» обрезающей функции м!у = (г). Притягивающий
А 3 п
член в (2) имеет вид: Vу =-м^ (г/)I Ще
п=1
в/ также зависят от типов атомов. Потенциал Еи — потенциал типа Леннарда-Джонса, который описывает
'в/г', где Ву и
дальнодеиствующее взаимодеиствие и представляется в виде:
УуА = у гу)С„ї“(г„) + (1 - ?(Гу ))С^(г„),
у'
А12
у у
(4)
где , ? и Су — функции обрезания; у/3 — классический потенциал Леннарда-Джонса. В формуле (2) Е1 описывает крутильный потенциал вокруг связей:
Е/к1 = (гу )Щк (ггк (г/1)у‘(юум X (5)
И 256
V =----------е,
405
„10
ЦЫ
cos
(Юу'й/2)"
1
10
-ЦЫ •
Численные значения в (2)-(5) для коэффициентов АШЕВО приведены в [7], радиус обрезания для части ИЕВО составлял 0.2 нм, для части LJ — 1.02 нм. Проблема выбора радиуса обрезания для части потенциала ИЕВО обсуждается ниже.
2.3. Постановка численных экспериментов
С помощью оригинального пакета были подготовлены молекулярные конфигурации пластин графена с размерами порядка 5 нм и необходимым углом кристаллической решетки. Затем атомы, расположенные на верхнем и нижнем краях (рис. 1), исключались из расчета в направлении растяжения. На атомы в этих зонах не влияли продольные силы, но сами они воздействовали на соседние атомы в расчетной области. Эти атомы могли перемещаться в поперечном направлении под действием соседних атомов, расположенных в расчетной области. Нижний край оставался неподвижным, а атомам в верхнем крае на каждом шаге задавалось постоянное малое смещение. Скорость растяжения верхнего края была порядка 10 м/с, что соответствует скорости деформации VI = 109 с-1, где V — скорость растяжения; I — характерный размер системы (в данной работе I — начальный продольный размер листа графе-на, который во всех экспериментах составлял 5 нм). Два
Рис. 1. Молекулярная конфигурация графена типа «зигзаг» перед растяжением. Стрелками указано направление растяжения
боковых края оставались свободными, либо на них задавались периодические граничные условия.
Перед каждым численным экспериментом система оптимизировалась и с помощью термостата плавно выводилась на нужную температуру. Характерная молекулярная конфигурация полностью оптимизированной системы представлена на рис. 1. Процесс растяжения продолжался до момента разрыва. Разрывом в данной работе называется момент, когда среднее напряжение в образце падает до нуля.
2.4. Расчет механических параметров графена
В расчетах определялись следующие параметры графена: модуль Юнга, критическое напряжение разрушения Е» и соответствующая деформация разрушения е». Деформация графена вдоль оси растяжения Оу рассчитывалась по формуле
е у =-
I -10
Iй
(6)
где I — текущиИ размер пластины в направлении растяжения, а 10 — начальный размер в том же направлении.
Напряжение в листе графена рассчитывалось следующим образом. Вычислялся симметричныИ тензор напряжении, приходящиИся на каждыИ атом, по формуле
(
то,о у +
1
1 (Ъ Г
\
• п=1
11* 1 у + Г2іГ2у )
(7)
где Уа — объем, занимаемый одним атомом углерода в графене. ПервыИ член в формуле (7) — это вклад кинетической энергии атома, второИ член — вклад энергии парного взаимодеИствия с Np соседями данного атома; г1 и г2 — координаты двух взаимодеИствующих атомов; Г и Г2 — соответственно силы парного взаимодеИствия, деИствующие на эти атомы. Затем рассчитывалась удельная продольная сила Г = /а на едини-
цу поперечного размера пластины по формуле
У -1 Syy, (8)
где V = abh — объем, занимаемый листом графена; а — поперечный размер; Ь — продольный размер; h — толщина графена; ак/N — поперечная площадь, приходящаяся на один атом углерода. С помощью формул (7),
(8) строились диаграммы растяжения — зависимость удельной силы от относительного удлинения.
Модуль Юнга определялся по формуле
1 аЕ1, = 1 аЕ
^ &еу к &еу ’ (9)
где Е,^ — полная продольная сила в образце; 51 = ah — площадь поперечного сечения листа графена. Сила F определялась по формуле (8), производная в формуле
(9) вычислялась при условии еу ^ 0. Толщина листа графена h = 0.335 нм выбиралась такой же, как в работе [5].
Критическая сила F» и критическая деформация разрушения е» в численных экспериментах на растяжение брались непосредственно перед разрушением гра-фена (момент разрушения t = I» определяется из условия F = 0). На диаграммах растяжения F(е) это соответствует максимальной силе, где здесь и ниже для сокращения будет опускаться индекс у компоненты деформации е = еу.
3. Результаты численных экспериментов
По описанной выше методике был произведен ряд численных экспериментов, в которых были выбраны параметры обрезания в потенциале АШБВО, исследовано влияние температуры, скорости деформации, угла между направлением растяжения и ориентации кристаллической решетки, граничных и поверхностных дефектов на разрушение графена.
3.1. Выбор радиуса обрезания в потенциале АШЕВО
Первоначально для потенциала АЖЕВО был предложен радиус обрезания потенциала в формуле (3) равный гс = 0.17 нм [7]. В работе [11] было отмечено, что в задачах, где происходит растяжение образца, необходимо увеличить радиус обрезания в формуле (3) до значения гс = 0.2 нм. Как отмечено в [12], в точке обрезания Гу = гс производная от потенциала дБ/дгу испытывает скачок, который приводит к нефизическому увеличению силы притяжения между атомами углерода Fу = —дБ/дгу . В последующих работах использовался радиус обрезания гс = 0.2 нм, однако в них отсутствует информация по влиянию ширины «размазывания» А обрезающей функции (3) на силу взаимодействия. На рис. 2 показана зависимость потенциальной энергии и силы взаимодействия между двумя атомами углерода, рассчитанная при одном значении радиуса обрезания гс = 0.2 нм, но для разных величин «размазывания»: А = = 0.001 и 0.03 нм. Видно, что при гс < г < гс + А для А = 0.03 нм возникает нефизическая сила притяжения между атомами углерода, а для А = 0.001 нм эта сила
1 2 1 3 1 4
1 \ : \ 1 1 1 \ ‘ "-Ч. ■ -"■* \
\ /
0.10 0.15 г 0.20 0.25
г, нм
Рис. 2. Зависимости потенциальной энергии Е (1, 2) и силы взаимодействия £ (3, 4) от расстояния г между атомами углерода, полученные для радиуса обрезания Гс = 0.2 нм при различной ширине «размазывания» А = 0.03 (1, 3) и 0.001 нм (2, 4). Е = 2.0-2.3 (1), 2.0-2.01 (2), £ = 2.0-2.3 (3), 2.0-2.01 (4)
отсутствует. На рис. 3 показаны зависимости удельной силы от деформации £(е), рассчитанные при растяжении графена для различных значений гсигс + А. Расчеты проводились при температуре графена 300 К, скорость растяжения графена — 100 м/с, боковые свободные границы — типа «зигзаг». Видно, что в диапазоне радиусов обрезания от 0.19 до 0.21 нм значения критической удельной силы F» и деформации е» существенно зависят от ширины «размазывания» А. Из сравнения результатов расчетов (рис. 3) с результатами рассмотренных выше работ [2-5] следует, что для потенциала АЖЕВО в формуле (3) нужно использовать радиус обрезания 0.19 нм < гс < 0.21 нм и ширину «размазывания» 0 < А < 0.001 нм.
3.2. Влияние скорости деформации и температуры на разрушение графена
На рис. 4 показаны диаграммы растяжения графена £(е), рассчитанные при различных скоростях деформа-
Рис. 3. Зависимость удельной силы от деформации £(е) при растяжении графена для различных радиусов обрезания Гс при 0 < А < 0.03 нм (а) и 0 < А < 0.001 нм (б)
Рис. 4. Зависимость £(е) при растяжении графена (периодические граничные условия на боковых границах) при различных скоростях растяжения
ции е = VЬ, где значения скоростей растяжения V указаны на рис. 4. Видно, что скорость деформации слабо влияет на критические параметры разрушения. Это связано с тем, что до момента разрушения графен деформируется как нелинейно-упругое тело. Небольшое увеличение критической силы АF*/ F» = 0.05 при увеличении скорости деформации в 103 раз связано с конечным временем разрушения графена. Из рис. 4 были получены механические параметры графена: модуль Юнга (9) Е = 0.96 ± 0.05 ТПа, критическая сила F» = = 34 Н/м, критическая деформация разрушения е» ~ 0.2. Экспериментальное значение модуля Юнга, критической силы и деформации равны соответственно [1] 1.02 ± 0.03 ТПа, 42 ± 4 Н/м и 0.25. Видно, что численные расчеты с помощью потенциала АГОБВО дают заниженное на 20 % значение критической силы разрушения, остальные рассчитанные параметры близки к экспериментальным значениям.
На рис. 5 приведены результаты расчетов деформации листа графена с границей типа «зигзаг» при различ-
ных температурах. Видно, что с увеличением температуры наблюдается монотонное уменьшение критической силы и деформации разрушения пластины графена. При температуре выше 2000 К лист графена теряет устойчивость и разрушается при сколь угодно малых нагрузках.
3.3. Влияние граничных и поверхностных дефектов на разрушение графена
В работах [2-5] было показано, что критические параметры разрушения F*, е» существенно уменьшаются при изменении типа боковой границы листа графена от «зигзага» к «креслу». В обоих этих предельных случаях лист графена состоит из шестиугольников и в нем отсутствуют дефекты. Представляет интерес рассмотреть растяжение листа графена при произвольной ориентации его границы по отношению к оси растяжения. С учетом симметрии решетки графена угол а = 0 ° будет соответствовать границе типа «зигзаг», угол а = 30° — границе типа «кресло». Для промежуточных углов граница листа графена будет содержать дефекты в виде недостроенных шестиугольников. На рис. 6, а показаны результаты расчетов £(е) для различных углов ориентации решетки относительно оси растяжения. Максимальные значения критических параметров достигаются для границы типа «зигзаг». Кругами, соединенными линией, на рис. 6, а показана зависимость £(е), рассчитанная в работе [5] при скорости деформации е = = 109 с-1 для потенциала А1ИБВО с помощью пакета LAMMPS [6]. Наблюдается хорошее совпадение между результатами нашей работы и работы [5] для границы типа «зигзаг». В отличие от [5], мы рассмотрели промежуточные типы границы листа графена между «креслом» и «зигзагом». Видно, что даже небольшой угол поворота на а = 3° относительно направления «зигзаг» приводит к значительному уменьшению критических параметров разрушения F*, е». Это связано с наличием дефектов типа незаполненных шестиугольников на границе, которые уменьшают энергию зарождения пор и
Рис. 5. Зависимости £(е) при растяжении графена (граничные условия на боковых границах типа «зигзаг») с различными температурами (а); зависимость критической удельной силы и деформации разрушения от температуры (б)
Рис. 6. Зависимость £(е) при растяжении графена с различными углами а (а); зависимость критической удельной силы F* и деформации разрушения е» от угла а (б); молекулярная конфигурация графена с углом а = 9° перед растяжением (в) и разрушением (г)
трещин на границе. В работе [13] показано, что краевая энергия монотонно уменьшается при увеличении угла от «зигзага» к «креслу». Это должно приводить к монотонному уменьшению критического напряжения и деформации разрушения графена. Данные результаты расчетов (рис. 6, б) показывают, что при увеличении угла величины Е», е» уменьшаются, однако их зависимость от угла не является монотонной функцией. Это отличие связано с тем, что в работе [13] рассматриваются хиральные углы с целыми периодами (т, и), когда граница состоит из правильных шестиугольников. В данном случае моделируется растяжение с произвольными углами, поэтому граница состоит из неполных шестиугольников (рис. 6, в, а = 9°), на которых зарождается трещина (рис. 6, г, а = 9°). Граничная энергия определяется не только числом шестиугольников на границе, но и степенью их заполнения. Для всех рассмотренных углов (за исключением а = 9°) критические значения лежат между значениями для границы типа «кресло» (а = 30°)
Е» = 27.2 Н/м, е» = 0.128 и границы типа «зигзаг» (угол а = 0°) Е» = 35 Н/м, е» = 0.2. Отметим, что в работе [5] получены близкие критические значения напряжения и деформации для границы типа «зигзаг» Е» = 36 Н/м, е» = 0.2 и для границы типа «кресло» — Е» = 30.2 Н/м, е» = 0.13.
Кроме рассмотренных выше граничных дефектов, в листе графена могут существовать поверхностные дефекты. На рис. 7 показаны рассчитанные зависимости £(е) для четырех типов дефектов D1-D4. Дефект D1 был образован путем перемещения двух атомов углерода из вершин двух шестиугольников в два соседних шестиугольника. В результате внутри отмеченного круга
Рис. 7. Зависимость £(е) при растяжении графена с различными типами дефектов
Рис. 8. Молекулярная конфигурация графена с дефектом DI (дефект отмечен кругом) перед растяжением (а); молекулярная конфигурация графена с дефектом DI перед разрушением (кругами отмечены вершины трещин) (б)
(рис. 8, а) образовались два семиугольника и два пятиугольника. Дефект D2 является порой, образованной удалением двух атомов углерода (рис. 9, а), дефект D3 является порой, образованной удалением четырех атомов углерода, дефект D4 (рис. 9, в) образован из двух пор D3. Образование дефектов производилось при нулевой температуре, затем производилось нагревание листа графена с дефектом до 300 К. Размеры листа графена
5x5 нм, граница — типа «зигзаг». При растяжении листа графена, содержащего один дефект D1, произошло значительное уменьшение критических параметров разрушения (рис. 7). Для листа графена с дефектами D2-D4 критические значения разрушения близки между собой и меньше, чем для дефекта D1. Такое поведение графена связано с тем, что наличие поры значительно сильнее уменьшает энергию зарождения трещины, чем
-2
-2
Рис. 9. Моле^лярная конфигурация графена с дефектом D2 (дефект отмечен кругом) перед растяжением (а); моле^лярная конфигурация графена с дефектом D2 перед разрушением (кругами отмечены вершины трещин) (б); молекулярная конфигурация графена с дефектом D4 (дефект отмечен двумя кругами) перед растяжением (в); моле^лярная конфигурация графена с дефектом D4 перед разрушением (кругами отмечены вершины трещин) (г)
Рис. 10. Молекулярная конфигурация графена без дефектов (кругами отмечены вершины трещин)
дефект D1. Разрушение графена с дефектами зарождается в окрестности дефектов и происходит путем распространения трещины в направлении перпендикулярном к оси растяжения (рис. 8, б, 9, б). Интересная особенность наблюдается при разрушении дефекта D4, образованного двумя порами. Поперечная трещина зарождается у одного дефекта, расположенного в центре листа графена. За счет упругой разгрузки на трещине происходит уменьшение напряжения в листе графена, поэтому второй дефект остается без изменений и не влияет на процесс разрушения. В случае бездефектного листа графена трещины зарождаются на границе и распространяются внутрь листа под углом 45° к оси растяжения (рис. 10).
4. Заключение
Приведенные выше результаты расчетов показывают, что метод молекулярной динамики с потенциалом А1ИБВО при соответствующем выборе параметров обрезания хорошо описывает поведение графена при одноосном растяжении в плоскости листа графена. Растяжение в графене происходит как в двухмерном нелинейно-упругом теле. Разрушение графена является хрупким, по этой причине критическое напряжение и деформация разрушения слабо зависят от скорости растяжения в широком диапазоне. Увеличение температуры приводит к монотонному уменьшению критического напряжения и деформации. Полученные критические значения напряжения и деформации для «бездефект-
ного» графена зависят от структуры свободной поверхности листа графена («кресло», «зигзаг») и согласуются с известными литературными данными. Очаги разрушения «бездефектного» графена зарождаются на границе и распространяются внутрь листа графена. Поверхностные и граничные дефекты в графене значительно уменьшают критическое напряжение и деформацию разрушения. Очаги разрушения «дефектного» графена зарождаются на дефектах и распространяются в поперечном направлении по отношению к оси растяжения. Поскольку реальные образцы графена содержат дефекты, полученные критические значения напряжения и деформации для «дефектного» графена могут быть полезны при исследовании графена.
Работа выполнена при частичной финансовой поддержке гранта РФФИ № 11-01-00357.
Литература
1. Lee C., Wei X., Kysar J.W., Hone J. Measurements of the elastic properties and intrinsic strength of monolayer graphene // Science. -2008.- V. 321. - P. 385-388.
2. Яновский Ю.Г., Никитин Е.А., Карнет Ю.Н., Никитин С.М. Кван-
тово-механическое исследование механизма деформации и разрушения графена // Физ. мезомех. - 2009. - Т. 12. - № 4. - С. 61-70.
3. Яновский Ю.Г., Никитина Е.А., Карнет Ю.Н., Никитин С.М. Моделирование деформации и разрушения графена: размерный эффект, влияние дефектов и модификации поверхности // Физ. мезомех. - 2010. - Т. 13. - № 5. - С. 139-147.
4. Liu F., Ming P., Li J. Ab initio calculation of ideal strength and phonon
instability of graphene under tension // Phys. Rev. B. - 2007. - V. 76.-P. 064120.
5. Zhao H., Min K., Aluru N.R. Size and chirality dependent elastic properties of graphene nanoribbons under uniaxial tension // Nano Letters. - 2009. - V. 9. - No. 8. - P. 3012-3015.
6. Plimpton S.J. Fast parallel algorithms for short-range molecular dynamics // J. Comput. Phys. - 1995. - V. 117. - P. 1-9.
7. Stuart S.J., Tutein A.B., Harrison J.A. Areactive potential for hydrocarbon with intermolecular interactions // J. Chem. Phys. - 2000. -V. 112. - P. 6472-6486.
8. Berendsen H.J.C., Postma J.P.M., van Gunsteren WF, DiNola A., Haak J.R. Molecular-dynamics with coupling to an external bath // J. Chem. Phys. - 1984. - V. 81. - No. 8. - P. 3684-3690.
9. Baskes M. Modified embedded-atom potentials for cubic materials and empurities // Phys. Rev. B. - 1992. - V. 46. - P. 2727-2742.
10. Tersoff J. Modeling solid-state chemistry: Interatomic potentials for multicomponent systems // Phys. Rev. B. - 1989. - V. 39. - No. 8. -P. 5566-5568.
11. Shenderova O.A., Brenner D.W., Omeltchenko A., Su X., Yang L.H. Atomistic modeling of the fracture of polycrystalline diamond // Phys. Rev. B. - 2000. - V. 61. - No. 6. - P. 3877-3888.
12. Belytscko T., Xiao S.P., Scatz G.C., Su X., Ruoff R.S. Atomistic simulations of nanotube fracture // Phys. Rev. B. - 2002. - V. 65. -P. 235430.
13. Gan C.K., Srolovitz D.J. First principles study of grapheme edge properties and flake shapes // Phys. Rev. B. - 2010. - V. 81. -P. 125440.
Поступила в редакцию 27.07.2011 г.
Сведения об авторах
Киселев Сергей Петрович, д.ф.-м.н., внс ИТПМ СО РАН, [email protected] Жиров Евгений Викторович, магистрант НГУ, [email protected]