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

УДК 004.94:539.62

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

А.Ю. Смолин12, Г.М. Еремина12, В.В. Сергеев13, Е.В. Шилько12, С.Г. Псахье123

1 Институт физики прочности и материаловедения СО РАН, Томск, 634055, Россия

2 Национальный исследовательский Томский государственный университет, Томск, 634050, Россия

3 Национальный исследовательский Томский политехнический университет, Томск, 634050, Россия

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

Ключевые слова: моделирование, метод подвижных клеточных автоматов, упрочняющие покрытия, наноиндентирование, склерометрия, трибоспектроскопия

Three-dimensional movable cellular automaton simulation of elastoplastic deformation and fracture of coatings in contact interaction with a rigid indenter

A.Yu. Smolin12, G.M. Eremina12, V.V. Sergeev13, E.V. Shilko12, and S.G. Psakhie123

1 Institute of Strength Physics and Materials Science, SB RAS, Tomsk, 634055, Russia 2 National Research Tomsk State University, Tomsk, 634050, Russia 3 National Research Tomsk Polytechnic University, Tomsk, 634050, Russia

In the work, the movable cellular automaton method is applied to consider 3D models of contact interaction such as nanoindentation, sclerometry (scratching), and tribospectroscopy. The system under study is a titanium substrate with a hardening coating. The substrate and coating materials are both described in the elastoplastic approximation. It is shown that the method chosen for numerical simulation is capable of correctly describing the contact interaction of elastoplastic materials under loads of different types with explicit allowance for fracture. The possibility to detect damages in material surface layers from friction force estimates is assessed.

Keywords: simulation, movable cellular automaton method, hardening coatings, nanoindentation, sclerometry, tribospectroscopy

1. Введение

Впервые задачу о контакте двух упругих тел с искривленными поверхностями решил Г. Герц в 1882 г. [1]. Позже с применением теории функций комплексной переменной были получены решения более сложных контактных задач, учитывающих такие факторы, как наличие неоднородностей и пустот в материале, вязко-упругих и пластических свойств, шероховатости поверхности (множественный контакт), их адгезии и т.д. [2-4]. Данная работа посвящена компьютерному моде-

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

Нанесение на поверхность материала специальных покрытий является эффективным способом повышения их функциональных свойств. Например, для имплан-

© Смолин А.Ю., Еремина Г.М., Сергеев В.В., Шилько Е.В., Псахье С.Г., 2014

татов используются многокомпонентные биоактивные наноструктурные покрытия на основе тугоплавких соединений ТЮ, (Т^Та)(С,К) с добавлением специальных элементов (Са, Zr, Si, О, Р), которые улучшают как трибологические, так и биоактивные свойства поверхности [5-8]. Следует отметить, что эти покрытия имеют наноструктурное состояние и малую толщину.

Для изучения механических свойств покрытий и пленок в основном используются такие методы, как наноиндентирование [9-12] и склерометрия [13-15]. Наноиндентирование представляет собой контролируемое внедрение сверхтвердого наконечника (индентора) в поверхность материала. В склерометрических испытаниях (царапание) после внедрения индентора осуществляется его движение вдоль покрытия, что приводит к его разрушению. Одним из методов неразрушающего контроля качества покрытия, позволяющим определить наличие в нем дефектов, является трибоспектроскопия [16], где в качестве измеряемого и анализируемого параметра используется сила трения при движении малого контртела вдоль поверхности.

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

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

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

который позволил бы моделировать на различных масштабах как процесс упругопластического деформирования, так и разрушение твердых тел. Как известно, лучшими возможностями для моделирования разрушения, в том числе зарождения и роста трещин, фрагментации и перемешивания вещества, обладают методы частиц, берущие свое начало из метода молекулярной динамики. Следует подчеркнуть, что из всех таких методов только метод подвижных клеточных автоматов (метод MCA) способен корректно описывать пластическую деформацию консолидированных тел [25], что и обусловило его выбор для построения модели.

2. Описание модели

Метод подвижных клеточных автоматов [25-29] является новым эффективным численным методом, основанным на концепции дискретных элементов (частиц), которая существенно отличается от концепции численных методов классической механики сплошных сред. Следует отметить, что развитие дискретных методов началось с метода молекулярной динамики, который был разработан для изучения материалов на атомном уровне. Однако возможности атомного описания на пространственных и временных масштабах, которые представляют интерес для инженерных приложений, существенно ограничены. Это привело к развитию методов, подобных методу молекулярной динамики для мезо- и макромасштабных уровней (методов частиц), в которых структурные элементы имеют конечный размер (в отличие от атомов, которые являются точечными массами) и взаимодействуют только с ближайшими соседями. Самым известным представителем этой группы методов является метод дискретных элементов [3034], который в настоящее время широко используется для изучения механического поведения гранулированных (сыпучих) и слабосвязанных сред, в частности реологических особенностей этих систем, их разрушения и перемешивания [32, 34]. В то же время до недавнего времени применимость метода дискретных элементов для изучения механического поведения консолидированных тел была ограничена главным образом хрупкими пористыми материалами [30-32, 34] в связи с недостаточным развитием математических моделей для вычисления сил взаимодействия дискретных элементов. В частности, большинство моделей дискретных элементов используют парные (двухчастичные) потенциалы (или силы) взаимодействия. Это упрощение приводит к ряду искусственных эффектов в поведении ансамбля частиц, таких как анизотропия, навязанное упаковкой значение модуля сдвига и коэффициента Пуассона и т.д. [30, 34]. Среди проблем, возникающих вследствие этих искусственных эффектов, важно отметить неадекватность моделирования накопления необратимых деформаций (пластичности) материалов.

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

Ц = 22 Ца1Г + Ц°.

3=1

Эта сила представлена как суперпозиция парных компонентов Цра1Г, зависящих от пространственного расположения автомата i по отношению к соседу j, и объемнозависящей компоненты обусловленной коллективными эффектами окружения.

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

г Н2К N

т = 2 ЦТ + Ц°,

з=1

I а20 N.

дт з=1

где К,, 9,, т{ и 31 — радиус-вектор, вектор вращения, масса и момент инерции '-го автомата соответственно; ЦРа1Г — парная сила механического взаимодействия автоматов i и j; — объемнозависящая сила, действующая на автомат i и обусловленная взаимодействием его соседей с другими автоматами. В последнем

уравнении М у = qij (п 3 х ЦРа1Г) + К 3, где qij — расстояние от центра /'-го автомата до точки его взаимодействия (контакта) с ]'-м автоматом; п 3 = (Я3 - Я,)Гу — единичный вектор ориентации пары и г3 — расстояние между центрами автоматов (рис. 1); К 3 — крутящий момент, обусловленный относительным вращением автоматов в паре (см. ниже).

Отметим, что автоматы пары могут представлять собой части различных тел или одного консолидированного тела. Поэтому их взаимодействие не всегда является контактным в том смысле, который имеется в виду в известных задачах механики [2-4]. Так, в случае автоматов одного консолидированного тела площадка контакта автоматов будет представлять собой просто элемент плоскости, разделяющей два конечных объема материала, аналогичный граням ячеек разностной сетки в численных методах континуальной механики. По этой причине здесь и далее мы будем слово «контакт» брать в кавычки. Кроме того, как показано на рис. 1, а, размер автомата характеризуется одним параметром di, но это не означает, что форма автомата представляет собой шар. Реальная форма автомата определяется областью его «контактов» с соседями. Например, если использовать в качестве начальной ГЦК-упаковку, то автоматы будут иметь форму ромбического додекаэдра, а если кубическую, то автоматы будут представлять собой кубики.

Для локально изотропных сред объемнозависящая компонента силы Ц может быть записана через давление Р3 в объеме соседнего автоматаj следующим образом:

Ц?=-А ЪР^ П3,

3=1

где SiJ — площадь контакта '-го автомата с j-м; А — некий материальный параметр.

С другой стороны, общая сила, действующая на автомат, может быть представлена в виде суммы нормальной ЦП и касательной (сдвиговой) ЦТ компонент:

F. =£ (Fj311 - APS nj) =

j =1

N

= E [(FT'" (hj) - APjSj )nj + F?^ (lfar)tj ] =

j=1

N.

= E (F" + Fj),

j=i

где ^?а1Г'и — нормальная, а рра1Г>т — касательная силы

и и

взаимодействия в паре, зависящие соответственно от межавтоматного перекрытия Ь- (рис. 1, а) и относительного тангенциального смещения 1-1еаг (рис. 1, б), рассчитанного с учетом вращения обоих автоматов [28, 29]. Отметим что, несмотря на то, что последнее выражение в уравнении (2) формально соответствует обычной записи силы взаимодействия в методе дискретных элементов [33-34], оно принципиально отличается от последнего вследствие многочастичного характера центрального взаимодействия автоматов.

Используя процедуру осреднения для тензора напряжений в частице, изложенную в [25, 34], выражение для компонент усредненного тензора напряжений в автомате i принимает вид:

1

J«p

Vi j=1

где а, X, Y, 2 лабораторной системы координат; Vi — текущий объем автомата г; п- а — а-компонента единичного вектора п -; Р- р — Р-компонента полной силы, действующей в точке «контакта» между автоматами

i и j. _

Давление Ц, или среднее напряжение а^ит, в объеме автомата может быть вычислено через компоненты тензора напряжений:

P =-oL

Знание компонент тензора напряжений позволяет вычислять все его инварианты в объеме автомата, в частности интенсивность напряжений:

= 1/V2 ((?„ - ^ )2 + tfyy - aZz )2 +

+ (^ -аХх)2 + 6[(аХу )2 + (а^)2 + (aixz)2])1/2. (5) Из уравнений (1)-(3) следует, что выражения для выи F?™'т определяют реологические

числения Fifair';

свойства модельной среды.

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

и j характеризуется соответствующими напряжениями % и т-

\ =ПуЯу ,

' г

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

Деформацию автомата i при его центральном взаимодействии с автоматом j можно характеризовать следующей безразмерной величиной (нормальной деформацией):

q.j- dil2 dj 2 '

В общем случае каждый из автоматов представляет собой различные материалы, и перекрытие в паре некоторым образом перераспределяется между г'-м иj-м автоматами:

АЬ- = Aqlj+ Aq- = А-, /2 + А^-, dJ /2,

где символ А означает приращение величины за временной шаг At численного интегрирования уравнений движения (1). Правило перераспределения деформации в паре тесно связано с выражением для вычисления силы взаимодействия автоматов. Это выражение для центрального взаимодействия похоже на соотношения Гука для диагональных компонент тензора напряжений:

А- 2С

где К — модуль объемной упругости; G — модуль сдвига материала автомата г; Ц — давление автомата г, которое может быть вычислено с помощью формул (3) и (4) на предыдущем шаге по времени или по схеме предиктор-корректор.

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

где г- = (R- - Rг), V, — трансляционная скорость автомата г'; ю- — мгновенная скорость вращения пары как целого. Если мы умножим обе части уравнения (10) векторно слева на г- и пренебрежем вращением вокруг оси пары (т.е. пусть ю- • г- = 0, т.к. вращение вокруг оси пары не производит деформацию сдвига), то мы получим следующую формулу:

Г- Х

г2 'у

Г-

'у

Кроме такого вращения пары как целого (определяется разницей в поступательных скоростях автоматов), каждый из автоматов вращается со своей собственной угловой скоростью ю, (рис. 1, в). Разница между этими вращениями обусловливает деформацию сдвига. Так, приращение деформации сдвига автоматов г и ]

за шаг Дt определится относительным тангенциальным

смещением в точке «контакта» пары Д131еаг, деленной на расстояние между автоматами:

shear

дг

ДУ а +дУ ц = —

(Чц («ц - Щ)х па + q а («а - « а)х пц)д

лу а

Выражение для вычисления силы тангенциального взаимодействия подвижных клеточных автоматов похоже на соотношения Гука для недиагональных компонент тензора напряжений:

Дх„. = 2С

и является чисто парным.

Разные скорости вращения каждого из автоматов также приводят к деформации относительного «изгиба» и «кручения» (последнее присутствует только в 3Э) в паре. Очевидно, что сопротивление относительному вращению в паре вызывает крутящий момент, величина которого пропорциональна разнице между поворотами автоматов:

К 3 +GJ )

Формулы (1)-(4), (6)-(9), (11)—(14) описывают механическое поведение линейно-упругого тела в рамках метода подвижных клеточных автоматов. Отметим, что соотношения (8), (9), (12) и (13) записаны в приращениях, т.е. в гипоупругой форме. В работе [35] показано, что эта модель дает результаты, соответствующие численному решению методом конечных разностей обычных уравнений механики сплошной среды для изотропной линейно-упругой среды. Это обстоятельство позволило использовать метод подвижных клеточных автоматов совместно с численными методами механики сплошных сред в рамках единого дискретно-континуального подхода. В работе [29] показано, что именно учет вращения автоматов совместно с многочастичной формой центрального взаимодействия позволяет методом подвижных клеточных автоматов адекватно описывать изотропный отклик материала на приложенную нагрузку.

Для описания упругопластического поведения в рамках метода подвижных клеточных автоматов предлагается использовать теорию пластического течения, а именно модель идеальной пластичности с критерием Мизеса. Для этого к методу подвижных клеточных автоматов [25] был адаптирован известный алгоритм Уил-кинса [36, 37]. Этот алгоритм состоит в решении упругой задачи на каждом временном шаге и последующем «сбросе» компонент девиатора тензора напряжений Дар = аар -1/3 акк5ар на поверхности текучести Мизеса в случае, когда интенсивность напряжений превышает заданную предельную величину (рис. 2):

дав = Дав М,

Рис. 2. Схема работы алгоритма Уилкинса

где М = ар1/ а1п1 , аы — интенсивность напряжений; ар1 — радиус круга текучести Мизеса.

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

I (ааа ) = (ааа - атеап )М, + атеап '

Ыр )' = 50р М,,

где а, в = X, Y, 2 и а ф в; (а^р)' — скорректированные компоненты осредненного тензора напряжений; агар — компоненты тензора напряжений, полученные в результате решения упругой задачи на текущем временном

шаге; М{ = 0^/ а[п1 — текущее значение коэффициента «сброса»; ар1 — текущее значение радиуса круга текучести Мизеса для автомата

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

К = (% - Р) М1 + Р,

{ т 3 = т 3М,,

где и т 3 представляют собой скорректированные значения удельных сил (напряжений). Легко показать, что подстановка этих соотношений в выражения (2) и (3) автоматически обеспечивает корректировку компонент осредненного тензора напряжений в автомате к кругу текучести [25].

Таким образом, реологические свойства материала автомата ' определяются заданием единой кривой упрочнения а|п1 =Ф(ё^) (здесь ёП — интенсивность осредненного тензора деформаций, компоненты которого могут быть вычислены аналогично [25, 34]), эта зависимость также называется функцией отклика автомата.

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

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

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

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

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

3. Результаты расчетов по наноиндентированию

Наноиндентирование — это процесс контролируемого внедрения сверхтвердого наконечника определенной формы (индентора) под действием нарастающей нагрузки в плоскую поверхность образца на глубину менее 100 нм, при этом в процессе нагружения постоянно измеряется сила Р, действующая на индентор, и глубина его погружения в материал h [8]. Иногда наноин-дентирование называют измерительным индентирова-нием. На основе анализа измеряемой Р-^диаграммы можно получать такие характеристики материала, как модуль упругости, упругое восстановление и нанотвер-дость. При таком анализе в основном используется методика Оливера-Фарра [10]. Однако, как показали, например, авторы работы [11], корректно определить модуль упругости и нанотвердость тонких покрытий и пленок по данным наноиндентирования с использованием стандартной методики Оливера-Фарра возможно только при условии совпадения этих характеристик у покрытия и подложки. Очевидно, что в большинстве практически важных приложений это условие не выполняется. Одним из способов решения этой проблемы может быть использование компьютерного моделирования, в рамках которого возможно получение достаточно точных зависимостей для индентирования покрытий на любых подложках.

При разработке модели, а также ее тестирования были использованы результаты натурных экспериментов, опубликованные в работах [5-8, 38-40]. Так, геометрические характеристики модельной системы «покрытие - переходный слой - подложка» определялись по фотографии поперечного среза реального образца [38]. Полученные по результатам натурного наноинден-тирования значения модуля Юнга и предела текучести учитывались при задании параметров модельных материалов.

Механическое поведение материалов и подложки, и покрытия описывалось в рамках упругопластического приближения. Упруго-идеальнопластическое тело характеризуется следующими параметрами: плотность р, модуль сдвига G, модуль объемной упругости К и предел текучести ау. На основе данных, приведенных в [39, 40], значения указанных параметров для подложки составили р = 4500 кг/м3, G = 41 ГПа, К = 100 ГПа и ау =2 ГПа, а для покрытия р = 4700 кг/м3, G = 76 ГПа, К = 167 ГПа и ау = 15 ГПа.

В экспериментах по наноиндентированию и склерометрии обычно используются алмазные инденторы. Модуль упругости алмаза (1141 ГПа) приблизительно в пять раз превышает модуль упругости для многокомпонентных биоактивных наноструктурных покрытий (200-250 ГПа). Поэтому в данных расчетах материал индентора принимался абсолютно жестким (неде-формируемым). Форма индентора соответствовала трехгранной пирамиде Берковича с углом между гранями и осью пирамиды равным 65.3°.

Рассматривалось индентирование двух типов образцов, каждый из которых имел форму параллелепипеда. Первый тип образцов состоял полностью из нанострук-турного титана, а второй — из нескольких слоев, нижний из которых представлял собой подложку из титана, а верхний — покрытие. Толщина такого покрытия составляет 350-2000 нм [5-8], при этом общая высота образцов варьировалась от 700 до 2500 нм. Между покрытием и титановой подложкой находился промежуточный слой толщиной порядка 50 нм. В модели свойства промежуточного слоя задавались равными средним величинам между значениями для покрытия и подложки. Размер автоматов, исходя из предварительных расчетов [41, 42], брался равным 10 нм.

Общий вид построенной модели представлен на рис. 3, а в виде ГЦК-упаковки сферических частиц. Процесс нагружения имитировался путем задания одинаковой скорости движения для всех автоматов инден-тора. Для того чтобы предотвратить смещение нагружаемого материала как целого, нижний слой автоматов подложки считался жестко закрепленным.

Анализ поля скоростей в начале индентирования, когда материал деформируется упруго, показывает, что смещения в материале локализуются под инденто-

Рис. 3. Общий вид модели (а), распределение скоростей автоматов (б, абсолютные значения в м/с), интенсивности напряжений (Па) по сечению образца при его упругой (в) и пластической (г) деформации

1.5-10

1.0-10^ —

0.75 и 0.50 -0.25 " 0.00

ром и направлены преимущественно вглубь материала (рис. 3, б). Полученные результаты согласуются с данными натурных экспериментов, которые показывают, что результаты индентирования тонкого покрытия при больших нагрузках в значительной степени определяются свойствами подложки [9, 11, 39, 40]. Распределение упругих напряжений в образце (рис. 3, в) качественно соответствует аналитическому решению для задачи Герца [4]. На рис. 3, г показано поле интенсивности напряжений при максимальной глубине проникновения, когда покрытие и подложка деформируются уже пластически. Видно, что максимальные напряжения находятся в непосредственной близости от индентора, однако из-за меньшего предела текучести у подложки она также деформируется пластически. Поля напряжений, представленные на рис. 3, г, находятся в хорошем качественном согласии с литературными данными [43, 44].

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

На рис. 4 представлены расчетная и экспериментальная [39] Р-^диаграммы индентирования титанового образца на глубину 650 нм. Их анализ позволяет сделать вывод о возможности нахождения свободных параметров модели для описания свойств материала

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

Для определения адекватности параметров модельного покрытия был сделан расчет с глубиной внедрения в образец с покрытием, равной 60 нм [42]. Результаты этого расчета хорошо соответствуют экспериментальным данным при максимальной нагрузке на индентор 2 мН [39, 40]. Для изучения влияния подложки на Р-^ диаграмму моделировалось погружение на большую глубину. В связи с этим использовались образцы с общей высотой 2.5 и 5.0 мкм и размером автоматов 70 и 140 нм соответственно, толщина покрытия составляла 1000 нм. Р-^диаграммы, полученные в результате моР, мЩ А

20

10

Рис. 4. Экспериментальная [39] (1) и расчетная (2) диаграммы индентирования титанового образца

Р, мН " 80  60  40  20

0 -

О

Рис. 5. Экспериментальные [39, 40] (1, 3) и расчетные (2, 4) диаграммы индентирования титанового образца с нанострук-турным покрытием

делирования и эксперимента [39, 40], показаны на рис. 5. Согласие результатов расчетов с имеющимися экспериментальными данными свидетельствуют о том, что предложенная модель хорошо описывает особенности деформационного поведения таких систем.

4. Результаты расчетов по склерометрии

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

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

Модель процесса склерометрии (рис. 6, а) в целом подобна модели для наноиндентирования и состоит из четырех блоков: основание, подложка, переходный слой, покрытие и индентор (контртело). Образец исследуемого материала, состоящий из трех материалов, имеет геометрическую форму параллелепипеда, высота покрытия 1.5 мкм, переходного слоя — 0.2 мкм и подложки из титана — 2.5 мкм. Основание представляет собой пластину со сторонами равными сторонам образца и высотой равной размеру автомата. Высота контртела в форме конуса составляет 0.2 мкм. Процесс нагружения имитировался заданием скоростей в вертикальном (ось 2) и горизонтальном (ось Y) направлениях всем автоматам контртела. То есть вместо силы, действующей на индентор по оси 2 в эксперименте, в расчетах задавалась постоянная вертикальная скорость равная -1 м/с до погружения индентора на заданную глубину, после чего вертикальная составляющая нагрузки задавалась равной нулю. В рассматриваемом случае глубина погружения индентора составляла расстояние до промежуточного слоя. Для движения контртела вдоль модельного образца всем его автоматам задавалась постоянная скорость в горизонтальном направлении по оси Y равная -5 м/с. Отметим, что в реальных экспериментах по склерометрии скорость движения индентора относительно материала, как правило, не превышает нескольких сантиметров в секунду.

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

Рис. 6. Общий вид модели для расчетов склерометрии в виде сетки межавтоматных связей

О 10 20 30 40 мкс

Рис. 7. Зависимость коэффициента трения

туры межавтоматных связей в исходном состоянии (рис. 6, а) и после прохождения индентора (рис. 6, б) позволяет увидеть, что вдоль пути прохождения индентора переходный слой и примыкающая к нему часть подложки деформировались пластически. При этом покрытие растрескивается и разрушается вплоть до переходного слоя, однако полного отслоения материала покрытия не происходит. Качественно полученные изображения соответствуют начальной стадии разрушения покрытия в натурном эксперименте по измерительному царапанию наноструктурных покрытий на титановой подложке, описанному в [8].

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

Полученную в результате проведения численного эксперимента зависимость (рис. 7) сравнивали с данными натурного эксперимента [8]. Средняя величина коэффициента трения, полученная в результате компьютерного моделирования (0.61), соответствует экспериментальному значению в момент разрушения покрытия (0.63, при этом в покрытии появляются большие трещины, но отслоения от подложки не происходит).

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

5. Моделирование трибоспектроскопии

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

Методика моделирования в рамках метода МСА процессов, имеющих место при трибоспектроскопии, описана в работах [16, 45, 46]. Их сущность сводится к моделированию движения контртела по поверхности образца. Вследствие искусственной шероховатости, вызванной дискретностью метода моделирования, между контртелом и поверхностью образца возникают силы трения, а от пятна контакта генерируются упругие волны, распространяющиеся по образцу. Вследствие динамического характера такого процесса регистрируемая сила трения скольжения меняется во времени. Преобразование Фурье данной силы для бездефектного образца позволяет выделить некоторые характерные частоты [45]. К ним относятся частота, соответствующая искусственной шероховатости поверхности образца (она связана с размером элементов модели — клеточных автоматов), и собственные частоты моделируемой системы.

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

Моделировалось керамическое покрытие на титановом образце (рис. 8). Толщина покрытия была выбрана равной Н= 60 нм, размер автомата d = 3 нм, длина образца Ь = 250 нм, ширина образца М= 125 нм. В рамках используемого приближения покрытие из керамики ZrO2 считалось однородным и изотропным.

Контртело имело форму усеченного конуса с диаметром основания 60 нм, его движение моделировалось заданием для автоматов верхней поверхности контртела постоянной скорости V = 5 м/с в направлении оси У (рис. 8). При этом нижняя поверхность образца (слой автоматов, моделирующий подложку) была неподвижной, а его боковые поверхности — свободными. При

Рис. 8. Образцы с различными видами пор: общий вид моделируемого образца со сквозными трещинами, расположенными на глубине 30 нм

движении контртела проводилась регистрация силы сопротивления его движению по поверхности (силы трения Р). Затем строилось преобразование Фурье этой функции [45].

Поврежденность поверхностного слоя моделировалась заданием наноскопических нарушений сплошности покрытия. Анализировались протяженные повреждения — нанопоры. Рассматривалось периодическое расположение нанопор с расстоянием Р = 75 нм различного вида: сквозные повреждения, несквозные расположенные посередине образца и несквозные нанопоры, расположенные сбоку образца (рис. 8). При этом длина всех пор составляла w = 48 нм, высота h = 3 нм, а ширина пор, расположенных сбоку и посередине, равна их дли-

s <

0.05 0.15 0.25

Частота, ГГц

Рис. 10. Спектры Фурье силы, действующей на контртело, в различное время расчета

не и составляла 48 нм. В работе изучалась возможность определения расстояния между такими дефектами.

Частота, соответствующая периоду следования нанопор, вычисляется по формуле V|P, где V — скорость движения контртела, а Р—период следования нанопор, и в нашем случае приблизительно равна 48 МГц. Данная частота хорошо видна на спектрах Фурье для составляющей силы сопротивления в направлении оси 2, независимо от глубины расположения и вида повреждений, но имеет тенденцию к уменьшению амплитуды при увеличении глубины их расположения от поверхности (рис. 9).

Частота, соответствующая размеру пор (48 нм), четко идентифицируется на спектрах горизонтальной составляющей силы сопротивления и равна 104 МГц. Следует отметить, что пик на этой частоте определяется достаточно четко на начальной стадии движения инден-тора (до прохождения второй поры, пунктирная линия на рис. 10), после чего амплитуда этого пика значительно уменьшается и на спектре появляется частота, характеризующая период расположения пор (сплошная линия на рис. 10).

Рис. 9. Спектры Фурье силы, действующей на контртело, для образцов с нанопорами, расположенными посередине образца на глубине от поверхности 6

Рис. 11. Спектры Фурье вертикальной составляющей силы для образцов с различными видами нанопор: сквозные (1), посередине (2) и боковые (3); глубина пор от поверхности 30 нм

д 100

<

25

ÜJ 005 ' 0Л5 ' Ö25

Частота, ГГц

Рис. 12. Спектры Фурье вертикальной составляющей силы для образцов с нанопорами, расположенными сбоку образца с толщиной 60

Величина амплитуды пиков также зависит от вида пор. Наибольшее значение она имеет при сквозных порах, меньшее — при расположенных посередине и минимальное — сбоку (рис. 11).

Кроме того, исследовалось влияние толщины покрытия, значения которой составляли 60, 70 и 90 нм, при условии одинаковой глубины расположения нанопор от исследуемой поверхности. Было установлено, что чем больше толщина покрытия, тем меньше амплитуда пика, соответствующего периоду следования нанопор. Для примера на рис. 12 представлены спектры Фурье для боковых трещин с различной толщиной покрытия.

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

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

6. Заключение

На основании результатов проведенных исследований можно сформулировать следующие выводы.

Представленные трехмерные модели процессов на-ноиндентирования и измерительного царапания, по-

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

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

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

Данные исследования выполнены в рамках Программы фундаментальных научных исследований государственных академий наук на 2013-2020 гг., работ по государственному контракту № 16.523.12.3002, скоординированному с проектом Евросоюза ViNaT, Contract No. 295322, FP7-NMP-2011-EU-Russia, NMP.2011.1.4-5.

Литература

1. Ландау Л.Д., Лифшиц Е.М. Теория упругости. - М.: Наука, 1987. -

246 с.

Landau L.D., Lifshitz E.M. Theory of Elasticity. - Oxford: Pergamon Press, 1986.

2. Развитие теории контактных задач в CCCP / Под ред. Л.А. Галина. - М.: Наука, 1976. - 493 с.

Development of the theory of contact problems in the USSR / Ed. by L.A. Galin. - Moscow: Nauka, 1976. - 493 p.

3. Горячева И.Г. Механика фрикционного взаимодействия. - М.: Наука, 2001. - 477 с.

Goryacheva I.G. Mechanics of Friction Interaction. - Moscow: Nauka, 2001. - 477 p.

4. Попов В.Л. Механика контактного взаимодействия и физика трения: от нанотрибологии до динамики землетрясений. - М.: Физ-матлит, 2013. - 350 с.

Popov V.L. Contact mechanics and friction. Physical Principles and Applications. - Berlin: Springer-Verlag, 2010. - 362 p.

5. Shtansky D. V., Kiryukhantsev-Korneev Ph. V., Bashkova I.A., Shevei-ko A.N., Levashov E.A. Multicomponent nanostructured films for various tribological applications // Int. J. Refract. Met. Hard Mater. -2010. - V. 28. - P. 32-39.

6. Shtansky D. V., Gloushankova N.A., Bashkova I.A., Petrzhik M.I., Sheveiko A.N. Multifunctional biocompatible nanostructured coatings for load-bearing implants // Surf. Coat. Tech. - 2006. - V. 201. -P. 4111-4118.

7. Shtansky D. V., Levashov E.A., Glushankova N.A. et al. Structure and properties of CaO- and ZrO2-doped TiCxNy coatings for biomedical applications // Surf. Coat. Tech. - 2004. - V. 182. - P. 101-111.

8. Левашов Е.А., Петржик М.И., Тюрина М.Я., Кирюханцев-Корне-ев Ф.В., Цыганков П.А., Рогачев А.С. Многослойные нанострук-турные тепловыделяющие покрытия. Получение и аттестация механических и трибологических свойств // Металлург. - 2010. -№ 9. - C. 66-74.

Levashov E.A., Petrzhik M.I., Tyurina M.Ya., Kiryukhantsev-Korneev F.V., Tsygankov P.A., Rogachev A.S. Multilayer nanostructured

heat-generating coatings. Preparation and certification of mechanical and tribological properties // Metallurgist. - 2011. - V. 54. - No. 910. - P. 623-634.

9. Головин И.Ю. Наноиндентирование и его возможности. - М.: Машиностроение, 2009. - 316 с.

Golovin I.Yu. Nanoindentation and its Possibilities. - Moscow: Mashi-nostroenie, 2009. - 316 p.

10. Oliver W.C., Pharr G.M. An improved technique for determining hardness and elastic modulus using load and displacement sensing indentation experiments // J. Mater. Res. - 1992. - No. 7. - P. 1564-1583.

11. Шугуров A.P., Панин A.B., Оскомов K.B. Особенности определения механических характеристик тонких пленок методом наноинден-тирования // ФТТ. - 2008. - Т. 5. - № 6. - С. 1007-1012. Shugurov A.R., Panin A.V., Oskomov K.V. Specific features of the determination of the mechanical characteristics of thin films by the nanoindentation technique // Phys. Sol. St. - 2008. - V 50. - No. 6. -P. 1050-1055.

12. Venkatesh T.A., van Vliet K.J., Giannakopoulos A.E., Suresh S. Determination of elasto-plastic properties by instrumented sharp indention: guidelines for property extraction // Scripta Mater. - 2000. -V. 42. - No. 9. - P. 833-839.

13. Flores S.E., Pontin M.G., ZokF.W. Scratching of elastic/plastic materials with hard spherical indenters // J. Appl. Mech. - 2008. - V 75. -P. 061021-1-7. - DOI: 10.1115/1.2966268.

14. Шугуров A.P., Панин A.B., ШестериковE.B. Исследование гальванических покрытий AuNi и AuCo методом склерометрии // Письма в ЖТФ. - 2011. - Т. 37. - № 5. - C. 64-71.

Shugurov A.R., Panin A.V., Shesterikov E.V Sclerometric study of galvanic AuNi and AuCo coatings // Tech. Phys. Lett. - 2011. - V 37. -No. 3. - P. 223-225.

15. Азаренков H.A., Литовченко C.B., Береснев B.M., Чишкала B.A., Гравнова Л..B, Огиенко СИ. Особенности разрушения силицид-ных покрытий при динамических механических нагрузках // Физическая инженерия поверхности. - 2012. - Т. 10. - № 4. - C. 411417.

Azarenkov N.A., Litovchenko S.V, Beresnev V.M., Chishkala V.A., Gravnova L.V., Ogienko S.I. Peculiarities of fracture of silicide coatings at mechanical dynamic loadings // Fiz. Inzhen. Pover. - 2012. -V. 10. - No. 4. - P. 411-417.

16. Псахье С.Г., Попов B.M., Шилько E.B., Смолин AЮ., Дмитри-евAИ. Изучение поведения и диагностика свойств поверхностного слоя твердого тела на основе спектрального анализа. Нано-трибоспектроскопия // Физ. мезомех. - 2009. - Т. 12. - № 4. -C. 27-42.

Psakhie S.G., Popov V.L., Shilko E.V, Smolin A.Yu., Dmitriev A.I. Spectral analysis of the behavior and properties of solid surface layers. Nanotribospectroscopy // Phys. Mesomech. - 2009. - V. 12. -No. 5-6. - P. 221-234.

17. Dao M., Chollacoop N., van Vliet K.J., Venkatesh T.A., Suresh S. Computational modeling of the forward and reverse problems in instrumented sharp indentation // Acta Mater. - 2001. - V. 49. - P. 38993918.

18. Bucaille J.L., Stauss S., Felder E., Michler J. Determination of plastic properties of metals by instrumented indentation using different sharp indenters // Acta Mater. - 2003. - V. 51. - P. 1663-1678.

19. Ogasawara N., Chiba N., Chen X. Measuring the plastic properties of bulk materials by single indentation test // Scripta Mater. - 2006. -V. 54. - P. 65-70.

20. Sreeranganathan A., Gokhale A., Tamirisakandala S. Determination of local constitutive properties of titanium alloy matrix in boron-modified titanium alloys using spherical indentation // Scripta Mater. -2008. - V 58. - No. 2. - P. 114-117.

21. Zimmerman J.A., Kelchner C.L., Klein P.A., Hamilton J.S., Foiles S.M. Surface step effects on nanoindentation // Phys. Rev. Lett. - 2001. -V. 87. - P. 165507-165511.

22. Saraev D., Miller R.E. Atomic-scale simulations of nanoindentation-induced plasticity in copper crystals with nanometer-sized nickel coatings // Acta Mater. - 2006. - V. 54. - P. 33-45.

23. Mei J., Li J., Ni Y., Wang H. Multiscale simulation of indentation, retraction and fracture processes of nanocontact // Nanoscale Res. Lett. - 2010. - V. 5. - P. 692-700.

24. Solar M., Meyer H., Gauthier C. Analysis of local properties during a scratch test on a polymeric surface using molecular dynamics simulations // Europ. Phys. J. E. - 2013. - V. 36. - No. 29. - P. 1-18.

25. Psakhie S.G., Horie Y., Shilko E.V., Smolin A.Yu., Dmitriev A.I., Astafurov S.V Development of discrete element approach to modeling heterogeneous elastic-plastic materials and media // Int. J. Terra-space Sci. Eng. - 2011. - V. 3. - No. 1. - P. 93-125.

26. Псахье С.Г., Остермайер Г.П., Дмитриев A.^, Шилько E.B., Смолин AЮ., Коростелев С.Ю. Метод подвижных клеточных автоматов как новое направление дискретной вычислительной механики. I. Теоретическое описание // Физ. мезомех. - 2000. -Т. 3. - № 2. - С. 5-13.

Psakhie S.G., Ostermeyer G.P., Dmitriev A.I., Shilko E.V., Smolin A.Yu., Korostelev S.Yu. Method of movable cellular automata as a new trend of discrete computational mechanics. I. Theoretical description // Phys. Mesomech. - 2000. - V. 3. - No. 2. - P. 5-12.

27. Попов B.M., Псахье С.Г. Теоретические основы моделирования упругопластических сред методом подвижных клеточных автоматов. I. Однородные среды // Физ. мезомех. - 2001. - Т. 4. - № 1.-С. 17-28.

Popov V.L., Psakhie S.G. Theoretical principles of modeling elasto-plastic media by movable cellular automata method. I. Homogeneous media // Phys. Mesomech. - 2001. - V. 4. - No. 1. - P. 15-25.

28. Псахье С.Г., Моисеенко Д.Д., Дмитриев A.^, Шилько E.B., Коростелев С.Ю., Смолин A.Ю., Дерюгин E.E., Кульков С.Н. О возможности компьютерного конструирования материалов с высокопористой и каркасной структурой на основе метода подвижных клеточных автоматов // Письма в ЖТФ. - 1998. - Т. 24. - №4.-С. 71-76.

Psakhie S.G., Moiseenko D.D., Dmitriev A.I., Shilko E.V., Korostelev S.Yu., Smolin A.Yu., Deryugin E.E., Kulkov S.N. A possible method of computer-aided design of materials with a highly porous matrix structure based on the method of moving cellular automata // Tech. Phys. Lett. - 1998. - Т. 24. - № 2. - С. 154-156.

29. Смолин AM., Роман H.B., Добрынин C.A., Псахье С.Г. О вращательном движении в методе подвижных клеточных автоматов // Физ. мезомех. - 2009. - Т. 12. - № 2. - С. 17-22.

Smolin A.Yu., Roman N.V., Dobrynin S.A., Psakhie S.G. On rotation in the movable cellular automaton method // Phys. Mesomech. -2009. - V. 12. - No. 3-4. - P. 124-129.

30. Cundall P.A., Strack O.D.L. A discrete numerical model for granular assemblies // Geotechnique. - 1979. - V. 29. - No. 1. - P. 47-65.

31. Jing L., Stephansson O. Fundamentals of Discrete Element Method for Rock Engineering: Theory and Applications. - Oxford: Elsevier, 2007. - 562 p.

32. Sibille L., NicotF., Donze F. V., Darve F. Material instability in granular assemblies from fundamentally different models // Int. J. Numer. Analyt. Meth. Geomech. - 2007. - V. 31. - No. 3. - P. 457-481.

33. Martin C.L., Bouvard D. Study of the cold compaction of composite powders by the discrete element method // Acta Mater. - 2003. -V. 51.- No. 2. - P. 373-386.

34. Potyondy D.O., Cundall P.A. A bonded-particle model for rock // Int. J. Rock Mech. Mining Sci. - 2004. - V 41. - No. 8. - P. 13291364.

35. Псахье С.Г., Смолин A.M., Стефанов Ю.П., Макаров П.B., Шилько E.B., Чертов M.A., Eвтушенко EH. Моделирование поведения сложных сред на основе комбинированного дискретно-континуального подхода // Физ. мезомех. - 2003. - Т. 6. - № 6. - С. 1121.

Psakhie S.G., Smolin A.Yu., Stefanov Yu.P., Makarov P. V., Shilko E.V., Chertov M.A., Evtushenko E.P. Simulation of behavior of complex media on the basis of a discrete-continuous approach // Phys. Mesomech. - 2003. - V. 6. - No. 5-6. - P. 47-56.

36. Wilkins M.L. Computer Simulation of Dynamic Phenomena. - Berlin: Springer-Verlag, 1999. - 246 p.

37. Уилкинс М.Л. Расчет упругопластических течений // Вычислительные методы в гидродинамике. - М.: Мир, 1967. - C. 212-263. Wilkins M.L. Calculation of Elastic-Plastic Flow // Methods in Computational Physics, Vol. 3 / Ed. by B. Alder, S. Fernbach, M. Rotenberg. - New York: Academic Press, 1964. - P. 211.

38. ЛевашовЕ.А., ШтанскийД.В., Кирюханцев-Корнеев Ф.В., Петр-жик М.И., Тюрина М.Я., Шевейко А.Н. Многофункциональные наноструктурные покрытия: получение, структура и обеспечение единства измерений механических и трибологических свойств // Деформация и разрушение. - 2009. - № 11. - C. 19-35. Levashov E.A., Shtanskii D.V., Kiryukhantsev-Korneev F.V., Petr-zhikM.I., Tyurina M.Ya., Sheveiko A.N. Multifunctional Nanostruc-tured Coatings: Formation, Structure, and the Uniformity of Measuring Their Mechanical and Tribological Properties // Russ. Metall. -2010. - No. 10. - P. 917-935.

39. Левашов Е.А., Петржик М.И., Кирюханцев-Корнеев Ф.В., Штанский Д.В., Прокошкин С.Д., Гундеров Д.В., Шевейко А.Н., Коро-тицкий А.В., Валиев Р.З. Cтруктура и механическое поведение при индентировании биосовместимых нано-структурированных титановых сплавов и покрытий // Металлург. - 2012. - № 5. -C.79-89.

Levashov E.A., Petrzhik M.I., Kiryukhantsev-Korneev F.V., Shtans-kiiD.V., Prokoshkin S.D., Gunderov D.V, Sheveiko A.N., Korotits-kiiA. V., Valiev R.Z. Structure and mechanical behavior during indentation of biocompatible nanostructured titanium alloys and coatings // Metallurgist. - 2012. - V. 56. - No. 5. - P. 395-407.

40. Levashov E.A., Petrzhik M.I., Shtansky D.V, Kiryukhantsev-Kor-neevF.V, Sheveiko A.N., Gunderov D.V., Prokoshkin S.D., Korotits-kiy A.V., Smolin A.Yu. Nanostructured titanium alloys and multicom-ponent bioactive films: Mechanical behavior at indentation // Mater. Sci. Eng. A. - 2013. - V. 570. - P. 51-62.

41. Смолин А.Ю., Аникеева Г.М., Шилько Е.В., Псахье С.Г. Моделирование деформации наноструктурных покрытий на титановой под-

ложке при наноиндентировании // Вестник ТГУ. Математика и механика. - 2013. - Т. 24. - № 4. - C. 111-125. Smolin A.Yu., Anikeeva G.M., Shilko E.V., Psakhie S.G. Simulation of deformation of nanostructural coatings on the titanium substrate at nanoindentation // Vest. TGU. Mat. Mekh. - 2013. - V. 24. - No. 4. -P. 111-125.

42. Psakhie S.G., Smolin A.Yu., Shilko E.V., Anikeeva G.M., Pogo-zhev Yu.S., Petrzhik M.I., Levashov E.A. Modeling nanoindentation of TiCCaPON coating on Ti substrate using movable cellular automaton method // Comput. Mater. Sci. - 2013. - V. 76. - P. 89-98.

43. Muliana A., Steward R., Haj-ali R.M., Saxena A. Artificial neural network and finite element modeling of nanoindentation tests // Metall. Mater. Transact. A. - 2002. - V. 33. - P. 1939-1948.

44. Feng Z.-Q., Zei M., Joli P. An elasto-plastic contact model applied to nanoindentation // Comput. Mater. Sci. - 2007. - V. 38. - P. 807-813.

45. Добрынин С.А., Колубаев Е.А., Смолин А.Ю., Дмитриев А.И., Псахье С.Г. Частотно-временной анализ акустических сигналов, генерируемых при трении // Письма в ЖТФ. - 2010. - Т. 36. -№ 13. - C. 47-53.

Dobrynin S.A., Kolubaev E.A., Smolin A.Yu., Dmitriev A.I., Psakhie S.G. Time-frequency analysis of acoustic signals in the audiofrequency range generated during Hadfield's steel friction // Tech. Phys. Lett. - 2010. - V. 36. - No. 7. - P. 606-609.

46. Сергеев В.В., Смолин А.Ю., Добрынин С.А., Коростелев С.Ю., Псахье С.Г. Изучение возможности идентификации наноскопичес-ких пор на основе анализа силы трения // Вестник ТГУ. Математика и механика. - 2010. - Т. 12. - № 4. - C. 116-122. Sergeev V.V., Smolin A.Yu., Dobrynin S.A., Korostelev S.Yu., Psakhie S.G. Possibilities of detection ofnanoscopic pores using the analysis of friction forces // Vest. TGU. Mat. Mekh. - 2010. - V. 12. -No. 4. - P. 116-122.

Поступила в редакцию 22.02.2014 г.

Сведения об авторах

Cмолин Алексей Юрьевич, д.ф.-м.н., проф., внс ИФПМ CO РАН, проф. ТГУ, asmolin@ispms.tsc.ru

Еремина Галина Максимовна, инж. ИФПМ CO РАН, асп. ТГУ, anikeeva@ispms.tsc.ru

Cергеев Валерий Викторович, мнс ИФПМ CO РАН, инж. ТПУ, svalera@ispms.tsc.ru

Шилько Евгений Викторович, д.ф.-м.н., доц., внс ИФПМ CO РАН, проф. ТГУ, shilko@ispms.tsc.ru

Псахье Cергей Григорьевич, д.ф.-м.н., чл.-корр. РАН, дир. ИФПМ CO РАН, проф. ТГУ, зав. каф. ТПУ, sp@ispms.tsc.ru