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

ВВЕДЕНИЕ

Трение играет центральную роль в разнообразных системах и явлениях которые на первый взгляд могут показаться не связанными друг с другом, но у которых при более внимательном осмотре обнаруживаются общее особенности, которые свойственны всем трибологическим процессам в областях технологии, геологии и биологии. Разработка износостойких поверхностей с низким трением и тонких смазочных пленок стала важным фактором в миниатюризации подвижных частей во многих технических приборах. Последние включают микроэлектромеханические системы (МЭМС), компьютерные запоминающие системы, миниатюрные двигатели и др. Старые, простые эмпирические законы трения не всегда работают в таких системах. Это происходит потому, что они характеризуются высоким отношением поверхности к объему, и для них более важна химия поверхностей, адгезия и структура (или шероховатость) последних. Обычные трибологические методы, что используются для макроскопических объектов, могут оказаться не эффективными на нанометровом масштабе, что требует новых подходов для контролирования изучения подобных систем. Особенно быстро развевается область трибологии, что изучает биосистемы, и также смазочные механизмы в суставах.

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

Первым шагом к решению данной задачи есть изучение трибологических свойств пленок воды. Так как реальные контакты макроскопических тел происходят на большом количестве микрои нанонеровностях и шероховатостей поверхностей, для пояснения макроскопического трения необходимо понять поведение одного наноконтакта. Экспериментально это стало возможным в 80-х годах ХХ столетия одновременно с изобретением сканирующей зондовой микроскопии и аппарата поверхностных сил (АСП). Исследования ультратонких (толщиной 1−6 молекулярных диаметров) пленок разнообразных жидкостей с использованием АПС определили в них ряд необычных свойств, обусловленными главным образом малыми размерами систем. Кроме того, было установлено, что пленка толщиной один молекулярный диаметр, сжатая между гидрофильными поверхностями слюды, имеет очень хорошие смазочные свойства.

С целью поиска теоретической модели, которая смогла бы дать возможность объяснить трибологические свойства тонкой пленки воды, в работе проведены компьютерные эксперименты с использованием метода молекулярной динамики. Исследовалось поведение молекул воды, сжатых между атомарно-гладкими абсолютно твердыми алмазными поверхностями. Пленка состоит с одного или двух слоев молекул, для которых использована модель TIP4P. В зависимости от приложенной к поверхностям нагрузке, силы сдвига и количества шаров жидкости исследовано равновесные и динамические свойства пленки. Полученные результаты указывают на способность модели отобразить общие свойства, принадлежащие ультратонким пленкам жидкости с простыми молекулами, сжатыми между атомарно-гладкими твердыми поверхностями.

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

1. СВОЙСТВА ТОНКИХ ПЛЕНОК ЖИДКОСТЕЙ, СЖАТЫХ МЕЖДУ ТВЕРДЫМИ ПОВЕРХНОСТЯМИ

1.1 История исследований и экспериментальные инструменты Трение является одной из самых старых проблем в физике с огромной практической значимостью. Оно издавна признавалось ограничивающим фактором для многочисленных практических применений и для многих технологических областей. Поэтому всегда существовал значительный интерес в понимании и контролировании этого процесса. Например, древние египтяне изобрели новые технологии передвижения камней, использовавшихся для построения пирамид; французский ученый Кулон изучал трение из-за необходимости в легком перемещении кораблей с суши на воду; исследования Джонсона автомобильных стеклоочистителей улучшили понимание контактной механики и поверхностных энергий. На сегодняшний день интерес и развитие трибологии — науки о трении — сфокусированы на микрои наномасштабные машины с подвижными частями, что продолжает требовать понимания принципов, лежащих в основе трения. Это мотивировало исследователей изучать данное явление на очень малых масштабах, определять его фундаментальные причины, и обусловило появление нового направления трибологии, названного нанотрибологией. Исследования трения на атомарном уровне также проложили путь для некоторых инноваций, например, самосмазывающихся поверхностей и стойких к износу материалов. Известно, что механизмы трения на наноуровне иногда полностью отличаются от механизмов, доминирующих на макромасштабах. Это имеет непосредственное значение для устройств с нанометровыми режимами работы, например, магнитных запоминающих дисков, размеры которых постепенно уменьшились на протяжении нескольких последних лет, и микроэлектромеханических систем, где трение на атомарном уровне, адгезия и износ являются доминирующими процессам.

Научные исследования трения на атомарном масштабе стали проводиться с конца 1980;х годов. Это произошло благодаря одновременному развитию новых экспериментальных инструментов для измерения трения на нанометровых расстояниях при небольших нагрузках, быстрому росту вычислительных мощностей компьютеров, а также дозреванию теоретических методов, необходимых для изучения процессов в материалах реалистическими способами. Например, аппарат поверхностных сил (АПС) и зондовая сканирующая микроскопия обеспечили исследователей новой информацией о трении и смазывании многих жидких и твердых систем с необыкновенно высокой разрешающей способностью. С помощью современной экспериментальной аппаратуры можно исследовать поверхности трения на атомарном уровне и, в конце концов, связать наблюдаемое поведение с явлениями, имеющими место на макроскопическом уровне.

Экспериментальное исследование граничных смазок проводится с помощью АПС, модифицированного для экспериментов по трению.

Схематически он представлен на рисунке 1.1. Основная часть АПС состоит из двух молекулярно-плоских поверхностей слюды каждая толщиной около 2 мкм, покрытых полупрозрачным 50 — 60 нм слоем серебра и приклеенных к цилиндрическим кремниевым дискам радиуса 10 мм. Оси цилиндров взаимно перпендикулярны, используется схема перекрещивающихся цилиндров. Поверхности слюды находятся в герметичной камере, заполненной исследуемой жидкостью или газом. Расстояние между поверхностями контролируется с точностью до 0.1 нм (для измерений используется оптическая интерференция полос равного хроматического порядка), силы — с точностью до 10 нН.

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

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

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

1.2 Основные трибологические понятия и законы Прежде всего, необходимо отличать два разных явления: статическое трение и кинетическое трение. Статическая сила трения определяется как минимальная сила, необходимая для начала скольжения. Ее значение зависит от атомной структуры границы раздела поверхностей и адгезионных взаимодействий. Чтобы начать скольжение, необходимо или разорвать межатомные связи, или инициировать пластическое течение на границе жидкость — поверхность. Очевидно, что этот процесс будет происходить сначала в некоторых «слабых» местах.

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

Как статическое, так и кинетическое трение являются очень важными для практических применений, и в различных ситуациях требуется как высокое, так и низкое их значение. Без статического трения мы не могли бы ходить, а автомобили передвигаться. Высокое статическое трение необходимо для стабильного удержания механических конструкций, соединенных гайками и болтами. Низкое статическое трение требуется в подвижных частях машин, например, автомобильных двигателях, и наилучшим решением было бы нулевое значение с жидкими смазками.

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

Выделяют два различных режима — гидродинамическое (жидкостное) трение и граничное смазывание (или граничное трение). В первом случае поверхности разделены толстой (приблизительно более 0.01 мм) жидкой пленкой смазки. Физическая задача в этом случае сводится к решению гидродинамического уравнения Навье — Стокса с соответствующими граничными условиями и геометрией контакта, а кинетическое трение определяется главным образом вязкостью жидкой смазки. Во втором случае поверхности разделены тонкой пленкой смазки (толщиной несколько атомарных диаметров). Случай «сухого трения» также относится к этому классу. Граничное трение, очевидно, является наиболее важным в микромашинах. Однако, даже в макромашинах, где, как правило, действует гидродинамическое трение, граничное смазывание также является важным в моменты остановки/запуска, когда смазка выдавливается из зоны контакта, и поверхности переходят в непосредственное соприкосновение.

Опишем трибологические законы и свойства граничной смазки, хорошо установленные экспериментально. Прежде всего, два основных закона Амонтонса, датированные 1699 годом, утверждают:

— сила трения прямо пропорциональна нагрузке. Этот закон также известен как закон трения Кулона;

— трибологический коэффициент трения не зависит от площади контакта и нагрузки и имеет значение порядка единицы.

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

Для граничного трения, когда поверхности разделены тонкой пленкой смазки толщиной несколько нанометров, экспериментально измеряемое трение намного превышает, например в 10², или даже в 10⁵ — 10⁷ раз гидродинамическое трение. Оно зависит в большей мере не от вязкости смазки, а от ее химического состава. Как правило, хороший смазочный материал — это вещество, которое адсорбируется на твердых поверхностях, поскольку это препятствует выдавливанию смазки из зоны контакта.

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

Рисунок 1.2 — Экспериментально полученные временные зависимости силы трения в режиме прилипания — скольжения для поверхностей слюды, разделенных слева — молекулами поверхностно-активных веществ (для различных скоростей скольжения), рисунок справа — простыми сферическими молекулами (для различного числа слоев молекул) Режим прилипания — скольжения является признаком наличия твердоподобного состояния смазки и считается одной из основных причин износа трущихся поверхностей. Рассмотрев общие трибологические понятия и законы, перейдем к более детальному описанию свойств ультратонких пленок жидкостей.

1.3 Свойства тонких пленок жидкостей, сжатых между твердыми поверхностями Поведение жидкостей, ограниченных в очень малых объемах, например, узких порах или тонких пленках, в общем случае полностью отличается от их объемного поведения. Как теоретически, так и экспериментально выяснено, что в пленках, толщина которых меньше десяти молекулярных диаметров, молекулы становятся более упорядоченными, и их подвижность резко уменьшается [1 — 4].

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

Что касается динамических характеристик, то установлено, что когда расстояние между поверхностями превышает приблизительно десять молекулярных диаметров, простая жидкость в промежутке между ними сохраняет свое ньютоновское поведение (т. е. наблюдается гидродинамическое трение). Однако «эффективная» вязкость более тонких пленок значительно возрастает. В экспериментах с простыми жидкостями, например, циклогексаном, найдено, что когда две молекулярно-гладкие поверхности слюды, между которыми находится только один или два слоя молекул жидкости, скользят друг по другу (течение Куэтта), «эффективная» вязкость превышает объемное значение на 5 — 7 порядков. Более того, выяснено, что вообще вся концепция ньютоновской вязкости перестает быть справедливой для таких тонких пленок: например, у них появляются пределы текучести, и сдвиговые напряжения больше не зависят от скорости сдвига.

Также во время скольжения, когда расстояние между пластинами приблизительно совпадает с минимумами энергии осциллирующих силовых кривых, поверхности разделены целым числом слоев молекул жидкости, и сдвиговые напряжения и другие трибологические свойства являются «квантованными» по числу слоев. Свойства «жидких» пленок толщиной только один — два молекулярных диаметра оказываются больше твердоподобными, чем жидкоподобными.

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

Для объяснения описанных свойств ультратонких пленок смазки широко применяется компьютерное моделирование. Рассмотрим некоторые результаты, полученные с его помощью.

1.4 Результаты компьютерных экспериментов Моделирования методами Монте Карло и МД, примененные для исследования сдвига в молекулярно-тонких пленках жидкостей, также указывают на их твердоподобные свойства [1, 4, 9, 18 — 21]. Результаты моделирования двух пластин, разделенных леннард — джонсовскими (ЛД) жидкостями, показывают, что для пленок, толщина которых менее 6 — 10 молекулярных диаметров, молекулы между стенками упорядочиваются в квазидискретные слои, при этом коэффициент диффузии уменьшается, а вязкость возрастает. Для атомарно-структурированных (т. е. шероховатых) поверхностей предсказано существование резкого фазового перехода жидкость — твердое тело в пленках толщиной менее 6 молекулярных диаметров. При этом молекулы становятся упорядоченными как перпендикулярно, так и параллельно пластинам. Однако горизонтальное упорядочение исчезает для неструктурированных, т. е. математически гладких стенок. Как и в экспериментах, были выявлены сильно выраженные квантованные пределы текучести для шероховатых поверхностей. Они также исчезали для неструктурированных пластин. Моделирования ясно указывают на большую важность атомарной структуры или «гранулярности» реальных поверхностей при изучении пленок, которые сами имеют такие же атомарные размеры.

В работе проведено исследование методом МД ультратонкой пленки ЛД жидкости между твердыми гладкими поверхностями. Наблюдался режим прилипания — скольжения, причину которого авторы усматривают в переходах пленки между твердоподобным и жидкоподобным состояниями. МД моделирование, проведенное Перссоном с сотрудниками, показывает механизм, по которому этот переход происходит. А именно, в случае скольжения диэлектрических кристаллических поверхностей твердоподобная смазка может находиться в «суперсмазочном» состоянии, в котором трение становится пренебрежимо малым. Однако из результатов моделирования ясно, что даже малая концентрация поверхностных дефектов будет переводить смазку из твердоподобного в жидкоподобное состояние. Когда скольжение происходит на металлических поверхностях при температурах, превышающих криогенную, имеется ненулевой электронный вклад в трение, и суперсмазочное состояние оказывается невозможным.

Если молекулы жидкости не являются сферическими, то им сложнее упорядочиваться и затвердевать. МД моделирования показывают, что сферические молекулы имеют большие критические скорости, чем разветвленные молекулы. Например, МД моделирования, проведенные Томпсоном и Роббинсом [18, 19] показывают, что когда молекулы разветвленные, не сила трения изменяется со скоростью скольжения, а изменяется количество времени, на протяжении которого различные части системы находятся в режиме прилипания или скольжения.

Таким образом, и реальные, и компьютерные эксперименты указывают на необычные свойства граничных смазок, которые кратко можно сформулировать следующим образом:

— молекулы жидкостей, сжатых между двумя атомарно-плоскими твердыми поверхностями, становятся более упорядоченными и стремятся образовать слоистые структуры. Плотность жидкости в ультратонкой пленке неоднородна и имеет осциллирующий профиль. Период осцилляций близок к величине диаметра молекул и отражает вынужденное упорядочение молекул в квазидискретные слои;

— подвижность молекул в ультратонких пленках значительно снижается по сравнению с объемными жидкостями. Это проявляется в уменьшении коэффициента диффузии и в увеличении сдвиговой вязкости и молекулярных времен релаксации;

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

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

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

1.5 Уникальные свойства воды Вода — необычное вещество. Ее свойства занимают особое место среди явлений, связанных с межмолекулярными силами. Два типа «специальных взаимодействий» — водородная связь и гидрофобный эффект — присущи именно молекулам воды.

Для жидкости с низкой молекулярной массы вода имеет неожиданно высокие температуры плавления и кипения и удельную теплоту испарения. Конечно, существует много других веществ с низкой молекулярной массой и высокими соответствующими температурами, но все они без исключения ионные кристаллы или металлы, атомы которых взаимодействуют друг с другом сильными кулоновскими или металлическими связями. Максимум плотности при 4 ^оС и необычное явление, что твердая вода (лед) легче жидкой, указывают на то, что молекулы льда предпочитают находиться на большем расстоянии друг от друга, чем в жидкости. Можно сделать вывод, что в жидком состоянии сохраняются сильные межмолекулярные связи, образующиеся у льда. Они должны сильно зависеть от взаимной ориентации молекул, поскольку вода выбирает тетраэдрическую координацию (четыре ближайших соседа на молекулу) вместо плотно упакованной (12 ближайших соседа на молекулу). Вода имеет и другие особенности, например, очень низкую сжимаемость и исключительные свойства растворимости как относительно растворимого, так и растворителя.

Перечисленные необычные свойства воды указывают на то, что некоторые очень сильные и зависимые от ориентации связи включены во взаимодействие между ее молекулами. Эти связи называются водородными. Они играют особенно важную роль в воде, поскольку каждый атом кислорода с двумя атомами водорода данной молекулы воды могут участвовать в четырех таких соединениях с другими молекулами — два включают ее собственные атомы водорода, и два — неспаренные электроны с атомами водорода других молекул. Необходимо отметить, что тетраэдрическая координация молекул воды является сердцем ее необычного поведения, более чем сами по себе водородные связи.

Вода также имеет необычные смазочные свойства и, как правило, дает очень низкое трение без изнашивания и режима прилипания — скольжения [4, 11]. Это проиллюстрировано на рисунке 1.3 для двух необработанных поверхностей, скользящих с мономолекулярным (толщиной 0.25 нм) слоем воды между ними.

Рисунок 1.3 — Зависимости силы трения от нагрузки для двух поверхностей слюды, скользящих друг по другу и одновременно погруженных в 1% водный раствор соли КCl (неадгезионные условия) Пленка воды является атомарно-тонкой с толщиной от 0.25 до 0.5 нм, межповерхностная сила очень мала, и коэффициент трения (до разрушения поверхности). После того как поверхности начинают разрушаться, коэффициент трения приобретает значение около 0.3

Эффективность пленки воды толщиной 0.25 нм для снижения силы трения была ассоциирована с «гидрофильностью» поверхностей слюды (которые «смачиваются» водой) и с существованием в водных растворах между такими поверхностями сильно отталкивающих короткодействующих сил гидратации, которые эффективно устраняют вклад в силу трения, контролируемый адгезией. Также интересным является то, что наличие пленки воды толщиной 0.25 нм между двумя поверхностями слюды достаточно для того, чтобы коэффициент трения упал до 0.01 — 0.02 — значения, отвечающего необычно низкому трению льда. Очевидно, один слой молекул воды может быть очень хорошей смазкой — намного лучшей, чем большинство других мономолекулярных жидких пленок.

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

2. ИССЛЕДУЕМАЯ СИСТЕМА

2.1 Основные принципы МД

2.1.1 Общие понятия Теоретические модели и моделирование, еще называемое компьютерным экспериментом, помогают в интерпретации экспериментальных данных и обеспечивают предсказание явлений, которые со временем могут быть подтверждены или опровергнуты экспериментами. Среди этих методов можно выделить аналитические модели и крупномасштабные моделирования методом молекулярной динамики (МД).

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

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

Однако эффективное использование МД моделирования для изучения трения на атомарном масштабе требует понимания многих деталей, не являющихся очевидными в этом простом анализе. Рассмотрим некоторые из них.

2.1.2 Межатомные потенциалы взаимодействия Как отмечалось, МД моделирования — это компьютерные программы, в которых ньютоновские уравнения движения интегрируются для отслеживания движения атомов во времени в ответ на приложенные силы:

(2.1а)

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

Преимущество этого метода в том, что он дает возможность следить за индивидуальным движением всех атомов в данной системе в режиме реального времени. Недостатком является то, что временные масштабы очень ограничены (от пикодо наносекунд). К тому же размеры системы, которые можно рассматривать на данный момент, ограничены 10⁶ — 10⁸ атомов, которое хотя и впечатляет, однако все еще далеко отстоит от реальных систем, содержащих 10²³ атомов и более. Следовательно, хотя МД моделирования используются для изучения трения на атомарном масштабе и сделали хорошую работу для обеспечения понимания этого процесса, однако они все еще ограничены пространственными и временными масштабами, которые значительно меньше экспериментальных значений.

Для атомистического моделирования материала необходимы математические выражения его потенциальной энергии (в (2.1б)). На сегодня существуют два подхода касательно выбора вида потенциала. В первом предполагается, что потенциальная энергия атомов может быть представлена как функция только их относительных атомных положений. Эти функции, называемые эмпирическими потенциалами, обычно основываются на упрощенных трактовках общих принципов квантовой механики и, как правило, содержат некоторое число свободных параметров. Последние выбираются таким образом, чтобы как можно лучше воспроизводился набор физических свойств интересуемой системы. Однако существует много трудностей, связанных с поиском пригодных функций потенциальной энергии. Например, параметры, используемые для подгонки потенциальной энергии, определяются из ограниченного набора известных свойств системы. Для металлов свойства, под которые подгоняется функция потенциальной энергии, могут включать постоянную решетки, когезионную энергию, постоянные упругости и энергию образования вакансий. Следствием этого является то, что другие свойства, включая и те, которые могут быть решающими для определения результата данного моделирования, определяются единственно предполагаемой функциональной формой. Предсказанные свойства для металла могут включать реконструкции поверхности, энергетику дефектов и реакции (эластическую и пластическую) на приложенную нагрузку. Форма потенциала, таким образом, является ключевой, если от моделирования требуется способность адекватно отражать физику системы.

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

Поскольку вычисленные силы основываются на квантовых принципах, вид межатомных взаимодействий может быть вполне переносимым между различными атомными средами. Однако это не обязательно гарантирует, что силы из вычислений полуэмрирической электронной структуры являются точными. Неудачно выбранные параметры и функциональные формы могут также давать неверные результаты. Недостатком этого подхода является то, что используемые потенциалы требуют значительно больших вычислительных затрат, чем используемые в эмпирических потенциальных функциях. Большие времена моделирования обуславливают меньшие доступные размеры и временные масштабы исследуемых систем, чем при использовании эмпирических потенциалов. Например, самые большие системы, которые на данный момент могут быть исследованы с помощью первых принципов и полуэмпирических МД моделирований состоят из нескольких сотен и нескольких тысяч атомов соответственно. Таким образом, данные методы еще не получили распространения в крупномасштабных МД моделированиях.

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

(2.2)

где индексы и являются метками атомов, — скалярное расстояние между ними, — предполагаемая функциональная форма для энергии. Некоторые традиционные формы для потенциала даются выражением

(2.3)

где параметр определяет минимум энергии для пар атомов. Две общие формы этого выражения представлены потенциалами Морсе () и Леннард-Джонса (ЛД) «12 — 6». Здесь и — произвольные параметры, используемые для подгонки потенциала под наблюдаемые свойства. Короткодействующая экспоненциальная форма функции Морсе обеспечивает разумное описание отталкивающих сил между атомными ядрами, в то время как член в ЛД потенциале описывает основной член дальнодействующих дисперсионных сил.

2.1.3 Термостаты Обычная МД отличается от большинства экспериментальных исследований тем, что в моделированиях фиксированными являются энергия и объем, а не температура и давление. В терминах статистической механики обычная МД дает средние по микроканоническому ансамблю (- количество молекул), в то время как эксперименты с постоянной температурой соответствуют каноническому ансамблю .

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

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

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

Более утонченной техникой для поддержания температуры системы является динамика Ланжевена. Этот метод изначально использовался для описания броуновского движения и получил широкое распространение в МД моделированиях. В этой технике дополнительные члены, соответствующие слагаемому, отвечающему трению, и случайной силе, добавляются в уравнения движения. Уравнение движения для атомов, подчиненных термостату Ланжевена, имеет вид:

(2.4)

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

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

Случайная сила в (2.4) обычно дается гауссовским распределением, в котором выбор ширины учитывает удовлетворение флуктуационно-диссипационной теоремы. Она определяется уравнением:

(2.5)

где функция — это случайная сила в (2.4), — масса частицы, — требуемая температура, — постоянная Больцмана, — время, и — коэффициент трения. Отметим, что случайные силы не зависят от сил на предыдущих шагах (что отмечено дельта-функцией), и ширина гауссовского распределения, из которого получается случайная сила, изменяется с температурой.

Ланжевеновский подход, представленный выше, не требует никакой обратной связи от текущей температуры системы. Случайные силы определяются единственно из выражения (2.5).

2.2 Общая схема компьютерного эксперимента

2.2.1 Моделирование воды С целью поиска модели трибологической системы, которая объясняла бы свойства ультратонкой пленки воды, были проведены компьютерные эксперименты, описываемые в данном разделе. Рассмотрим их общую схему.

Моделирования проводились для плоской геометрии Куэтта, напоминающей экспериментальные системы, исследуемые с помощью АПС. Рассматривалась ультратонкая пленка воды, ограниченная двумя твердыми стенками с периодическими граничными условиями в плоскости пластин, примененными как к жидкости, так и к поверхностям. Изучались атомарно-плоские и атомарно-шероховатые пластины. В первом случае каждая стенка состоит из 1152 атомов углерода, составляющих две (001) кристаллические поверхности с алмазной решеткой. Во втором шероховатость (периодический рельеф) обеспечивалась удалением 144 атомов, находящихся непосредственно на поверхности пластин. Все значения параметров моделирования поверхностей с различным рельефом одинаковы, за исключением того, что в случае шероховатых поверхностей введено слабое взаимодействие поверхностей друг с другом. Хотя изучение трения алмазных поверхностей, разделенных ультратонкой пленкой воды, еще не проводилось экспериментально, однако алмаз прозрачен для света и может быть использован в исследованиях с АПС. Для упрощения задачи пластины считаются абсолютно жесткими, и модель не включает упругости поверхностей. Но, учитывая то, что алмаз является одним из самых твердых материалов, в данной работе решено проверить это приближение. Исследовались пленки воды толщиной в один и два молекулярных диаметра. Один слой молекул в начале моделирований содержит молекул. Максимальное количество молекул, использовавшееся в моделированиях, составляло .

Рабочая программа содержит код, выполняющий вычисления, и интерфейс, позволяющий непосредственно наблюдать анимацию трехмерного изображения системы в режиме реального времени, и написана в среде Microsoft Visual Studio 2005 с использованием компилятора Microsoft Visual C++ 8.0.

На рисунке 2.2 представлены начальные конфигурации системы для пленок толщиной в один и два молекулярных слоя между гладкими поверхностями. На рисунке 2.3 -вид «сверху» (в отрицательном направлении оси z) начальной конфигурации молекул воды и кристаллической поверхности.

Рисунок 2.2 — Конфигурации исследуемой системы для пленки воды из одного и двух слоев в начале моделирования, сжатой между атомарно плоскими поверхностями Атомы углерода, кислорода и водорода изображены соответственно синими, красными и зелеными шариками Рисунок 2.3 — Расположение молекул воды (слева) и вид алмазной поверхности (справа) в плоскости ху в начале моделирований Для визуализации атомов использовались эмпирические значения их ковалентных радиусов, которые для углерода, кислорода и водорода соответственно составляют 77, 73 і 37 пм (1 пм = 10^-12 м), постоянная решетки алмаза равна 356.68 пм. Молекулы воды в начале моделирования размещались в вершинах кубической решетки с постоянной решетки, отвечающей значению 1048 кг/м³ плотности воды, являющимся близким к объемному (при нормальных условиях) значению 1000 кг/м³. При этом молекулярная ориентация принимает случайные значения, угловая и линейная скорости каждой молекулы имеют фиксированные величины, зависимые от температуры, и случайно выбранное направление. Пластины и пленка воды несоизмеримы, поскольку постоянная решетки алмаза не равна расстоянию между молекулами воды. Начальное расстояние (щель) между поверхностями для одного слоя составляет 0.3104 нм, для двух слоев — 0.6207 нм. Горизонтальные размеры ячейки моделирования вдоль осей х и у одинаковы и составляют 42.806 .

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

Температура составляла 25^оС, а нагрузка моделировалась приложением постоянной силы к каждому атому плоскостей вдоль оси z в направлении молекул воды. Сдвиг также имитировался приложением постоянной горизонтальной силы к каждому атому верхней поверхности вдоль оси х, соответствующей направлению [010], и в данной работе использовался алгоритм постоянной сдвигающей силы.

2.2.2 Моделирование аргона Моделирования аргона проводились только для атомарно-гладких поверхностей. Параметры алмазных пластин для аргона полностью аналогичны параметрам при моделировании воды. Также исследовались пленки в начальный момент содержащие два и один молекулярных слоя, в каждом из слоев по 196 молекул. Начальные конфигурации системы представлены на рисунке 2.4.

Для визуализации молекулы аргона использовалось расчетное значение ее радиуса, равное 71 пм. Так же, как и для воды молекулы аргона в начале моделирования размещались в вершинах кубической решетки с постоянной решетки, отвечающей значению 942 кг/м³ плотности жидкого аргона. При этом скорость каждой молекулы имеет фиксированные величины, зависимые от температуры, и случайно выбранное направление. Температура составляет 120 К. Внешняя нагрузка и сдвиг имитировались абсолютно так же, как и при моделировании воды.

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

2.3 Модели молекул исследуемых жидкостей

2.3.1 Молекула воды Для объяснения свойств воды необходима картина распределения заряда в ее молекуле. Были предложены разнообразные модели, например, ST2, TIP3P и др., но до сих пор еще не существует единой модели, которая была бы способной удовлетворительно учесть свойства воды во всех трех фазах (лед, жидкость и пар). Основным требованием к любой модели воды является то, что она должна воспроизводить тетраэдрическую структуру, обусловленную водородными связями.

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

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

Рисунок 2.6 — Схематическое изображение TIP4P молекулы воды (слева), вид молекулы при визуализации в программе (справа) Расстояния и угол, необходимые для полного определения координат сайтов, имеют значения:, ,. Энергия взаимодействия между двумя молекулами и состоит из двойной суммы по всем сайтам обеих молекул. Члены с индексами и учитывают кулоновское взаимодействие между электрическими зарядами, связанными с сайтами, а также вклады ЛД типа:

(2.6)

где,. Соответствующая сила определяется выражением:

(2.7)

Заряды имеют значения:, ,, где Кл. Как часть молекулярного проектирования отрицательный заряд из сайта О был смещен на небольшую величину в сайт М, введенный именно с этой целью. Параметры ЛД части потенциала, действующей только между О сайтами, имеют значения: (ккал/моль), (ккал/моль), и 1 ккал/моль = 4184 Дж/моль.

Перейдем к безразмерным МД единицам измерения, удобным для решения задачи. Определим единицу измерения расстояния как имеющую значение, удовлетворяющее равенству, откуда. Единица измерения энергии. За единицу измерения массы выбираем массу молекулы воды г. Физические и безразмерные единицы измерения тогда связаны следующим образом:, ккал/моль, и единица измерения времени равна с. Также определим безразмерные единицы измерения заряда, в которых. Поскольку единица измерения энергии соответствует К, то типичная температура 298К равна 3.8 в безразмерном виде. Использовавшийся временной шаг имел значение, в размерных единицах это составляет с.

Приведем безразмерные и размерные значения сил, использовавшихся в моделированиях. Значение силы на каждый атом изменялось от 2 (в системе СИ это 6.838 пН) до 50 (0.171 нН), что соответствует значениям давления на каждую поверхность соответственно от 0.43 ГПа до 10.751 ГПа. Учитывая, что порог пластичности алмаза порядка 100 ГПа [1], для рассматриваемых значений давления можно считать приемлемым принятое в модели приближение абсолютной жесткости стенок. Горизонтальная сдвигающая сила на каждый атом верхней пластины изменялась от 0.5 (1.71 пН) до 10 (34.2 нН), или полная горизонтальная сила на поверхность изменялась от 1.97нН до 394 нН.

Для взаимодействия вода — алмаз учитывались взаимодействия атомов поверхностей только с сайтом О молекулы воды. Потенциал имеет ЛД вид

(2.8)

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

2.3.2 Молекула аргона Молекулы аргона моделировались упругими шарами, взаимодействие для которых имеет ЛД вид

(2.9)

где первое равенство выполняется при, а второе — при; (безразмерные единицы измерения те же, что и для воды). Параметры в (2.9) имеют значения, , масса атома аргона 2.24. Взаимодействие молекул аргона с атомами поверхностей считалось таким же, как и для воды. Значения силы сдвига изменялись в тех же пределах, что и для воды. Значение силы на каждый атом изменялось от 20 (в системе СИ это 68.38 пН) до 300 (1.026 нН), что соответствует значениям давления на каждую поверхность соответственно от 4.3 ГПа до 64.506 ГПа.

2.4 Уравнения движения В моделированиях использовались классические уравнения движения. Вследствие абсолютной жесткости поверхностей рассматривалось движение их центра масс. Для верхней стенки уравнения имеют вид:

(2.10)

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

Уравнения поступательного движения (без учета слагаемого, отвечающего за поддержание постоянной температуры) для центра массй молекулы имеют вид:

(2.11)

где — сила, действующая со стороны всех атомов углерода, находящихся от данной молекулы на расстоянии, меньшем, и определяется аналогично как и сила для поверхностей; - сила, действующая наю молекулу со стороны других молекул, находящихся от данной не далее расстояния; определяется по формуле (2.7).

В модели также учитывается вращательное движение молекул воды. Для описания вращения используются гамильтоновы кватернионы. Угловая скорость связана с производной по времени от компонент, ,, кватерниона соотношением

(2.12)

где (2.13)

Уравнения Эйлера, описывающие вращение твердого тела (без слагаемого, отвечающего за постоянство температуры), имеют вид:

(2.14)

и аналогично для остальных компонент. В (2.14) — момент сил, действующих на тело, — главные моменты инерции молекулы,, ,, массы и ассоциируются с сайтами О и Н,. Уравнение движения через кватернионы имеет вид

(2.15)

Можно исключить из правой части (2.15), используя уравнения Эйлера (2.14). После замены компонент вектора линейными комбинациями величин из (2.12) получаем набор уравнений движения, выраженных полностью через кватернионы и их производные.

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

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

2.5 Контроль температуры Одним из наиболее распространенных способов контролирования температуры при моделировании граничной смазки является применение термостата Ланжевена [1, 18, 19, 21]. В данной работе используются неравновесные уравнения движения, или, другими словами, на систему накладывается механическое ограничение.

Ограничение основывается на поступательной и вращательной кинетической энергиях. Для каждой молекулы в поступательные уравнения движения (2.11) вводится слагаемое, содержащее множитель Лагранжа, а также в каждое из уравнений Эйлера (2.14) добавляется член общего вида. Поскольку кинетическая энергия равна

(2.16)

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

(2.17)

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