Научная статья на тему 'МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ФИЗИЧЕСКИХ ПРОЦЕССОВ В МЕТАЛЛИЧЕСКИХ КРИСТАЛЛИЧЕСКИХ РЕШЕТКАХ'

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ФИЗИЧЕСКИХ ПРОЦЕССОВ В МЕТАЛЛИЧЕСКИХ КРИСТАЛЛИЧЕСКИХ РЕШЕТКАХ Текст научной статьи по специальности «Математика»

CC BY
234
32
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ / ПОТЕНЦИАЛ ЛЕННАРДА-ДЖОНСА / ПАКЕТ ПРОГРАММ MATLAB / ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ / АТОМ / КРИСТАЛЛИЧЕСКАЯ РЕШЕТКА / МЕТАЛЛЫ И УПОРЯДОЧЕННЫЕ СПЛАВЫ / ЭНЕРГИЯ ВЗАИМОДЕЙСТВИЯ / АМПЛИТУДНО-ЧАСТОТНАЯ ХАРАКТЕРИСТИКА / РАЗНОСТНАЯ АППРОКСИМАЦИЯ

Аннотация научной статьи по математике, автор научной работы — Семёнова Мария Николаевна, Семёнов Александр Сергеевич, Бебихов Юрий Владимирович, Якушев Илья Анатольевич

Колебания кристаллических решеток определяют такие важные свойства материалов, как теплопроводность, теплоемкость, тепловое расширение и многие другие, поэтому их изучение является актуальной и важной задачей. Наряду с экспериментальными методами изучения нелинейной динамики кристаллической решетки широко используются эффективные методы компьютерного моделирования, такие как первопринципное моделирование и метод молекулярной динамики. Реже используется метод математического моделирования, так как погрешность расчетов при его применении может достигать 10%. Вместе с тем он является наименее затратным в вычислительном плане. В настоящей работе представлен процесс и описаны результаты математического моделирования физических процессов в металлах и упорядоченных сплавах при помощи потенциала Леннарда-Джонса в пакете программ MatLab, хорошо зарекомендовавшего себя для решения задач технических вычислений. В теоретической части описаны кристаллические решетки, дифференциальные уравнения движения для моделирования, их начальные и граничные условия, а также разностная аппроксимация. В качестве метода моделирования выбран принцип молекулярно-динамического моделирования при помощи одного из парных потенциалов. В практической части представлены вычислительный алгоритм, линеаризация числа операций, термодинамические расчеты, описаны скоростная схема Верле и потенциал межатомного взаимодействия, а также приведена поэтапная разработка модели в среде моделирования. Получены следующие результаты: построен график трехмерного распределения атомов по расчетной ячейке, доказывающий возможность перемещения до пяти межатомных расстояний; произведена оценка амплитудно-частотной характеристики методом Уэлча с относительной среднеквадратической ошибкой, не превышающей 30%; получена графическая зависимость энергий связи между настоящей моделью и стандартными данными для ячейки ГПУ металла с погрешностью немногим более 3%; произведено вычисление оптимальной модели для кусочно-линейной аппроксимации, и построена её трехмерная интерполяция. Все проведенные исследования показывают хорошую степень применимости математического моделирования к задачам изучения динамических процессов в физике кристаллов.

i Надоели баннеры? Вы всегда можете отключить рекламу.

Похожие темы научных работ по математике , автор научной работы — Семёнова Мария Николаевна, Семёнов Александр Сергеевич, Бебихов Юрий Владимирович, Якушев Илья Анатольевич

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

MATHEMATICAL MODELING OF PHYSICAL PROCESSES IN METAL CRYSTAL LATTICES

Vibrations of crystal lattices determine such important properties of materials as thermal conductivity, heat capacity, thermal expansion, and many others; therefore, their study is an urgent and important task. Along with experimental methods of studying the nonlinear dynamics of a crystal lattice, effective methods of computer modeling, such as ab initio modeling and the method of molecular dynamics, are widely used. The method of mathematical modeling is used less often, since the calculation error during its application can reach 10%. At the same time, it is the least computationally expensive method. This paper presents the process and describes the results of mathematical modeling of physical processes in metals and ordered alloys using the Lennard-Jones potential in the MatLab software package, which is well-proven for solving technical computing problems. The theoretical part describes crystal lattices, differential equations of motion for modeling, their initial and boundary conditions, as well as difference approximation. The principle of molecular dynamics modeling using one of the paired potentials was chosen as a modeling method. In the practical part, a computational algorithm, linearization of the number of operations, thermodynamic calculations are presented, the Verlet velocity scheme and the potential of interatomic interaction are described, as well as a step-by-step development of a model in a modeling environment. The following results were obtained: a graph of the three-dimensional distribution of atoms over the computational cell was constructed, which proves the possibility of moving up to five interatomic distances; the estimation of the amplitude-frequency characteristic by the Welch method with the relative root-mean-square error not exceeding 30% was made; a graphical dependence of the binding energies between the present model and the standard data for a metal hcp cell was obtained with an error of slightly more than 3%; the calculation of the optimal model for piecewise-linear approximation was made and its three-dimensional interpolation was constructed. All conducted studies show a good degree of applicability of mathematical modeling to problems of studying dynamic processes in crystal physics.

Текст научной работы на тему «МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ФИЗИЧЕСКИХ ПРОЦЕССОВ В МЕТАЛЛИЧЕСКИХ КРИСТАЛЛИЧЕСКИХ РЕШЕТКАХ»

ФИЗИКО-МАТЕМАТИЧЕСКИЕ НАУКИ

УДК 531.37

М. Н. Семёнова, А. С. Семёнов, Ю. В. Бебихов, И. А. Якушев

Математическое моделирование физических процессов в металлических кристаллических решетках

СВФУ им. М.К. Аммосова, Политехнический институт (филиал), г. Мирный, Россия

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

СЕМЁНОВА Мария Николаевна - старший преподаватель кафедры ФиПМ, Политехнический институт (филиал) СВФУ им. М.К. Аммосова.

E-mail: mariya_semyonova86@mail.ru

SEMENOVA Mariya Nikolaevna - Senior lecturer of the Department of fundamental and applied mathematics, Polytechnic Institute (branch), M.K. Ammosov North-Eastern Federal University.

Семёнов Александр Сергеевич - к. ф.-м. н., доцент, заведующий кафедрой ЭиАПП, Политехнический институт (филиал) СВФУ им. М.К. Аммосова.

E-mail: sash-alex@yandex.ru

SEMENOV Alexander Sergeevich - Candidate of Physics and Mathematics, Associate Professor, Head of the Department of electrification and automation of industrial production, Polytechnic Institute (branch), M.K. Ammosov North-Eastern Federal University.

БЕБИХОВ Юрий Владимирович - к. ф.-м. н., доцент кафедры ЭиАПП, Политехнический институт (филиал) СВФУ им. М.К. Аммосова.

E-mail: bebikhov.yura@mail.ru

BEBIKHOV Yuriy Vladimirovich - Candidate of Physics and Mathematics, Associate Professor of the Department of electrification and automation of industrial production, Polytechnic Institute (branch), M.K. Ammosov North Eastern Federal University.

ЯКУшЕВ Илья Анатольевич - к. ф.-м. н., доцент кафедры ФиПМ, Политехнический институт (филиал) СВФУ им. М.К. Аммосова.

E-mail: yakushevilya@mail.ru

YAKUSHEV Ilya Anatolievich - Candidate of Physics and Mathematics, Associate Professor of the Department of fundamental and applied mathematics, Polytechnic Institute (branch), M.K. Ammosov NorthEastern Federal University.

результаты математического моделирования физических процессов в металлах и упорядоченных сплавах при помощи потенциала Леннарда-Джонса в пакете программ MatLab, хорошо зарекомендовавшего себя для решения задач технических вычислений. В теоретической части описаны кристаллические решетки, дифференциальные уравнения движения для моделирования, их начальные и граничные условия, а также разностная аппроксимация. В качестве метода моделирования выбран принцип молекулярно-динамического моделирования при помощи одного из парных потенциалов. В практической части представлены вычислительный алгоритм, линеаризация числа операций, термодинамические расчеты, описаны скоростная схема Верле и потенциал межатомного взаимодействия, а также приведена поэтапная разработка модели в среде моделирования. Получены следующие результаты: построен график трехмерного распределения атомов по расчетной ячейке, доказывающий возможность перемещения до пяти межатомных расстояний; произведена оценка амплитудно-частотной характеристики методом Уэлча с относительной среднеквадратической ошибкой, не превышающей 30%; получена графическая зависимость энергий связи между настоящей моделью и стандартными данными для ячейки ГПУ металла с погрешностью немногим более 3%; произведено вычисление оптимальной модели для кусочно-линейной аппроксимации, и построена её трехмерная интерполяция. Все проведенные исследования показывают хорошую степень применимости математического моделирования к задачам изучения динамических процессов в физике кристаллов.

Ключевые слова: математическое моделирование, потенциал Леннарда-Джонса, пакет программ MatLab, дифференциальные уравнения, атом, кристаллическая решетка, металлы и упорядоченные сплавы, энергия взаимодействия, амплитудно-частотная характеристика, разностная аппроксимация.

DOI 10.25587/SVFU.2021.84.4.002

Благодарности

Семёнов А. С. (моделирование динамических процессов в металлах) выражает благодарность Российскому научному фонду, грант № 21-12-00275.

M. N. Semenova, A. S. Semenov, Yu. V. Bebikhov, I. A. Yakushev

Mathematical modeling of physical processes in metal crystal lattices

M.K. Ammosov North Eastern Federal University, Polytechnic Institute (branch), Mirny, Russia

Abstract. Vibrations of crystal lattices determine such important properties of materials as thermal conductivity, heat capacity, thermal expansion, and many others; therefore, their study is an urgent and important task. Along with experimental methods of studying the nonlinear dynamics of a crystal lattice, effective methods of computer modeling, such as ab initio modeling and the method of molecular dynamics, are widely used. The method of mathematical modeling is used less often, since the calculation error during its application can reach 10%. At the same time, it is the least computationally expensive method. This paper presents the process and describes the results of mathematical modeling of physical processes in metals and ordered alloys using the Lennard-Jones potential in the MatLab software package, which is well-proven for solving technical computing problems. The theoretical part describes crystal lattices, differential equations of motion for modeling, their initial and boundary conditions, as well as difference approximation. The principle of molecular dynamics modeling using one of the paired potentials was chosen as a modeling method. In the practical part, a computational algorithm, linearization of the number of operations, thermodynamic calculations are presented, the Verlet velocity scheme and the potential of interatomic interaction are described, as well as a step-by-step development of a model in a modeling environment. The following results were obtained: a graph of the three-dimensional distribution of atoms over the computational cell was constructed, which proves the possibility of moving up to five interatomic distances; the estimation of the amplitude-frequency characteristic by the Welch method

with the relative root-mean-square error not exceeding 30% was made; a graphical dependence of the binding energies between the present model and the standard data for a metal hcp cell was obtained with an error of slightly more than 3%; the calculation of the optimal model for piecewise-linear approximation was made and its three-dimensional interpolation was constructed. All conducted studies show a good degree of applicability of mathematical modeling to problems of studying dynamic processes in crystal physics.

Keywords: mathematical modeling, Lennard-Jones potential, MatLab software package, differential equations, atom, crystal lattice, metals and ordered alloys, interaction energy, amplitude-frequency characteristic, difference approximation.

Введение

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

Кристаллические тела (все металлы и металлические сплавы) характеризуются упорядоченным расположением атомов. В металлах и металлических сплавах атомы находятся в узлах пространственных кристаллических решеток. В процессе кристаллизации металлов и сплавов могут образовываться кристаллические решетки разного типа [2]. Большинству металлов свойственно образование высокосимметричных решеток с плотной упаковкой атомов [3]. К основным типам кристаллических решеток металлов относятся: объемно-центрированная кубическая (ОЦК), гранецентрированная кубическая (ГЦК), гексагональная плотноупакованная (ГПУ). Элементарные ячейки таких кристаллических решеток показаны на рис. 1.

Наряду с экспериментальными методами изучения нелинейной динамики кристаллической решетки широко используются эффективные методы компьютерного моделирования, такие как первопринципное моделирование и метод молекулярной динамики [4]. Реже используется метод математического моделирования, так как погрешность расчетов при его применении может достигать 10%. Вместе с тем он является наименее затратным в вычислительном плане. В настоящее время математическое моделирование рассматривается как мощное средство теоретических исследований нелинейных проблем в физике [5].

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

Рис. 1. Элементарные ячейки объемно-центрированной кубической (а), гранецентрированной кубической (б) и гексагональной плотноупакованной (в) кристаллических решеток

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

В настоящей работе авторами впервые будет применен пакет программ МаЛаЬ для моделирования перемещения атомов в трехмерной кристаллической решетке при помощи потенциала взаимодействия Леннарда-Джонса. Также новизной будет являться определение зависимости энергии связи в исследуемой модели и её сравнение с данными, полученными авторами ранее для ячейки ГПУ металла Ве, подтверждающими высокую точность моделирования с погрешностью немногим > 3%.

Материалы и методы

Математическая модель состоит из системы дифференциальных уравнений, их разностного аналога, потенциала межатомного взаимодействия, специфически определяемых начальных и граничных условий [6].

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

du

т.

dt

= F + Fe

dri dt

(1)

= и

где mi, Ц, Г - масса, скорость и радиус-вектор соответственно; Fi - сила взаимодействия с остальными частицами; Fi - сила взаимодействия с внешними полями, i = 1... N ; N - количество точечных частиц.

Зная координаты и скорости всех N частиц в начальный момент времени t = 0: ( )=0 , можем произвести интегрирование системы уравнений (1).

Поскольку в начальный момент времени моделируемая среда представляет собой кристалл, подготовку начальных условий будем производить в три этапа: 1) расчет приблизительного значения постоянной решетки аТ при температуре Т; 2) задание

координаты частиц; 3) задание скорости частиц в соответствии с распределением Максвелла.

При молекулярно-динамическом моделировании (МДМ) необходимо рассчитать траектории огромного количества частиц, обеспечив как вычислительную эффективность, так и требования точности. Для этого применим алгоритм Верле в скоростной форме [8], формулируемый из условий известности всех величин на момент времени tk с переходом к следующему моменту времени +1 = + At по формулам:

-к+-

= V

1 -М 2

ак м

-к+-

-к+1

. , 2м

grad (ик+1)

(2)

т

+ -

ак м

1+

2

где верхний индекс к - номер шага по времени; и - потенциальная энергия взаимодействия системы.

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

Ключевым элементом МДМ является выбор потенциала межатомных и межмолекулярных взаимодействий частиц. Имеются два основных подхода: первопринципный [9, 10] и метод полуэмпирических потенциалов [11, 12]. Первый основан на развитии идей квантовой механики, в том числе на решении уравнения Шредингера для конденсированной среды. Второй оценивает способность того или иного потенциала хорошо описывать соединения с различным атомным окружением.

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

((г ) = 4е

^12

V г )

^6

V г )

(3)

где а - значение межатомного расстояния, при котором = 0, е - глубина

потенциальной ямы, расположенной на расстоянии а

. В

е

этом выражении слагаемое ~ г"6 доминирует на больших расстояниях и соответствует дисперсионному диполь-дипольному притяжению. Слагаемое ~ г~12 моделирует сильное отталкивание между парой атомов за счет обменного взаимодействия, если они оказываются очень близко друг к другу.

Характерный вид потенциала Леннарда-Джонса показан на рис. 2, его минимум лежит в точке г . = оЩ2 .

тт *

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

2

2

Рис. 2. Характерный вид потенциала Леннарда-Джонса

свойств одномерных нелинейных колебательных мод в треугольной решетке с Леннард-Джонсовскими взаимодействиями [14].

Результаты

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

В вычислительном плане наиболее затратным мог оказаться первый этап, так как для вычисления потенциалов необходимо выполнить в общем случае N арифметических операций, где N - число составляющих систему частиц. Однако выбранный нами потенциал Леннарда-Джонса является короткодействующим. Это означает, что существует некоторое критическое расстояние гсг, такое, что значениями потенциала для г > гсг можно пренебречь. Благодаря этому был построен алгоритм расчета с линейно растущим числом частиц, который раз в М шагов по времени пересчитывает количество «соседей» во взаимодействующей среде. Соседями при этом считаем пару взаимодействующих атомов, находящуюся друг от друга на расстоянии не более /тах = гсг + А/. Таким образом, линеаризованный алгоритм и очевидность скоростной схемы Верле позволили нам производить ряд вычислений параллельно.

Вычисление основных термодинамических величин (давления и температуры) при МДМ производили, исходя из быстрого установления локального термодинамического равновесия, время которого может составлять единицы пикосекунд. Температуру со степенью свободы а = X, у, z получили, усредняя хаотическую составляющую кинетической энергии по ансамблю частиц с установившимся равновесным распределением. После чего общую температуру определили как среднюю по степеням свободы (4). Давление рассчитали аналогично путем усреднения по объему и по времени (5).

тх+Т + т

Т -У-£ , (4)

3

р + р + р

р _ хх уу гг (5)

" 3 '

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

С помощью пакета программ Ма1ЬаЬ, хорошо зарекомендовавшего себя для решения задач технических вычислений, было смоделировано поведение (хаотическое движение) п-го количества атомов в кристаллической решетке методом молекулярной динамики. Ранее авторы для подобных целей использовали непосредственно программы МДМ, к примеру, свободно распространяемый пакет LAMMPS [15, 16]. Для реализации расчетов выбрали скоростную схему Верле. В общем виде формулы для траектории (6), скорости (7), ускорения (8), скорости с новым ускорением (9) выглядят следующим образом:

г ( + М ) = г () + vAt + аМ2/2, (6)

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

у ( + Дt /2 ) = V + aAt /2, (7)

а () = -grad (и (г ())))2, (8)

V ( + At) = V (t + At/2) + aAt/2. (9)

Здесь в формуле ускорения (8) и - потенциал Леннарда-Джонса (3), который был выбран из множества других парных потенциалов, как оптимально подходящий для условий нашей задачи.

Моделирование производилось в трехмерном (3D) пространстве, поэтому градиент в формуле (8) принимался по трем направлениям. Преобразование потенциала и (г (?) в и (х (?), у (7), г () осуществлялось по формуле:

гуо = XI + yj + zk. (10)

В качестве исходных данных для 3D моделирования были выбраны 125 атомов (Ы) с размерами кристаллической решетки (ячейки) L х W х Н = 5 х 5 х 5, структура которой показана на рис. 3. К основным исходным данным для моделирования также относятся: Тк - время окончания расчетов; dt - шаг времени; г - расстояние между центрами частиц, е - глубина потенциальной «ямы», а - расстояние, на котором энергия взаимодействия становится равной нулю. Параметры е и а являются характеристиками атомов соответствующего вещества, которые можно найти через коэффициент Джоуля-Томсона, либо сравнивая экспериментальное значение коэффициента вязкости со значением, полученным из формулы для потенциальной энергии. Таблицы значений параметров е и а для некоторых веществ есть в книге [17]. Для ускорения расчётов будем обрывать потенциал Леннарда-Джонса на расстоянии гс = 2.5 • ст . Такой выбор обусловлен тем, что на этом расстоянии значение энергии взаимодействия составляет лишь -0,0163 от глубины ямы е.

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

Рис. 4. Выгрузка кода программы из окна редактора пакета программ Ма1ЬаЬ

ускорения; вычисление следующего цикла; повторная проверка граничных условий; вычисление скорости с новым ускорением. Вычислительный процесс может повторяться до бесконечности. Программа была разработана и сохранена в т-файле пакета программ Ма1ЬаЬ [18]. На рис. 4 представлен код программы с исходными данными и двумя циклами для моделирования в 3D пространстве.

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

Рис. 5. Трехмерное распределение атомов (а) и двумерный график разброса данных (б)

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

Приведем снимки с графопостроителя Plot и анализ-фильтров с вычислительными операциями в трехмерном пространстве, а именно: трехмерное распределение атомов при хаотическом перемещении внутри расчетной ячейки и двумерный график разброса данных в векторах x и y с отображением предельных распределений в виде одномерных гистограмм на горизонтальной и вертикальной осях графика рассеяния (рис. 5); сравнение энергий связи и амплитудно-частотные характеристики потенциала

Рис. 7. Результаты аппроксимации методом линейной интерполяции

взаимодействия, включая оценку спектральной плотности мощности по методу Уэлча (рис. 6); непараметрическую аппроксимацию с использованием сглаживающих сплайнов и линейным методом интерполяции (рис. 7).

Обсуждение

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

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

На рис. 6 (а) представлена амплитудно-частотная характеристика при нормализованной частоте (известна как частота Найквиста). При такой обработке непрерывная временная переменная t, выраженная в секундах, заменяется дискретной целочисленной переменной п с единицами отсчетов. В нашем случае временная переменная в секундах была нормализована (разделена) на интервал выборки Т (секунды/выборка), благодаря чему время приняло удобные целочисленные значения в моменты выборки. Максимальная величина на частотной области достигает 50 дБ, минимальная - немного не доходит до -10 дБ. Так как частота колебаний для каждого атома имеет разные значения, применим для оценки спектральной плотности мощности метод Уэлча. В результате моделирования 125 атомов (сэмплов) получили очень сильно «изрезанный» график (рис. 6 (в)), что обуславли-

вается недостаточно длительным интервалом наблюдения и приводит к высокой дисперсии оценки. Относительная среднеквадратическая ошибка (коэффициент утечки) мощности не превышает 30%.

При сравнении энергий связи, полученных в результате моделирования, со стандартными данными для ГПУ металлов мы увидели большую абсолютную и относительную погрешности при моделировании дальних взаимодействий. На рис. 6 (б) исследованы взаимодействия близлежащих атомов. По осям абсцисс и ординат откладываются квантили расчетного значения и смоделированного результата энергии связи для ГПУ металла. Так для ГПУ металла Be, который авторы изучали в работе [19], погрешность моделирования составила немногим более 3%. Таким образом, мы подтвердили, что потенциал Леннарда-Джонса хорошо применим для парных близлежащих взаимодействий. В результате моделирования получили положительное значение потенциала, выборочное среднее значение которого определили по формуле:

- X"и,

U = = 1.0978,

n

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

На рис. 7 (в) представлены результаты аппроксимации кривой из 125 точек, состоящей в определении функции z = f (x, у). Была найдена оптимальная общая модель с уравнением вида:

f (x,y) = а + р-sin(/л-п-x■ y) + y-e~(y) ,

где коэффициенты в векторной форме с 95% вероятности равны:

a =(2,376 2,961);в =(-0,6042 0,1607) = (-1,795 0,6895);^ = (0,9121 1,03)) = (-6,557 3,131).

На рис. 7 (а) показана линейная интерполяция, которая является частным случаем интерполяционной формулы Лагранжа и интерполяционной формулы Ньютона. Она применяется для приближенного представления данных в виде кусочно-линейной функции. На рис. 7 (б) показана её трехмерная проекция, вычисленная из плоскости с нормированными средним и стандартным значениями для x: 2,546 и 1,453; для y: 2,542 и 1,464.

Все проведенные исследования показывают хорошую степень применимости математического моделирования к задачам изучения динамических процессов в физике кристаллов. Разработанная программа математического моделирования была апробирована в учебном процессе в курсах физики и информатики [20, 21], а также получено свидетельство о государственной регистрации данной программы для ЭВМ [22].

Заключение

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

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

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

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

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

Основными результатами моделирования являются: трехмерное распределение атомов по расчетной ячейке, доказывающее возможность перемещения до пяти межатомных расстояний; оценка амплитудно-частотной характеристики методом Уэлча с относительной среднеквадратической ошибкой, не превышающей 30%; получение графической зависимости энергий связи между настоящей моделью и стандартными данными для ячейки ГПУ металла с погрешностью немногим более 3%; вычисление оптимальной модели (уравнений) для кусочно-линейной аппроксимации и построение её трехмерной интерполяции.

Л и т е р а т у ра

1. Борн, М. Динамическая теория кристаллических решеток / М. Борн, Х. Кунь. - Москва : Издательство иностранной литературы, 1958. - 488 с. - Текст : непосредственный.

2. Греков Ф.Ф. Структурная кристаллография / Ф. Ф. Греков, Г. Б. Рябенко, Ю. П. Смирнов.

- Ленинград : Издательство ЛГПИ, 1988. - 80 с. - Текст : непосредственный.

3. Гуляев, А. П. Металловедение / А. П. Гуляев. - изд. 6-ое. - Москва : Металлургия, 1986.

- 544 с. - Текст : непосредственный.

4. Хеерман, Д. В. Методы компьютерного эксперимента в теоретической физике / Д. В. Херман.

- Москва : Наука, 1990. - 176 с. - Текст : непосредственный.

5. Якушев, И. А. Математическое моделирование сложных технических систем в среде MATLAB / И. А. Якушев, М. Н. Семенова, Ю. В. Бебихов и др. - Москва : Издательство «Спутник +», 2019.

- 126 с. - Текст : непосредственный.

6. Мажукин, В. И. Температурная зависимость кинетической скорости плавления и кристаллизации алюминия / В. И. Мажукин, А. В. Шапранов, М. М. Демин и др. // Краткие сообщения по физике ФИАН. - 2016. - № 9. - С. 37-43. - Текст : непосредственный.

7. Маркидонов, А. В. Расчет термодинамических характеристик системы Fе-P методом молекулярной динамики / А. В. Маркидонов, Д. А. Лубяной, В. В. Коваленко и др. - Текст : непосредственный // Известия высших учебных заведений. Черная металлургия. - 2019. - Т. 62.

- № 9. - С. 725-731.

8. Verlet, L. Computer Experiments on Classical Fluids. I. Thermodynamical Properties of Lennard-Jones Molecules / L. Verlet // Phys. Rev. - 1967. - V. 159. - Pp. 98-103.

9. Ракитин, М. С. Первопринципное и термодинамическое моделирование влияния упругих напряжений на энергию растворения водорода в альфа-железе / М. С. Ракитин, А. А. Мирзоев, Д. А. Мирзаев. - Текст : непосредственный // Известия высших учебных заведений. Физика. - 2017.

- Т. 60. - № 12. - С. 76-82.

10. Боровик, А. М. Моделирование электрических характеристик графенового полевого транзистора на основе данных расчетов из первых принципов / А. М. Боровик, М. С. Баранова, Д. Ч. Гвоздовский. - Текст : непосредственный // Инфокоммуникационные и радиоэлектронные технологии. - 2020. - Т. 3. - № 1. - С. 63-74.

11. Zalizniak, V. E. Mathematical model of conducting nanopore for molecular dynamics simulations / V. E. Zalizniak, O. A. Zolotov, O. P. Zolotova // Siberian Journal of Science and Technology. - 2018.

- V. 19. - Is. 4. - Pp. 677-682.

12. Sanchez-Badillo, J. Potential of mean force calculations for an s„2 fluorination reaction in five

' N

different imidazolium ionic liquid solvents using quantum mechanics/molecular mechanics molecular dynamics simulations / J. Sanchez-Badillo, M. Gallo, R. A. Guirado-Lopez et al. // Journal of Physical Chemistry B. - 2020. - V. 124. - Is. 21. - Pp. 4338-4357.

13. Jones, J. E. On the Determination of Molecular Fields / J. E. Jones // Proceedings of the Royal Society of London. - 1924. - V. 106. - Is. 738. - Pp. 463-477.

14. Sunagatova, I. R. Properties of one-dimensional nonlinear vibrational modes in triangular lattice with Lennard-Jones interactions / I. R. Sungatova, A. M. Subkhangulova, M. N. Semenova et al. // IOP Conference Series: Materials Science and Engineering. - 2020. - V. 1008. - No. 012073.

15. Krylova K.A. Spherically localized discrete breathers in bcc metals V and Nb / K. A. Krylova, I. P. Lobzenko, A. S. Semenov et al. // Computational Materials Science. - 2020. - V. 180. - No. 109695.

16. Babicheva R. I. Discrete breathers in a triangular P-Fermi-Pasta-Ulam-Tsingou lattice / R. I. Babicheva, A. S. Semenov, E. G. Soboleva et al. // Physical Review E. - 2021. - V. 103. - Is. 5.

- No. 052202.

17. Цянь Сюэ-Сень. Физическая механика / Цянь Сюэ-Сень. - Москва : Мир, 1965. - 544 с.

- Текст : непосредственный.

18. Bebikhov, Y. V. The application of mathematical simulation for solution of linear algebraic and ordinary differential equations in electrical engineering / Y. V. Bebikhov, A. S. Semenov, I. A. Yakushev et al. // IOP Conference Series : Materials Science and Engineering. - 2019. - V. 643. - No. 012067.

19. Bachurina, O. V. Properties of Moving Discrete Breathers in Beryllium / O. V. Bachurina, R. T. Murzaev, A. S. Semenov et al. // Physics of the Solid State. - 2018. - V. 60. - Is. 5. - Pp. 989-994.

20. Татаринов, П. С. Сборник лабораторных работ по разделу «Механика» курса «Физика» / П. С. Татаринов, В. Д. Яковлев, Д. Ч. Ким и др. - Москва : Издательство «Спутник +», 2021.

- 175 с. - Текст : непосредственный.

21. Васильева, А. В. Практикум по работе с пакетами прикладных программ на ЭВМ / А. В. Васильева, М. Н. Семенова, И. А. Якушев. - Москва : Издательство «Спутник +», 2020. - 144 с. - Текст : непосредственный.

22. Семёнов, А. С. Программа математического моделирования физических процессов в металлах и упорядоченных сплавах / А. С. Семенов, М. Н. Семенова, И. А. Якушева и др. // Свидетельство № 2021615775 от 13.04.2021, заявка № 2021614436 от 02.04.2021. - Текст : непосредственный.

R e f e r e n c e s

1. Born, M. Dinamicheskaya teoriya kristallicheskih reshetok / M. Born, H. Kun'. - Moskva : Izdatel'stvo inostrannoj literatury, 1958. - 488 s. - Tekst : neposredstvennyj.

2. Grekov F.F. Strukturnaya kristallografiya / F. F. Grekov, G. B. Ryabenko, Yu. P. Smirnov. - Leningrad : Izdatel'stvo LGPI, 1988. - 80 s. - Tekst : neposredstvennyj.

3. Gulyaev, A. P. Metallovedenie / A. P. Gulyaev. - izd. 6-oe. - Moskva : Metallurgiya, 1986.

- 544 s. - Tekst : neposredstvennyj.

4. Heerman, D. V. Metody komp'yuternogo eksperimenta v teoreticheskoj fizike / D. V. Herman. - Moskva : Nauka, 1990. - 176 s. - Tekst : neposredstvennyj.

5. Yakushev, I. A. Matematicheskoe modelirovanie slozhnyh tekhnicheskih sistem v srede MATLAB / I. A. Yakushev, M. N. Semenova, Yu. V. Bebihov i dr. - Moskva : Izdatel'stvo «Sputnik +», 2019.

- 126 s. - Tekst : neposredstvennyj.

6. Mazhukin, V. I. Temperaturnaya zavisimost' kineticheskoj skorosti plavleniya i kristallizacii alyuminiya / V. I. Mazhukin, A. V. Shapranov, M. M. Demin i dr. // Kratkie soobshcheniya po fizike FIAN.

- 2016. - № 9. - S. 37-43. - Tekst : neposredstvennyj.

7. Markidonov, A. V. Raschet termodinamicheskih harakteristik sistemy Fe-P metodom molekulyarnoj dinamiki / A. V. Markidonov, D. A. Lubyanoj, V. V. Kovalenko i dr. - Tekst : neposredstvennyj // Izvestiya vysshih uchebnyh zavedenij. Chernaya metallurgiya. - 2019. - T. 62. - № 9. - S. 725-731.

8. Verlet, L. Computer Experiments on Classical Fluids. I. Thermodynamical Properties of Lennard-Jones Molecules / L. Verlet // Phys. Rev. - 1967. - V. 159. - Pp. 98-103.

9. Rakitin, M. S. Pervoprincipnoe i termodinamicheskoe modelirovanie vliyaniya uprugih napryazhenij na energiyu rastvoreniya vodoroda v al'fa-zheleze / M. S. Rakitin, A. A. Mirzoev, D. A. Mirzaev. - Tekst : neposredstvennyj // Izvestiya vysshih uchebnyh zavedenij. Fizika. - 2017. - T. 60. - № 12. - S. 76-82.

10. Borovik, A. M. Modelirovanie elektricheskih harakteristik grafenovogo polevogo tranzistora na osnove dannyh raschetov iz pervyh principov / A. M. Borovik, M. S. Baranova, D. Ch. Gvozdovskij.

- Tekst : neposredstvennyj // Infokommunikacionnye i radioelektronnye tekhnologii. - 2020. - T. 3. - № 1.

- S. 63-74.

11. Zalizniak, V. E. Mathematical model of conducting nanopore for molecular dynamics simulations / V. E. Zalizniak, O. A. Zolotov, O. P. Zolotova // Siberian Journal of Science and Technology. - 2018.

- V. 19. - Is. 4. - Pp. 677-682.

12. Sanchez-Badillo, J. Potential of mean force calculations for an sN2 fluorination reaction in five different imidazolium ionic liquid solvents using quantum mechanics/molecular mechanics molecular dynamics simulations / J. Sanchez-Badillo, M. Gallo, R. A. Guirado-Lopez et al. // Journal of Physical Chemistry B. - 2020. - V. 124. - Is. 21. - Pp. 4338-4357.

13. Jones, J. E. On the Determination of Molecular Fields / J. E. Jones // Proceedings of the Royal Society of London. - 1924. - V. 106. - Is. 738. - Pp. 463-477.

14. Sunagatova, I. R. Properties of one-dimensional nonlinear vibrational modes in triangular lattice with Lennard-Jones interactions / I. R. Sungatova, A. M. Subkhangulova, M. N. Semenova et al. // IOP Conference Series: Materials Science and Engineering. - 2020. - V. 1008. - No. 012073.

15. Krylova K.A. Spherically localized discrete breathers in bcc metals V and Nb / K. A. Krylova, I. P. Lobzenko, A. S. Semenov et al. // Computational Materials Science. - 2020. - V. 180. - No. 109695.

16. Babicheva R. I. Discrete breathers in a triangular P-Fermi-Pasta-Ulam-Tsingou lattice / R. I. Babicheva, A. S. Semenov, E. G. Soboleva et al. // Physical Review E. - 2021. - V. 103. - Is. 5. - No. 052202.

17. Cyan' Syue-Sen'. Fizicheskaya mekhanika / Cyan' Syue-Sen'. - Moskva : Mir, 1965. - 544 s.

- Tekst : neposredstvennyj.

18. Bebikhov, Y. V. The application of mathematical simulation for solution of linear algebraic and ordinary differential equations in electrical engineering / Y. V. Bebikhov, A. S. Semenov, I. A. Yakushev et al. // IOP Conference Series : Materials Science and Engineering. - 2019. - V. 643. - No. 012067.

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

19. Bachurina, O. V. Properties of Moving Discrete Breathers in Beryllium / O. V. Bachurina, R. T. Murzaev, A. S. Semenov et al. // Physics of the Solid State. - 2018. - V. 60. - Is. 5. - Pp. 989-994.

20. Tatarinov, P. S. Sbornik laboratornyh rabot po razdelu «Mekhanika» kursa «Fizika» / P. S. Tatarinov, V. D. Yakovlev, D. Ch. Kim i dr. - Moskva : Izdatel'stvo «Sputnik +», 2021. - 175 s. - Tekst : neposredstvennyj.

21. Vasil'eva, A. V. Praktikum po rabote s paketami prikladnyh programm na EVM / A. V. Vasil'eva, M. N. Semenova, I. A. Yakushev. - Moskva : Izdatel'stvo «Sputnik +», 2020. - 144 s. - Tekst : neposredstvennyj.

22. Semyonov, A. S. Programma matematicheskogo modelirovaniya fizicheskih processov v metallah i uporyadochennyh splavah / A. S. Semenov, M. N. Semenova, I. A. Yakusheva i dr. // Svidetel'stvo № 2021615775 ot 13.04.2021, zayavka № 2021614436 ot 02.04.2021. - Tekst : neposredstvennyj.

^iMir^ir

i Надоели баннеры? Вы всегда можете отключить рекламу.