CC BY
43
где угловые скобки означают усреднение по всем возможным перестановкам г и г и интегрирование проводится по 3(М — 1) координатам. В силу полноты одночастичных решений уравнения Шредингера, для связанных состояний многочастичное решение можно представить как сумму произведений различных ортогональных одночастичных решений с соответствующими статистике коэффициентами:
Ф = X] ^ П X (гз)
/ (Х (гз ^
2
фк (Г]дУз = 1.
Тогда для плотности зарядов электронного облака получается выражение вида
р(г) = Я Х а2ф2(т). (2)
г
Самый простой вид функции распределения зарядов — сферически симметричный. Как известно [14], решение уравнение Шредингера для водородоподобного атома в сферических координатах имеет вид
фп^т(г,в,р) = О^гЩ^О^У-™, (3)
где У1,т(6,ф) —сферические гармоники, (г) — полиномы, которые выражаются через обобщённые полиномы Лаггера, п,1,т — квантовые числа и постоянные тп находятся из условия существования нетривиальных решений при заданных граничных условиях. Подставляя волновые функции (3) в (2) и усредняя полученное выражение по углам, в силу ортонормированности сферических гармоник получим:
р(т) = <Э Х а*Р?(г)е-^.
г
Здесь учтено, что сумма квадратов полиномов Ип\г) с положительными коэффициентами будет неотрицательным полиномом, который приближаем в виде квадрата некоего полинома Рг(г).
При расчете энергии взаимодействия двух атомов физически интересны достаточно большие расстояния между ними, где основной вклад дают члены с малыми а^, то есть внешние оболочки атома. Поскольку в дальнейшем мы намерены находить коэффициенты из сравнения расчетных данных и результатов эксперимента, полином следует брать достаточно небольшой степени, чтобы ограничиться разумным количеством подгоночных параметров. Тем более, что степень полинома определяет число его нулей, а после сложения квадратов полиномов Ип\г) различных степеней результирующий полином не будет иметь действительных корней при г > 0. Тогда в интересующей нас весьма ограниченной области положительных действительных г, искомый полином может быть с достаточной точностью приближен полиномом много меньшей степени.
Таким образом, вышеизложенное дает основания полагать, что уже функция вида
р(г) « (1 + /Зг)2е-аг
может быть неплохой аппроксимацией плотности заряда электронного облака. С учетом нормировки, плотность зарядов электронного облака будет иметь вид
После интегрирования уравнения Пуассона получим создаваемый такой плотностью потенциал:
Ф(г)
Q (2 (а$2(г2а2 + 2га + 2) 2а@(а + 2р )(га + 1) 2
+ а + -
8^о Vг \ + 1202 + а2) (6а@ + 12$2 + а2)
где £0 — диэлектрическая проницаемость вакуума.
Рассмотрим электростатическую энергию взаимодействия двух атомов. Она состоит из энергии взаимодействия двух ядер Wnn, энергии взаимодействия их электронных оболочек и энергии взаимодействия ядер с электронными оболочками другого атома Шпз:
У (г) = №пп + Ш88 + Шп8, (4)
№пп (г)
ЯгЯ2
4п£оГ
, (г) = ^!Ф2 (г) + (г)
Wяя (г)
(г) = 2тг ^ ^ (г') р2 ('г2 + (г')2 - 2г (г') сов(0)^ ^ яп(г1)2 Аг'+
ОО / Г'Ж
Ф2 (г')рг (\/г2 + (г')2 - 2г (г') сов(0) ) ) вт(0)с$) (г1)2 ¿г
2к /
./о Ч./0
Здесь г — расстояние между центрами электронных оболочек атомов, где находятся ядра с зарядами Q1 и Q2. Отметим, что в общем случае разных атомов мы имеем две пары подгоночных параметров а1, /31 и а2, /32. С помощью интегрирования можно получить аналитическое выражение для но оно получается очень громоздким даже для одинаковых атомов и поэтому здесь не приводится.
Очевидно, что при перекрытии электронных облаков квантово-механические эффекты обменного взаимодействия играют существенную роль. Как было уже указано выше, точное аналитическое или полуаналитическое решение такой квантово-механической задачи для общего случая не представляется возможным. Поэтому будем учитывать эти эффекты с помощью подгоночного параметра 7 и тогда потенциал парного взаимодействия примет вид
р (г) = 7 (Шпп + Ш88 + Шп8).
Подставляя аналитическое выражение для в (4) и преобразуя получающееся выражение, потенциал парного взаимодействия для одинаковых атомов можно записать в виде
(р (г) = £ ■ ехр (—аг) ^^ ап (аг)г
п= — 1
где параметры ап зависят от а и ¡3 :
ап (а,Р) =
^2 е = "I->
4ъ£о
а-1 = 1, Рп+1 )
п
0,..., 6.
Функции рг,
Ро (а,р
Р1 (а,Р Р2 (а,Р Рз (а,Р Р4 (а,Р Р5 (а,Р
ро (а,Р)
возникающие после преобразования выражения (4), имеют следующий вид:
(а2 + 6ар + 12$2)2,
0.3125а4 + 5.125а3@ + 31.5а2@2 + 88.875а@3 + 98.1562@4, —0.1875а4 — 0.875а3 @ + 1.5а2/32 + 16.845ар3 + 72@4, —0.02а4 — 0.5а3 @ — 2.44а2@2 — 3.875ар3 + 2.125@4, —0.04а3 @ — 0.4375а2 @2 — 1.5а$3 — 1.87504, —0.15а@3 — 0.029а2 @2 — 0.22504,
г
Рб (а,р) = —0.008«^3 — 0.01875,54, рт (а,р) = — 0.00089^4.
Из квантово-механических расчётов [7] следует, что функция Р(р), описывающая энергию вложенного атома, должна удовлетворять следующим свойствам:
(1) стремиться к нулю, когда электронная плотность стремится к нулю:
Р (0) = 0, (5)
(2) иметь отрицательное значение, так как она описывает притяжение атомов:
Р (р) < 0, (6)
(3) иметь отрицательный наклон:
дР
Ж(р)<0. (7)
(4) иметь положительную кривизну:
д2 Р
( ) > 0, (8)
др2
для характерных электронных плотностей р, наблюдаемых в металлах. Функция Р(р) задавалась в следующем виде
4 / \ п
^(рНЕМ^"1 ' (9)
п=0 '
где ре — равновесная электронная плотность, а способ определения коэффициентов сп будет представлен ниже.
3. Параметризация
Таким образом необходимо подобрать только два параметра а и ¡3 для того чтобы определить потенциал взаимодействия между атомами металла одного сорта. Рассмотрим идеальный, однородный кристалл при нулевой температуре, тогда все атомы будут эквивалентны. Поэтому мы можем связать начало координат с одним из атомов и переписать (1) в виде
1 2
Еп = Р(Ре) + -^2^(гт), (Ю)
т
Ре = (Гт) '
т
где суммирование производится по всем окружающим атомам, которые находятся на расстояниях гт от начала координат. Эти расстояния определяются типом и параметрами решётки для каждого конкретного металла. Основные свойства твёрдого тела могут быть вычислены из уравнения (10). Условие равновесия решётки имеет вид
\ X Тт% ^ + ^ ^ X = ^ ^
Энергия связи на атом определяется выражением
1 2
т
Энергия образования вакансии приближённо может быть вычислена следующим образом:
= 0 ^ (Гт)) ~ F (Ре)]
m m
+ ^ р - р2 •
т т т
Также можно вычислить упругие постоянные твёрдого тела. Приближённые значения трёх независимых упругих постоянных вычисляются следующим образом:
_ _ ___т
-11
М2?г£(<гД(Гт) г,„<fr<Гт))+ % w <--) - (->) + f ы (Ç ££ (г->) "}
д£ ^ у xjyj (_ J_dp \ + d^F fy x^dp \ /y yidp dp e rm V dr2 m rm dr m 1 dp2 e \ rm dr m I \ ^ rm dr m
m m m m
44 К rl V dr2 [Гт) rm dr [Гт)) +
m m
2
m
f ы £ ^ @ Ы - В Ц + 0 <*> (?:fï (-»)2} ■
о vre/ / , о l Ï о V m/ Ï v m/ I 1 о о \ie/ \ / , Ï
dp ' rm V «r2 rm dr ) dp2 y ' rm dr
где ^ — атомный объём. Тогда приближённое значение модуля всестороннего сжатия равно
= — ¡-Т(г2^(г )-г ^(г )1 + 9У„ т dr2 [т) т йг [т) 1
te (*£<'->-'-24
m
Е (гф Ы - rmfr (rm)) + Ц (Ре) fe Г-Тг {Гт}) f ■
m m
Процедура параметризации состоит из трёх частей:
(1) для различных пар значений (а, [3) вычисляются значения 7, ре, F(ре), F'(pe), F"(pe) при которых выполняется условие (11) и
(Ес - Е^)2 + (Evf - Е^У + (В - В(а))2 = 0,
(2) определяется множество пар значений (а,@) при которых выполняется условие ^ > 0 и условия (6)-(8) для р = pe,
(3) из полученного множества выбирается пара оптимальных значений (а, /3) при которых величина
{си — с^^ + ^С12 — + {с44 — éll^j
принимает минимальное значение.
Коэффициенты сп функции вложенной энергии (9) определяются из условия (5), вычисленных значений ре, F(ре), F' (ре), F' '(ре) и условия F"(2ре) = 6. Величина 6 подбирается так, чтобы зависимость энергии кристалла Е от параметра решётки а* наилучшим образом приближала уравнение состояния Rose et al [18] при 0.8а < а* < 1.6а.
При параметризации учитывалось влияние только ближайших атомов, при этом радиус усечения задавался равным 2а для ГЦК и ОЦК металлов, и 2с для ГПУ металлов. Далее приводятся результаты параметризации для трёх металлов с различными кристаллическими структурами. Значения физических величин, используемых для нахождения параметров а и /3 для этих металлов, приведены в Таблице 1.
Таблица 1
Величина А1(ГЦК) Fe(OU,K) Mg(rny)
а, А 4,05 [19] 2,87 [19] 3,176 [19]
с, А - - 5,145 [19]
Ес, эВ 3,34 [19] 4,28 [19] 1,519 [19]
Evf, эВ 0,64 [20] 1,6 [21] 0,58 [22]
В, эВ/А3 0,451 [23] 1,052 [23] 0,2306 [23]
Си, эВ/А3 0.6688 [23] 1.4766 [23] 0.3969 [23]
cía, эВ/А3 0.3813 [23] 0.8486 [23] 0.1644 [23]
с44, эВ/А3 0.1813 [23] 0.7545 [23] 0.1150 [23]
Q, е 13 26 12
Вычисленные параметры парного потенциала взаимодействия приведены в Таблице 2. Распределение электронной плотности для рассмотренных атомов показано на рис. 1. Зависимость энергии парного взаимодействия двух атомов от расстояния показана на рис. 2. Вычисленные коэффициенты функций Р(р) приведены в Таблице 3. Зависимость энергии Р от безразмерной плотности р* = р/ре показана на рис. 3.
ТАБЛИЦА 2. Параметры парного потенциала взаимодействия
Параметр А1(ГЦК) Fe(OU,K) Mg(rny)
а, 1/А 1.8206 2.1025 1.4798
/3, 1/А -2.452 -2.7959 -3.1749
£, ЭВ 4.59254 13.21523 4.77938
a_i 1 1 1
а0 0.725 0.72551 0.71066
ai 0.9248 0.93268 0.72737
а2 0.05880 0.05959 0.03888
а3 -0.01504 -0.01506 -0.01427
а4 -1.977305Е-03 -1.982755Е-03 -1.818449Е-03
а5 -1.954224Е-04 -1.963332Е-04 -1.692550Е-04
а6 -1.357733Е-05 -1.372200Е-05 -1.002825Е-05
Поверхностная энергия является одним из основных свойств поверхностей металлов. Она может быть вычислена следующим образом:
Е3 = ^(Е1(М)-ЕЬ(М)),
РИС. 2. Зависимость энергии парного взаимодействия от расстояния
Р ; эВ
РИС. 3.
Таблица 3. Параметры функции Р (р)
Параметр А1(ГЦК) Ре(ОЦК) М§(ГПУ)
ре, е/А3 0.67022 1.81908 0.40755
с0, эВ -3.337719 -4.27535 -1.51658
си эВ -0.53775 -1.59524 -0.57753
с2, эВ 1.01860 1.075809Е-03 5.626435Е-04
с3, эВ -0.73678 -0.44412 1.115650Е-01
С-4, ЭВ 1.04460 2.23491 1.05005
где Еъ^) —энергия N атомов в неограниченном теле, Е1(М) —энергия N атомов в бесконечном слое и Б - площадь поверхности слоя содержащего N атомов. Толщина слоя задавалась равной 7а для А1 и Бе, и 7с для М^.
В Таблице 4 приведено сравнение экспериментальных и вычисленных значений основных параметров металлов. Первая строка соответствует экспериментальным значениям, вторая — вычисленным значениям. Число в круглых скобках означает индекс Миллера поверхности, а символ «*» означает поликристаллическую поверхность.
4. Заключение
В работе предложен потенциал межатомного взаимодействия для чистых металлов. Для определения параметров взаимодействия в конкретном металле требуется подобрать всего два параметра распределения электронной плотности. Это значительно облегчает обобщение предложенного потенциала на случай взаимодействия между атомами различных металлов. В качестве примера приводятся параметры потенциала взаимодействия для
Зависимость энергии Б от безразмерной электронной плотности
ТАБЛИЦА 4. Вычисленные значения физических параметров
Величина А1(ГЦК) Бе(ОЦК) Mg(rny)
Ec, эВ 3,34 4,28 1,519
3,34 4,28 1,519
El, эВ 0,64 1,6 0,58
0,64 1,6 0,58
В, эВ/А3 0,451 1,052 0,2306
0,451 1,052 0,2306
Си, эВ/А3 0.6688 1.4766 0.3969
0.6688 1.5542 0.4000
с12, эВ/А3 0.3421 0.8486 0.1644
0.3421 0.8009 0.1624
с44 эВ/А3 0.1813 0.7545 0.1150
0.1816 0.6036 0.0940
Es, Дж/м2 1.28(100)[24], 1.16(*)[20] 2.48(*) [20] 0.785(0001) [25]
1.24(100) 3.05(100) 0.929(0001)
Al, Fe и Mg, которые представляют три основные кристаллические структуры, наблюдаемые в металлах. Наблюдается хорошее соответствие между вычисленными и экспериментальными значениями различных физических параметров металлов. Из процедуры параметризации следует, что вычисленные значения энергии связи, энергии образования вакансии и модуля объёмного сжатия равны экспериментальным значениям этих параметров. Вычисленные значения упругих постоянных с11 и с12 незначительно отличаются от их экспериментальных значений. Менее точно вычисляются значения упругой постоянной с44 для Fe и Mg. Для выбранных металлов предложенный потенциал взаимодействия позволяет довольно точно оценить поверхностную энергию.
Точность описания свойств металлов с помощью предложенного потенциала ограничивается погрешностями самого метода вложенного атома и также зависит от того насколько точно предложенное распределение электронной плотности описывает распределение электронной плотности реального атома.
Приведённая процедура параметризации позволяет определить параметры потенциала, которые дают достаточно точные результаты при моделировании явлений, когда средняя электронная плотность системы атомов близка к равновесной. Применение потенциала с этими параметрами для моделирования явлений, когда р/ре мало (например, малые кластеры, процессы напыления) даёт менее точные результаты. В этом случае для определения параметров потенциала следует применять процедуру параметризации, основанную на других условиях (например, имеются обширные данные по длине и энергии связи димеров металлов).
Литература
[1] Daw M.S. and Baskes M.I. Semiempirical, quantum mechanical calculation of hydrogen embrittlement in metals // Phys. Rev. Letters. — 1983. — 50(17). — 1285.
[2] Daw M.S. and Baskes M.I. Embedded-atom method: Derivation and application to impurities, surfaces, and other defects in metals // Phys. Rev. B. — 1983. — 29. — 6443.
[3] Finnis M. W. and Sinclair J. E. A simple empirical N-body potential for transition metals // Phil. Mag. A. — 1984.— 50(1).— 45.
[4] Sutton A. P. and Chen J. Long-range Finnis-Sinclair potentials // Phil. Mag. Letters. — 1990. — 61(3). — 139.
[5] Ercolessi F., Tosatti E. and Parrinelo M. Au (100) surface reconstruction // Phys. Rev. Letters. — 1986. — 57(6).— 719.
[6] Ercolessi F., Parrinelo M. and Tosatti E. Simulation of gold in the glue model // Phil. Mag. A. — 1988. — 58(1). — 213.
[7] Foiles S.M., Baskes M.I. and Daw M.S. Embedded-atom-method functions for the fcc metals Cu, Ag, Au, Ni, Pd, Pt, and their alloys // Phys Rev B. — 1986. — 33(12). — 7983.
[8] Johnson R.A. Analytic nearest-neighbor model for fcc metals // Phys Rev B. — 1988. — 37(8). — 3924.
[9] Doyama M., Kogure Y., Embedded atom potentials in fcc and bcc metals // Comput. Mater. Science. — 1999. — 14(1-4), 80.
[10] Dai X.D., Kong Y., Li L.H. and Lin B.X. Extended Finnis-Sinclair potential for bcc and fee metals and alloys // J. Phys.: Condensed Metter. — 2006. — 18. — 4527.
[11] Karolewski M.A. Tight-binding potentials for sputtering simulations with fcc and bcc metals // Radiation Effects and Defects in Solids. — 2001. — 153(3). — 239.
[12] Cleri F. and Rosato V. Tight-binding potentials for transition metals and alloys // Phys. Rev. B. — 1993. — 48(1). -22.
[13] Shunling Chen, Jincheng Xu, Hongsheng Zhang, A new scheme of many-body potentials for hcp metals // Comput. Mater. Science. — 2004. — 29. — 428.
[14] Ландау Л.Д., Лифшиц Е.М., Теоретическая физика. Т. 3. Квантовая механика. — М.: Наука, 1989.
[15] Фок В.А., Начала квантовой механики. — М.: Наука, 1976.
[16] Jones R.O. and Gunnarsson O. The density functional formalism, its application and prospects // Rev. Mod. Phys. — 1989. — 61(3). — 689-746.
[17] Баранов М.А., Дятлова И.В., Бумажникова K.H., Бразовский В.Е., Электростатическое приближение и его применение к описанию устойчивости многокомпонентных металлических соединений // Горизонты Образования. Научно-Образовательный Журнал АлтГТУ. — 2003. — 5. — 26-44.
[18] Rose J.H., Smith J.R., Guinea F., Ferrante J. Universal features of the equation of state of metals // Phys. Rev. B. — 1984. — 29. — 2963-2969.
[19] Kittel C. Introduction to solid state physics. — New York : Wiley, 1996.
[20] de Boer F.R., Boom R., Mattens W.C.M., Miedema A.R., Niessen A.K. Cohesion in metals. — Amsterdam: North Holland, 1988.
[21] Puska M. and Nieminen R.M. In: Density functional methods in chemistry and materials science. — New York: Wiley, 1997.
[22] Doyama M. and Koehler J.S. The relation between the formation energy of a vacancy and the nearest neighbor interactions in pure and liquid metals // Acta Metall. — 1976. — 24(9). — 871-879.
[23] Simmons R.O. and Wang H. Single Crystal Elastic Constants and Calculated Aggregate Properties: A Handbook. — Cambridge MA: MIT Press, 1977.
[24] Jacobs P.W.M., Zhukovskii Yu.F., Mastrikov Yu., Shunin Yu.N. Bulk and surface properties of metallic aluminium: DFT simulations // Computer Modelling & New Technologies. — 2002. — 6(1). — 7-28.
[25] Tyson W.R. and Miller W.A. Surface free energies of solid metals: Estimation from liquid surface tension measurements // Surface Science. — 1977. — 62(1). — 267.
