Использование приближения линейной связи при расчете энергий образования точечных дефектов ионных кристаллов Текст научной статьи по специальности «Физика»

УДК 548.4:548.31

С. И. Гращенков, В. Л. Вейсман ИСПОЛЬЗОВАНИЕ ПРИБЛИЖЕНИЯ ЛИНЕЙНОЙ СВЯЗИ

при расчете энергий образования точечных дефектов ионных кристаллов

В работе предлагается новый метод вычисления энергий образования точечных дефектов в ионных кристаллах. Этот метод использует процедуру Мотта-Литтлтона и приближение линейной связи для энергии взаимодействия ионов области II с ионами области I и дефектом. Предложенный метод достаточно прост и эффективен для быстрого анализа различных комплексов точечных дефектов, т. к. в полученном базовом уравнении для энергии образования дефекта в явном виде отсутствует энергия взаимодействия ионов области II с ионами области I и энергия взаимодействия ионов области II друг с другом. Разработанный метод не содержит (в отличие от общеизвестных) жёстких ограничений на конфигурацию области I, что может оказаться полезным для моделирования сильно протяжённых дефектов. При этом точность вычисления энергетических параметров дефектов можно повысить простым увеличением области I, что делает данный метод особенно актуальным в связи с постоянным ростом мощности компьютеров.

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

ВВЕДЕНИЕ

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

В настоящее время наиболее интересные данные в теоретических исследованиях точечных дефектов в ионных кристаллах получены методами атомного моделирования с помощью ЭВМ. Одним из основных методов моделирования характеристик дефектов в кристаллах является вариационный метод (метод молекулярной статики) [1-5]. Этот метод позволяет провести расчёт энергий образования и миграции структурных дефектов (т. е. вакансий, ионов внедрения, замещения и их комплексов) в ионных кристаллах при Т = 0 К. В случае ионных кристаллов вариационный метод, как правило, используется совместно с хорошо известной процедурой Мотта-Литтлтона [6]. В этих моделях [7-11] кристалл разбивается на две области. В первой (I), ближней к дефекту области локальные поля определяются микроскопически в рамках одного из известных методов расчёта: оболочечной модели иона, модели точеч-но-поляризуемых ионов и др. Силы взаимодействия между ионами определяются потенциалами, выраженными в аналитической форме. Область II — континуум. В

этой области кристалла расчёт ведется приближенными методами. Основную сложность в численных расчётах представляет учёт энергии взаимодействия ионов области II с ионами области I. В настоящее время наиболее эффективным для учёта этой энергии является подход, основанный на методе Эвальда [8-10]. Однако следует заметить, что использование метода Эвальда приводит к уравнениям, аналитический вид которых весьма сложен и требует много времени при создании соответствующих компьютерных программ. Это не позволяет быстро оценить пригодность метода молекулярной статики для исследования различных сложных систем дефект — кристалл.

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

теоретическая часть

Энергию образования дефекта можно представить в виде [7, 8]:

где Е0 — энергия идеального кристалла, Е' — энергия кристалла с дефектами, Е1 (х)

— изменение энергии взаимодействия ионов кристаллической решетки в области I,

Е2 (х,Е) — изменение энергии взаимодействия ионов области I с ионами области

II, Е3 (Е) — изменение энергии взаимодействия ионов в области II, W0 — энергия

дефекта в недеформированной кристаллической решетке, W1 (х , х') — изменение

энергии взаимодействия ионов области I с дефектом, W2 (х', Е) — изменение энергии взаимодействия дефекта с ионами области II.

Параметр х характеризует смещение ионов области I и их поляризацию, е характеризует смещение и поляризацию ионов области II (для оболочечной модели х

и Е — смещения ядер и оболочек ионов из узлов кристаллической решетки), х'

— координаты дефекта. Основную трудность при вычислении энергии образования

дефекта по (1) составляет определение слагаемых Е2 (х,Е) и Е3 (Е), т. к. область II рассматривается в континуальном приближении.

Для преодоления этих трудностей в ряде работ [2, 7-10] Е2 и Е3 переопределяют таким образом, что для переопределенной энергии Е3 становится верным гармоническое приближение по Е , при условии достаточной удаленности области II от дефекта. Это позволяет уравнение (1) представить в виде:

здесь е2 (х ,Е) — переопределение Е2 (х ,Е) Е — смещение, соответствующее равновесной конфигурации ионов области II. Кроме того, в работах [8-10] использу-

Е = Е' - Е0 = Е1 (х) + Е2 (х , Е) + Е3 (Е) + W0 + Щ (х , х') + W2 (х', Е),

(1)

~ Е + Щ0 + Щ1(х ,х) +

(2)

ется метод Эвальда, который позволяет при вычислениях Е2(х ,Е) и дЕ2 (х ,Е) учи-

тывать только ионы области НА, непосредственно прилегающие к области I.

В настоящей работе мы не используем указанное выше переопределение

Е2 (х ,Е) и Е3 (Е), а представляем Е3 (Е) в виде:

Е3 =

4=0-4 +1 4-4 -4

где 4 — матрица силовых постоянных [1].

Запишем условие равновесия для ионов области II:

дЕ2 (* ,4)

дЕз (4)

4=4

54

дW2 (* ' ,4)

4=4

д4

4=4

.= 0

(3)

(4)

С другой стороны, условие равновесия для бездефектного кристалла дает:

дЕ2 (* ,4)

54

4=0 *=0

дЕз (4)

д4

4=0

= 0

(5)

С учётом (3, 4, 5) уравнение (1) приводится к виду:

Е = Е1( *) + Е2 (* ,4) -

1 дЕ2 (* ,4)

2 д4

+ W0 + №1 (* ,*') + W2 (*' ,4) --

*=0 -4 4=0

1 дW2 (*' ,4)

1 дЕ2 (*,4)

2

д4

4=4

•4 +

(6)

2 д4

4=4

Отметим, что уравнение (4) не содержит внутренних противоречий, несмотря на явную зависимость от х только первого члена, т. к. мы, как и авторы работ [8-10], считаем, что ионы области II адиабатически следуют движению ионов области I, и поэтому Е являются функциями х [1].

Воспользуемся линейным приближением по параметрам, характеризующим

смещение и поляризацию ионов (приближение линейной связи) для Е2 (х,Е) и

Щ2 (х', Е). В этом случае

дх

А • * *=0

4=0

дЕ2 (* ,4)

д4

*=0-4

4=0

_дит2 (*' ,4)

д4

Кроме того, в приближении линейной связи можно считать, что

дЕ2 (* ,4) дЕ2 (* ,4)

д4

д4

*=0 4=0

и

дW2 (*' ,4)_ дW2 (*' ,4)

д4

д4

4=0

(7)

(8)

(9)

С другой стороны, для идеальной решетки

5Е2 (* ,4)

е2 =

дЕ1( х)| + дЕ2 (х ,Е)

дх х—0 дх

Е=0

С учетом уравнений (7-11) уравнение (6) приводится к виду

0 = 0 (11) х—0

Е — Е1( х) -^ дх

х + Щ0 + Щ1( х, х' ) +1Щ2 (х' ,Е) (12)

х—0

Для нахождения равновесной конфигурации системы дефект-кристалл в области I используется условие минимума энергии образования дефекта, получаемое из (12):

дЕ1 (х)

дх

дЕ1 (х)

х—X дх

х—0 дх

— 0

х—X

(13)

где ~ — параметры, характеризующие смещение и поляризацию ионов области I, соответствующие равновесной конфигурации кристалла с дефектом.

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

при вычислении энергии Щ2 взаимодействия дефекта с кристаллом. Понятно, что

данное замечание не относится к расчёту Щ0. Если же рассматривается нейтральный дефект (например, бивакансия), то весь расчёт по нашему методу проводится, естественно, без учёта координат ионов II области. При расчёте же по формуле (2) во всех случаях необходимо учитывать положение ионов II области

из-за присутствия в уравнении Е'2 (х ,Е).

РЕЗУЛЬТАТЫ РАСчЕТОВ И Их ОБСУЖДЕНИЕ

Используя полученные выражения (12, 13), мы провели расчёт энергий образования катионной и анионной вакансий для щелочно-галоидных кристаллов структуры №С1 (Таблица 1). Расчёты выполнены с использованием оболочечной модели иона Дика-Оверхаузера [12] и потенциалов межионного взаимодействия Сангстера-Атву-да [13]. Область I выбиралась в виде куба и содержала 729 ионов. Согласно оболо-чечной модели каждый ион моделируется двумя материальными точками, представляющими собой ядро и оболочку. При взаимодействии ядер между собой и ядер с оболочками эти потенциалы учитывают кулоновское взаимодействие, а для оболочек потенциал парного взаимодействия разделяется на кулоновскую и короткодействующую части. Ядро и оболочка каждого иона связаны между собой упругой силой.

Вакансия моделировалась наложением на ион в центре куба «антииона», т. е. иона с такими же по модулю, но противоположными по знаку потенциалами взаимодействия с другими ионами. При расчёте Щ0 энергия электростатического взаимодействия вакансии с бесконечным бездефектным кристаллом рассчитывалась через постоянную Маделунга, а энергия короткодействующего взаимодействия находилась непосредственным суммированием. Для суммирования выбиралась область таких размеров, чтобы энергией короткодействующего взаимодействия между вакансией и ионами, не входящими в область, можно было пренебречь. В нашем расчёте мы пре-

небрегали короткодействующим взаимодействием между ионами, расстояние между которыми было больше 7а, где а — межионное расстояние.

Равновесная конфигурация системы в области I находится итерационным способом. Исходя из условий минимума энергии (13), величина —

VE (y)-VEx (у)

У=У0 +VWl (y ,х ')

рассматривается как сила, действующая

на ядро или оболочку. Здесь у — радиус-вектор соответственно данного ядра или

оболочки, у0 — радиус-вектор, соответствующий положению иона или оболочки в недеформируемом кристалле. Итерации прекращаются, когда силы становятся достаточно малыми или же когда энергия рассчитана с нужной точностью и дальнейшие итерации не приводят к существенному изменению результата. Энергии Е1 и W1 находятся непосредственным суммированием по области I. Для расчёта энергии W2 область II разбивается на две подобласти. Первая прилегает непосредственно к

области I и для нее соответствующая часть энергии W2 находится непосредственным суммированием. Оставшаяся подобласть представляет собой континуум и оставшаяся часть W2 находится интегрированием.

При расчете W2 мы использовали смещения в форме Мотта-Литтлтона [6], при этом параметры Мотта-Литтлтона вычислялись по формулам, приведённым в работе [14]. В Таблице 1 также представлены результаты расчётов энергий образования вакансий, выполненных по формуле (2) с использованием метода Эвальда и тех же потенциалов межионного взаимодействия (программа HADES). Как видно из таблицы, результаты расчетов, выполненные обоими методами, хорошо согласуются друг с другом. Хорошее согласие наших расчётов с результатами вычислений по формуле (2) подтверждает справедливость возможности использования приближения линейной

связи при расчете энергии E2 (х, £). Сама же возможность использования линейной

связи для Е2 обусловлена тем, что при достаточно больших размерах области I смещения ионов малы по сравнению с расстояниями между ионами I и II области. При этом точность вычислений повышается простым увеличением размеров области I.

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

решетки (I и II области), т. е. для расчёта энергетических слагаемых W1 (х ,х') и

W2 (х ' ,<Е). Расчёты, выполненные авторами [16] для ряда щелочно-галоидных кристаллов, показали, что такое использование приближения линейной связи не эффективно и может приводить к большим погрешностям. Это связано с тем, что при расчете в приближении линейной связи энергии взаимодействия дефекта с ионами решетки смещения ионов должны быть малы по сравнению с их расстояниями до дефекта, что в общем случае не выполняется. Как отмечают и сами авторы работы, удовлетворительные результаты получаются только для некоторых систем (кристалл-дефект), для которых смещения ионов меньше 10 % от постоянной решетки.

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

Естественно, полученные значения энергий образования катионной (Ек) и анионной (Еа) вакансий позволяют рассчитать энергии дефектов по Шоттки (Е) для щелочно-галоидных кристаллов.

Es = Ек + Еа + Е0 (14)

где Е0 — энергия связи, приходящаяся на пару ионов. С другой стороны,

Е =- 2(^к +W0a ) (15)

где W0k0a — энергии, соответственно, катионной и анионной вакансий в не-деформируемой кристаллической решетке.

В Таблице 2 представлены результаты соответствующих теоретических расчетов и для сравнения приведены значения Едля ряда ЩГК, полученные экспериментально в нашей лаборатории двумя, принципиально различными методами исследования термодинамических параметров точечных дефектов.

Таблица 1

Энергия образования одиночных вакансий для ионных кристаллов

структуры NaC1 ( 5 = EHADES Е 100 % )

Е HADES

Кристалл: дефект ЕэВ наш расчёт Е эВ HADES, данные по [15] 5 Е эВ данные по [16] 5

иР : Ук 7,43 7,46 0,4 — —

ир : г; 6,27 6,29 0,3 — —

: У- 5,53 5,57 0,7 — —

ЫС1 : У: 5,00 5,04 0,8 — —

и1Бг : У- 5,06 5,11 1,0 — —

ихВг : У+ 4,93 4,95 0,4 — —

-У Н-5 3 4,47 4,48 0,2 4,63 3,4

н-5 3 4,59 4,62 0,6 5,24 13,4

: Ук 6,32 6,37 0,8 — —

NaF : У а+ 6,28 6,28 0 — —

Продолжение Таблицы 1

№С1 : У- 5,15 5,19 0,8 — —

№С1 : ¥а+ 5,26 5.31 0,9 — —

№Вг : Ук 4,81 4,85 0,8 4,88 0,6

NaBr: У+ 5,19 5,22 0,6 5,47 4,8

NaJ : Ук 4,34 4,38 0,9 — —

NaJ : У; 5,00 5,00 0 — —

^ : Ук 5,54 5,60 1,1 — —

^ : у: 5,28 5,32 0,8 — —

КС1 : Ук 4,81 4,84 0,6 — —

КС1 : У а+ 4,93 4,97 0,8 — —

КВг : Ук 4,58 4,61 0,7 4,66 1,1

КВг : У+ 4,94 4,97 0,6 5,06 1,8

К е-Н У 4,23 4,26 0,7 — —

КК : У+ 4,83 4,87 0,8 — —

RbF : Ук 5,34 5,41 1,3 — —

RbF : У+ 4,80 4,84 0,8 — —

RbCl : Ук 4,69 4,72 0,6 — —

RbCl : Уа+ 4,68 4,72 0,8 — —

RbBr: Ук 4,50 4,53 0,7 4,59 1,3

RbBr: У+ 4,73 4,76 0,6 4,81 1,1

RbJ : Ук 4,19 4,22 0,7 — —

ИВК : У+ 4,68 4,71 0,6 — —

CsF : Ук 5,14 5,21 1,3 — —

О* : у+ 4,11 4,19 1,9 — —

Таблица 2

Энергии связи и образования пар Шоттки для ЩГК структуры ЫаС1

Крис- талл эВ 0а, наш расчет ^ эВ наш расчет Е№ эВ расчет по [16] ^ эВ данные по [13] Е, эВ s, наш расчет Е, эВ s, экспер. по [17-21]

LiF 10.882 11.427 -11.155 -11.165 2.55 —

LiCl 8.676 8.945 -8.811 -8.822 1.72 —

LiBr 8.371 8.347 -8.359 -8.375 1.63 —

LiI 7.913 7.694 -7.804 -7.824 1.26 —

NaF 9.565 9.493 -9.529 -9.538 3.07 —

ШС1 8.087 7.927 -8.007 -8.015 2.40 2.20

ШВг 7.854 7.480 -7.667 -7.676 2.33 1.78

7.483 6.937 -7.210 -7.221 2.13 —

О 8.205 8.294 -8.250 -8.258 2.57 —

КС1 7.297 7.137 -7.217 -7.225 2.52 2.44

КЬг 7.155 6.808 -6.982 -6.991 2.54 2.39

Ю 6.896 6.370 -6.633 -6.639 2.43 —

RbF 7.703 7.981 -7.842 -7.850 2.30 —

RbC1 6.955 6.886 -6.921 -6.926 2.45 —

RbBr 6.853 6.599 -6.726 -6.735 2.49 —

RbI 6.635 6.195 -6.415 -6.422 2.46 —

CsF 7.016 7.598 -7.307 -7.320 1.94 —

ВЫВОДЫ

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

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

Литература

1. Стоунхэм А. М. Теория дефектов в твердых телах. Т. 1. М.: Мир, 1978. 569 с.

2. Lidiard A. B., Norgett M. J. In: Computational solid state Physics / Eds. Herman F., Dalton N. W. and Koehler T. R. New York: Plenum, 1972. P. 385-412.

3. Computer simulation of solids. V. 166 / Eds. C. R. A. Catlow, W. C. Mackrodt. Berlin: Spring-er-Verlag, 1982. 320 p.

4. Дефекты в кристаллах и их моделирование на ЭВМ. Л.: Наука, 1980. 210 с.

5. Gale J. D. GULP - a computer based program for the symmetry adapted simulation of solids // J. Chem. Soc. Faraday Trans. 1997. V 93. P. 629-637.

6. Mott N. F., Littleton M. J. Conduction in polar crystals. I. Electrolytic conduction in solid salts // Trans. Faraday Soc. 1938. V. 34. P. 485-489.

7. Boswarva I. M., Lidiard A. B. The energy of formation of Schottky defects in ionic crystals // Phil. Mag. 1967. V 16. P. 805-826.

8. Norgett M. J. // Report AERE - 7650, 1974; AERE Harwell Report. M2880, 1977.

9. Колмогоров Ю. Н., Вараксин А. Н. Программа MOLSTAT для расчета характеристик структурных дефектов в ионных кристаллах // Журн. структур. химии. 1991. Т. 32. С. 162-165.

10. Gavartin J. L., Catlow C. R. A., Shluger A. L., Varaksin A. N., Kolmogorov Yu. N. Modelling Simul. // Matter. Sci. Eng. 1992. V 1. P. 29-38.

11. Lidiard A. B. Atomistic calculations of defects in ionic solid - their development and their significance // J. Phys.: Condens. Matter. 1993. V. 5. P. 137-148.

12. Dick B. G., Overhauser A. W. Theory of the dielectric constants of alkali halide crystals // Phys. Rev. 1958. V. 112. P. 90-103.

13. Sangster M. J. L., Atwood R. M. Interionic potentials for alkali halides: II. Completely crystal independent specification of Born - Mayer potentials // J. Phys. C.: Solid State Phys. 1978. V 11. No. 8. P. 1541-1555.

14. Sangster M. J. L., Schroder U., Atwood R. M. Interionic potentials for alkali halides: I. Crystal independent shell parameters and fitted Born - Mayer potentials // J. Phys. C.: Solid State Phys. 1978. V. 11. No. 8. P. 1523-1540.

15. Rowell D. K., Sangster M. J. L. Calculations of intrinsic defect energies in the alkali halides // J. Phys. C: Solid State Phys. 1981. V 14, No. 21. P. 2909-2921.

16. Колмогоров Ю. Н., Вараксин А. Н. Определение энергий образования точечных дефектов в ионных кристаллах в приближении линейной связи // ФТТ. 1991. Т. 33. № 2. C. 597-598.

17. Першиц Я. Н., Павлов Е. В. Определение параметров проводимости в галогенидах калия // ФТТ. 1968. Т. 10. № 5. C. 1418-1424.

18. Першиц Я. Н., Вейсман В. Л. Определение параметров проводимости в галогенидах натрия // ФТТ. 1970. Т. 12. № 11. C. 3175- 3178.

19. Першиц Я. Н., Каленникова Т. А. Определение энергии образования дефектов по Шот-тки в кристаллах KCl и NaCl путем электолитического введения одновалентных ионов // ФТТ. 1981. Т. 23. № 9. С. 2553-2559.

20. Першиц Я. Н., Вейсман В. Л. Определение энергетических параметров точечных дефектов в кристалле NaBr путем введения одновалентных и двухвалентных примесей // ФТТ. 1983. Т. 25. № 5. С. 1379-1385.

21. Вейсман В. Л., Першиц Я. Н. Исследование дефектов в щелочногалоидных кристаллах на основе изучения процессов электропереноса // Тез. докладов Международной научнотехнической конференции «Диэлектрики - 97». Санкт-Петербург, 1997. Т. 1. C. 28 - 30.

22. Физика электролитов / Ред. Дж. Хладик. М.: Мир, 1978. 555 с.

S. Grashchenkov, V. Veisman

CALCULATION OF FORMATION ENERGIES OF pOINT DEFECTS IN IONIC CRYSTALS BASED ON THE LINEAR BOND approximation

A novel method of the calculation ofpoint defects formation energies in ionic crystals has been suggested. This method uses the Mott - Littleton procedure and the linear approximation for the interaction energy of ions of the region II with those located in the region I and with the point defect. The suggested method is efficient and simple enough for rapid analysis of different point defect complexes due to the absence of the interaction energy between ions located in the region II as well as between those located in the regions II and I, as shown in our basic equation for point defect formation energy. In contrast with the well-known methods, our approach does not strictly limits the configuration of the region I, which may be useful for long-range defects modeling. This method seems to be especially actual taking into consideration computer technique progress: in order to enhance the accuracy of point defect energy parameters calculations one can simply increase the region I.

Key words: alkali halide crystals, point defects, formation energy, migration energy, atomic modeling, variation method, molecular statics.