Численное моделирование гетероперехода

ЧЕБЫШЕВСКИЙ СБОРНИК Том 11 Выпуск 2

УДК 511.9.

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ГЕТЕРОПЕРЕХОДОВ1

А. В. Ермолов, Н. М. Добровольский

Аннотация

В данной работе рассмотрены некоторые вопросы численного моделирования гетеропереходов с применением теоретико-числовых методов приближенного анализа.

Введение.

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

Среди всего многообразия полупроводниковых материалов следует особо выделить структуры на основе редких земель. Они обладают рядом свойств, которые определяются особенностями электронной структуры атомов Бт, Ей и т.д., а именно наличием недозаполненной 4/-оболочки, Значительное место среди таких материалов занимают халькогениды редкоземельных элементов, часть из которых является магнитными полупроводниками (одновременно проявляют как полупроводниковые, так и магнитные свойства). Магнитные свойства соединений редкоземельных элементов обусловлены прежде всего обмен/

редких земель. Магнитное упорядочение сказывается на характере движения носителей заряда, а концентрация носителей, в свою очередь, существенно влияет на магнитное упорядочение [3].

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

1Работа выполнена по гранту РФФИ 08-01-00790

их электронной и магнитной подсистем, позволяет целенаправленно управлять свойствами материала, и следовательно гетероструктуры на их основе, воздействием на образец внешнего магнитного поля.

Проблема транспорта епин-поляризованных носителей в гетероструктурах привлекает большое внимание в связи с интенсивным развитием спинтроники. Создание спиновых инжекторов и детекторов является одним из главных направлений исследований данной отрасли электроники. Наиболее эффективным методом как детектирования так и создания спиновой поляризации носителей является резонансное туннелирование через двухбарьерную гетероструктуру

[4], При этом рассеяние спина на гетерограницах будет минимально, если материалы гетеропар при одинаковой кристаллической структуре будут обладать практически совпадающими постоянными решеток. Данным условиям удовлетворяют гетероетуктуры на основе моносульфидов европия, самария и свинца кристаллизующиеся в структурном типе МаС1 с необычайно близкими по величине постоянными решеток, составляющими соответственно 5.968, 5.967 и 5.940 А, В качестве ферромагнитного материала предполагается использовать моносульфид европия, обладающем кристаллической структурой - типа МаС1 (группа кристаллической симметрии Гт3т - 0^), и 110 РЯДУ магнитных параметров превосходящем материалы данного класса,

БтБ

тероетруктуре к ЕпБ. Данные материалы, при одинаковом типе кристаллической структуры обладают практически совпадающими постоянными решеток, различие в которых составляет менее 0,05 %, Следует отметить, что внесение Бт ЕпБ

Еп БтБ

ферромагнитное состояние. Данный факт объясняется переходом иона самария в присутствии двухвалентного иона европия в магнитоактивное состояние /6 ^ /5 + 1е- Также следует отметить, что топкие пленки феромагнитного материала зачастую обладают повышенной точкой Кюри, по сравнению с объемными материалами. Это объясняется поверхностным магнетизмом, обусловленным дополнительным косвенным обменом через таммовские уровни, В случае гетероструктур стабилизации ферромагнитного состояния материала следует ожидать благодаря влиянию состояний локализованных на границах гетеропереходов.

Всё это позволяет предположить, что гетероструктуры ЕпБ-БтБ являются весьма перспективными в свете как повышения точки Кюри системы, так и получения епин-поляризованного потока электронов.

Моделирование гетероперехода

Описываемые материалы обладают кристаллической структурой типа МаС1, т,е, катионы К (Еп, Бт, РЬ) расположены в узлах Я гранецентрирован-ной кубической решетки и окружены шестью анионами А (ионы серы Б), находящимися на расстоянии ё, = 0, 5а (а - постоянная решетки) в направлениях, задаваемых 6 векторами октаэдра Т (г = 1,..., 6). Будем считать анион ста = Т

вторым атомом в элементарной ячейке.

Заметим, что ионное состояние таких соединений будет Ме2+Б2-, Электронная конфигурация нейтральных атомов:

1бБ ^ 1в22в22р63в23р4;

63Еп ^ оболочка 54Хе + 4/76з2;

62Бт ^ оболочка 54Хе + 4/656з2;

Образование ионной связи соответствует переходу 2 электронов с атома металла на атомы халькогенида, что соответствует электронной структуре ионов: 16Б2- ^ 1в22в22р63в23р6;

63Еп2+ ^^^очка 54Хе + 4/76з0;

62Бт2+ ^ оболочка 54Хе + 4/66з0;

Таким образом, минимально базисную систему можно ограничить в-орбита-лями валентных оболочек всех катионов и анионов и тремя рХ-, ру-, р^-орбита-лями валентной оболочки каждого аниона. Однако реальный кристалл не будет полностью ионным, но типы симметрии занятых орбиталей будут соответствовать ионной структуре. Заполненные зоны образуются из з- и р-епин-орбиталей аниона, способных вместить восемь внешних электронов на молекулу, но имеется сильная гибридизация орбиталей аниона с р и ^-орбиталями атомов металла в валентной зоне, как следствие перенос значительной величины заряда на изначально пустые 5^-орбитали этого металла [2], При этом согласно теории кристаллического поля лигандов в октаэдрическом окружении ^-состояния металла испытывают расщепление на дуплет и триплет, перекрытие оболочек наиболее значительно для состояний дуплета - симметрия которых х2 — у2, г2, как следствие именно эти состояния вносят определяющий вклад в ковалентную составляющую связи. Валентные зоны, образуемые из связывающих орбиталей, и зоны проводимости, образуемые из антисвязывающих орбиталей, являются комбинацией рподобных орбиталей отрицательного нона и ^-подобных орбиталей иона металла, результат этого наличие сильной ковалентной связи наряду с ионной. Состояния лигандов также претерпят смешения вследствие воздействия поля окружающих ионов металла.

Рассмотрим состояния на гетерогранице, В простейшем случае достаточно

Еп Бт Б

образом, достаточно расширить систему базисных функций, учетом орбитали дуплета каждого катиона и аниона, обеспечивающей максимальное перекрытие. Тогда 11 различных базисных з-, р-, и /-орбитадей (3-для аниона Б и по 4

Еп Бт

ф1(х) = (х); ф4,8(х) = ^а (х);

ф2 (х) = ^2 (х) ; ф5,9 (х) = (х);

фз(х) = (х); ф6,ю (х) = ^2(х);

ф7,11(х) = ^а (х);

где

х/у/х (Х

1^Лх2/

1~5~х2 — у2, г2 — г2

V 4п г2

Др/К

радиальные части орбиталей Лёвдина,

Из данного базиса конструируется волновая функция

Использование орбиталей Левдина позволило систему неортоганальных орбиталей, находящихся на разных атомах, преобразовать в ортоганальную с сохранением атомной симметрии, и как следствие из гамильтониана системы выпадают интегралы перекрытия,

ЕпБ — БтБ

будет сверхобмен посредством лигандов (механизм Крамерса-Андерсона), Задача соответственно имеет вид 3-х центров и остановимся на 8-ми электронах. Распределяя всеми возможными способами спиновые функции по наборам орбит основного и возбужденного состояний, получаем N базисных функций. Так как операторы энергии И, квадрата полного момента Б2 и проекции на ось г полного момента Бх попарно коммутируют друг с другом, а также матричные

И

лю, то матрица энергии разобьется на квадратные блоки, соответствующие всем мультиплетноетям системы.

Рассмотрим матрицу Р (Бх) всех возможных распределений электронов по спинам для заданного Бх,

Например

а а а ... в

а а а ... в

а а а ... в

в в а ... а

в в в ... а

в в в ... а

здесь Р(0)^ определяет направление спина г-того электрона в к-том наборе распределений.

Также рассмотрим матрицу К(б^.) всех возможных распределений электронов по орбитам левого катиона /Еи1, рЕи1, ^Еи1, зЕи1, правого катиона /5т2, Р^т2, ^5т2, в5т2 и центрального аниона ^для заданного

к

3 Уе«1 3 3

Ре«1 Рём1 Ре«1 Рём1 •••

Ре«1 РЄм1 Ре«1 ^Ем1 •••

Рз Рз •••

Рз Рз Рз •••

Рзт2 Рзт2 Рзт2 Рзт2 •••

Рзт2 Рзт2 Рзт2 Рзт2 •••

Узт2 Узт2 Узт2 Узт2 •••

Тогда детерминант Слетера, для 1-1 ого распределения по спину, ]-'юго распределения по орбиталям и для заданного можно записать в виде:

из*

/е«1

/е«1

Рз

Рз

Р3т2

Рзт2

/е«1

Ре«1

Рз

р (Б* )|,

/е«1

где справа указана лишь сокращенная форма записи и не подразумевается умножение.

Введем обычным образом спиновые операторы Б+а = 0, Б+в = а, Б-а =

в, Б-в = 0 [1], Результат действия оператора повышающего мультиплетноеть

системы 5+ на Ф0 - полный набор волновых функций с мультиплетноетыо 8^=0,

будет равен Б+Ф0 =01 Ф1, где элементы матрицы 0, равны 1, если путем замены

Ф0

Ф1

типлетность системы оператора Б- па Ф1 - полный набор волновых функций с

мультиплетноетыо Б*=1, будет равен Б Ф1

10 Ф0

(01)Т Ф

Ф можно представить в виде Ф

и

1,1

и

1,р

и

и

и

,,1

,,р

Важно отметить, что оператор квадрата полного спина для состояний с проекцией спина Б*=0 можно представить как произведение спиновых операторов

понижающих и повышающих суммарный спин системы

52Ф0 = (5-5+ + 52 + 5,) Ф0 = 5-5+Ф0 = (01) (01)Т Ф0 = ММ0Ф0

Собственные функции V5,52 состояния с мультиплетностью Б, Б, будем искать как С5,52 • Ф5г, Тогда собственные функции триплетного состояния (8=1, 8,2=0) приобретут вид V1,0 = С1,0 • Ф0, для которой с одной стороны 5^1,0 = 5(5 + 1^1>0 = ^ 1>^ а с другой стороны 5^ 1>0 = (5-5+ + ^ + 5,) V 1>0 = в-в+V1’0, Тогда для ьой функции триплета имеем:

1,0 1,0 0

V = С • Ф = 2_^ Ск«к

к

5V/’0 = 5^С'0 • Ф0 = С1’0 • ММ0 • Ф0 = 2Сг1,0 • Ф0 В силу линейной независимости ик составляющих вектор Ф0 получаем

С1>0 • ММ0 = 2С1>0,

С 1,0

му числу 2 = 5(5 + 1). Ранг матрицы |ММ0 — 5(5 + 1)| определяет количество

С 1,0

1,0

соответственно V

По аналогии собственные функции синглетного состояния

V0,0 = С0,0 • Ф°, а для ^ой функции имеем:

С0,0 • ММ0 = 0 • С^0,0 Остается только найти матричные элементы оператора Н вида,

г,0*

3 А / 7 2 ^2 ^2 ^2

ЕЕЕЕ(^.- е е е

_г=1 .7 = 1 А=1 В=1

2<т ^Аг ^АВ

V г>0^т,

что позволит в свою очередь определить энергию синглетного Е и триплетного Е состояний и соответственно энергию эффективного обменного взаимодействия атомов металла </эфф = \ (Е3 — Е4),

Методика численного расчета

Гетеропереходы в этих структурах практически не изучены, хотя от них зависят все основные электрические, магнитные и оптические характеристики гетеросистемы. Решение подобного рода задач возможно в рамках приближения сильной связи, базирующемся на использовании метода ЛКАО исследуемых многоатомных систем, В расчетах из первых принципов в данном методе требуется вычисление многоцентровых кулоновеких и обменных интегралов, подразумевающими привлечение теоретико-числовых методов (например - метод Коробова), Теоретическое исследование подобного рода систем наиболее

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

При использовании неэмпирических методов расчета основные затраты времени ЭВМ уходят на вычисление интегралов межэлектронного взаимодействия. При увеличении размеров молекулы число этих интегралов растет пропорционально п4, где п - размер базиса молекулярной орбитали.

Расчет двух-, трех- и четырехцентровых интегралов представляет сегодня одну из наиболее сложных проблем, встречающихся в квантовой теории молекул, В виду большого объема вычислений необходимо этот процесс запрограммировать для использования на ЭВМ,

Обменный и кулоновекий интегралы обладают рядом специфических особенностей, и попытки вычислить их аналитически приводят к достаточно громоздким выражениям. Однако, наряду с аналитическим подходом, получило распространение прямое численное интегрирование. Преимущества очевидны: нет нужды выполнять предварительные, весьма трудоемкие аналитические выкладки для получения окончательных формул, пригодных к программированию; существенно расширяется класс волновых функций, интегралы от которых могут быть найдены одним и тем же методом. Не менее очевидны и недостатки прямого численного интегрирования молекулярных интегралов: до сих пор не достигнута приемлемая скорость счета; трудно контролировать точность вычислений. Эти обе важные стороны метода в значительной степени определяют метод численного интегрирования.

Исторически для приближенного расчета многоцентровых интегралов применялись следующие два метода: метод Барнетта и Коулеона и метод гауссовых функций, или метод Бойза, Эти методы, хотя они и совершенно различны, имеют одну общую черту: в основе каждого из них лежит довольно простое преобразование, которое позволяет записать функцию, локализованную около данной точки, через функции, локализованные в другой точке, В методе Барнетта и Коулеона все волновые функции, входящие в многоцентровые интегралы, разлагаются в ряд около одного из атомов, принятого за центр. После этого проблема вычисления этих интегралов становится формально схожей с проблемой вычисления соответствующих интегралов для атома.

Метод Барнетта и Коулеона основан на представлении атомных орбиталей в виде линейных комбинаций членов вида, г“ ехр(—аг), умноженных на сферические функции углов, где а и а — постоянные. Известно, что, составляя суперпозицию достаточного числа функций такого типа, можно получить очень хорошее приближение к правильным атомным орбиталям. Тогда проблема сводится к тому, чтобы начать с атомной орбитали, локализованной у атома В и

переписать ее через функции, локализованные у атома Л, Как только это будет сделано, дальнейший ход вычисления ясен, хотя практически этапы расчета могут быть очень трудоемкими,

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

Метод Бойза проще метода Барнетта и Коулеона, поскольку он сводится к вычислению интегралов для набора гауссовых функций. Однако он имеет и недостаток: для того чтобы выразить приемлемым образом атомные орбитали как линейные комбинации гауссовых функций, последних требуется гораздо больше, чем функций вида, г“ехр(-аг). Это означает, что нужно скомбинировать большее их число, чтобы получить приемлемое значение интегралов,

В настоящей работе для расчета кулоновеких и обменных интегралов, которые встречаются в расчетах энергетического спектра гетероструктур использовался метод параллелепипедальных сеток Коробова,

В данном методе для вычисления интегралов произвольной кратности применяются квадратурные формулы с параллелепипедальными сетками:

где а1;а2,.., ,ав - специальным образом выбранные целые числа (оптимальные коэффициенты). Для вычисления оптимальных коэффициентов мы использовали алгоритмы из работы [5], а для организации вычислений руководствовались подходами из работы [5],

Для вычисления кратных интегралов по квадратурной формуле

использовалась частичная периодизация подинтегральной функции, так как по некоторым переменным она была выполнена.

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

Прежде всего мы должны выразить прямоугольные координаты х», у* и г»

через эллиптические Аі, ^ И фі.

К сов Рі [~^~2 ~ ~ КдШ'Рі [~^~2 ~ “ ~ ХіЦіЯ

Хг = -------1у(\2 - !) • = ---^V/^2 ~~ ^ ^ ~ ^

Соответственно расстояния от ядер о и Ь до і- того электрона в эллиптической системе координат будут иметь вил:

Гаі = !^\Лі + 5 Ты = "^"(Аг — Ці),

а элемент объема

Д3

ІУі = —

и расстояние между электронами -

Г12 = I [(А| - 1) • (1 - (4) + (А? - 1) • (1 - ц\) + (А2//2 - А1/Х1)2-

______________________________________.. I

2 сов((/?2 - ¥?1)^/(А| - 1) • (1 - /X2) • (А2 - 1) • (1 - /X2) 2.

Для расчета интегралов воспользуемся квадратурной формулой (1) и возможностями среды Ма11н-а<1. Однако перед тем, как применять квадратурную формулу, область интегрирования подлежащих вычислению интегралов необходимо свести к единичному многомерному кубу, В нашем случае это реализуется путем следующей замены переменных:

Аі = ± /її = 2х3 - 1 <~рі = 2тгх5 А2 = -3- = 2х4 - 1 ір2 = 2ттх6

Результаты расчета

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

СПИСОК ЦИТИРОВАННОЙ ЛИТЕРАТУРЫ

[1] Бердышев А, А, Введение в квантовую теорию магнетизма, Екатеринбург, изд-во Уральского университета, 1992 г., 276 С,

[2] Губанов В,А,, Жуков В,П., Литинский А,О, Полуэмпиричеекие методы молекулярных орбиталей в квантовой химии, — М.: Наука, 1976,

[3] Головнев Ю.Ф, Наноразмерные ферромагнитные гетеросистемы, — Тула: Изд-во ТГПУ, 2007. - 262 с.

[4] Головнев Ю.Ф., Ермолов А.В. Влияние ферромагнитного порядка на квантоворазмерные эффекты в гетероструктурах на основе ферромагнитного полупроводника ЕиЭ / / Химия твердого тела и современные микро- и нанотехнологии: Материалы VII Международной научной конференции. — Кисловодск, 2007. — С. 261 — 263.

[5] Бочарова (Добровольская) Л. П. Алгоритмы поиска оптимальных коэффициентов // Чебышевекий сборник 2007 Т. 8. Вып. 1(21). Ту-ла, Из-во ТГПУ им. Л. Н. Толстого. С. 4 — 109.

[6] Добровольская Л. П., Добровольский Н. М., Симонов А. С. О погрешности приближенного интегрирования по модифицированным сеткам // Чебышевекий сборник 2008 Т. 9. Вып. 1(25). Тула, Из-во ТГПУ им. Л. Н. Толстого. С. 185 - 223.

ТГПУ им. Л. Н. Тол сто го Поступило 27.11.2010