РАЗДЕЛ 2.МЕТОДИКА КОМПЬЮТЕРНОГО МОДЕЛИРОВАНИЯ ДЕФОРМАЦИИ И РАЗРУШЕНИЯ В НАНООБЪЁМЕ
2.1 Общее описание модели
В силу малости размеров нанокристаллов, исследование прочностных свойств этих объектов связано со значительными методическими трудностями. Поэтому, в качестве основного инструмента использовался оптимальный для таких размеров метод молекулярной динамики [65, 66].
Такой вычислительный эксперимент предоставляет возможность дополнить экспериментальные исследования, и даже заменить их, когда данные недоступны для непосредственного измерения или явления очень сложны. Выполняемая на компьютере программа выступает одновременно и как прибор, и как исследуемая система. Свойства такой системы определяются экспериментатором (с учетом ограниченных возможностей компьютера). Такой прибор для вычислительного эксперимента обеспечивает полностью контролируемую среду. В отличие от физического эксперимента измерение физических величин не вызывает затруднений и возмущений, что особенно ценно при исследовании явлений на наноуроуровне.
2.1.1 Метод молекулярной динамики
Метод молекулярной динамики основан на представлении материала совокупностью взаимодействующих частиц (материальных точек или твердых тел), для которых записываются классические уравнения динамики. Взаимодействие частиц описывается посредством потенциалов взаимодействия, основным свойством которых является отталкивание при сближении и притяжение при удалении. Перед началом моделирования задается некоторое начальное распределение частиц в пространстве (исходная структура материала) и начальное распределение скоростей частиц (механическое и тепловое движение системы в исходном состоянии). Далее задача сводится к решению задачи Коши для системы обыкновенных дифференциальных уравнений.
Для простых металлов можно считать атомы материальными точками, таким образом имеем систему из 6N дифференциальных уравнений первого порядка [66-67]:
(2.1)
здесь: i - номер атома; N - количество атомов в системе; x - радиус вектор; p - импульс; m - масса частицы; f - сила; V - потенциал межатомного взаимодействия.
2.1.2 Выбор программного обеспечения
Для проведения расчётов выбор был остановлен на программном пакете XMD [68], разработанном Джоном Рифкиным в Центре Моделирования материалов Коннектикутского Университета специально для моделирования металлов и керамики методом молекулярной динамики. Программа распространяется в исходных кодах, благодаря чему есть возможность производить в них изменения и добавлять собственный программный код.
Также, по сравнению с другими распространяющимися программами, данный пакет наименее зависит от программной и аппаратной платформы, что позволило максимально использовать имеющиеся вычислительные мощности.
2.1.3 Алгоритм Гира
Для численного решения уравнений движения XMD использует предикторно-корректорный алгоритм Гира [69-70]. При использовании такого алгоритма хранятся значения координаты, а также пяти её производных для предикторного и корректорного этапов.
Наиболее удобно рассматривать алгоритм Гира в масштабированных шагом времени скоростях, ускорениях и т.д. Так r0 - координата частицы, то масштабированная скорость - r1=?t(dr0/dt); ускорение - r2=1/2 ?t2(d2r0/dt2); r3=1/6 ?t3(d3r0/dt3) и так далее.
На первом шаге итерации вычисляются предварительные (предикторные) значения координаты и производных:
(2.2)
где rp - предикторные значения; такая форма значительно облегчает применение метода, здесь используется матрица Паскаля. Таким образом, когда нам известны предварительные значения координат, мы можем рассчитать силу, которая действует на атом в момент времени (t+?t) и найти точное значение ускорения , а масштабированная по шагу времени производная:
(2.3)
Тогда, учитывая эту ошибку , откорректируем найденные ранее значения:
(2.4)
Сi - матрица корректорных коэффициентов, она имеет вид:
.
Таким образом, использование предикторно-корректорного алгоритма позволяет увеличить сходимость и уменьшить погрешность расчётов, которая возникает только за счет пренебрежения производными координаты шестого и больше порядков.
2.1.4 Параметры термостата
Для поддержания постоянной температуры в модельном объёме, применялся метод масштабирования скорости атомов на каждом шаге времени:
(2.5)
значение выбрано равным 33. Для системы материальных точек, каковой являются простые металлы, такой алгоритм является адекватным[67]. Выбранное значение меньше единицы, принятой в молекулярной статике. Это связано с тем, что использование предикторно-корректорного алгоритма позволяет избежать дрейфа температуры, связанного с погрешностью расчётов (рис. 2.1).
Рис. 2.1. Зависимость температуры системы из 10000 атомов железа от продолжительности расчётов при отсутствии внешнего воздействия без применения алгоритма термостата.
С другой стороны, такое значение более плавно релаксирует температуру, адекватно моделирует процесс остывания и позволяет избежать нежелательного эффекта "замораживания" скоррелированого движения атомов (образования и движения дефектов, разрушения) обусловленного тем, что температура рассчитывается на основании кинетической энергии поступательного движения атомов в лабораторной системе координат (рис. 2.2). Следует заметить, что такой выбор приводит к значительному времени первоначальной релаксации температуры - 2·10-11 с, а поскольку при всех расчётах величина шага времени =10-15с, то это составляло 2·104 итераций.
Рис. 2.2. Зависимость температуры системы из 5795 атомов железа от продолжительности расчётов при одноосном растяжении. Первый пик соответствует процессу первоначальной релаксации, второй -