Для розрахунку енергетичної зонної структури кристалів останнім часом набув поширення метод апріорних атомних псевдопотенціалів (ПП). У загальних рисах цей підхід ґрунтується на самоузгодженому пошуку ПП у наближенні функціонала локальної спінової густини. Стартова точка цієї процедури базується на релятивіському рівнянні Дірака для хвильової функції Gl(r) і Fl(r):

dFl(r)/dr - (g/r).Fl(r) + a.[El - V(r)]Gl(r) = 0 (1a)

dGl(r)/dr + (g/r).Gl(r) - a[(2/a2) + El - V(r)].Fl(r) = 0, (1b)

де l = a-1=137.07 – обернене значення константи надтонкої структури; g - ненульове ціле число. Розв’язки рівняння (1) визначають густину заряду:

r(r) = å [|Gl(r)|2 + | Fl(r)|2]. (2)

El < EF

Використовуючи процедуру нелінійної інтерполяції, визначимо псевдопотенціал:

Vps(l)(r) = Vост(r) + Vioн(l)(r)+ Vсо(l)(r), (3)

де Vост(r) – остовний потенціал, V(l)ioн(r) визначає іонну корекцію псевдопотенціалу і V(l)со(r) – спін-орбітальна корекція. Кожний з вищенаведених доданків може бути виражений в аналітичній формі і тому відповідні матричні елементи можна точно обчислити (замість числового інтегрування):

Vîńň(r) = (-Zv.e2/r).Ci.erf[ai1/2.r]; (4)

Vłîí(r) = (Ai+r2.Ai+3).exp(-gi.r2). (5)

Інтерполяційні коефіціенти Ai, ai, CI, gi визначаються як розв`язки самоузгодженого рівняння Дірака – Хартрі – Фока – Слетера для конкретних атомів з відповідними орбітальними числами з подальшою нелінійною інтерполяцією (4, 5).

Повний псевдопотенціал є сумою нелокальних псевдопотенціалів, які перекриваються і розміщені в точках tp,q. Взаємодія електрона з остовами, які описуються періодичними іонними ПП, визначається оператором:

Vps(r, r') = å Vq,s(r - Rp - tq,s, t - Rp - tq,s), (6)

p,q,s

де Rp – вектор прямої ґратки, який визначає розташування елементарної комірки; tq,s – вектор, який визначає розташування s-го іона сорту q в елементарній комірці. При обчисленні матричних елементів секулярного рівняння в базисі плоских хвиль необхідно перейти від r-простору до оберненого G-простору за допомогою Фур’є перетворення:

<k+Gi|Vps(r, r’)|k+Gj>=N-1exp[-iq(Rp + tq,s)].V-1.

ò dr.dr’.exp[-i(k+Gj).r].Vq(r, r').exp[i(k+Gj).r'], (7)

де q = Gi - Gj і V – об’єм першої зони Бріллюена. Враховуючи, що qj.ti = 2.p.dij (ti, qj – основні вектори прямої й оберненої ґратки) знаходимо, що при довільних значеннях p exp(-iqRp) = 1, а тому сума по p просто дає множник N i остання рівність у локальному наближенні набуває вигляду:

<k+Gi|Vps(r, r’)|k+Gj>= å exp[-iq.tq,s)].V-1.ò d3r.exp[-i(q.r].Vq(r) (8)

q, s

Зауважимо, що тут інтегрування по r здійснюється в основній сфері кристала.

Розглянемо локальну частину псевдопотенціала. Форм-фактор потенціалу іона (4) дорівнює:

Vq = V-1. ò Vîńň(r).exp(-iqr).d3r.

Для обчислення цього виразу використаємо загальну процедуру. Запишемо розклад плоских хвиль за Релеєм:

exp(iqr) = (2l+1).il.jl(|q.r|.Pl(cosQq^r), (9)

де jl(x) – сферичні функції Бесселя, Pl(cosq^r) – поліноми Лежандра l-го порядку і q^r – кут між q і r. Через сферичну симетрію s-орбіталі основний вклад в інтеграл у рівнянні (8) надає лише член ряду з l = 0. Враховуючи, що P0(cosq^r) = 1, одержимо:

¥ 2

V(q) = 4p/V ò r2.{-Zve2/r) å Ci.erf[ai1/2.r]+(Ali+r2.Ali+3).exp(-gli.r2)j0(qr)}d3r. (10)

0 i=1

Форм-фактор нелокальної частини псевдопотенціалу (перші два доданки 4) може бути визначений із:

< k + Gi|Vqíë (r, r’)| k + Gj > = V-1 ò d3r.d3r’.exp[-i(k + Gi).r].{(Ail + r2.Ail+3).

exp(-gil.r2) - (A(0)l + r2.A(0)l+3).exp(-g(0)l r2).Pl}.exp[-i(k + Gj).r]. (11)

Проекційний оператор l-ої компоненти моменту має вигляд:

Pl = Ylm(Qi,j).Y*lm(Q‘i,j‘), (12)

де Ylm(Qi,j) – сферичні гармонічні. За визначенням, дія проекційного оператора на довільну функцію виражається так:

Pl exp(ikr) = Ylm(Qi,j).{ò dQ'.sinQ'.dj'.Y*lm(Q‘,j‘).exp(ikr)}, (13)

де Q, Q', j, j' – полярний і азимутальний кути векторів k i r. Використовуючи (13), рівняння (11) можна записати:

l ¥ 2p p

< k + Gi|Vqíë (r, r’)| k + Gj > = V-1 å ò r2dr ò dj ò sinQ dQ´

m=-l 0 0 0

3

exp[-i(k + Gi).r].Ylm(Q,j){å (Ail + r2.Ail+3).exp(-gil.r2) - (A(0)l + r2.A(0)l+3).

l=1

2p p

exp(-g0)l r2)}.ò dj' ò sinQ'dQ Y*lm(Q‘,j‘).exp[i(k + Gj).r] (11a)

0 0

Враховуючи умову ортонормованості сферичних функцій:

2p 2p

òdj ò sinQ.dQ. Ylm(Q, j).Y*l'm'(Q, j) = dll' .dmm',

0 0

а також теорему додавання сферичних гармонік:

l

å Ylm(Q, j).Y*lm(Q, j) = [-(2l+1)/4p].Pl(cosQGi^Gj),

m=-l

одержимо:

< k+Gi|Víë(r, r’)| k+Gj>= 4p/V.(2l+1).Pl(cosQGi^Gj).Tl, (14)

де

cosQGi^Gj = GiGj/|Gi |.|Gj|, P0 = 1, P1 = cosQGi^Gj, P2 = (3.cosQGi^Gj - 1)/2;

а оператор Tl подається таким виразом:

Tl = [All' M1l +All'+3.Mll' - A0l' M0l - A0l'+3.M0l']; (15)

M1l = d3r.r2..jl(Gl.r).exp(-gil r2). jl(Gj r); (16)

M2(l) = d3r.r4.jl(Gir) exp(-gil r2).jl(Gj r). (17)

При проведенні самоузгоджених розрахунків враховувався екрануючий потенціал Хартрі-Слетера. Цей потенціал можна подати як звичайний кулонівський потенціал у формі, що визначається рівнянням Пуасона:

2.VH(r) = -4.p.e2.r(r), (18)

де r(r) - електронна густина. Фур`є образ кулонівського потенціалу дорівнює :

VH(Gi - Gj) = [4.p.e2./|Gi - Gj|].r(|Gi - Gj|) (19)

де

¥

r(Gi - Gj) = V-1 ò dr.r(r).exp[-i(Gi - Gj).r] (20)

0

Крім заданого потенціалу іона Vlіон(r) та потенціалу Хартрі VH(r), необхідно також враховувати потенціал, пов`язаний з обміном та кореляцією електронів, які беруть участь в екрануванні. У літературі описані різні обмінно-кореляційні потенціали та наведено результати розрахунків, одержаних при їх застосуванні [1]. Найбільш поширеним є обмінно-кореляційний потенціал Слетера без урахування спінової поляризації електронних станів. Останні залежать від сорту атома і враховуються множником:

VOK(r) = -3e2.b.[(3/8.p).r(r)]1/3 (21)

Параметер b змінюється в межах 0.5 - 1. Збільшення параметра посилює обмінно-кореляційний потенціал біля атома у порівнянні з міжатомною сферою. З цього випливає, що вибір обмінного потенціалу лише визначає розташування енергетичного спектра на шкалі енергій. Після перетворення Фур’є одержимо:

N

VOK(Gl - Gj) = -3e2.b.(3/8.p)1/3.N-1 å r(rl)1/3.exp[-i(Gl - Gj).rl].

l = 1

Отже, матричні елементи ермітової матриці, яка визначає гамільтоніан, записуються наступним чином :

<k+Gi|H| k+Gj>=( k + Gi)2.dGi^Gj + [v(q) (|Gi - Gj|) +

+< k + Gi|V(l)łîí| k+Gj>.r(q).(Gi - Gj) + VH(Gi - Gj) + VOK(Gi - Gj) (22)

Структурний фактор іона типу q є:

r(q)(|Gi - Gj|) = 1/2 å exp[i.(|Gi - Gj|.tp,q)].

p

Знаючи вираз гамільтоніана для старту процедури обчислення, необхідно вибрати певну початкову густину валентних електронів r(r), яку подають як суперпозицію електронних густин нейтральних атомів:

r(r) = å Cq,n,l < n; q, l; n; q, l >,

n,q,l

де Cq,n,l – фактор заповнення орбіталі n, q, l атома n. Хвильова функція кожної орбіталі задається рівнянням:

Yq(r) = Rqnl(r).Ulm(Q, j)/r.

Аналітичні вирази та коефіцієнти для функції Rqnl(r) наведені в роботах [2, 3].

Обчислення Фур’є перетворень для екрануючого кулонівського потенціалу VH(r) та обмінно-кореляційного потенціалу VOK[r(r)], виходячи з функції Rqnl(r), дуже трудомісткі і, як показала практика, не завжди доцільні. Тому при проведенні обчислень замість суми локального потенціалу Vlq(q), кулонівського екрануючого потенціалу VH(q) та обмінно-кореляційного потенціалу VOK(q) використовувалася незалежна модель [4]:

Vq(q) = V(q)/e(q),

де діелектрична проникність e(q) визначається рівнянням:

e(q) = 1 - [8.p.e2/Vŕň.q2].(1 - f(q)].c(q).

Щоб отримати повну густину валентних електронів, треба додати r(r) за всіма валентними зонами і за всіма різними станами k у першій зоні Бріллюена. Для визначення точного значення електронної густини здійснювалося додавання за багатьма точками ЗБ. Одержане значення використовується для обчислення Фур’є-образів кулонівського і обмінно-кореляційного потенціалів (19 - 21). Ітераційний процес самоузгодження зводився до розв’язування рівняння (7) на кожному етапі ітерації. Це дозволяло визначити густину валентних електронів, яка використовувалася на наступному етапі для обчислення екрануючого потенціалу. Показником ступеня самоузгодження є максимальна різниця власних значень гамільтоніана En,k на двох послідовних етапах ітерації. Точність розрахунків енергії була менша, ніж 0.2 eV, залежно від часу комп’ютерного розрахунку. Орбіталі йоду 5s-, 5px-, 5py-, 5pz і 5s-Cd вважалися базовими для секулярних рівнянь. Орбіталі 5s-J вводилися в першому порядку теорії збурень. Ефекти екранування приймалися до уваги для досягнення хорошого збігу між обчисленою енергетичною щілиною та її експериментальним значенням. Додаткове прискорення процедури самоузгодження досягалося внаслідок проведення варіації за параметром b. Зміна параметра веде до зміни відповідного енергетичного терма і менше впливає на дисперсію.