ЖЭТФ, 2021, том 160, вып. 2 (8), стр. 249-274

ЭНЕРГИИ ОБРАЗОВАНИЯ ВАКАНСИЙ КУБИЧЕСКОЙ ФАЗЫ

МАГНЕТИТА В РАМКАХ DFT+U

М. И. Шутиковаa,b*, В. В. Стегайлов^a,b

^a Объединенный институт высоких температур Российской академии наук

125412, Москва, Россия

^b Московский физико-технический институт (Национальный исследовательский университет)

141701, Долгопрудный, Московская обл., Россия

Поступила в редакцию 22 марта 2021 г.,

после переработки 22 марта 2021 г.

Принята к публикации 20 мая 2021 г.

В рамках теории функционала спиновой плотности с учетом сильных электронных корреляций (метод

DFT+U) вычислены свойства катионных вакансий кубической фазы магнитного оксида Fe3O4 (маг-

нетита). Проанализировано влияние ключевых параметров модели (значения эффективной поправки

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

ной плотности) на описание электронной структуры кубической фазы магнетита. Предложена процедура

поиска решения с наименьшей полной энергией, решающая свойственную для DFT+U проблему ло-

кальных минимумов. Охарактеризованы найденные варианты зарядового и орбитального упорядочения

основного состояния электронной системы. Описаны различные варианты оптимизации структуры мо-

делей катионных вакансий и определены соответствующие энергии формирования.

DOI: 10.31857/S0044451021080101

которой требуется рассмотреть несколько дефект-

ных конфигураций и учесть эффект релаксации ло-

кального окружения дефекта.

1. ВВЕДЕНИЕ

Метод DFT+U был применен для описания низ-

Магнетит (Fe₃O₄) — магнитный оксид железа

котемпературной фазы магнетита [8,12-22] и его ку-

с сильными корреляциями в электронной подси-

бической фазы [4-7, 9-11, 16, 18]. Он требует обос-

стеме; их присутствие сильно усложняет и обого-

нованного выбора значения полуэмпирической по-

щает физику его электронного строения [1, 2]. По-

правки. Один из способов ее верификации заклю-

этому для вычисления свойств магнетита требует-

чается в сравнении вычисленного значения ширины

ся теоретический метод, способный учитывать эф-

запрещенной зоны с ее экспериментальным значе-

фекты электронных корреляций. Стандартные при-

нием [5, 18, 24]. Однако существующие в настоящее

ближения обменно-корреляционных функционалов

время экспериментальные данные не позволяют сде-

в методе теории функционала плотности (density

лать однозначного утверждения о том, существует

functional theory, DFT) сделать этого не позволяют.

ли в кубической фазе магнетита запрещенная зона

или нет. Во Введении к работе [5] приводятся ссылки

В последнее время в работах по ab-initio расче-

там свойств магнетита для учета сильных электрон-

на экспериментальные данные, интерпретация кото-

рых позволяет говорить о ее существовании.

ных корреляций в основном используется два под-

хода: 1) полуэмпирическая поправка Хаббарда (ме-

Методу DFT+U присуща проблема поиска реше-

тод DFT+U [3]) [4-22] и 2) гибридные функционалы

ния, обладающего наименьшей полной энергией сре-

[5, 9, 23]. Расчеты с гибридными функционалами —

ди множества локальных минимумов. Данная про-

чрезвычайно ресурсоемкие, потому в настоящее вре-

блема обсуждается в ряде работ, посвященных ис-

мя практически невозможно использовать этот ме-

следованию соединений других d- и f-металлов с

тод для получения модели дефектов магнетита, в

недостроенными оболочками [25-28]. Основным от-

личием подхода DFT+U от стандартных приближе-

^* E-mail: shutikova.mi@phystech.edu

ний теории функционала плотности является его за-

249

ЖЭТФ, вып. 2 (8)

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Таблица 1. Литературные данные. Энергии обра-

зования EFVA и EFVB (эВ) вакансий в A- и B-под-

решетках магнетита и их разница ΔAB = EVA -

- EFVB (эВ)

Метод

E^FV

Δ_AB

−

1.10

[4]

DFT+U

3.45

2.30

1.15

[6]

−

1.01

[7]

−

1.27

[10]

DFT

2.69

0.83

1.86

[29]

висимость от орбитальных степеней свободы систе-

мы. Эвристические методы поиска решений с наи-

меньшей полной энергией в рамках DFT+U были

предложены в работах [26, 27]. Можно назвать их

основой контроль за получаемыми в DFT+U чис-

Рис. 1. Кристаллическая структура кубической фазы маг-

лами заполнения локализованных орбиталей. Попы-

нетита (для визуализации используется VESTA [31]). Ионы

ток их применения для исследования кубической

железа и кислорода представлены крупными и маленьки-

ми шарами соответственно. Катионы A-подрешетки отме-

фазы магнетита в литературе найдено не было.

чены темно-синим цветом, двух- и трехвалентные катио-

Свойства дефектов магнетита в рамках DFT

ны B-подрешетки — голубым и золотистым соответствен-

и DFT+U получены в небольшом числе работ

но (показан один из вариантов распределения Fe²⁺ и Fe³⁺

[4, 6, 7, 10, 11, 29]. В работе [29] расчет свойств кати-

по узлам B-подрешетки)

онных дефектов кубической фазы магнетита про-

водился без учета сильных электронных корреля-

составляют 3.45 эВ и 2.30 эВ соответственно, в то

ций. В работах [10,11] для рассмотренных в рамках

время как без поправки для тех же дефектов были

DFT+U катионных и анионных дефектов были по-

получены значения 2.69 эВ и 0.83 эВ [29] (табл. 1).

дробно описаны свойства локального окружения, но

Для тетраэдрического и октаэдрического междоуз-

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

лий в работе [6] были получены энергии образова-

(DFT+U)-работах [4, 7, 10] приводится только раз-

ния 1.25 эВ и 0.81 эВ соответственно, но без уче-

ница в энергиях образования вакансий A- и B-под-

та поправки в работе [29] их значения составляют

решеток¹⁾ магнетита (табл. 1) и отмечается, что об-

5.40 эВ и 3.26 эВ. Концентрации дефектов в равно-

разование B-вакансии более вероятно, чем образо-

весных условиях экспоненциально зависят от энер-

вание A-вакансии. В недавней работе по исследова-

гий образования, поэтому такие разногласия в тео-

нию кислородных вакансий в перовските SrTiO₃ [28]

ретических предсказаниях дают огромные неточно-

обсуждаются дополнительные приближения (метод

сти в определении концентраций дефектов.

DFT+U+V), которые могут потребоваться для ис-

Противоречивые литературные данные опреде-

следования точечных дефектов в оксидах; однако в

лили цель настоящей работы — разобраться в при-

перечисленных выше работах по исследованию де-

чинах возникновения таких отличий и получить

фектов магнетита такой подход не применялся.

энергии образования катионных вакансий кубиче-

Значения энергий образования дефектов, полу-

ской фазы Fe₃O₄ в рамках DFT+U. Однако с пер-

ченные в рамках DFT+U в работе [6], очень отли-

вой попытки сделать это не удалось: были получе-

чаются от значений, полученных ранее в работе [29]

ны отрицательные энергии образования дефектов.

без поправки Хаббарда. Действительно, энергии об-

При этом вычисленные ранее в других работах зна-

разования A- и B-вакансий, полученные в работе [6],

чения не были отрицательными (табл. 1). Поэтому

мотивацией работы стал поиск таких приближений

¹⁾ В современных работах по исследованию свойств магне-

и параметров расчета, которые позволили бы полу-

тита [5,10,16] для обозначения его тетраэдрической и октаэд-

рической катионных подрешеток используются буквы A и B

чить адекватные свойства дефектов кубической фа-

соответственно.

зы в DFT+U. Пристальное внимание было уделено

250

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

существующим в литературе ab-initio моделям без-

ра низкотемпературной фазы при 90 K описывается

дефектной кубической фазы магнетита.

с помощью не менее 168 координат (замороженных

фононных мод с амплитудами менее 0.24Å).

В связи c тем, что основным состоянием магне-

1.1. Свойства бездефектной кубической

тита является состояние не с кубической, а с низ-

фазы магнетита в рамках DFT+U

котемпературной симметрией, вычисление свойств

бездефектной кубической фазы в рамках DFT+U

Известно, что в нормальных условиях магнетит

остается возможным только с условием сохране-

имеет структуру обращенной шпинели (Fd3m). Эле-

ния ее симметрии. В противном случае, при опти-

ментарная ячейка этой структуры (рис. 1) содержит

мизации геометрии без сохранения симметрии, бу-

56 атомов (24 катиона железа и 32 аниона кисло-

дет найдено решение с искажением решетки. Та-

рода согласно ионной модели). Анионы кислорода

кое искажение решетки с потерей кубической сим-

образуют гранецентрированную подрешетку, иска-

метрии при оптимизации координат атомов в рам-

женную внедренными в ее междоузлия катионами

ках DFT+U показано в работе [18]. Искажение да-

железа. Одну восьмую имеющихся тетраэдрических

ют также вычисления с гибридным функционалом

междоузлий (позиции 8 a в нотации Вайкоффа) за-

B3LYP [23]. При оптимизации геометрии в методе

нимают трехвалентные катионы железа Fe3+A, обра-

DFT без учета сильных электронных корреляций

зуя тетраэдрическую или A-подрешетку. Одну вто-

такого искажения не наблюдается [18,23]; в литера-

рую имеющихся октаэдрических междоузлий (пози-

туре это объясняется тем, что стандартные обменно-

ции 16 d) занимают в произвольном порядке двух-

корреляционные функционалы DFT не дают точно-

и трехвалентные катионы железа, Fe2+B и Fe3+B, об-

го описания низкотемпературной фазы магнетита.

разуя октаэдрическую или B-подрешетку. Формулу

магнетита в этих обозначениях можно записать как

Только в очень давней работе [38] по исследова-

(Fe³⁺)_A(Fe²⁺Fe³⁺)_B O2-4 [30]. До температуры Нееля

нию кубической фазы магнетита оптимизация гео-

(850 К) магнетит находится в состоянии ферримаг-

метрии не проводилась совсем; расчеты выполня-

нитного упорядочения: магнитные моменты катио-

лись с фиксированной ионной конфигурацией, соот-

нов A- и B-подрешеток направлены противополож-

ветствующей экспериментальным значениям кисло-

ным образом и нескомпенсированы.

родной координаты и постоянной решетки. В насто-

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

ящее время в рамках DFT+U для кубической фа-

требуются два параметра: постоянная решетки a₀ и

зы магнетита вычисляются оптимизированные па-

кислородная координата x. Экспериментальные зна-

раметры a₀ и x [5, 6, 9, 10, 16, 18]. При этом опти-

чения этих параметров приведены в табл. 2. Кисло-

мизированные в рамках DFT+U значения постоян-

родная координата показывает степень искажения

ной решетки оказываются примерно на 1-2 % выше

анионной подрешетки из-за размещения в ней кати-

экспериментального значения (табл. 2), и чем боль-

онов железа. В работе [32] отмечается, что ее зна-

ше используемое значение поправки Хаббарда, тем

чение остается практически постоянным вплоть до

больше отклонение от экспериментальных данных

600 K. Если бы кислородная подрешетка была иде-

[5]. Оптимизированные значения кислородной коор-

альной гранецентрированной (x = 0.2500), то углы

динаты также могут отличаться от эксперименталь-

между положительными направлениями трех осей

ного значения (табл. 2). В некоторых работах по-

всех октаэдров Fe_BO₆ составляли бы 90^◦. При на-

дробного описания того, какие именно параметры

блюдаемом для магнетита экспериментальном зна-

оптимизировались (объем, векторы решетки, коор-

чении (x = 0.2549) эти углы — меньше 90^◦, и окта-

динаты атомов), а какие были фиксированы, не при-

эдры Fe_BO₆ искаженны одинаковым образом.

водится. Зачастую указывается только оптимизиро-

Кубическая фаза магнетита ограничена снизу

ванное значение постоянной решетки, но не значе-

низкотемпературной фазой, в которую она пере-

ние кислородной координаты. Не указывается и то,

ходит при температуре Вервея (120 ± 5 K) [33].

какая начальная ионная конфигурация использова-

Считается, что низкотемпературная фаза магнети-

лась при запуске расчетов с оптимизацией геомет-

та имеет не кубическую, а моноклинную симметрию

рии. В работах [5, 16] проводилась оптимизация и

[34-37]. Элементарная ячейка низкотемпературной

постоянной решетки, и кислородной координаты; в

√

фазы имеет размеры (

2a₀,

2a₀, 2a₀), где a₀ — по-

этих работах и в работе [9] указывается, что опти-

стоянная решетки кубической фазы, и содержит 224

мизация проводилась с сохранением симметрии ку-

атома. Согласно уточненным данным [37], структу-

бической фазы.

251

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Таблица 2. Литературные данные. Свойства кубической фазы магнетита, полученные в рамках DFT+U, гибридных

функционалов, DFT (без поправки Хаббарда) и экспериментальные данные. Обозначения: a0 (Å) — постоянная

решетки, x — кислородная координата, μ (μB /f.u.) — магнитный момент формульной единицы Fe3O4, μ_Fe3+ ,

μ_Fe2+

, μ_Fe3+ , μ_O2- (μB) — магнитные моменты частиц A-, B- и O-подрешеток, Eg (мэВ) — ширина запрещенной

зоны, U_eff (эВ) — значение эффективной полуэмпирической поправки Хаббарда, Nat — число атомов в расчетной

сверхячейке

a₀

μ_Fe3+

μ_Fe2+

μ_Fe3+

μ_O2-

E_g

Метод U_eff N_at

Год

8.446

4.0

DFT+U

3.80

2006

[18]

8.446

0.2549

-3.97

+3.87

+0.04

DFT+U

3.20^a

2007

[16]

8.271

4.0

-3.68

+3.62

DFT+U

3.94

2013

[10]

8.488

4.0

-4.09

+(3.96 ÷ 3.98)

+0.03

DFT+U

4.00

2014

[9]

8.477

0.2548

-3.93

+3.64

+4.01

DFT+U

3.00

8.491

0.2548

-3.96

+3.62

+4.06

177

DFT+U

3.50

2017

[5]

8.505

0.2549

-3.99

+3.64

+4.09

474

DFT+U

4.00

8.488^b

-4.06

+(3.49 ÷ 4.37)

+0.03

HSE 25 %

2014

[9]

8.389

0.2547

-4.21

+3.79

+4.27

421

HSE 25 %

2017

[5]

8.488^b

-3.92

(+3.82 ÷ 3.85)

+0.06

HSE 15 %

2014

[9]

8.395

0.2548

-4.10

+3.92

+4.03

HSE 15 %

2017

[5]

8.370

4.0

-3.63

+3.60

+0.11

DFT

2003

[29]

8.377

0.2545

-3.40

+3.52

+0.08

DFT

2007

[16]

8.387

4.0

-3.47

+(3.55 ÷ 3.58)

+0.08

DFT

2014

[9]

8.380

0.2546

-3.58

+3.63

DFT

2017

[5]

8.394

эксп.

XRD^d

1951

[43]

8.394

эксп.

XRD

1953

[44]

8.398

эксп.

XRD

1977

[34]

8.394

0.2549

эксп.

XRD

1981

[45]

8.396

0.2549

эксп.

XRD

1996

[32]

−

100*

эксп.

PES^e

1997

[41]

−

70*

эксп.

STS^f

2020

[39]

Примечание.^a Вместо U_eff в работе использовались значения U = 4.0 эВ, J = 0.8 эВ.

^b Расчет с гибридными функционалами проводился для геометрии, оптимизированной в DFT+U.

^c Запрещенная зона есть, но ее ширина не указана в статье.

^d Рентгеноструктурный анализ.

^e Фотоэмиссионная спектроскопия.

^f Сканирующая туннельная спектроскопия.

^∗ В обоих случаях речь идет о наблюдаемой в эксперименте псевдощели.

252

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

Экспериментальные данные [39-41] указывают

пользование малых сверхячеек внести в модель ар-

на то, что ниже температуры перехода Вервея маг-

тефакты, связанные с их размером.

нетит является полупроводником. Для ширины его

В большинстве работ, выполненных для кубиче-

запрещенной зоны в этой области температур бы-

ской фазы магнетита в DFT+U [9,10,16,18], как вид-

ли получены оценки: 150 мэВ [41], 100 мэВ [40] и

но из (табл. 2), не предсказывается разница магнит-

320 мэВ [39]. Предполагается, что в низкотемпера-

ных моментов двух- и трехвалентных B-катионов и

турной фазе магнетита кроме ферримагнитного мо-

запрещенная зона. Катионы октаэдрической подре-

жет наблюдаться также зарядовое, орбитальное и

шетки в этих моделях представляются одинаковыми

тримеронное упорядочения [42]. Однако остается от-

частицами, подобно тому, какими они получаются

крытым вопрос о том, какие из черт низкотемпера-

в модели DFT, где сильные электронные корреля-

турной фазы и в какой мере характерны для куби-

ции не учитываются. Однако в работе [5] в рамках

ческой фазы магнетита.

DFT+U с очень близкими к используемым в осталь-

ных работах значениями U_eff разница между μ_Fe3+

Обоснование выбранного значения полуэмпири-

и μ_Fe2+ была получена и составила 0.44 μB2). Разни-

ческой поправки для учета сильных электронных

корреляций является важной деталью вычислений

цу Fe2+B и Fe3+B для кубической фазы предсказывают

в методе DFT+U. Значения U_eff = U - J, исполь-

также расчеты с гибридным функционалом HSE06

зуемые в работах по моделированию бездефектной

[5,9]. При этом в зависимости от доли точного обме-

кубической фазы магнетита, находятся в основном

на в функционале Хартри-Фока и в зависимости от

в интервале от 3.2 до 4.0 эВ (табл. 2). В некоторых

выбранного значения U_eff в DFT+U в этих подхо-

работах [9, 10] при выборе поправки U_eff ссылают-

дах может быть также предсказана ненулевая ши-

ся на значения, полученные ранее в других работах.

рина запрещенной зоны в кубической фазе магне-

Авторы [5,18] выбирали значение U_eff исходя из ши-

тита. При условии, что ширина запрещенной зоны

рины запрещенной зоны (E_g). В работе [18] значе-

не будет слишком большой, модель кубической фа-

ние E_g, вычисленное после оптимизации геометрии

зы с запрещенной зоной не вызовет противоречий

(без сохранения кубической симметрии), сравнива-

с наблюдаемым на эксперименте порядком скачка

лось с экспериментальным значением ширины за-

проводимости при переходе Вервея [41].

прещенной зоны для низкотемпературной фазы. В

Авторы [5] указывают, что причина, по которой

работе [5] значение E_g для кубической фазы маг-

отличия B-катионов и ширина запрещенной зоны

нетита, полученное после оптимизации геометрии с

появились в модели кубической фазы, — это пред-

сохранением кубической симметрии, сравнивалось

положение о том, что волновые функции и элек-

со значением, вычисленным с помощью гибридных

тронная плотность необязательно должны удовле-

функционалов — метода, не требующего полуэмпи-

творять дополнительным требованиям симметрии,

рической поправки. По мнению авторов [5] значение

накладываемым на них симметрией решетки куби-

U_eff = 3.5 эВ дает более точные предсказания для

ческой фазы. В работе [5] показано, что такие не

ширины запрещенной зоны и магнитных моментов

обладающие дополнительной симметрией решения

катионов кубической фазы магнетита.

находятся ниже по энергии, чем симметричные.

Ввести такое предположение кажется естествен-

Интересно отметить, что авторы [18] при вери-

ным при расчетах свойств точечных дефектов. По-

фикации U_eff по ширине запрещенной зоны низко-

этому вызывают интерес данные по расчету дефек-

температурной фазы проводили оптимизацию гео-

тов в работе [10], где, несмотря на значение поправ-

метрии в сверхячейках с 56 атомами, в то время

ки, при котором в работе [5] уже была получена за-

как элементарная ячейка низкотемпературной фа-

прещенная зона и разница в магнитных моментах

зы гораздо больше. Возникает вопрос, насколько ве-

B-катионов, ничего такого получено не было. Неко-

лики погрешности верификации U_eff при модели-

торые детали вычислений в работе [10] не упомина-

ровании низкотемпературной фазы в сверхячейках,

ются, поэтому установить точнее причины отличий

меньших ее элементарной ячейки, особенно при уче-

в предсказаниях работ [5,10] не представляется воз-

те ее сложной геометрии [37]. В работе [5] при вери-

фикации U_eff по E_g и магнитным моментам катио-

нов кубической фазы использовались некубические

²⁾ Точнее, наверное, сказать, что отличия в вычисленных

значениях магнитных моментов (и зарядов) B-катионов поз-

сверхячейки с 14 атомами, меньшего размера, чем

воляют разделить их на две группы, интерпретировать как

элементарная ячейка кубической фазы с 56 атома-

двух- и трехвалентные и обозначить Fe2+B и Fe3+B, соответ-

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

ственно, поскольку разница составляет менее 1 μ_B .

253

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

можным. В работе [9] в DFT+U с также уже доста-

метрам вычислений, после чего в рамках получен-

точно большим значением поправки U_eff ни отли-

ного набора допущений вычисляются свойства изо-

чий в магнитных моментах B-катионов, ни E_g полу-

лированных катионных вакансий магнетита.

чено не было (в этой работе рассматривался толь-

ко бездефектный магнетит, поэтому можно предпо-

ложить, что ее авторы при оптимизации геометрии

2. ДЕТАЛИ ВЫЧИСЛЕНИЙ

накладывали на волновые функции и плотность те

требования симметрии, от которых авторы [5] отка-

Все вычисления в работе проводились в рамках

зались).

теории функционала спиновой плотности с учетом

К числу приближений, о которых не упомина-

сильных электронных корреляций (метод DFT+U).

ется в деталях вычислений большинства работ по

Расчеты выполнялись в коде VASP [46, 47] с бази-

исследованию кубической фазы магнетита в рам-

сом плоских волн и приближением псевдопотенциа-

ках DFT+U (явные упоминания о них встречаются,

ла. PAW-псевдопотенциалы, полученные в прибли-

пожалуй, только в работах без поправки Хаббарда

жении обменно-корреляционного функционала PBE

[29, 38]), относятся начальные приближения к вол-

и включающие восемь и шесть валентных состояний

новым функциям и электронной (спиновой) плот-

для железа и кислорода соответственно, были взяты

ности, необходимые для запуска самосогласованного

из библиотеки этого кода.

расчета в методе теории функционала плотности. В

Начальные приближения к волновым

работе [18] вскользь отмечается, что некоторые ре-

функциям и электронной плотности. На-

шения для кубической фазы, обладающие наимень-

чальные приближения к волновым функциям и

шей полной энергией, были получены с использо-

электронной плотности задавались двумя спо-

ванием волновых функций и плотности, вычислен-

собами. В первом случае задавались начальные

ных ранее для геометрии с более низкой симметри-

магнитные моменты катионов и анионов, а на-

ей. В работе [7] для магнетита приводятся разницы

чальные приближения к волновым функциям и

в полных энергиях для состояний с ферро- и ферри-

электронной плотности строились в соответствии с

магнитным упорядочением. Такие данные позволя-

ними. Для этого случая рассматривались различ-

ют предположить, что задание произвольно выбран-

ные схемы начального магнитного упорядочения.

ных начальных приближений к волновым функци-

Во втором случае в качестве начальных приближе-

ям и электронной плотности вовсе не гарантирует

ний к волновым функциям и плотности при запуске

нахождения решения с наименьшей полной энерги-

нового расчета с плавным изменением одного из

ей, но может приводить к локальному минимуму,

параметров расчета брались волновые функции

причем скорее всего к тому, ближе к которому они

и плотность решений, полученных ранее. Допол-

оказались. Иначе не совсем ясно, как авторам [7]

нительных, определяемых симметрией решетки

удалось найти решение с ферромагнитным упоря-

кубической фазы, требований симметрии на волно-

дочением, не характерным для магнетита.

вые функции и плотность в расчетах дефектной и

Проведенный литературный обзор показал, что

бездефектной конфигураций не налагалось.

использование метода DFT+U для вычисления

Учет сильных электронных корреляций.

свойств кубической фазы магнетита и ее дефектов

Учет сильных электронных корреляций осуществ-

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

лялся по схеме, предложенной в работе [48]. Боль-

получение отрицательных энергий образования

шая часть представленных результатов была полу-

дефектов связано с тем, что для бездефектной

чена с использованием U_eff = 3.5 эВ. Дополнитель-

кубической фазы в рамках DFT+U еще не было до-

но обсуждаются результаты, полученные при U_eff =

стигнуто состояние с наименьшей полной энергией.

= 0 и U_eff = 5.0 эВ. Выбор значения U_eff для рас-

При этом, поскольку расчеты проводились без на-

чета свойств дефектов магнетита проводился в мо-

ложения дополнительных требований симметрии на

дели кубической фазы, в предположении, что в ней

волновые функции и плотность, причина того, что

существует запрещенная зона, а порядок ее шири-

основное состояние не было достигнуто, кроется,

ны близок к наблюдаемому экспериментально выше

по-видимому, в существовании множества локаль-

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

ных минимумов системы и трудности нахождения

в магнетите (табл. 2).

основного состояния среди них.

Параметры сходимости и оптимизация гео-

В данной работе в рамках DFT+U исследуется

метрии. Свойства бездефектного кристалла вычис-

чувствительность свойств кубической фазы к пара-

лялись с использованием кубической сверхячейки с

254

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

56 атомами, соответствующей элементарной ячей-

го ранее, но внешнее давление в ОЦК-сверхячейке

ке кубической фазы магнетита (рис. 1), а также в

железа составляет 100 кбар. В работе рассматрива-

сверхячейке с 14 атомами, которая была использо-

лись два типа вакансий: вакансии A- и B-подреше-

вана ранее для расчета свойств кубической фазы в

ток магнетита.

работе [5]. Свойства дефектов вычислялись в куби-

ческих сверхячейках, содержащих 56 атомов и соот-

ветствующий дефект. Базисный набор плоских волн

3. РЕЗУЛЬТАТЫ

во всех случаях был ограничен энергией 550 эВ. Для

суммирования по зоне Бриллюэна использовались

3.1. Возможность получения отрицательных

Γ-центрированные сетки с 8×8×8 и 6 × 6 × 6 k-точ-

энергий образования

ками для сверхячеек с 14 и 56 атомами соответствен-

но. Порог сходимости при поиске самосогласованно-

Перед изложением основных результатов пред-

го решения был задан равным 10-6 эВ, что обеспе-

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

чило сходимость полной энергии и внешнего давле-

получению отрицательных энергий образования де-

ния лучше 1.0 мэВ и 0.05 кбар соответственно. По-

фектов. Расчеты проводились следующим образом.

рог оптимизации геометрии во всех случаях опре-

В кубической сверхячейке с 56 атомами с фик-

делялся по силам на атомы и составлял 10-2 эВ/Å.

сированным значением постоянной решетки (a =

Оптимизация геометрии дефектных конфигураций

= 8.490Å) и фиксированным значением кислород-

проводилась при фиксированном равновесном зна-

ной координаты (x = 0.2549) создавалась вакан-

чении постоянной решетки, как это было сделано

сия, после чего все координаты атомов (и только

ранее для дефектов магнетита в работах [10, 29].

они) оптимизировались для вычисления значения

Энергии образования вакансий. Энергии об-

E_V в формуле (1). В бездефектной сверхячейке так-

разования вакансий вычислялись по аналогии с тем,

же проводилась оптимизация координат атомов (и

как это было сделано для дефектов магнетита ранее

только их) для получения значения E₀ в формуле

в работе без поправки Хаббарда [29]; при этом была

(1). Вычисленные значения E_V , E₀, E_Fe и E^FV для

использована формула, приведенная в работе [24],

вакансий A- и B-подрешеток представлены в табл. 3.

В этих расчетах начальные приближения к вол-

E^FV = E_V - E₀ + E_Fe,

(1)

новым функциям и электронной плотности задава-

где E^FV - энергия образования вакансии, E_V , E₀ —

лись с помощью начальных магнитных моментов

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

катионов и анионов. Начальные магнитные момен-

фекта соответственно, E_Fe — химический потенциал

ты катионов A- и B-подрешеток были заданы рав-

железа в магнетите. Как и в работе [29], в данной ра-

ными -3.8 μ_B и +3.8 μ_B соответственно, а началь-

боте за химический потенциал железа в магнетите

ные магнитные моменты анионов кислорода были

(E_Fe) была принята энергия, приходящаяся на атом

заданы равными нулю. Такие значения близки к су-

в ферромагнитном ОЦК-железе. Однако в настоя-

ществующим литературным данным (например, ре-

щей работе вычисления свойств дефектов проводи-

зультаты (DFT+U)-расчетов [9,10,16,18] в табл. 2),

лись с учетом полуэмпирической поправки Хаббар-

что, как предполагалось, могло бы способствовать

да. Поэтому для того чтобы можно было воспользо-

скорейшему нахождению минимума. В дефектных

ваться формулой (1), химический потенциал железа

сверхячейках задавались те же начальные значения

вычислялся также с учетом поправки. При U_eff =

магнитных моментов, что и в бездефектном слу-

= 3.5 эВ полученные значения оптимизированной

чае, за исключением магнитного момента удаленно-

постоянной решетки и магнитного момента железа

го катиона. Все вычисления проводились при U_eff =

равны 2.936Å и 2.87 μ_B соответственно, что несколь-

= 3.5 эВ.

ко выше по сравнению с данными без поправки,

2.830Å и 2.19 μ_B [49], и экспериментальными данны-

Таблица 3. Возможность получения отрицательных

ми, 2.870Å и 2.22 μ_B [49]. Полученное в этих усло-

энергий образования дефектов магнетита

виях значение E_Fe = -5.57 эВ выше, чем в расче-

тах без учета поправки (-8.24 эВ). Также был про-

defect

E_V

E₀

E_Fe

E^FV

веден расчет с использованием экспериментальных

V_A

-380.532

-386.513

-5.570

0.41

значений постоянной решетки и магнитного момен-

та ОЦК железа в рамках DFT+U. В этом случае

V_B

-381.976

-1.03

химический потенциал на 0.03 эВ выше полученно-

255

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Как видно из табл. 3, расчет с такими прибли-

жениями приводит к тому, что энергия образования

B-вакансии оказывается отрицательной. Хотя раз-

ница между энергиями образования вакансий A- и

B-подрешеток составляет Δ_AB = 1.44 эВ, что в це-

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

в других работах (табл. 1). Разница между пол-

ным магнитным моментом дефектной и бездефект-

ной сверхячеек составляет +8 μ_B и -2 μ_B для A- и

B-вакансии соответственно, что также согласуется с

литературными данными [10, 29].

Анализ вычисленных магнитных моментов ка-

тионов показал, что в дефектных сверхячейках су-

ществуют две группы катионов B-подрешетки, в то

время как в бездефектной сверхячейке все катионы

октаэдрической подрешетки имеют одинаковые зна-

чения магнитных моментов. В плотностях электрон-

ных состояний, полученных для случаев с A- и B-ва-

кансиями, существуют запрещенные зоны, однако в

бездефектном случае запрещенной зоны нет. Такие

расчетные данные и результаты работы [5] заста-

вили в первую очередь пересмотреть приближения,

закладываемые в модель бездефектной кубической

фазы магнетита.

3.2. Свойства бездефектной кубической

фазы

3.2.1. Множество локальных минимумов

системы

Расчеты в сверхячейке с

атомами.

Сверхячейки с 14 атомами использовались ранее

Рис. 2. Полная энергия E и внешнее давление p в зависи-

мости от объема сверхячейки V с 14 атомами. Точки —

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

вычисленные значения, кривые — аппроксимация урав-

в рамках DFT+U в работах [5, 18]. Они — мень-

√

нением состояния Винета. В легенде p(V ) указаны вы-

шего размера (a₀/

2, a₀/

2), чем ячейка

численные значения равновесной постоянной решетки a0:

кубической фазы с 56 атомами, и некубические

p(a0) = 0 для соответствующих случаев

(рис. 13). Уменьшение числа атомов в расчетной

сверхячейке существенно снижает ресурсоемкость

вычислений, что позволяет рассмотреть большее

ной координаты (x = 0.2549). После этого вычислен-

число интересных случаев. В то же время отличия

ные зависимости полной энергии E и внешнего дав-

расчетной ячейки от ячейки кубической фазы с

ления p от объема сверхячейки V аппроксимирова-

атомами могут стать источником расчетных

лись уравнением состояния Винета [50] (рис. 2).

артефактов при ее использовании в периодических

граничных условиях.

Начальные приближения к волновым функциям

Вычисления проводились при U_eff = 3.5 эВ. По-

и электронной плотности задавались с помощью на-

иск равновесной постоянной решетки (a₀: p(a₀) = 0)

чальных магнитных моментов анионов и катионов.

осуществлялся следующим образом. Для значений

Начальные магнитные моменты анионов во всех

постоянной решетки a из интервала 8.460÷8.510Å с

случаях задавались равными нулю. Для катионов

шагом 0.005Å проводились самосогласованные рас-

рассматривалось два варианта задания начальных

четы фиксированных ионных конфигураций с дан-

магнитных моментов, соответствующих (1) ферро-

ным значением постоянной решетки a и фиксиро-

магнитному (обозн. ferrO) и (2) ферримагнитному

ванным экспериментальным значением кислород-

упорядочению (обозн. ferrI). Ферромагнитное на-

256

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

Таблица 4. Свойства решений, полученных для слу-

чальное приближение дает решение также обладаю-

чая с постоянной решетки a = 8.475Å: разница в

щее ферромагнитным упорядочением, не характер-

значениях полной энергии данного решения и реше-

ным для магнетита. Из рис. 2 видно, что для всех

ния с наименьшей полной энергией среди найден-

рассмотренных значений постоянной решетки это

ных ΔE (мэВ/f.u.) и магнитные моменты катионов

ферромагнитное решение выше по энергии, чем ре-

шения, полученные с ферримагнитным начальным

μ (μB)

приближением. Это говорит о том, что с использо-

ванием ферромагнитного начального приближения

Начальные

ΔE

найти основное состояние системы невозможно.

приближения

Fe3+A Fe2+B Fe3+B

Для случая (2), с ферримагнитным начальным

ferrO, 555

+646

+4.19

+4.12

приближением (ferrI), дающим решения с более глу-

бокими минимумами, рассматривалась зависимость

ferrI, -888

+278

-4.00

+3.94

результатов от конкретных значений начальных

ferrI, -555

+235

-4.05

+3.93

+3.95

магнитных моментов. В обозначениях легенды E(V )

на рис. 2 три цифры рядом с обозначением заданно-

wfu5, -878

+30.4

-4.02

+3.72

+4.10

го начального магнитного упорядочения (ferrO или

ferrI, -878

+29.7

-4.04

+3.85

+4.00

ferrI) — это значения начальных магнитных момен-

тов катионов μ(Fe3+A), μ(Fe2+B), μ(Fe3+B) в магнето-

ferrI, -656

0.0

-4.03

+3.73

+4.09

нах Бора (μ_B )³⁾. Для сравнения на рис. 2 пред-

ставлены также результаты расчетов, выполненных

при U_eff = 3.5 эВ с использованием в качестве на-

ницей B-катионов (wfu5, -878, ferrI, -878 и ferrI,

чальных приближений решений, найденных при той

-656). Последние, обладающие более низкими пол-

же постоянной решетки при U_eff = 5.0 эВ (обозн.

ными энергиями, отличаются между собой не более

wfu5, -878). Из рис. 2 отчетливо видно, что различ-

чем на 0.03 эВ/f.u., что на порядок меньше указан-

ные начальные приближения к волновым функциям

ных выше отличий. Из рис. 2 видно, что только одно

и плотности позволяют обнаружить только некото-

из этих решений (ferrI, -878) близко к p = 0. Реше-

рые локальные минимумы системы. При этом от ис-

ние, обладающее наименьшей полной энергией при

пользуемых начальных приближений зависит также

a = 8.475Å, полученное с начальным приближением

вычисленное значение равновесной постоянной ре-

ferrI, -656, характеризуется разницей в магнитных

шетки.

моментах двух- и трехвалентных катионов, состав-

Характер полученных локальных минимумов

ляющей 0.36 μ_B. Почти такой же разницей характе-

рассмотрен подробнее для случая с постоянной ре-

ризуется и решение, полученное с wfu5, -878, с чуть

шетки a = 8.475Å. Свойства решений, полученных с

более высокой полной энергией. Отличия этих реше-

разными начальными приближениями при этом зна-

ний видны на рис. 3, где показаны плотности элек-

чении постоянной решетки, представлены в табл. 4.

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

Значения вычисленных магнитных моментов кати-

Действительно, случай ferrI, -656 дает запрещен-

онов в табл. 4 показывают, что решения, найден-

ную зону и для состояний со спином вверх, и для

ные с использованием ферримагнитного начально-

состояний со спином вниз, а случай wfu5, -878 —

го приближения без отличий B-катионов (случаи

только для состояний со спином вверх.

ferrI, -888 и ferrI, -555), приводят к решению, так-

же обладающему ферримагнитным упорядочением

Из рис. 2 видно, что не все решения, полученные

без отличий B-катионов (магнитные моменты ка-

при одном и том же начальном приближении, мож-

тионов Fe2+B и Fe3+B практически одинаковы). При

но аппроксимировать уравнением состояния (см.,

этом из табл. 4 видно, что ферромагнитное (ferrO,

например, случай ferrI, -555: при различных посто-

555) и ферримагнитные без отличий B-катионов ре-

янных решетки это начальное приближение приво-

шения на 0.65 эВ/f.u. и 0.25 эВ/f.u. соответственно

дит к решениям, не укладывающимся на одну кри-

выше по энергии, чем решения, полученные с фер-

вую в рамках общей для всех точности). Так проис-

римагнитными начальными приближениями с раз-

ходит и в случае ferrI, -656. На рис. 4 представлены

свойства решений, полученных с ferrI, -656 при раз-

ных значениях постоянной решетки (это увеличен-

³⁾ Так, запись ferrI,-878 обозначает, что при запуске рас-

чета было задано ферримагнитное начальное приближение с

ное изображение зависимости E(V ) с ferrI, -656 c

μ(Fe3+A) = -8.0 μ_B , μ(Fe2+B) = 7.0 μ_B , μ(Fe3+B) = 8.0 μ_B .

рис. 2). Видно, что одно и то же начальное прибли-

257

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Рис. 4. Свойства решений, полученных с начальным при-

ближением ferrI, -656: зависимость E(V ) — в верхней ча-

Рис. 3. Плотность электронных состояний со спином вверх

сти рисунка, плотность электронных состояний для точек,

(в положительном направлении оси ординат) и со спи-

укладывающихся на верхнюю и нижнюю кривые соответ-

ном вниз (с обратным знаком в отрицательном направ-

ственно — внизу. Плотность электронных состояний для

лении оси ординат) для случаев, представленных в табл. 4

«выпавшей» точки (это случай ferrI, -656 с a = 8.475Å)

(a = 8.475Å). Красной вертикальной линией отмечен уро-

показана на рис. 3 внизу

вень Ферми

ются полуметаллическим состоянием (верхняя кри-

жение ferrI, -656 может давать решения, которые

вая), решения с меньшей полной энергией демон-

распадаются на три случая. Два из них можно ап-

стрируют псевдощелевое поведение (нижнаяя кри-

проксимировать двумя различными кривыми E(V )

вая), решение для выпавшей точки отличается за-

с разницей в 17 мэВ/f.u. между значениями полной

прещенной зоной (это тот самый случай ferrI, -656

энергии в минимумах (это верхняя и нижняя кри-

при a = 8.475Å в табл. 4 и самый нижный график

вые на рис. 4), а третий случай — «выпавшая» точ-

на рис. 3).

ка с a = 8.475Å. То, что решения, укладывающие на

каждую из кривых, относятся к одному типу, вид-

На рис. 5 представлены числа заполнения лока-

но по плавно меняющейся с изменением постоянной

лизованных орбиталей d-состояний для каждого из

решетки плотности электронных состояний в рам-

шести катионов сверхячейки с 14 атомами для слу-

ках каждой серии (рис. 4, внизу). При этом реше-

чаев с a = 8.470Å, a = 8.480Å и a = 8.475Å с

ния с чуть большей полной энергией характеризу-

рис. 4. Видно, что все эти решения действительно

258

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

Рис. 5. Степень заселенности d-орбиталей со спином вверх

(в положительном направлении оси ординат, выделена зе-

леным цветом) и со спином вниз (с обратным знаком

в отрицательном направлении оси ординат, выделена си-

ним цветом) для катионов A- и B-подрешеток с номерами

N = 1 - 2 и N = 3 - 6 соответственно в сверхячейке с 14

Рис. 6. Степень заселенности d-орбиталей для обсужда-

атомами для случаев рис. 4

емых случаев. Обозначения тождественны обозначениям

рис. 5

обладают ферримагнитным упорядочением: у кати-

онов A-подрешетки (с номерами N = 1, 2) заполне-

мощью метода работы [26] с шагом δU_eff = 0.25 эВ

ны все состояния со спином вниз, а у катионов B-

при U_eff = 3.5 эВ. Это решение на 0.25 эВ/f.u. выше

подрешетки (с номерами N = 3-6) — все состояния

решения, полученного при начальном приближении

со спином вверх. Для случая с a = 8.470Å c верх-

ferrI, -878 (третья картинка на рис. 6), что совпада-

ней кривой на рис. 4 числа заполнения для одной

ет с разницей между полными энергиями решений с

из локализованных орбиталей у каждого B-катиона

и без отличий B-катионов, представленными ранее.

составляют 0.5. Для случая с a = 8.480Å c ниж-

Из рис. 6 видно, что вычисления без учета сильных

ней кривой наблюдается состояние с заселенностью

электронных корреляций не дают локализованных

более 0.8 у одного катиона B-подрешетки (N = 4).

состояний у B-катионов со спином вниз. Их исполь-

Наконец, в случае с выпавшей точкой a = 8.475Å,

зование в качестве волновых функций и плотности

имеющей наименьшую полную энергию, локализо-

при Ueff = 0.25 эВ, и далее использование реше-

ванные состояния наблюдаются у двух B-катионов

ний от U_eff = 0.25 эВ при U_eff = 0.50 эВ, и далее,

(это половина общего числа B-катионов), однако их

вплоть до U_eff = 3.5 эВ, дает решения также без

числа заполнения немного отличаются (на 0.1).

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

В работе [26] был предложен метод поиска реше-

ми начальными приближениями при U_eff = 3.5 эВ

ний с наименьшей полной энергией в DFT+U. Его

можно найти решения с более низкими значениями

алгоритм заключается в последовательном исполь-

полной энергии. Эти решения отличаются локали-

зовании решений, полученных сначала без поправ-

зованными состояниями со спином вниз у половины

ки, а затем с увеличением U_eff в качестве началь-

B-катионов. Однако числа заполнения этих локали-

ных приближений в расчетах с плавным увеличени-

зованных орбиталей не равны в точности единице,

ем U_eff . В работе была сделана попытка его при-

а несколько ниже (0.65). Интересно отметить, что

менения в сверхячейке с 14 атомами с a = 8.490Å

при тех же условиях при U_eff = 5.0 эВ получаются

и x = 0.2549. На первой картинке рис. 6 показаны

близкие к единице числа заполнения (рис. 6, самая

числа заполнения для случая с U_eff = 0, а на вто-

нижняя картинка). Поэтому, вероятно, их использо-

рой картинке сверху — решение, полученное с по-

вание в качестве начальных приближений для рас-

259

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Таблица 5. Вычисленные значения полной энергии

чета с U_eff = 3.5 эВ дает решения из области с более

E (эВ/f.u.), внешнего давления p (кбар) и разни-

низкими полными энергиями (wfu5, -878 на рис. 2).

цы между магнитными моментами Fe3+B и Fe2+B,

Однако, как будет обсуждаться далее, ширина за-

прещенной зоны, полученная при U_eff = 5.0 эВ, на

δμ (μB), для случая c a = 8.490Å, x = 0.2549 при

U_eff = 3.5 эВ

порядок выше, чем при U_eff = 3.5 эВ (в то время

как при U_eff = 3.5 эВ ширина запрещенной зоны

ближе к ее экспериментальному порядку).

δμ

Подводя промежуточный итог, можно отметить,

Без отличий B-катионов

что рассматриваемая система обладает множеством

-555 (isym2)

-48.313

-20.7

0.05

локальных минимумов, найти некоторые из них поз-

воляет задание различных начальных приближений

−888 (isym2)

-48.313

-20.7

0.05

к волновым функциям и плотности. При этом ре-

−555

-48.597

-7.9

0.37

шениями с наименьшими полными энергиями явля-

ются решения с отличиями октаэдрических катио-

−888

не сошлась

нов. При вычислении энергий образования дефек-

−555 (14)

-48.561

-7.8

0.30

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

0.02-0.03 эВ/f.u. для величины E₀ (и, возможно, для

−888 (14)

-48.313

-20.0

0.00

величины E_V ) между решениями, имеющими отли-

С отличиями B-катионов

чия октаэдрических катионов, по-видимому, вносит

Начальная схема Вервея (mv)

вклад в неточность, с которой эти энергии образова-

ния могут быть известны (порядка 0.16-0.24 эВ при

-545

-48.571

-7.6

0.36

исследовании дефектов в сверхячейках с 8 f.u. (56

−878

-48.564

-10.5

0.17

атомами)). Она кажется не столь существенной по

сравнению с разницей в 0.25 эВ/f.u. между решени-

−545 (14)

-48.549

-9.7

0.15

ями с и без отличий B-катионов, дающей 2.0 эВ. По-

−878 (14)

-48.564

-10.2

0.15

следняя может быть основной причиной появления

Случайная начальная схема (m2)

отрицательных энергий образования дефектов.

Расчеты в сверхячейке с 56 атомами. Далее

-545

-48.605

-6.6

0.44

в статье, кроме особо оговоренных случаев, будут

−878

-48.605

-7.2

0.46

представлены результаты, полученные в сверхячей-

ках с 56 атомами. В табл. 5 представлены результа-

ты расчетов в сверхячейке с 56 атомами при фикси-

рованной ионной конфигурации с a = 8.490Å и x =

волновые функции и электронную плотность были

= 0.2549 и значении U_eff = 3.5 эВ. Рядом для срав-

наложены дополнительные требования симметрии,

нения представлены результаты, полученные ранее

определяемые симметрией решетки. Как ранее от-

при аналогичных условиях в сверхячейке с 14 ато-

мечалось в работе [5], в таком случае решения име-

мами (обозн. (14)).

ют более высокую полную энергию, чем решения, не

При запуске расчетов в сверхячейке с 56 атома-

обладающие симметрией. Разница в полных энерги-

ми, как и в случае с 14 атомами, задавались различ-

ях решений, полученных в случаях -555 (и обла-

ные начальные приближения. Их задание осуществ-

дающих отличиями B-катионов) и -555 (isym2) (и

лялось с помощью начальных магнитных моментов

не обладающих такими отличиями), в данной рабо-

катионов. При этом ферромагнитная начальная схе-

те составляет 0.28 эВ/f.u., что сравнимо с разницей,

ма уже не рассматривалась, только ферримагнит-

полученной при вычислениях в сверхячейке с 14 ато-

ная. Обозначения аналогичны обозначениям в слу-

мами (0.25 эВ/f.u.). По представленным на рис. 7

чае с 14 атомами. Для начального приближения с

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

ферримагнитным упорядочением рассматривались

но, что для решений, обладающих симметрией (см.

случаи с отличиями и без отличий B-катионов. В

-888 (isym2)), локализованных состояний со спином

случае с 56 атомами начальные приближения без от-

вниз у B-катионов не наблюдается. Это выражается

личий B-катионов могут приводить (-555), а могут

также в том, что разница вычисленных магнитных

и не приводить (-555 (isym2), -888 (isym2)) к ре-

моментов катионов в табл. 5 близка к нулю. Видно,

шениям с отличиями B-катионов. Это объясняется

что в случае с 14 атомами такой локальный мини-

тем, что в случаях (-555 (isym2), -888 (isym2)) на

мум также был найден (-888 (14)).

260

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Проведенный ранее анализ решения, обладаю-

и (2) случайная начальная схема зарядового упоря-

щего наименьшей полной энергией среди найденных

дочения (обозн. m2), показанная на рис. 1.

в сверхячейке с 14 атомами, показал, что оно ха-

Из табл. 5 видно, что решения с наиболее низ-

рактеризуются локализованными состояниями у по-

кими полными энергиями были получены с исполь-

ловины B-катионов. Показалось естественным пред-

зованием начальной случайной схемы m2. Если ана-

положить, что этим будут отличаться и решения с

лизировать представленные для случаев табл. 5 чис-

наименьшими полными энергиями в сверхячейке с

ла заполнения локализованных орбиталей (рис. 7),

56 атомами. В таком случае в качестве начальных

можно заметить, что решения с наименьшей пол-

приближений естественно было бы задать прибли-

ной энергией (-545 (m2) и -878 (m2)) отличают-

жения, учитывающие разницу B-катионов (полови-

ся от других решений следующими особенностями:

ну B-катионов с большим начальным моментом, по-

1) у них наблюдаются локализованные d-состояния

ловину — с меньшим). Однако при задании началь-

со спином вниз у B-катионов (в отличие от случая

ных приближений с отличиями B-катионов в сверхя-

-888 (isym2), в котором их нет), 2) числа их за-

чейке с 56 атомами возникают существенные труд-

полнения стремятся к наибольшим возможным зна-

ности. Дело в том, что в сверхячейке с 14 атомами

чениям (в отличие от случая -878 (mv)), 3) кати-

четыре B-катиона. Их можно распределить по двум

онов с такими орбиталями — половина от общего

группам с меньшим и большим начальными момен-

числа B-узлов (в отличие от случая с -545 (mv)),

тами шестью способами. При этом все шесть спо-

4) и в рамках этой половины числа заполнения име-

собов приводят к шести решениям с неотличимой в

ют практически одинаковые значения (в отличие от

рамках погрешности полной энергии. Решения отли-

случая -555). Из табл. 5 видно также, что решения

чаются только схемами распределения двух- и трех-

-545 (m2) и -878 (m2) не отличаются полной энер-

валентных катионов по узлам B-подрешетки. Одна-

гией (в рамках погрешности). Однако на рис. 7 их

ко в сверхячейке с 14 атомами все эти схемы, при

отличия видны. Дело в том, что в случае -545 (m2)

условии нахождения минимума с четким разделе-

решение характеризуется схемой зарядового упоря-

нием на две группы, эквивалентны и представляют

дочения Вервея (рис. 11), а в случае -878 (m2) —

схему, известную как схема Вервея.

нет. По-видимому, для решения со схемой Вервея

Понятие «схема Вервея» возникает при обсужде-

характерны одинаковые значения чисел заполнения

нии зарядового упорядочения в магнетите, т. е. рас-

локализованных орбиталей у половины B-катионов

пределения двух- и трехвалентных катионов по уз-

(в отличие от решения с -878 (m2), где среди них

лам B-подрешетки (или локализации заряда). В схе-

есть небольшие отличия). Интересно отметить, что

ме Вервея, рассматриваемой, например, в элемен-

к схеме Вервея не привели и случаи с начальной

тарной ячейке кубической фазы магнетита, плоскос-

схемой Вервея (-545 (mv) и -878 (mv)). При этом

ти, перпендикулярные оси [001] и содержащие це-

характерные отличия в полных энергиях между схе-

почки катионов Fe2+B, чередуются вдоль этой оси с

мами с отличиями B-катионов составляют порядка

плоскостями, также перпендикулярными [001] и со-

0.03 эВ/f.u., что согласуется с данными, полученны-

держащими цепочки катионов Fe3+B (см. решение со

ми в сверхячейках с 14 атомами. Ясно, однако, что

схемой Вервея на рис. 11).

решений типа тех, которые в сверхячейке с 56 ато-

Трудность в сверхячейке с 56 атомами заключа-

мами дали случаи -545 (mv) и -878 (mv), в сверхя-

ется в том, что в ней не четыре, а шестнадцать узлов

чейке с 14 атомами из-за ее геометрии получить

B-подрешетки. Распределить по ним двух- и трехва-

нельзя. Для случая с 56 атомами исследовалась так-

лентные катионы можно большим числом способов,

же начальная схема -545 (m1), которая отличалась

причем далеко не все из них эквивалентны схеме

от схемы -545 (m2) перестановкой двух началь-

Вервея. Такие начальные приближения, отличаю-

ных магнитных моментов катионов. Полученное в

щиеся не только начальным магнитным упорядоче-

этом случае решение на 4 мэВ выше, чем в случае

нием и значением начальных магнитных моментов

-545 (m2) и не характеризуется ни зарядовым, ни

катионов, а еще и схемой начального распределения

орбитальным упорядочением.

катионов разной валентности по B-узлам, могут да-

Результаты оптимизации геометрии сверхячей-

вать новые локальные минимумы в сверхячейках с

ки с сохранением симметрии кубической фазы для

56 атомами. Поскольку вычисления в сверхячейках

отдельных случаев представлены в табл. 6. При

с 56 атомами более ресурсоемкие, в данной работе

этой процедуре варьироваться мог только объем

было рассмотрено два случая: (1) начальное зарядо-

сверхячейки. В качестве начальной ионной конфи-

вое упорядочение в виде схемы Вервея (обозн. mv)

гурации использовалась описанная выше, с a

262

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

Таблица 6. Вычисленные значения полной энер-

гии E (эВ/f.u.), равновесной постоянной решетки

a0 (^Å) и разницы между магнитными моментами

Fe3+B и Fe2+B, δμ (μB ), для равновесного случая c

x = 0.2549

a₀

δμ

Без отличий B-катионов

-888 (isym2)

-48.319

8.456

0.01

С отличиями B-катионов

Случайная начальная схема (m2)

-545

-48.606

8.479

0.43

−878

-48.606

8.478

0.44

= 8.490Å, x = 0.2549. Из сравнения данных в табл. 6

и табл. 5 видно, что свойства равновесных реше-

ний очень близки к тем, которые были получе-

Рис. 8. Полная и спроектированные на d-состояния отдель-

ны ранее для неравновесной конфигурации с a =

ных катионов A- и B-подрешетки плотности электронных

= 8.490Å, x = 0.2549 при тех же начальных при-

состояний для равновесного случая при U_eff = 0 в сверхя-

ближениях. Полученные при начальных приближе-

чейке с 56 атомами

ниях -545 (m2) и -878 (m2) решения не отличают-

ся в рамках погрешности, а приближение без отли-

Таблица 7. Вычисленные значения магнитных мо-

чий B-катионов ferrI, -888 дает другую оптимизи-

ментов катионов μ (μB ) и равновесной постоянной

рованную постоянную решетки и более высокую (на

решетки a0^Å для случаев с разными Ueff (эВ)

0.29 эВ/f.u.) полную энергию. Равновесная постоян-

ная решетки в случае -545 (m1) составляет 8.477Å;

полная энергия для этого случая -48.598 эВ, что

U_eff

a₀

на 7 мэВ выше, чем для равновесного случая -545

Fe3+A

Fe2+B

Fe3+B

(m2), который дал схему Вервея.

0.0

-3.48

3.56

8.393

3.5

-4.02

3.68

4.10

8.479

3.2.2. Подбор Ueff

5.0

-4.16

3.72

4.24

8.487

U_eff = 0. В этом предельном случае, без учета

сильных электронных корреляций, при x = 0.2549

вычисленное значение равновесной постоянной ре-

шетки составляет 8.393Å, что в точности соответ-

ком окружении — t2g, что согласуется с теоретичес-

ствует экспериментальному значению для кубичес-

кими предсказаниями [51].

кой фазы магнетита (табл. 2). При этом ферримаг-

U_eff = 5.0 эВ. На рис. 9 показан другой предель-

нитные начальные схемы с и без отличий B-катио-

ный случай, с достаточно большим значением по-

нов дают одно и то же ферримагнитное решение без

правки, U_eff = 5.0 эВ. Видна предсказываемая мо-

отличий B-катионов. Начальное приближение, не

делью запрещенная зона, ширина которой, однако,

включающее отличия B-катионов, сходится к этому

почти на порядок превосходит наблюдаемую в экс-

решению быстрее. Модель с U_eff = 0 не предсказы-

периментах с кубической фазой ширину псевдоще-

вает запрещенной зоны для кубической фазы маг-

ли. При U_eff = 5.0 эВ в сверхячейке с 14 атомами

нетита. Уровень Ферми (рис. 8) пересекает состоя-

чувствительности к конкретным значениям началь-

ния B-катионов со спином вниз. На рис. 8 отчетливо

ных магнитных моментов для ферримагнитной на-

видно расщепление d-состояний в кристаллическом

чальной схемы с отличиями B-катионов не наблю-

поле на t2g (d_xy, d_xz, d_yz) и e_g (d_z2_-r2 , d_x2_-y2 ) состоя-

дается, а схема без отличий B-катионов дает бо-

ния, причем в тетраэдрическом окружении ниже по

лее высокую полную энергию. В отличие от слу-

энергии оказываются e_g состояния, а в октаэдричес-

чая с U_eff = 0, где отличий B-катионов не наблю-

263

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Рис. 9. Полная и спроектированные на d-состояния отдель-

Рис. 10. Полная и спроектированные на d-состояния от-

ных катионов A- и B-подрешетки плотности электронных

дельных катионов A- и B-подрешетки плотности электрон-

состояний для равновесного случая при U_eff = 5.0 эВ в

ных состояний для равновесного случая при U_eff = 3.5 эВ

сверхячейке с 56 атомами

в сверхячейке с 56 атомами

дается, при U_eff = 5.0 эВ у половины B-катионов

ями Fe2+B. Распределение спиновой плотности для

наблюдаются локализованные состояния со спином

состояний из интервала энергий, соответствующих

вниз. На плотности состояний B-катионов с мень-

этому пику, показано на рис. 11. Из рис. 11 видно,

шим магнитным моментом (и потому интерпрети-

что это решение характеризуется также орбиталь-

рованных и обозначенных как Fe2+B) вблизи уров-

ным упорядочением. При этом характер располо-

ня Ферми можно заметить характерный ассоцииро-

жения орбиталей кажется соответствующим иска-

ванный с ними пик (на рис. 9 показан случай для

жению решетки (x = 0.2549). Интересно отметить,

одного из восьми Fe2+B), которого нет у Fe3+B. Раз-

что при x = 0.2500 для схемы Вервея наблюдается

личия в магнитных моментах катионов Fe2+B и Fe3+B

орбитальное упорядочение без искажений, и орби-

при U_eff = 5.0 эВ выше, чем при U_eff = 3.5 эВ и со-

тали Fe2+B выстраиваются в идеальные цепочки; од-

ставляют 0.52 μ_B (табл. 7). Однако при U_eff = 5.0 эВ

нако, как уже упоминалось, случай с x = 0.2500 —

случайно заданная начальная схема -545 (m2), дав-

выше по энергии. В случае сверхячейки с 14 атома-

шая при U_eff = 3.5 эВ решение со схемой Вервея,

ми орбитали также выстраиваются в идеальные це-

схему Вервея не дала. Вычисленные значения рав-

почки (рис. 13); но в этом случае такой эффект, по-

новесной постоянной решетки и магнитных момен-

видимому, может быть обусловлен размерным эф-

тов катионов при U_eff = 5.0 эВ оказываются выше,

фектом сверхячейки, включающим, возмножно, ар-

чем при U_eff = 0 и U_eff = 3.5 эВ (табл. 7).

тефактное взаимодействие атомов сверхячейки с их

U_eff = 3.5 эВ. В равновесном случае при U_eff =

изображениями в ее образах. Парциальная плот-

= 3.5 эВ, полученном c -545 (m2), ширина за-

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

прещенной зоны на порядок меньше, чем в случае

случая -878 (m2) показана на рис. 12; можно заме-

Ueff = 5.0 эВ (рис. 10). В плотности состояний со

тить, что двух- и трехвалентные катионы в цепочках

спином вниз, как и в случае с U_eff = 5.0 эВ, ниже

чередуются, образуя отличный от -545 (m2) харак-

уровня Ферми есть пик, связанный с t2g состояни-

тер упорядочений.

264

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

Рис. 11. Зарядовое и орбитальное упорядочения для рав-

Рис. 12. Парциальная плотность локализованных состоя-

новесного случая -545 (m2) при U_eff = 3.5 эВ в сверхя-

ний для равновесного случая -878 (m2) при U_eff = 3.5 эВ

чейке с 56 атомами. На нижней картинке показан вид

в сверхячейке с 56 атомами

сверху

U_eff = 3.2 ÷ 3.7 эВ. Для проверки того, насколь-

и ширина запрещенной зоны кубической фазы E_g

ко сильно равновесная постоянная решетки и ши-

увеличиваются (рис. 14). Экспериментальный поря-

рина запрещенной зоны могут изменится при вари-

док наблюдаемой для кубической фазы псевдощели

ации значения U_eff в небольшом интервале, были

(100 мэВ [41], 70 мэВ [39]) может быть получен при

найдены равновесные решения при U_eff от 3.2 эВ

U_eff в диапазоне от 3.55 эВ до 3.65 эВ. Судя по ха-

до 3.7 эВ. Для этого, с использованием волновых

рактеру зависимости E_g(U_eff ), можно также пред-

функций и плотности равновесного решения, полу-

положить, что при U_eff менее 3.0 эВ модель уже

ченного при U_eff = 3.5 эВ c -545 (m2), в качестве на-

не будет предсказывать запрещенную зону для ку-

чальных приближений, проводился расчет сначала

бической фазы магнетита. Стоит отметить, что не

с U_eff = 3.6 эВ, а затем, с использованием получен-

представленные здесь расчеты в сверхячейке с 14

ного при этом U_eff решения, — для U_eff = 3.7 эВ, и

атомами при U_eff менее 3.0 эВ при отсутствии за-

аналогично — в сторону уменьшения U_eff . Все полу-

прещенной зоны дают небольшие (видимо, недоста-

ченные таким образом решения также обладали схе-

точные для открытия запрещенной зоны) отличия

мой Вервея и орбитальным упорядочением. С уве-

B-катионов. В данной работе решено было оставить

личением U_eff равновесная постоянная решетки a₀

значение U_eff = 3.5 эВ, поскольку характер полу-

265

ЖЭТФ, вып. 2 (8)

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

ченных с ним решений не сильно отличается от слу-

чаев 3.55 эВ до 3.65 эВ. Значение ширины запре-

щенной зоны модели при U_eff = 3.5 эВ составляет

60 мэВ.

3.2.3. Релаксация атомов

При получении свойств дефектов кубической

фазы магнетита планировалось учесть эффект ре-

лаксации их локального окружения. В работах

[10, 29] для этого в сверхячейке с дефектом опти-

мизировались координаты атомов при фиксирован-

ной равновесной постоянной решетки. Принимая во

внимание наблюдаемый в магнетите переход Вервея

и то, что модель DFT+U способна описывать низ-

котемпературное состояние магнетита, можно пред-

положить, что при оптимизации координат атомов

в DFT+U даже система без дефектов перейдет в

состояние с более низкой полной энергией, чем по-

лученное равновесное состояние кубической фазы.

Это произойдет потому, что наложенные на ионную

конфигурацию условия симметрии кубической фа-

зы будут сняты.

Действительно, полученные в случае -545 (m2)

Рис. 13. Зарядовое и орбитальное упорядочения для рав-

данные показали, что оптимизация координат ато-

новесного случая при Ueff = 3.5 эВ в сверхячейке с 14

мов в сверхячейке с 56 атомами из равновесного со-

атомами. Ребра расчетной сверхячейки показаны черны-

стояния кубической фазы с фиксированным значе-

ми сплошными линиями. Видна более четкая ориентация

нием постоянной решетки приводит к понижению

орбиталей по сравнению со случаем 56 атомов, что, воз-

полной энергии на 106 эВ/f.u. (0.85 эВ) при повыше-

можно, является проявлением размерного эффекта

нии внешнего давления на 2.3 кбар. При этом проис-

ходит искажение решетки с понижением ее симмет-

рии, а ширина запрещенной зоны увеличивается от

60 мэВ до релаксации до 588 мэВ после (рис. 20). Ха-

рактер зарядового упорядочения в случае -545 (m2)

при этом не изменился (рис. 18), но характер орби-

тального упорядочения стал другим (рис. 15). По

расположению орбиталей (рис. 15) можно предпо-

ложить, что при релаксации атомов они разверну-

лись так, чтобы область их перекрытия в цепоч-

ке Fe2+B была меньше по сравнению с тем, что бы-

ло до релаксации (рис. 11); это могло способство-

вать увеличению ширины запрещенной зоны. После

релаксации атомов разница между зарядами двух-

и трехвалентных B-катионов немного уменьшается

(рис. 19). Релаксация атомов в случаях -878 (m2) и

-545 (m1) дает понижение полной энергии на 112 и

129 мэВ/f.u. соответственно; запрещенная зона при

Рис. 14. Зависимость равновесной постоянной решетки

релаксации увеличивается на порядок; зарядового

a0 и ширины запрещенной зоны кубической фазы Eg от

значения U_eff в диапазоне 3.2-3.7 эВ. Эксперименталь-

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

ный порядок ширины запрещенной зоны достигается при

конфигурациях нет.

U_eff = 3.646 эВ

Симметрия решетки после релаксации атомов не

контролировалась. Не исключено, что для релак-

266

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

сированной конфигурации мог быть найден также

только локальный минимум. Кажется существен-

ным только то, что при релаксации атомов найдена

конфигурация с меньшей полной энергией, чем до

релаксации. Возможно, в задаче вычисления энер-

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

тельный эффект релаксации атомов, связанный не

с релаксацией дефектного окружения, а со стрем-

лением структуры кубической фазы перейти в низ-

комемпературную, — еще один фактор, способству-

ющий получению неадекватных энергий образова-

ния дефектов кубической фазы магнетита. При ис-

следовании дефектов в структуре низкотемператур-

ной фазы такого эффекта, казалось бы, возникнуть

не должно. Однако вычисления со сверхячейками,

отвечающими структуре низкотемпературной фазы,

включающими 224 атома, — довольно ресурсоемки,

и получение основного состояния в них при данном

U_eff с учетом чувствительности результатов к пара-

метрам расчета — непростая задача. Другой путь —

создание дефектной конфигурации после релакса-

ции атомов в структуре кубической фазы. Однако,

как будет показано ниже, в этом случае дефектная

конфигурация может придти к минимуму с более

высокой полной энергией, чем в случае, когда де-

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

мов.

3.3. Свойства вакансий магнетита

3.3.1. Энергии образования

Дефекты создавались в узлах, выделенных си-

ним (V_A) и голубым (V_B ) полиэдрами на рис. 1. Рас-

сматривалось два случая. В первом из них катион

удалялся до, а во втором — после релаксации ато-

мов в бездефектной конфигурации. Затем в том же

предположении о начальных магнитным моментах,

что и для бездефектного массива (-545 (m2)), за

исключением магнитного момента удаленного кати-

она, проводился расчет с и без (обозн. как стати-

ческий расчет) оптимизации координат атомов де-

фектной сверхячейки. В табл. 8 и табл. 9 приведе-

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

и бездефектных конфигураций, а также вычислен-

ные по формуле (1) значения энергий образования

дефектов.

Рис. 15. Плотность электронных состояний и орбитальное

Интересно отметить, что в рамках подхода с

упорядочение после релаксации атомов. При сравнении с

уточненными параметрами расчета отрицательных

рис. 11 (внизу) видна переориентация орбиталей после ре-

энергий образования дефектов нет. Эффект релак-

лаксации и связанное с ним изменение спроектированной

сации локального окружения дефекта оказывается

плотности состояний; характер зарядового упорядочения

важен. Действительно, в статических расчетах, без

соответствует показанному на рис. 11 (вверху)

оптимизации координат атомов дефектной конфи-

267

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Таблица 8. Энергии образования A- и B-вакансий

ных электронных корреляций полученные значе-

магнетита при U_eff = 3.5 эВ. Начальная конфигу-

ния энергий образования A- и B-вакансий состави-

рация без релаксации атомов

ли 2.89 эВ и 1.02 эВ соответственно, что согласует-

ся с данными работы [29], полученными также без

Defect

E_V

E₀

E_Fe

E^FV

поправки Хаббарда (табл. 1). При U_eff = 5.0 эВ вы-

Статический расчет

численные значения энергий образования A- и B-ва-

кансий равны 2.27 эВ и 1.46 эВ соответственно, что

V_A

-378.558

-388.847

-5.570

4.72

очень похоже на данные, полученные в данной ра-

V_B

-379.734

3.54

боте при U_eff = 3.5 эВ (табл. 8).

Расчет с релаксацией атомов

Полученные при U_eff = 3.5 эВ с использовани-

в сверхячейке с дефектом

ем начальных приближений -545 (m2) и -878 (m2)

энергии образования дефектов отличаются не более

V_A

-380.779

-388.847

-5.570

2.50

чем на 0.03 эВ/f.u. и 0.09 эВ/f.u. для расчетов с

V_B

-382.017

1.26

и без оптимизации координат атомов соответствен-

V_B*

-381.836

1.44

но. Если задавать начальные магнитные моменты

катионов такими, какие были получены в резуль-

^∗Был удален Fe3+B.

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

(за исключением магнитного момента удаленного

Таблица 9. Энергия образования B-вакансии магне-

катиона), то большинство задач не сходится. По-

тита при U_eff = 3.5 эВ. Начальная конфигурация с

видимому, такие малые значения ((3.72 - 4.12) μ_B)

релаксацией атомов

начальных магнитных моментов не подходят для

алгоритма расчета.

Defect

E_V

E₀

E_Fe

E^FV

Можно было бы предположить, что при нали-

Статический расчет

чии зарядового упорядочения в модели кубической

фазы могло бы иметь значение то, какой катион

V_B

-380.205

-389.689

-5.570

3.91

был удален при создании вакансии в октаэдриче-

Расчет с релаксацией атомов

ской подрешетке, Fe2+B или Fe3+B. В данной работе

в сверхячейке с дефектом

были сделаны попытки провести расчет дефектных

конфигураций с начальными магнитными момента-

V_B

-381.708

-389.689

-5.570

2.41

ми, соответствующими случаям, когда удален двух-

или трехвалентный B-катион. Различие в вычислен-

ных полных энергиях между этими случаями при

гурации, энергии образования вакансий оказывают-

U_eff = 3.5 эВ составляет 23 мэВ/f.u. (табл. 8, случай

ся завышенными примерно в два раза. При создании

VB *), что сравнимо с неточностью, которую можно

дефекта после релаксации атомов в бездефектной

получить при задании различных начальных маг-

конфигурации энергия образования B-вакансии ока-

нитных моментов в бездефектной сверхячейке. По-

зывается выше (E^F

= 2.41 эВ), чем в том случае,

этому из полученных данных, вероятно, нельзя сде-

когда вакансия создавалась в бездефектной конфи-

лать вывод о разнице в энергиях образования между

гурации до релаксации атомов (E^F

= 1.26 эВ), по-

случаями, когда удален Fe2+B или Fe3+B.

скольку найденное состояние дефектной конфигу-

При U_eff = 3.5 эВ разница в энергиях образова-

рации в этом случае обладает более высокой полной

ния вакансий в тетраэдрической и октаэдрической

энергией. Этот результат можно сравнить только с

подрешетках составляет 1.16 ± 0.20 эВ, что согласу-

результатом E^F

= 2.30 эВ работы [6] (табл. 1), в

ется с данными работ, выполненных ранее в рамках

которой проводились полная оптимизация решетки

DFT+U (табл. 1).

в сверхячейках с 56 атомами при U = 4.5 эВ и J =

= 0.89 эВ.

3.3.2. Локальное окружение дефекта

В случае, когда для создания вакансии атом уда-

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

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

без релаксации атомов, а затем проводилась релак-

A- и B-вакансий были получены значения 2.50 эВ и

сация атомов в сверхячейке с дефектом, в данной

1.26 эВ соответственно (табл. 8), рассматривались

работе были получены энергии образования так-

изменения локального окружения дефекта. Опти-

же при U_eff = 0.0 и U_eff = 5.0 эВ. Без учета силь-

мизированные конфигурации с A- и B-вакансиями

268

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

Рис. 16. Локальное окружение A-вакансии магнетита после

Рис. 17. Локальное окружение B-вакансии магнетита после

релаксации в случае, когда дефект создавался в равновес-

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

релаксации атомов в ней

показаны на рис. 16 и рис. 17 соответственно. Вид-

но, что в дефектных сверхячейках зарядовое и орби-

рошо согласуется с данными [10]. Это укорочение

тальное упорядочение нарушаются, а большинство

может говорить об усилении связи Fe_B - O₁. Объе-

двухвалентных катионов располагается как можно

мы вакансионных пор тетраэдрической и октаэдри-

дальше от центра вакансии.

ческой вакансии составляют 4.96Å3 и 15.44Å3, что

Полученная при U_eff = 3.5 эВ равновесная ку-

на 57 % и 15 % больше, чем объемы не занятых в без-

бическая фаза характеризуется длинами связей

дефектной решетке тетра- и октаэдрических междо-

= 1.908Å и l_Fe

= 2.079Å при a0 =

узлий. Внешнее напряжение в сверхячейке с А-ва-

lFeA-O

-O

= 8.479Å и x = 0.2549 (экспериментальные данные

кансией возрастает на 7.36 кбар, а в сверхячейке с

дают lFeA-O = 1.889Å, lFe_B -O = 2.058Å при a0 =

B-вакансией — снижается на 10.61 кбар.

= 8.394Å и x = 0.2549 [52]). Образование вакансий

Изменения полного магнитного момента сверхя-

в A- и B-подрешетках характеризуется в основном

чейки в случаях с тетраэдрической и октаэдриче-

укорочением длины связи Fe_B - O в первой коорди-

ской вакансией составляют +8 μ_B и -2 μ_B соответ-

национной сфере вакансии, причем анионы отдаля-

ственно, что согласуется с данными работ, выпол-

ются, а катионы становятся ближе к центру вакан-

ненных с учетом [10] и без учета [29] сильных элек-

сии. Такое укорочение lFeB -O1 составляет 6.4 % и

тронных корреляций. По вычисленным в работах

4.8 % для A- и B-вакансии соответственно, что хо-

[10, 29] значениям магнитных моментов катионов

269

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

Рис. 18. Степень заселенности d-орбиталей со спином вверх (в положительном направлении оси ординат) и со спином

вниз (с обратным знаком в отрицательном направлении оси ординат) для бездефектных и дефектных случаев

видно, что в этих работах разницы двух- и трехва-

лять μVA = 5 · 4μ_B + 11 · 5μ_B - 7 · 5μ_B = 40μ_B, что

лентных B-катионов найдено не было ни в дефект-

дает изменение 8 μ_B по сравнению с бездефектным

ном, ни в бездефектном случаях. В данной работе

случаем (μ = 32 μ_B). Аналогичным образом объяс-

эти отличия были получены (рис. 18, рис. 19). Из

няется разница в случае вакансии B-подрешетки.

рис. 18 видно, что в случаях с А- и B-вакансиями

Полученные данные согласуются с данными работы

число двухвалентных катионов в B-подрешетке сни-

[7], в которой отмечается, что в случае B-вакансии

жается с восьми до пяти штук. Иcходя из получен-

число трехвалентных B-катионов в сверхячейке уве-

ного числа двух- и трехвалентных катионов A- и

личивается на два, что по мнению авторов соот-

B-подрешеток, с помощью ионной модели магнети-

ветствует удалению B-катиона и перераспределению

та можно получить, что в случае A-вакансии пол-

двух электронов так, чтобы ближайшими соседя-

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

ми B-вакансии оказались Fe3+B, что восстановило бы

270

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

Рис. 19. Число электронов в сфере с радиусом Вигнера - Зейтца N^WZe , его d-составляющие N^WZe, эффективные заря-

ды ионов Z^WZ (Z = 8 - N^WZe) и магнитные моменты μ^WZ для катионов с номерами N. Использованный в расчетах

псевдопотенциал включал восемь валентных состояний для железа. А-катионы имеют номера 1-8, B-катионы — 16-24.

Заряд удаленного катиона обозначен равным 0

формальный заряд кислорода O2- у шести ближай-

дефектами ширина запрещенной зоны увеличивает-

ших к вакансии анионов.

ся (рис. 20).

В данной работе, как и в работе [7], ближайши-

4. ОБСУЖДЕНИЕ

ми соседями B-вакансии оказываются только трех-

валентные катионы. В случае A-вакансии большин-

Таким образом, в настоящей работе в рамках ме-

ство ближайших к ней B-катионов — также трех-

тода DFT+U были рассмотрены свойства кубичес-

валентные. При этом в области вакансионной поры

кой фазы магнетита и ее точечных дефектов. Полу-

снижается электронная плотность. Разницы в значе-

ченные данные позволяют сделать следующие вы-

ниях суммарного заряда между бездефектным и де-

воды.

фектным случаями, вычисленные при интегрирова-

1. Самосогласованное решение для равновес-

нии по сферам Вигнера-Зейтца, центрированных на

ной ионной конфигурации, отвечающей симметрии

атомах, составляют 6.8 e и 7.1 e для A- и B-вакансии

кубической фазы магнетита, может характеризо-

соответственно. Эти значения близки к вычислен-

ваться зарядовым упорядочением, отвечающим схе-

ными значениями зарядов катионов в бездефектной

ме Вервея, орбитальным упорядочением, а также

сверхячейке (рис. 19). Сравнение плотности элект-

иметь запрещенную зону, порядок которой соответ-

ронных состояний, вычисленных для дефектных и

ствует наблюдаемому экспериментально выше тем-

бездефектных случаев показывает, что в случаях с

пературы Вервея порядку псевдощели.

271

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

котемпературной симметрией. Он вносит дополни-

тельную неопределенность при интерпретации дан-

ных о релаксации локального окружения дефекта.

4. Создание дефекта в нерелаксированной струк-

туре кубической фазы позволило найти более глу-

бокий минимум для дефектной конфигурации, чем

создание дефекта в отрелаксированной структуре,

среди рассмотренных случаев. При этом эффект ре-

лаксации локального окружения дефекта кубиче-

ской фазы оказывается существенным: без релакса-

ции атомов в дефектной конфигурации энергии об-

разования дефектов кубической фазы оказываются

примерно в два раза выше, чем с релаксацией.

5. При образовании катионной вакансии в ку-

бической фазе магнетита, по-видимому, происходит

укорочение и усиление связей ближайших к вакан-

сии анионов кислорода с оставшимися катионами

Рис. 20. Изменение электронной плотности для случаев с

железа, что позволяет частично восстановить заряд

релаксацией атомов и дефектами

иона O2-. В общих чертах это согласуется с ионной

моделью магнетита. При этом в кубической сверхя-

Наличие этих черт в статической DFT+U модели

чейке с дефектом зарядовое и орбитальное упорядо-

кубической фазы магнетита, не включающей, кро-

чения полностью нарушены.

ме структуры, никаких температурных эффектов,

6. Полученные значения энергий образования ва-

может указывать на то, что вклад в наблюдаемую

кансий кубической фазы составляют 2.50 и 1.26 эВ

выше температуры Вервея проводимость магнетита

для вакансий A- и B-подрешеток соответственно.

вносят в основном динамические факторы, которые

Эти данные отличаются от энергий образования,

необходимо учесть для возможности сравнения рас-

представленных в DFT+U-работе [6]. В работе [6]

четных и экспериментальных данных.

указывается, что при проведении расчетов прово-

2. При вычислении свойств магнетита в малых

дилась полная оптимизация геометрии сверхячей-

сверхячейках, включающих 14 атомов, могут по-

ки, однако замечаний о том, сохранялась или нет

являться нежелательные эффекты, связанные с ее

симметрия кубической фазы при расчете бездефект-

размером и несоответствием элементарной ячейке

ной сверхячейки, не приводится. В данной работе

структуры шпинели; в то же время, наличие мно-

показано, что при оптимизации геометрии кубиче-

жества степеней свободы в больших сверхячейках

ской фазы без сохранения симметрии ее структура

не позволяет однозначно утверждать, что найден-

искажается. Поэтому можно предположить, что в

ный минимум является состоянием с наименьшей

работе [6] исследовались дефекты структуры, более

полной энергией для данной структуры. Вариация

близкой к структуре низкотемпературной, а не ку-

начальных приближений с последующим анализом

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

полученных решений выявила, что решения, обла-

в данной работе получены более надежные данные

дающие локализованными состояниями у половины

о дефектах кубической фазы, чем в работе [6].

B-катионов, с одинаковыми и как можно больши-

ми числами заполнения этих состояний, отличаются

5. ЗАКЛЮЧЕНИЕ

наименьшими полными энергиями среди найденных

для кубической фазы.

Выполненные расчеты и анализ полученных ре-

3. Релаксация атомов в сверхячейке кубической

зультатов позволяют сделать следующие выводы.

фазы при фиксированной постоянной решетки при-

1. Определены схемы зарядового и орбитального

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

водит к понижению полной энергии системы и уве-

личению запрещенной зоны. При этом решетка ис-

та с 14 и 56 атомами в элементарной ячейке. Пока-

зано наличие размерного эффекта: при сохранении

кажается с понижением симметрии. Такой эффект,

возникающий при снятии ограничений симметрии с

схемы Вервея в больших системах происходит час-

тичная разориентация орбиталей, снижающая пол-

позиций ионов, является отражением общего стрем-

ления системы перейти в основное состояние с низ-

ную энергию.

272

ЖЭТФ, том 160, вып. 2 (8), 2021

Энергии образования вакансий кубической фазы магнетита...

2. Проведено определение U_eff для модели с 56

J. Noh, O.I. Osman, S.G. Aziz, P. Winget, and

атомами на основе экспериментальных данных по

J. L. Brédas, Sci. Technol. Adv. Mater. 15, 044202

щели.

(2014).

3. Показано, что релаксация положений атомов

10.

R. Arras, B. Warot-Fonrose, and L. Calmels, J. Phys.

в кубической фазе в модели с 56 атомами приводит

Condens. Matter 25, 256002 (2013).

к изменению орбитального упорядочения при сохра-

11.

R. Arras, L. Calmels, and B. Warot-Fonrose, Appl.

нении зарядового упорядочения и к росту ширины

Phys. Lett. 100, 032403 (2012).

щели примерно на порядок. Данные изменения, по-

видимому, соответствуют стремлению системы пе-

12.

M. S. Senn, I. Loa, J. P. Wright, and J. P. Attfield,

рейти в низкотемпературную фазу.

Phys. Rev. B 85, 125119 (2012).

4. Показано, что удаление атома в расчетах энер-

13.

T. Fukushima, K. Yamauchi, and S. Picozzi, J. Phys.

гии образования вакансий должно проводиться из

Soc. Jpn. 80, 014709 (2011).

кубической фазы (ниже E_V ).

5. Проанализирована локальная структура ре-

14.

F. Zhou and G. Ceder, Phys. Rev. B 81, 205113

шетки после релаксации вакансий Fe_A/Fe_B.

(2010).

6. Получены энергии образования вакансий

15.

K. Yamauchi, T. Fukushima, and S. Picozzi, Phys.

E^F

для системы 56 атомов, существенно от-

VA /VB

Rev. B 79, 212404 (2009).

личающиеся от данных работы [6] по абсолютным

значениям. Подробный анализ сделанных прибли-

16.

P. Piekarz, K. Parlinski, and A. M. Oles, Phys. Rev.

жений позволяет считать результаты настоящей ра-

B 76, 165124 (2007).

боты более надежными, чем результаты работы [6].

17.

I. Leonov, A. N. Yaresko, V. N. Antonov, and

V. I. Anisimov, Phys. Rev. B 74, 165117 (2006).

Финансирование.

Работа

выполнена

18.

H. P. Pinto and S. D. Elliott, J. Phys. Condens.

при финансовой поддержке РФФИ (гранты

Matter 18, 10427 (2006).

№№ 20-21-00159, 20-32-90081).

19.

H. T. Jeng, G. Y. Guo, and D. J. Huang, Phys. Rev.

B 74, 195115 (2006).

ЛИТЕРАТУРА

20.

H. T. Jeng, G. Y. Guo, and D. J. Huang, Phys. Rev.

Lett. 93, 156403 (2004).

1. К. И. Кугель, Д. И. Хомский, УФН 136, 621 (1982).

21.

I. Leonov, A. N. Yaresko, V. N. Antonov, M. A. Koro-

2. С. В. Стрельцов, Д. И. Хомский, УФН 187, 1205

tin, and V. I. Anisimov, Phys. Rev. Lett. 93, 146404

(2017).

(2004).

3. V. I. Anisimov, J. Zaanen, and O. K. Andersen, Phys.

22.

V. I. Anisimov, I. S. Elfimov, N. Hamada, and K. Te-

Rev. B 44, 943 (1991).

rakura, Phys. Rev. B 54, 4387 (1996).

4. C. N. Lininger, C. A. Cama, K. J. Takeuchi,

23.

A. D. Rowan, C. H. Patterson, and L. V. Gasparov,

A. C. Marschilok, E. S. Takeuchi, A. C. West, and

Phys. Rev. B 79, 205103 (2009).

M. S. Hybertsen, Chem. Mater. 30, 7922 (2018).

24.

C. Freysoldt, B. Grabowski, T. Hickel, J. Neugebauer,

5. H. Liu and C. Di Valentin, J. Phys. Chem. C 121,

G. Kresse, A. Janotti, and C. G. Van de Walle, Rev.

25736 (2017).

Mod. Phys. 86, 253 (2014).

6. D. Li, B. Qu, H. Y. He, Y. G. Zhang, Y. Xu,

25.

M. Cococcioni and S. de Gironcoli, Phys. Rev. B 71,

B. C. Pan, and R. Zhou, Phys. Chem. Chem. Phys.

035105 (2005).

18, 7789 (2016).

26.

B. Meredig, A. Thompson, H. A. Hansen, C. Wolver-

7. C. L. Muhich, V. J. Aston, R. M. Trottier, A. W. Wei-

ton, and A. van de Walle, Phys. Rev. B 82, 195128

mer, and C. B. Musgrave, Chem. Mater. 28, 214

(2010).

(2016).

27.

J. P. Allen and G. W. Watson, Phys. Chem. Chem.

8. Y. Meng, X. W. Liu, C. F. Huo, W. P. Guo,

Phys. 16, 21016 (2014).

D. B. Cao, Q. Peng, A. Dearden, X. Gonze, Y. Yang,

J. Wang et al., J. Chem. Theory Comput. 12, 5132

28.

C. Ricca, I. Timrov, M. Cococcioni, N. Marzari, and

(2016).

U. Aschauer, Phys. Rev. Res. 2, 023313 (2020).

273

М. И. Шутикова, В. В. Стегайлов

ЖЭТФ, том 160, вып. 2 (8), 2021

29. S. Hendy, B. Walker, N. Laycock, and M. Ryan, Phys.

41.

J. H. Park, L. H. Tjeng, J. W. Allen, P. Metcalf, and

Rev. B 67, 085407 (2003).

C. T. Chen, Phys. Rev. B 55, 12813 (1997).

30. S. Zhang and J. Northrup, Phys. Rev. Lett. 67, 2339

42.

M. S. Senn, J. P. Wright, J. Cumby, and J. P. Attfield,

(1991).

Phys. Rev. B 92, 024104 (2015).

31. K. Momma and F. Izumi, J. Appl. Crystallogr. 44,

43.

N. C. Tombs and H. P. Rooksby, Acta Cryst. 4, 474

1272 (2011).

(1951).

32. H. Okudera, K. Kihara, and T. Matsumoto, Acta

44.

S. C. Abrahams and B. A. Calhoun, Acta Cryst. 6,

Crystallogr., Sect. B: Struct. Sci. 52, 450 (1996).

105 (1953).

33. E. J. W. Verwey, Nature 144, 327 (1939).

45.

M. E. Fleet, Acta Crystallogr., Sect. B: Struct.

Crystallogr. Cryst. Chem. 37, 917 (1981).

34. J. Yoshida and S. Iida, J. Phys. Soc. Jpn. 42, 230

(1977)

46.

G. Kresse and J. Hafner, Phys. Rev. B 47, 558 (1993)

35. M. Iizumi, T. F. Koetzle, G. Shirane, S. Chikazumi,

47.

G. Kresse and J. Furthmüller, Phys. Rev. B 54, 11169

M. Matsui, and S. Todo, Acta Crystallogr., Sect. B:

(1996).

Struct. Crystallogr. Cryst. Chem. 38, 2121 (1982).

48.

S. L. Dudarev, G. A. Botton, S. Y. Savrasov,

36. J. P. Wright, J. P. Attfield, and P. G. Radaelli, Phys.

C. J. Humphreys, and A. P. Sutton, Phys. Rev. B 57,

Rev. B 66, 214422 (2002).

1505 (1998).

37. M. S. Senn, J. P. Wright, and J. P. Attfield, Nature

49.

L. Zhang, M.

Šob, Z. Wu, Y. Zhang, and G. H. Lu,

481, 173 (2012).

J. Phys. Condens. Matter 26, 086002 (2014).

38. A. Yanase and K. Siratori, J. Phys. Soc. Jpn. 53, 312

50.

P. Vinet, J. Ferrante, J. R. Smith, and J. H. Rose,

(1984).

J. Phys. C: Solid State Phys. 19, L467 (1986).

39. A. Banerjee and A. J. Pal, J. Phys. Condens. Matter

51.

K. I. Kugel’ and D. I. Khomskii, Sov. Phys. Usp. 25,

32, 055701 (2020).

231 (1982).

40. L. V. Gasparov, D. B. Tanner, D. B. Romero, H. Ber-

52.

K. Glazyrin, C. McCammon, L. Dubrovinsky,

ger, G. Margaritondo, and L. Forró, Phys. Rev. B 62,

M. Merlini, K. Schollenbruch, A. Woodland, and

7939 (2000).

M. Hanfland, Am. Mineral. 97, 128 (2012).

274