Конструктивная математика

Акимов О.Е.

29. Закон дисперсии в рамках kp-модели

Закон дисперсии E(k), рассчитанный в приближении слабой связи, ведется через матричные элементы псевдопотенциала

< k + K' | V | k + K >,

входящие в матрицу гамильтониана, обычная размерность которой равна 100 × 100 плоским волнам. ЛКАО гамильтониан, составленный по модели сильной связи, имел у нас матрицу размерности 8 × 8. Сейчас мы познакомимся с моделью, которая дает матрицы гамильтонианов средней размерности 15 × 15 или 27 × 27, хотя потенциально они могут составляться очень небольшими, порядка 3 × 3, или, наоборот, большими, порядка 50 × 50; ограничений здесь никаких нет. Этот так называемый kp-метод лежит посередине между методами сильной и слабой связи. Его вычислительный потенциал огромен, но с точки зрения физики, он представляется достаточно искусственной конструкцией. По причине его простоты и эффективности им удобно пользоваться при анализе закона дисперсии. Впрочем, существует еще около десятка методов расчета с большим и средним числом матричных элементов, о которых мы умолчали.

Надо сказать, что положительной стороной теории конденсированного вещества является то, что она терпимо относится к самым различным конструктивным подходам. В ее арсенале мирно сосуществует множество модельных конструкций, которые не претендуют на истину в последней инстанции. Однако в теоретической физике есть области, где господствует консервативный формализм, когда ни одна буква в формулах не может быть изменена или иначе проинтерпретирована. Конструктивно настроенным исследователям там делать нечего, они должны всячески избегать подобных «заболоченных» мест науки. Но обратимся к нашему модельному kp-гамильтониану.

Свое название kp-метод получил от матричных элементов, которые представляют собой произведение волнового вектора k на величину импульса p. В этом методе, как и в предыдущих двух, действует так называемое адиабатическое приближение, которое позволяет пренебречь влиянием ядер на движение рассматриваемого электрона и, таким образом, позволяет отделить электронную систему кристалла от фононной. Здесь, как и в двух предыдущих случаях, так называемое зонное приближение предполагает, что кристаллический потенциал создан всей совокупностью атомных потенциалов, когда небольшие температурные отклонения отдельного атома или его дефектное расположение не может повлиять на строгую периодичность волновой функции. За счет взаимной компенсации случайных явлений, происходящих в кристаллической решетке, решение уравнения Шредингера ищется как для идеального монокристалла.

Если в исходное уравнение Шредингера подставить волновую функцию ψ(r + R), то для собственных значений энергии E_n(k) получим уравнение:

H_ku_nk = [p²/2m + ћkp/m + V(r)]u_nk = E_n(k)u_nk.

Слагаемое ћkp/m в гамильтониане H_k рассматривается как потенциал возмущения к основному потенциалу V(r). Решение последнего уравнения ищется в виде разложения u_nk по некоторой полной системе функций:

u_nk = Σ_m a_nm(k – k₀)u_mk0,

в которой точка k достаточно близко отстоит от точки k₀, так что разложение в ряд возможно.

Выделим в гамильтониане H_k невозмущенную часть H_k0, получим

[H_k0 + ћ(k – k₀)p/m + ћ²(k² – k₀² )p/2m]u_nk = E_n(k)u_nk.

где H_k0 = p²/2m + ћk₀p/m + ћ²k₀² /2m + V(r). Подстановка разложенной функции u_nk в последнее выражение, умножение обеих частей уравнения на сопряженную функцию u_nk и интегрирование по элементарной кристаллической ячейке дает систему уравнений в виде следующей суммы

Σ_m {[E_n(k₀) + ћ²(k² – k₀² )p/2m]δ_nm + ћ(k – k₀)P_nm(k₀)/m}a_nm = E_n(k)a_nm,

здесь P_nm(k₀) = < u_nk0 | p | u_mk0 > — матричный элемент оператора импульса, δ_nm — символ Кронекера. Данная система уравнений разрешима, если ее определитель равен нулю, т.е.

det{[E_n(k₀) – E_n(k) + ћ² (k² – k₀² )p/2m]δ_nm + ћ(k – k₀)P_nm(k₀)/m} = 0.

Чтобы найти закон дисперсии E_n(k) по всему k-пространству в пределах первой ЗБ, нужно знать энергии E_n(k₀) и матричные элементы импульса P_nm(k₀). Энергии E_n(k₀) берутся из эксперимента, а под них затем подгоняют матричные элементы P_nm(k₀). Для составления матрицы гамильтониана H_k, в соответствии с последним уравнением необходимо знать еще структуру кристаллического материала, в качестве которого мы рассматриваем полупроводники Ge, Si, GaAs и т.д.

Симметричная матрица kp-гамильтониана (приводится ее верхняя треугольная часть) размерностью 15 × 15 с подогнанными матричными элементами оператора импульса p_i = < u_n | p | u_m > для трех первых векторов обратной решетки: {000}, {111} и {200} представлены табл. 31. Ниже диагонали матрицы приведены значения энергии в симметричных точках ЗБ. Экспериментальные точки, по которым велась подгонка параметров импульса, выделены жирным шрифтом.

Таблица 31

В табл. 32 приведены матричные элементы импульса для Ge, Si и α-Sn по исследованиям других авторов, в частности, [53, с. 59]. Симплекс-методом подгонки параметров, который использовался нами во всех предыдущих случаях, можно добиться аналогичного результата меньшим числом параметров. В четвертой строке табл. 32 приведены параметры, найденные по симметричному направлению Δ, где известно максимальное число экспериментальных точек. Квадратичная ошибка для этого направления принимается за ноль (δ = 0). В пятой строке приведены десять оптимально найденных параметров, дающих минимальную ошибку (δ = 0,000721), вычисленную по всем известным экспериментальным точкам по Δ и L направлениям. Во всех последующих строках табл. 32 число параметров уменьшается. Например, в строке Ge(9) приведены девять параметров для германия. Два наиболее близких параметра, а именно: p₆ и p₉, принимались за один, и далее варьирование производилось уже по девяти параметрам.

Таблица 32

Таким образом, число параметров уменьшилось в два раза при ошибке δ ≈ 0,01. Тем самым была продемонстрирована неоднозначность параметров, поскольку вполне точную зонную структуру можно получить при совершенно различных значениях импульса. На рис. 47 энергетические зоны германия Ge(10) изображены сплошными линиями, Ge(5) (с ошибкой δ = 0,0075) — точечными линиями. Белыми кружками на этом рисунке обозначены экспериментальные точки.

Рис. 47

Число параметров для полупроводников, лишенных центра симметрии, увеличивается еще на шесть дополнительных. Приведем численные значения матричных элементов для GaAs с указанием их места расположения в гамильтониане:

< Г₁₅ | p₁₁ | Г '₂₅ > = 0.121,     < Г '₂ | p₁₂ | Г_1u > = – 0.224,
< Г '₂ | p₁₃ | Г₁ > = 0.324,     < Г₁₅ | p₁₄ | Г '_25u > = 0.156,
< Г '_2u | p₁₅ | Г₁ > = 0.001,     < Г '_2u | p₁₆ | Г_1u > = – 0.209,

Симметрия в центре ЗБ в точке Г в принятой системе обозначений будет следующая: Г '₂₅ для уровней Г₄, Г₈, Г₁₃; Г '_25u для Г₅, Г₁₁, Г₁₅; Г₁₅ для Г₃, Г₁₀, Г₁₄; Г '₁₂ для Г₉, Г₁₂; Г₁ для Г₁; Г_1u для Г₂; Г '₂ для Г₆; Г '_2u для Г₇ (u – верхний). На диагонали матрицы стоят значения энергии H_i = Г_i + k_x² + k_y² + k_z², где k_x, k_y, k_z — координаты волнового вектора k = 2π/r₀ (r₀ = 5,292·10^–9 — боровский радиус). В матрице фигурируют коэффициенты: , , , . Энергия E(k) измеряется в ридбергах (напомним, 1 Ry = 13,6 eV).

На матрице kp-гамильтониана разберем более детально, какие матричные элементы должны быть в принципе отличными от нуля, а какие равны нулю. Данный анализ распространяется на матрицы двух предыдущих приближений. Используя свойства линейного преобразования базиса можно факторизовать матрицу гамильтониана Г(15 × 15) по симметричным направлениям. Для направлений Λ и Σ будем иметь следующее разложение исходной матрицы на прямую сумму матриц меньшей размерности:

T_ΛГ(15 × 15)T_Λ^–1 = Λ₁(7 × 7) + 2Λ₃(4 × 4),

T_ΣГ(15 × 15)T_Σ^{– 1} = Σ₁(5 × 5) + Σ₂(3 × 3) + Σ₃(6 × 6) + Σ₄(1 × 1).

Матрица гамильтониана Г(15 × 15) за счет добавления вектора {220}, которое дает 12 новых энергетических уровней, может быть расширена до размерности Г(27 × 27). Найдем совокупность неприводимых представлений, на которые распадается уровень с вектором {220} в кристаллическом пространстве кубической симметрии O_h, которую имеет алмаз. Задачу разложения приводимого представления на сумму неприводимых удобно начать с примера 6-лучевой звезды {100} ПК решетки (рис. 43а).

Под действием преобразований группы симметрии O_h вектора {100} переходят друг в друга. Например, поворот на угол p вокруг оси k_x (C₄²) приводит к тому, что вектор (001), обозначенный через 1, переходит в вектор 4 (001), а этот, в свою очередь, займет положение 1. Векторы, обозначенные как 5 и 2, также меняются местами, а векторы 6 и 3 остаются на месте. Таким образом, это преобразование может быть записано матрицей, которая переставляет строки в векторе-столбце, а ее удобно представить подстановкой:

(14)(25) = = .

Характер этой матрицы равен h(C₄²) = 2. Характеры всех матриц-представителей 10 классов группы O_h для вектора {100} выписаны ниже:

h(E) = h[(1)] = 6,	h(I) = h[(14)(25)(36)] = 0,
h(C4) = h[(2653)] = 2,	h(IC4) = h[(1245)(36)] = 0,
h(C2) = h[(12)(36)(45)] = 0,	h(IC2) = h[(26)(35)] = 2,
h(C3) = h[(132)(456)] = 0,	h(IC3) = h[(123456)] = 0,
h(C42) = h[(14)(25)] = 2,	h(IC42) = h[(36)] = 4.

Аналогичным образом находим характеры для векторов {220}. Все результаты этих подсчетов занесены в последние строки таблицы характеров группы O_h (табл. 33)

Таблица 33

Для определения числа неприводимых представлений, содержащихся в приводимом представлении, необходимо воспользоваться соотношением:

где .

Здесь E_α — неприводимые представления группы O_h; E_β — приводимое представление группы векторов либо Г{100}, либо Г{220}; m_α — число, указывающее сколько раз представление E_α содержится в E_β; n = 48 — порядок группы O_h; h_β — характер приводимого представления (в нашем случае h_β равно либо h{100}, либо h{220}); h_α — характеры неприводимых представлений группы O_h. Характеры этой группы, группы октаэдра, в которой существует центр инверсии, приведены в табл. 33.

Найдем, сколько раз неприводимое представление Г₁ содержится в Г_{100}:

(6·1 + 8·0 + 3·2 + 6·2 + 6·0 + 1·0 + 8·0 + 3·4 + 6·0 + 6·2) = 1.

В приводимое представление Г_{220} неприводимое представление Г₁ тоже входит 1 раз. Но, например, представление Г '₁₅ уже не входит в Г_{220}:

[12·1·3 + 0·8·0 + 0·3·(–1) + 0·6·1 + 2·6·(–1) + 0·1·3 + 0·8·0 + 4·3·(–1) + 0·6·1 + 2·6·(–1)] = 0.

Окончательные результаты расчета дают следующие разложения:

Г_{100} = Г₁ + Г₁₂ + Г₁₅,     Г_{220} = Г₁ + Г₁₂ + Г₁₅ + Г₂₅ + Г '₂₅.

В табл. 34 приведены разложения и для других звезд обратной решетки.

Таблица 34

	Г1	Г2	Г12	Г '15	Г '25	Г '1	Г '2	Г '12	Г15	Г25
{000}	1	0	0	0	0	0	0	0	0	0
{100}	1	0	1	0	0	0	0	0	1	0
{111}	1	0	0	0	1	0	1	0	1	0
{200}	1	0	0	0	1	0	1	1	0	0
{220}	1	0	1	0	1	0	0	0	1	1
{311}	1	0	1	1	2	0	1	1	2	1
{222}	1	0	0	0	1	0	1	0	1	0
{400}	1	0	1	0	0	0	0	0	1	0

Пусть {φ_i} есть некоторое пространство функций, преобразующихся по приводимому представлению Г_β. Можно выбрать из {φ_i} такие линейные комбинации, которые бы преобразовывались только по соответствующим неприводимым представлениям Г_α. С этой целью из неприводимых представлений Г_α строим проекционный оператор, по стандартной методике [54], [55]:

,

где n — порядок группы O_h, n_α — размерность матриц неприводимого представления Г_α, P(R) — конкретный матричный оператор группы O_h, R — пробегает все элементы группы O_h, которой описывается алмазная решетка полупроводников Ge и Si.

Итак, нам нужен набор функций и набор матриц группы O_h. Можно было бы воспользоваться ортогональными матрицами группы октаэдра O [4, c. 234] однако, во-первых, группа O_h не совсем одно и то же, что группа O; во-вторых, в формуле для проекционного оператора ничего не говорится об ортогональности, т.е. матрицы могут быть общего вида, не обязательно ортонормированными на единицу; в-третьих, нам не мешало бы лучше усвоить на практике процедуру составления матриц представления по заданным базисным функциям. Набор симметризованных функций приведен в табл. 35.

Таблица 35

Гα	Состояния	Симметризованные базисные функции
Г1	s	1
Г15	px,  py,  pz	x, y, z
Г '25	d1,  d2,  d3,	xy, yz, zx
Г12	d3,  d4	3z2 – r2,    x2 – y2
Г25	f1,  f2,  f3,	z(x2 – y2),   x(y2 – z2),   y(z2 – x2)
Г '2	f4	xyz
Г '15	g1,  g2,  g3,	xy(x2 – y2),   yz(y2 – z2),   zx(z2 – x2)
Г '12	h1,  h2	xyz(3z2 – r2),   xyz(x2 – y2)
Г2	i	z4(x2 – y2) + x4(y2 – z2) + y4(z2 – x2)
Г '1	j	xyz[z4(x2 – y2) + x4(y2 – z2) + y4(z2 – x2)]

Теперь начинаем конструировать матрицы неприводимых представлений группы O_h. Матрицами всех одномерных представлений являются сами характеры этих представлений. В отношении трехмерных представлений Г₁₅ и Г '₂₅, отвечающие p- и d-орбиталям свободных атомов, т.е. не находящихся в кристаллической решетке (рис. 42), также никаких проблем не возникает. В частности, для преобразования, входящего в класс С₃, двумя конкретными матрицами представления являются

Г₁₅ (С₃) = , Г '₂₅ (С₃) = .

Несколько сложнее отыскать матрицы двумерного представления Г₁₂. Базисными функциями этого неприводимого представления являются {2z² – x² – y², x² – y²}. Преобразование С₃, по сути, производит циклическую перестановку координат x → y → z → x. Значит, нужно найти такие числа a, b, c и d, чтобы выполнялось условие:

.

Выписав из таблицы характеров (табл. 33) характер (–1), соответствующий представлению Г₁₂, для восьми (четыре прямых и четыре инверсных — 8C₃) поворотов куба вокруг его пространственных диагоналей, присваиваем a = d = –½ (хотя кому-то могут больше понравиться другие значения, например: a = –¾, d = –¼). Далее ищем по двум уравнениям два оставшихся неизвестных b и c. В итоге получим 

Г₁₂ (С₃) = .

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

Итак, все матрицы неприводимых представлений Г₁, Г₁₂ и Г₁₅, на которые распадается приводимое представление Г_{100}, найдены. Тогда, используя проекционный оператор и 6-лучевой базис {1, 2, 3, 4, 5, 6} состояния Г_{100}, найдем новый, штрихованный базис {1', 2', 3', 4', 5', 6'} для прямой суммы неприводимых представлений:

Г₁ : 1' = 1 + 2 + 3 + 4 + 5 + 6,
Г₁₂ : {2' = – 1 – 2 + 2·3 – 4 – 5 + 2·6,   3' = 1 – 2 + 4 – 5},
Г₁₅ :  {4' = 1 – 4,   5' = 2 – 5,  6' = 3 – 6}.

С помощью трансформационной матрицы Т, а лучше — построчного преобразования векторов, произведем разложение матрицы приводимого представления Г_{100} на прямую сумму матриц неприводимых представлений Г₁, Г₁₂ и Г₁₅. В качестве примера возьмем конкретную матрицу, записанную в базисе {1, 2, 3, 4, 5, 6}, и преобразуем ее в квазидиагональный вид посредством базиса {1', 2', 3', 4', 5', 6'}:

→ = .

Проверяем:

1' = 1 + 6 + 2 + 4 + 3 + 5 = (1) · 1',
2' = – 1 – 6 – 4 – 3 + 2 · 2 + 2 · 5 = (–½) · 2' + (–3/2) · 3',
3' = 1 – 6 + 4 – 3 = (–½) · 2' + (½) · 3',
4' = 1 – 4 = (1) · 4',
5' = 3 + 6 = (–1) · 6',
6' = 2 – 5 = (1) · 5'.

Поиск 12 базисных функций для неприводимых представлений Г₁, Г₁₂, Г₁₅, Г₂₅ и Г '₂₅ для интересующего нас состояния Г_{220} ничем не отличается, разве что своей громоздкостью. После того, как найдены неприводимые представления и базисные функции, можно приступать к составлению матрицы гамильтониана Г(27 × 27). Сейчас нам понадобятся правила отбора матричных элементов.

Матричные элементы только тогда отличны от нуля, когда прямое произведение представлений, стоящих в матричном элементе, содержит единичное представление. Для любого матричного элемента центральным множителем является импульс, преобразующийся по представлению радиус-вектора (x, y, z), т.е. по представлению Г₁₅. Кроме того, можно показать [55], что достаточно выполнить условие

Г_α × Г₁₅ = Г_β + Σ_δ Г_δ,

для матричного элемента импульса < Г_α | p | Г_β >, чтобы выполнялось и условие

Г_α × Г₁₅ × Г_β = Г₁ + Σ_δ Г_δ.

Раскладывать произведение представлений на прямую сумму мы уже умеем, поэтому сразу приведем результат: 

Г '₂₅ × Г₁₅ = Г₁₅ + Г₂₅ + Г '₁₂ + Г '₂.

Отсюда автоматически следует, что

< Г '₂₅ | p(Г₁₅) | Г '₁₂ > ≠ 0;     но < Г '₂₅ | p(Г₁₅) | Г₁₂ > = 0.

Проведенные вычисления показали, что общее число матричных элементов, отличных от нуля, в гамильтониане Г(27 × 27) равно 27 (числа совпали случайно). Для поиска числовых значений матричных элементов импульса матрицу Г(27 × 27) необходимо факторизовать по симметричным направлениям, поскольку одновременное варьирование по 27 свободным параметрам даст слишком большую неоднозначность. Под словом «факторизация» в данном случае понимается разложение одной большой матрицы на прямую сумму нескольких матриц меньшей размерности, т.е. сведение исходной матрицы к квазидиагональному виду. Однако технически эта операция выполняется намного проще, чем вышеописанная процедура перехода от нештрихованного базиса для Г_{100} и Г_{200} к штрихованному, где они тоже принимают квазидиагональный вид. Так, чтобы матрица Г(27 × 27) распалась на совокупность Δ-матриц, принимают k_x = k_y = 0, k_z — переменная; тогда получают 

2·Δ₁(6 × 6),     Δ₂(2 × 2),     2·Δ₅(6 × 6),     Δ'₁(1 × 1).

Но и эта простая процедура попадает под операцию подобия; для двух других симметричных направлений будем иметь следующие прямые суммы:

T_ΛГ(27 × 27)T_Λ^–1 = Λ₁(10 × 10) + Λ₂(1 × 1) + 2·Λ₃(8 × 8),

T_ΣГ(27 × 27)T_Σ^–1 = Σ₁(10 × 10) + Σ₂(5 × 5) + Σ₃(8 × 8) + Σ₄(4 × 4).

В табл. 36 верхняя часть является симметричной матрицей Λ₁(10 × 10), нижняя часть — матрица Λ₃(8 × 8), матрица Λ₂(1 × 1) = H₉; в табл. 37 верхняя часть — матрица Σ₁(10 × 10), нижняя — матрица Σ₃(8 × 8); в табл. 38 верхняя часть — матрица Σ₂(5 × 5), нижняя — матрица Σ₄(4 × 4); коэффициенты во всех таблицах равны: с₁ = 1.1547, с₂ = –0.5774, с₃ = 0.7071, с₄ = 0.5, с₅ = 0.866. Далее методом варьирования находят значения всех 27 параметров, опираясь на известные 10. В табл. 39 приведена совместимость центральной точки Г ЗБ с ее направлениями Δ, Λ и Σ.

Таблица 36

Таблица 37

Таблица 38

Таблица 39