CC BY
18
УДК 539.3 : 621.7
ИССЛЕДОВАНИЕ ДЕФОРМИРОВАННОГО СОСТОЯНИЯ а-ЖЕЛЕЗА МЕТОДОМ МОЛЕКУЛЯРНОЙ ДИНАМИКИ
НАГОРНЫХ И.Л., БУРНЫШЕВ И.Н., БЕСОГОНОВ В.В.
Институт прикладной механики УрО РАН, 426067, г. Ижевск, ул.Т.Барамзиной, 34
АННОТАЦИЯ. В ходе работы методом молекулярной динамики с использованием ЕАМ потенциалов межатомного взаимодействия получены зависимости механических напряжений от степени деформации идеальных кристаллов и кристаллов, содержащих трехмерный дефект, при комнатной температуре. Деформация проведена в различных кристаллографических направлениях с использованием различных потенциалов и программ. Показано, что объемный дефект кристаллического строения может служить зародышем образования двойника, и, как следствие, приводить к пластической деформации при меньших (чем для бездефектных кристаллов) механических напряжениях.
КЛЮЧЕВЫЕ СЛОВА: моделирование, молекулярная динамика, деформация, предел прочности, плоскость скольжения, двойникование.
ВВЕДЕНИЕ
В мире науки всегда наблюдался повышенный интерес к железу, как основному конструкционному материалу. Однако в поведениях систем на основе Fe остаётся много вопросов, в частности, изменение его свойств при деформировании. На настоящий момент нами сделана попытка изучения возникающих механических напряжений в атомных плоскостях при деформации систем на основе ОЦК-железа. Целью настоящей работы явилось получение диаграмм деформации кристаллов железа при использовании различных программ и потенциалов межатомного взаимодействия. В последующих работах будет произведена оценка изменения полученных напряжений при введении в систему различных дефектов кристаллического строения, в том числе примесей.
Метод молекулярной динамики (ММД) представляет собой мощный инструмент при исследовании свойств деформированных кристаллов. В настоящей работе для описания межатомного взаимодействия использован метод погруженного атома (the embedded-atom method, EAM), предложенный [1], [2]. Суть метода погруженного атома состоит в том, что каждый атом системы рассматривается как примесь, погруженная в ансамбль, состоящий из всех оставшихся атомов. Энергия, необходимая для погружения, принимается зависящей от электронной плотности в точке погружения.
Основное преимущество потенциалов в приближении метода погруженного атома по сравнению с потенциалами парного взаимодействия заключается в учете многочастичного взаимодействия, которое становится важным при наличии дефектов кристаллического строения, особенно при наличии химически активных примесей, таких как, например, водород.
Основные недостатки ЕАМ состоят в том, что, во-первых, электронная плотность в точке представляется как суперпозиция электронных плотностей каждого атома, во-вторых, взаимодействие предполагается сферически симметричным, что для переходных металлов, таких как Fe, является серьезным допущением.
ВЫЧИСЛИТЕЛЬНЫЙ ЭКСПЕРИМЕНТ
Расчеты проведены с использованием программных пакетов метода молекулярной динамики «LAMMPS» (21 May 2008) [3], [4] и «MDP», разработанного авторами данной статьи. Компоненты тензора напряжений рассчитывались по теореме вириала по формуле:
х(Г1 - гОР , Р=х, у, Z,
где Уа - средний объем, приходящийся на одну частицу; к - постоянная Больцмана; Тхуг - абсолютная температура, при вычислении которой учитываются только компоненты
\Х,у,2
скорости вдоль выбранного направления; N - число частиц в ансамбле; ¥у , (- о)
х, у 2
- проекции сил парного взаимодействия и расстояний между частицами на выбранное направление х, у или z прямоугольной декартовой системы координат.
Потенциальные функции взяты из [5] и [6]. С использованием данных [5] Веном и др. [7] сформирован набор потенциальных функций для системы РеОЦК-Н, который будет использован нами в дальнейшем при исследовании деформированных систем железо-водород. Источник [6] представляет собой самую последнюю доступную разработку (насколько нам известно) ЕАМ потенциалов для железа. Полуэмпирические потенциалы [6] удовлетворительно подогнаны, согласно авторам [6], под свойства низкотемпературной ОЦК модификации железа (а-железо) и более высокотемпературной ГЦК модификации (у-железо). Это важно, поскольку в деформированных системах на основе Fe могут наблюдаться одновременно фазы а-Ре и у-Ре [8].
Для проведения расчетов по деформированию вдоль направления <0\0> сформированы следующие системы:
\. Идеальный кристалл ОЦК-железа размером \5*\5х\5 в числах параметра решетки вдоль каждого направления трехмерной прямоугольной декартовой системы координат (рис. \, а). Граничные условия: периодические по всем направлениям. Число частиц: 6750.
2. Кристалл ОЦК-железа размером \5*\5*\5, в котором отсутствуют части плоскостей (0\0), то есть имеется объемный дефект кристаллического строения (рис. \, б). Граничные условия: периодические по всем направлениям. Число частиц: 6390. При соблюдении периодических граничных условий дефект имеет смысл протяженной поры.
а) б)
а) - идеальный кристалл б) - кристалл, содержащий Зх-мерный дефект
Рис.1. Вид сформированных кристаллов для вычислительного эксперимента.
Проекция на плоскость (001). Стрелками показаны направления растягивающих деформаций
Все расчеты проведены в приближении МРТ-ансамбля, в котором число частиц N механическое напряжение вдоль направлений, отличных от направления деформации (для процесса деформирования), и температура остаются постоянными. Параметр решетки: 0,286 нм.
ИССЛЕДОВАНИЕ ДЕФОРМИРОВАННОГО СОСТОЯНИЯ а-ЖЕЛЕЗА МЕТОДОМ МОЛЕКУЛЯРНОЙ _____________________________________ДИНАМИКИ_______________________________________
Для интегрирования использовалась схема Верле с шагом по времени 0,\ фс. Для поддержания необходимых температуры и давления в программе «МИР» использованы термостат и баростат Берендсена, в программе «LAMMPS» - термостат Нуза-Хувера и баростат Берендсена. Термодинамические параметры систем записывались через каждые 50 молекулярно-динамических (МД) шагов расчета. Все расчеты проведены по следующей схеме:
\. Релаксация (оптимизация) систем в течение \0000 МД-шагов методом молекулярной статики (на каждом МД-шаге производится обнуление скоростей частиц) для получения структуры с минимальной энергией при Т ~ 0 К.
2. Нагрев системы до Т = 300 К, затем выдержка в течение 30000 МД-шагов. Выдержка необходима для достижения системой состояния равновесия. Считается, что система достигает состояния равновесия, если ее параметры релаксируют к средним значениям. В качестве таких параметров использовались параметр решетки, объем, давление, температура, полная энергия. Результаты проведенных вычислительных экспериментов показывают, что выдержка в течение 30000 МД-шагов является достаточной для достижения системой состояния равновесия.
3. Фиксация размеров кристалла вдоль направления <0\0> (на рис. \ показано стрелками). Выдержка в течение 30000 МД-шагов. Фиксация необходима для оценки равновесного значения механического напряжения на плоскости (0 \ 0) недеформированного кристалла. Заметим, что линейные размеры моделируемой системы вдоль направлений <\00> и <00\> изменялись таким образом, что механическое напряжение на плоскости (\00) и (00\) составляло ~ 0 МПа.
4. Деформация кристаллов на \ % вдоль направления <0\0> с последующей фиксацией размера в течение 30000 МД-шагов. Деформация проводилась до \2 %. Выдержка необходима для релаксации.
Деформация происходила путем изменения линейного размера расчетной области в направлении <0\0> и соответствующего масштабирования координат атомов.
Аналогичные эксперименты проведены с направлением растяжений вдоль направлений <\ \ \>, <\2\>, <\0\>.
ОБСУЖДЕНИЕ РЕЗУЛЬТАТОВ
На рис. 2, 3 приведены полученные зависимости механического напряжения на плоскости (0 \ 0) (перпендикулярные направлению деформации) от степени деформации для идеального кристалла и кристалла, содержащего объемный дефект (пору) с использованием различных потенциалов и программ. Приведенные зависимости можно интерпретировать как диаграммы деформации с обратным знаком, то есть механические напряжения, которые необходимо приложить для того, чтобы, во-первых, для случая линейной зависимости, когда выполняется закон Гука (меньше 8 %), деформировать кристалл до соответствующего линейного размера, во-вторых, удержать систему в деформированном неизменном состоянии.
Примем, что предел прочности соответствует минимуму напряжения на диаграммах на рис. 2, 3. Усредненный по направлениям предел прочности для идеального кристалла при использовании потенциалов [5] составляет « \8 ГПа, при использовании потенциалов [6] « \4 ГПа. Теоретические оценки предела прочности бездефектного кристалла и экспериментальные результаты, полученные на нитевидных кристаллах железа, показывают значения (\3,5-\5) ГПа [9], что находится в согласии с результатами, полученными моделированием в настоящей работе. Данный факт отчасти свидетельствует о корректном поведении модели, сформированной для моделирования процесса деформирования.
pressure, MPa pressure, MPa
0
-2000
-4000
-6000
-8000
0
-1000
-2000
-3000
-4000
-5000
-6000
-7000
Simulation results
Simulation results
deformation, %
a)
Simulation results
deformation, %
deformation, %
б)
в)
а) - потенциалы [6], программа МИР; б) - потенциалы [6], программа ЬАММР8; в) - потенциалы [5], программа ЬАММР8
Рис. 2. Механические напряжения вдоль направления деформации <010>, возникающие при деформирования идеального кристалла
Simulation results
Simulation results
deformation, %
а)
deformation, %
б)
а) - потенциалы [6], программа ЬАММР8; б) - потенциалы [5], программа ЬАММР8
Рис. 3. Механические напряжения вдоль направления деформации <010>, возникающие при деформировании кристалла, содержащего пору
В процессе расчета проведены наблюдения за изменением структур систем. При степени деформации, превышающей линейную область (более 8 %) (см. рис. 2, 3) обнаружено образование двойниковых структур. Уменьшение напряжения по абсолютной величине на каждой из диаграмм напряжение-деформация сопровождается выделением двойников, что является признаком пластического характера деформирования. При использовании потенциалов [6] наблюдается появление двойников в плоскостях легкого скольжения (семейства плоскостей {\\0}, {2\\}) с границей двойников в направлениях <\\\> (рис. 4, а), что согласуется с представлениями о плоскостях и направлениях легкого скольжения в кристаллах с оцк структурой [\0]. При использовании потенциалов [5] наблюдается слабовыраженная двойниковая структура преимущественно в плоскостях {\00}.
а) б)
а) - плоскость (110). Деформация 10 %. Образования двойниковой прослойки в плоскости легкого скольжения (\ Ю) с границей двойника в направлении <\ \ \> (плоскость двойникования (\ Ш));
б) - плоскость (0Ю). Деформация 9 %.
Данные получены при использовании потенциалов [6], программы “ЬАММР8”
Рис. 4. Примеры образования двойниковых прослоек.
Наличие дефекта приводит к наступлению более ранней пластической деформации
При моделировании деформированного состояния систем, содержащих объемный дефект кристаллического строения, получены следующие результаты. При наличии объемного дефекта пластическая деформация наступает при деформациях, меньших, чем для случая идеальных кристаллов. Причина такого поведения заключается в том, что пора представляет собой зародыш образования двойника (рис. 4, б). Очевидно также снижение значений механических напряжений в области до \0 % по сравнению с бездефектными кристаллами.
Аналогичный эксперимент проведен для идеального кристалла ОЦК-железа, содержащего отличное число атомов (2000) с использованием потенциалов [6]. Результаты полностью совпадают с результатами, приведенными на рис. 2, что свидетельствует о независимости исхода эксперимента от числа частиц в идеальном кристалле.
При растяжении систем вдоль направлений <\\\>, <\20>, <\0\> качественный вид диаграмм деформации остается таким же, что и для случая деформирования вдоль направления <0\0>.
ЗАКЛЮЧЕНИЕ
Полученные в работе величины механических напряжений совпадают с теоретически предсказанными для идеальных кристаллов. Необходимо отметить, что, вероятно, при определенной степени деформации важным компонентом энергии кристалла становится энергия групп атомов, зависящая от углов между связями. В этом случае сферическая симметрия взаимодействия, лежащая в основе модели, становится некорректной. Использованные потенциалы не учитывают угловую зависимость энергий. Оценка данного вклада в энергию кристалла не проводилась.
СПИСОК ЛИТЕРАТУРЫ
1. Daw M.S., Baskes M.I. Semiempirical, quantum mechanical calculation of hydrogen embrittlement in metals // Phys. Rev. Lett. 1983. V. 50. P. 1285.
2. Daw M.S., Baskes M.I. Embedded-atom method: Derivation and application to impurities, surfaces, and other defects in metals // Phys. Rev. 1984. № 29 B. P. 6443.
3. Plimpton S.J. Fast parallel algorithms for short-range molecular dynamics // J. Comp. Phys. 1995. №1. P. 117.
4. URL: LAMMPS WWW Site: http://lammps.sandia.gov/
5. Johnson R.A., Oh D.J. Analytic embedded atom method model for bcc metals // J. Mater. Res. 1989. № 4 (5). P. 1195.
6. Chamati H., Papanicolaou N.I., Mishin Y., Papaconstantopoulos D.A. Embedded-atom potential for Fe and its application to self-diffusion on Fe(100) // Surface Science. 2006. № 600(9). P. 1793.
7. Wen M., Xu X-J., Fukuyama S. et.al.Yokogawa. Embedded-atom-method functions for the body-centered-cubic iron and hydrogen // J. Mater. Res. 2001. № 16. P. 3496.
8. Шабашов В.А., Коршунов Л.Г., Мукосеев А.Г. и др. Фазовые превращения в стали У13 при сильной холодной деформации // Проблемы нанокристаллических материалов. Екатеринбург : Изд-во УрО РАН, 2002. C. 111.
9. Сатель Э.А. Проблемы развития технологии машиностроения. М. : Машиностроение, 1968. 592 с.
10. Новиков И.И. Дефекты кристаллического строения металлов. М. : Металлургия, 1975. 208 с.
INVESTIGATION OF DEFORMED STATE OF ALPHA-IRON BY THE MOLECULAR DYNAMIC METHOD
Nagomykh I.L., Bumyshev I.N., Besogonov V.V.
Institute of Applied Mechanics Ural Branch of the Russian Academy of Sciences, Izhevsk, Russia
SUMMARY. We have tried to study properties of BCC-iron crystals in deformed states by the molecular dynamic method using different EAM potentials. In this work we have received relations between stresses and rates of deformation for perfect crystals and the crystals containing three-dimensional defect under the room temperature. Deformations have been applied in different crystallographic directions. It has been shown that the crystalline bulk defect appears be twin nucleation as well as leads to plastic deformation under smaller (than for perfect crystals) stresses.
KEYWORDS: simulation, molecular dynamics, deformation, tensile strength, slip plane, twining.
Нагорных Иван Леонидович, аспирант ИПМ УрО РАН
Бурнышев Иван Николаевич, кандидат технических наук, ведущий научный сотрудник ИПМ УрО РАН, e-mail inburn@mail.ru
Бесогонов Валерий Валентинович, кандидат технических наук, старший научный сотрудник ИПМ УрО РАН, тел. 8-3412-21-66-11
