Плоские волны и сетки
Метод плоских волн подходит для описания периодических кристаллов, где они обеспечивают интуитивное понимание, а также простые алгоритмы для практических расчетов. Методы, основанные на сетках в реальном пространстве, являются наиболее подходящими для конечных систем. Современные алгоритмы используют методы плоских волн и сеток с быстрым преобразованием Фурье [24].
1.7.1 Уравнение Шредингера для системы независимых частиц в базисе
плоских волн
Собственные состояния для любой независимой частицы в уравнениях Кона-Шэма, подобных уравнению Шредингера, в которых каждый электрон движется в эффективном потенциале Veff (г), отвечают требованиям уравнения на собственные значения:
Heff (Г)(г) = [~ — V2 + Vefr (r)](//i (г) = (г) • (1.74)
2гц,
В твердом теле (или любом другом состоянии конденсированной среды) удобно, чтобы состояния были нормализованы и подчинялись периодическим граничным условиям в большом объеме Q, что позволит перейти к бесконечности. Пользуясь тем, что любая периодическая функция может быть разложена в полный набор компонент Фурье, собственная функция может быть записана как:
^i(r) = 2У q X _exp(iq ? г) = qx I q >, (1.75)
q q
где C; n - коэффициенты разложения волновой функции по базису ортонормированных волн I q >:
< Ч'1 Ч >= И dr exp(-iq’-r)exp(-iq • г) = <>q q.
(1.76)
Подставляя (1.75) в (1.74), умножая слева на (1.77) Х< q’l Heff I q > ci>q = s, X< q'l q > ci q = . q q Матричный элемент оператора кинетической энергии (последнее выражение в единицах Хартри): Для кристалла потенциал Veff (г) периодический и может быть представлен в виде суммы компонент Фурье: Veff (г) = Xveff (Gm) exp(iGm • г), (1.79) m где Gm - вектор обратной решетки и: Veff (G) = J- fveff (r)exp(-iG • r)dr , (1.80) Qcell n „ “cell где Qcell - объем основной ячейки. Таким образом, матричные элементы потенциала: < Q I I Q >= (Gm )^q-q’,Gm ’ (1-81) m Отличаются от нуля тогда и только тогда, когда q и q' отличаются некоторым вектором обратной решетки Gm. Определим q = k+Gm и q'=k+Gm. (которые отличаются на вектор обратной решетки Gm„ = Gm -Gm>), тогда уравнение Шредингера для любого заданного к можно записать в виде матричного уравнения: XHm,m.(k)ci>m.(k) = ^(k^,m(k) , (1.82) m' где: Hm,m'(k)- й2 2Ше I к + Gm I2 <5m m. + Veff (Gm - Gm.) (1.83) Здесь обозначены собственные значения и собственные функции i = l,2,... для дискретного набора решений матричных уравнений для заданного к. Уравнения Теорема Блоха: Каждая собственная функция в уравнении Шредингера (1.82) для заданного к задается выражением (1.75) с суммой по q ограниченной q = k + Gm, которая может быть записана как: Пк<г) = Xci,m(k)х -7^exp(i(k + Gm) • г) m Т2 j , (1.84) = exp(ik-r)^==uik(r) л/Ncdi где Q=NccUQccU и: ui,k<r>= т= Xci,m(k)exP(iGm-r) . (1.85) / ^cell m которая имеет периодичность кристалла. Теорема Блоха устанавливает, что любой собственный вектор является произведением exp(ik-r) и периодической функции. Волновая функция k (г) должна быть ортонормирована по объему ?1, а функция ui k (г) ортонормирована в одной примитивной ячейке, то есть: 1 Jdru*i,k(r)ui.k(r) = X ^cell Ceii m Зоны собственных значений. Так как уравнение Шредингера (1.82) определяется для каждого к по отдельности, тогда каждое состояние может быть помечено волновым вектором к, а собственные значения и собственные векторы для каждого к независимы, если они не отличаются на вектор обратной решетки. В пределе большого объема ?1 точки к станут плотным континуумом, а собственные значения ^(к) становятся непрерывными зонами. Для каждого к существует дискретный набор собственных состояний обозначенных i = l,2,..., которые могут быть найдены из Гамильтониана (1.83) в базисе дискретных компонент Фурье k + Gm, m = l,2,... [24]. 1.7.3 Расчет электронной плотности: сетки Одной из наиболее важных операций является вычисление плотности электронов п. Общая форма для кристалла в теории независимых частиц может быть записана в виде: с1-87) Nk tr п,.к(г)=|^ж(г)|2> которая является средней по к точкам с i, которая обозначает зоны в каждой точке k, a f (?• к) - функция Ферми. В базисе плоских волн, используя (1.84): ni,k<r)=5 2Скп>(к)сьт'(к)ехР(«(Ст'-Ст)-Г)> (i-88) т,т' n,,k(G)= ’ С1-89) т где т" определяет вектор G, для которого Gm- = Gm +G. Несмотря на простоту выражения (1.89), это не самый эффективный способ для расчета плотности п(г) или n(G). Проблема в том, что нахождение всех компонент Фурье, используя (1.89), включает в себя двойную сумму, которая требует NG операций, где NG - это количество векторов G, необходимых для описания плотности. Для больших систем расчёты становятся трудоемкими. С другой стороны, если состояния Блоха известны на сетке NR точек в реальном пространстве, то плотность можно найти просто как квадрат NR операций. Можно использовать быстрое преобразование Фурье (FFT), которое позволяет преобразовать из одного пространства в другое N в N операций, где N = NG = NR. Блок-схема на рисунке 1.7 иллюстрирует алгоритм и общие черты для всех таких операций. На рисунке 1.7 {G} и {R} обозначает наборы N векторов G и N точек сетки R. Результат дается как в реальном, так и в обратном пространстве, что необходимо для расчета обменно-корреляционного и Хартри термов. Плотность п(г) необходима для нахождения обменно-корреляционных плотности ?хс(г) и потенциала Vxc(r). Обратное преобразование может быть использовано для нахождения плотности n(G). Рисунок 1.7 - Вычисление плотности с использованием преобразования Фурье и сеток Идеальные примером применения метода плоских волн являются кристаллы с небольшими примитивными ячейками и атомами, с точно представленными псевдопотенциалами. Если ячейка в реальном пространстве мала, то относительно 42 небольшой набор плоских волн k + Gm является эффективной основой. Операции просты, легко вычислить полную энергию, плотность электронов, силы и т.д. [24]. 1.7.4 Методы в реальном пространстве. Метод конечных разностей. При расчете плотности многие операции легче выполнять в реальном пространстве. Например, если ^/(г) явно представлена на сетке, тогда можно найти n(r) = ^1 ^/(r) I2 без необходимости быстрого преобразования Фурье (FFT) в i соответствии с требованиями метода плоских волн. Потенциал Хартри можно найти с помощью FFT или алгоритмов множественных сеток в реальном пространстве. При решения уравнений Пуассона и Шредингера методы в реальном пространстве являются предпочтительными для конечных систем, где волновые функции равны нулю вне некоторой границы и кулоновский потенциал не подчиняется периодическим граничным условиям. В методе конечных разностей (МКР) оператор Лапласа кинетической энергии оценивается из значений функции на множестве точек сетки. Для равномерной трехмерной сетки сточками (xt,yj,zk) в приближении m порядка: д2ц/ Эх2 m = XCm^(xi+mh,yj,zk) + O(h2m+2), (1.90) —m где h - шаг сетки, am- положительное целое число. На рисунке 1.8А показано, что Лапласиан в центральной точке вычисляется через значения функции в точках на пересекающихся осях. Размер 25 точек уменьшается согласно величине Ст. Расчеты с использованием алгоритмов конечных разницей были применены ко многим задачам, в том числе для кластеров и других конечных систем. Многосеточный метод, основанный на форме, изображенной на рисунке 1.8Б, для Лапласиана был применен ко многим периодическим и непериодическим системам. Вариант на рисунке 1.8Б является предпочтительнее варианта на рисунке 1.8А, особенно для конечных систем, где более расширенный «крест» приводит к более сложным краевым эффектам. Размеры точки обозначают их вес схематично: А) ортогональных «крест» с 25 очками; Б) более компактный куб из 27 точек. Конечные системы, такие как кластеры и молекулы могут быть изучены, используя плоские волны или сетки. Сетки являются очевидным выбором для локализованных функций, так как нужно учитывать только ту часть сетки, где функция отлична от нуля. В случае плоских волн важно построить суперячейку, в которой локализована система [1, 3, 4, 6, 24]. А Б Рисунок 1.8- Два примера трафаретов для расчета Лапласиана методом конечных разностей Как отмечалось в пункте 1.33, собственные значения Кона-Шэма не имеют физического смысла, разность между собственными значениями не соответствует действительным энергиям возбуждения. Тем не менее, собственные значения Кона-Шэма были использованы для обсуждения спектров твердых тел, атомов и молекул. Это приводит к проблеме запрещенной зоны: в приближении LDA предсказания запрещенной зоны, как правило, на 30 - 50 процентов (и более) ниже ширины запрещенной зоны, полученной из эксперимента. Например, для кремния ширина запрещенной зоны, полученная методом функционала плотности, составляет 2.6 эВ, тогда как экспериментальное значение 3.4 эВ. Несмотря на эти трудности, в первом приближении исследование электронных спектров в пределах статической теории функционала плотности полезно. Можно записать уравнение для поляризуемости системы независимых частиц /{) (приближение случайных фаз (RPA)), как сумму по независимым переходам: |^clexp(iq-r)ly/v}| / = 2Е v,c (1-91) w-(fc -sv) + irj где i// - собственные одноэлектронные волновые функции, s - собственные значения, v и с обозначает валентные зоны и зоны проводимости. Диэлектрическая проницаемость (ь> - Кулоновский потенциал): s = l-u/ (1.92) мнимая часть которой описывает спектр поглощения: Abs = Im{^} = Im{l - }, (1.93) а мнимая часть обратной 8 1 описывает спектр характеристических потерь энергии электронами (EELS): (1.94) Спектр поглощения более явно _ 2 (1.95) Abs = Im{f}cc2jH lexp(iq-r)l^v)| J(w-(ffc-?v)). v,c Несмотря на простоту, это приближение может быть использовано для качественного, а иногда и количественного описания электронных спектров [32].
Электронные спектры в теории функционала плотности
