Микроэлектроника, 2020, T. 49, № 5, стр. 323-333

Моделирование резистивного переключения в мемристорных структурах на основе оксидов переходных металлов

О. О. Пермякова ab*, А. Е. Рогожин a**

a Физико-технологический институт им. К.А. Валиева Российской академии наук
117218 Москва, Нахимовский проспект, 34, корп. 1, Россия

b Московский физико-технический институт (государственный университет)
141701 Долгопрудный, Институтский переулок, 9, Россия

* E-mail: o.permyakova@phystech.edu
** E-mail: rogozhin@ftian.ru

Поступила в редакцию 30.01.2020
После доработки 31.01.2020
Принята к публикации 19.02.2020

Полный текст (PDF)

Аннотация

Рассмотрены основные подходы к моделированию резистивного переключения в системах со слоями оксидов переходных металлов. Кратко изложены алгоритмы и результаты ab-initio моделирования таких систем, в частности значений энергетических барьеров для генерации, рекомбинации и миграции дефектов в пленках HfOx и Ta2O5. Подробно описаны подходы к моделированию с помощью метода конечных элементов и метода Монте-Карло и его результаты, особое внимание уделено изменению характера переключения в мемристорной структуре при изменении кинетического барьера для диффузии ионов кислорода, который определяется границей слоев, участвующих в резистивном переключений.

Ключевые слова: мемристор, резистивное переключение, моделирование Монте-Карло

ВВЕДЕНИЕ

В последнее время возрос интерес к энергонезависимой памяти с произвольным доступом на основе резистивного переключения (ReRAM). Высокая скорость переключения такой памяти позволяет рассматривать ее в качестве альтернативы динамической памяти DRAM. Еще больший интерес вызывает возможность создания нейроморфных систем на основе мемристоров [14].

Основным преимуществом ReRAM является ее способность к масштабированию – минимальный размер ячейки такой памяти составляет 4f  2 (f – размер элемента ячейки) [5, 6]. К резистивной памяти предъявляют следующие технологические требования: рабочее напряжение (<1 В), потребление энергии (~10 пДж/переключение), время переключения (<10 нс), количество циклов переключений (>109 циклов), время сохранения состояния (>10 лет), площадь ячейки (<576 нм2). Отдельные требования достигнуты для различных структур металл–изолятор–металл. Наилучшие результаты получены для оксидной памяти (OxRAM) на основе оксидов переходных металлов HfOx или TaOx [7].

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

Основные параметры мемристора: RLRS(RON) – сопротивление в состоянии с низким сопротивлением, RHRS(ROFF) — сопротивление в состоянии с высоким сопротивлением, VSET — напряжение переключения в состояние с низким сопротивлением, VRESET — напряжение переключения в состоянии с высоким сопротивлением, VFORMING – напряжение электроформовки, необходимое для первичного формирования филамента в диэлектрике. Соответственно можно выделить три основных процесса (рис. 1): процессы переключения между состояниями (SET и RESET) и электроформовка [5, 6].

Рис. 1.

Типичная ВАХ биполярного резистивного переключения с указанием напряжений переключения. Черным цветов выделен этап формовки.

Электрофизические характеристики структур с резистивным переключением определяются параметрами используемых материалов, а также распределением и свойствами дефектов в них. Для расчета параметров материалов и анализа свойств дефектов можно воспользоваться моделированием из первых принципов (ab-initio) [810].

Параметры, полученные при ab-initio моделировании, используют для кинетических моделей Монте-Карло и метода конечных элементов [1114]. С помощью этих подходов можно, оценивая изменение распределения дефектов, моделировать процесс формирования филаментов и рассчитывать электрические характеристики приборов.

Расчет характеристик цепей, состоящих из нескольких мемристоров, производиться с помощью компактных моделей [15]. При этом предполагается, что параметры мемристоров и статистические отклонения этих параметров известны.

В данной работе рассмотрены основные подходы к моделированию процессов переключения и электроформовки мемристорных структур на основе оксидов гафния HfOx и тантала TaOx. Кратко рассмотрены особенности ab-initio моделирования этих материалов и компактного моделирования цепей из мемристоров на основе HfOx и TaOx.

1. AB-INITIO МОДЕЛИРОВАНИЕ

На первом этапе ab-initio моделирования, то есть метода молекулярной динамики и теории функционала плотности, задается структура вещества, например, кристаллическая структура. Затем задается взаимодействие между атомами решетки с помощью набора псевдопотенциалов, называемого силовым полем (force field), рассчитанного из теории функционала плотности, либо выбранного с учетом экспериментальных результатов [5]. Силовое поле представляет собой функцию потенциальной энергии молекулы от координат ядер атомов, оно позволяет учесть электрон-электронное и электрон-ядерное взаимодействия в веществе. В подавляющем большинстве работ по ab-initio расчету параметров материалов мемристорных структур для построения модели используется прикладной пакет VASP (Vienna ab-initio simulation package).

В результате ab-initio моделирования можно получить значения энергии генерации и рекомбинации дефектов и энергии миграции (энергии активации) [8, 16], определить диффузионную подвижность дефектов внутри диэлектрика, значения энергии релаксации и термической ионизации, а также рассчитать электронную структуру в высокоомном и низкоомном состояниях. Верифицировать результаты ab-initio моделирования можно, сравнив их с экспериментальными данными о деградации подзатворного диэлектрика во времени [17], переменном изменением температуры и характеристикой произвольного телеграфного шума [18].

С помощью ab-initio моделирования были получены значения энергий активации вакансий кислорода на границе Hf/HfO2, а также в объеме HfO2 для моноклинной структуры m-HfO2 [8]. Рассчитанное значение энергии образования кислородных вакансий на границе металл-оксид оказалось ниже на 1.1 эВ по сравнению с энергией образования вакансий в объеме оксида. Полученное значение начальной энергии миграции кислородных вакансий на границе Hf/m-HfO2 составляет 1.30 эВ, в объеме оксида – 1.76 эВ, а вблизи границ зерен в объеме оксида – 0.87 эВ. Минимальное значение энергии необходимое для преодоления диффузионного барьера вдоль границ зерен в m-HfO2 зависит от направления движения [16]. Так при движении дважды положительных вакансий в объеме зерна минимальное значение энергии составляет 0.71 эВ, значение энергии при движении к границе составляет 1.14 эВ и при движении к центру зерна 0.42 эВ. Необходимо отметить, что ввиду кулоновского отталкивания наибольшая концентрация вакансий вдоль границ зерен ограничена значением 2.5 × 1014 см–2.

Напряжение формовки резистивного переключение в слое оксида сильно зависит от энергии генерации вакансий в отсутствии электрического поля. С помощью легирующих примесей можно изменить энергию образования вакансий в отсутствии поля в оксидном слое HfOx [17]. Так среди примесей Al, Ti и Si, наибольший эффект на снижение энергии генерации, а значит и на напряжения формовки оказывает Al, затем Ti и Si. Из-за того, что кубическая фаза HfO2 наиболее плотная снижение энергии генерации в этом случае оказывается больше по сравнению с моноклинной и ромбической фазой.

Образование дефектов на границе материалов и переход ионов кислорода в кислородный обменный слой, например Hf для HfO2, или переход ионов кислорода в резервуар, т.е. активный электрод, является одной из гипотез образования филамента при биполярном переключении. В мемристорных структурах на основе оксидов HfO2 и Ta2O5 влияние кислородного обменного слоя вблизи верхнего электрода (Hf, Ti, Ta) и влияние нижнего электрода (TiN, Ru) на время сохранения состояния и повторяемость ВАХ в зависимости от числа циклов определяются кинетическими барьерами для диффузии кислорода между слоями вещества [19]. Эти барьеры также определяют скорость переключения и энергию необходимую для резистивного переключения.

Результаты ab-initio расчетов говорят о том, что кислородные вакансии, формирующие проводящий канал, имеют заряд 0 или +1, тогда как мобильные изолированные вакансии внутри объема диэлектрика имеют заряд +2 из-за сильной зависимости когезионной энергии вакансии от ее заряда [20].

Рассматривалась возможность влияния катионов металла в оксиде на механизм образования канала [21, 22]. Например, отношение коэффициентов диффузии ионов тантала и кислорода в TaOx достигает 0.25–0.5, а отношение коэффициентов диффузии ионов гафния и кислорода в HfOx не превышает 0.013–0.1 [22].

Проводилось моделирование и более сложных систем, например, структур TiN/HfO2/Ti/TiN и TiN/Ta2O5/TaOx/TiN [23]. Показано, что, несмотря на более высокое отношение сопротивлений для структуры с оксидом тантала (на порядок больше), число циклов переключения в слое Ta2O5/TaOx составляет только 104 по сравнению со структурой на основе оксида гафния, число циклов переключения в которой более 108. Авторы объясняют данный эффект большей энергией формирования вакансий в стеке HfO2/Ti.

С помощью ab-initio моделирования рассчитаны энергии образования дефектов внутри оксида тантала [9, 10, 24]. Значения энергии формирования вакансий в значительной степени зависят от конфигурации λ фазы Ta2O5 [9]. Так для конфигурации 2f-a энергии активации нейтральных вакансий 5.74 эВ и дважды положительных вакансий 0.69 эВ являются наибольшими, тогда как энергия активации для положительных вакансий 2.89 эВ, что является средним значением по сравнению с энергиями активаций для других конфигураций рассмотренных в работе. В тонком промежуточном слое Ta2O3 слое миграция вакансий может происходить с меньшими энергиями активаций, кроме того мобильные вакансии с зарядом +2 имеют меньшую энергию активации по сравнению с неподвижными, образующими филамент, нейтральными и положительно заряженными вакансиями [10]. Результаты сравнения параметров дефектов (кислородных вакансий и междоузельных ионов Ta) в Ta2O5 говорят о том, что с уменьшением химического потенциала атомов кислорода вероятность формирования ионов Ta становиться сравнимой с вероятностью формирования кислородных вакансий [24].

2. МОДЕЛЬ МОНТЕ-КАРЛО И МЕТОД КОНЕЧНЫХ ЭЛЕМЕНТОВ

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

Рис. 2.

Типичная блок-схема метода Монте-Карло для моделирования процессов переключения в ReRAM структурах.

Модель полного резистивного переключения, построенную с помощью метода Монте-Карло и метода конечных элементов удобно разделить на следующие модули [21], которые также выделены в блок-схеме (рис. 2):

• Начальная структура: представляет собой начальное распределение вакансий/ионов внутри диэлектрика;

• Кинетика дефектов: описывает механизмы формирования, движения и рекомбинации дефектов внутри моделируемой структуры;

• Модель электронной проводимости: необходима для получения вольтамперных характеристик (ВАХ) моделируемой структуры;

• Модель распространения тепла;

• Стохастическая компонента: алгоритм МК, применяемый к кинетической модели, в котором определяются частоты или времена релаксации диффузии, образования и рекомбинации дефектов для создания статистической модели.

Процесс резистивного переключения включает в себя несколько этапов, а именно, процесс формовки, процесс переключения в состояние с низким сопротивлением (SET) и процесс переключения в состояние с высоким сопротивлением (RESET). Модель может описывать только один или несколько этих процессов.

Начальная структура

В зависимости от моделируемых процессов задается начальное распределение вакансий. Например, для модели, описывающей все три процесса, задается небольшое количество вакансий в объеме диэлектрика, моделирующие границы зерен [11]. В отсутствии процесса формовки необходимо начальное распределение вакансий характерное для обратимого пробоя диэлектрика. Например, начальным состоянием ячейки TiN/HfO2/Hf/TiN с низким сопротивлением может быть филамент, имеющий ассиметричную форму песочных часов [25]. Количество вакансий, образующих такой проводящий канал, обычно калибруют с помощью экспериментальной ВАХ.

Кинетическая модель и стохастическая компонента

Кинетические барьеры для диффузии ионов О определяют скорость резистивного переключения [19], энергию необходимую для переключения и способность к сохранению состояния. Правильно задав кинетические барьеры и условия на границах между слоями, можно добиться, чтобы резистивное переключение происходило только в одной области филамента. Такой режим переключения, характерный для ассиметричных структур, может обеспечить более стабильные параметры процессов переключения. Для создания асимметрии филамента один из электродов должен быть активным, т. е. обладать высоким сродством с кислородом, а другой инертным к оксидации и обладать высокой работой выхода. Асимметрия приводит к кинетической модели, в которой вакансии образуются вблизи одного из электродов и затем диффундируют в сторону другого электрода (рис. 3) [11]. Впрочем, генерация вакансий может происходить как на границе с активным электродом, так и в объеме диэлектрика [12], но энергетический барьер необходимый для генерации дефекта в объеме диэлектрика значительно выше. В работе [13] также предполагается один активный электрод (резервуар), но считается, что генерация вакансий может происходить во всем объеме диэлектрика, тогда как их рекомбинация с ионами кислорода происходит только на границе с резервуаром.

Рис. 3.

Движение дефектов в диэлектрике, рассматриваемое при моделировании: зеленый круг – вакансия кислорода, сплошной круг – ион кислорода, КОС – кислородный обменный слой; 1 – генерация вакансии на границе с резервуаром; 2 – рекомбинация вакансии на границе с резервуаром; 3 – генерация вакансии в объеме диэлектрика; 4 – рекомбинация вакансии и иона кислорода в объеме диэлектрика; 5 – диффузия кислородных вакансий в объеме диэлектрика. Диффузия ионов кислорода рассматривается как последовательность 3 и 4.

Существуют и другие кинетические модели образования вакансий. Например, оба электрода могут считаются активными [25, 26], тогда генерация дефектов происходит у каждого из них. Параметры переключения структуры TiN/HfO2/Hf/TiN зависят от изменения количества вакансий образующих филамент, а также емкости (количество электронов, которые может поглотить электрод) верхнего и нижнего резервуаров [25]. Наиболее стабильные параметры переключения при низком токе (высокое сопротивление) происходят при асимметрии емкостей резервуаров, при этом филамент приобретает коническую форму. С другой стороны, при равных емкостях процесс переключения может становиться нестабильным, т.е. увеличивается дисперсия напряжений VSET и VRESET.

В модели [27] рассматривается разрыв в середине сформированного канала униполярной структуры. Процессы генерации и рекомбинации вакансий рассматривались только в пределах разрыва.

Кинетическая модель может описывать не только механизм образования и движения вакансий [11, 13, 2629], но и движение ионов, например, ионов кислорода [12, 3034]. При этом движение ионов может рассматриваться как во всем объеме структуры [3033], так и только в кислородном обменном слое [34]. В этом случае образование вакансии связано с перемещением иона кислорода со своего положения в кристаллической решетке в междоузлие, а движение вакансии – с переходом иона кислорода в соседнюю вакансию или в кислородный обменный слой. Однако, возможен подход, при котором модель предполагает, что вакансии неподвижны, но могут появляться и исчезать, а ионы кислорода могут двигаться в объеме диэлектрика или переходить в кислородный обменный слой [30, 32, 34].

Формирование цепочки вакансий – не единственный возможный механизм образования проводящего канала [35]. Для TaOx-ReRAM, основным механизмом образования проводящего канала является не только образование большого количества вакансий в диэлектрике Ta2O5, но и переход в проводящую фазу TaO2.

Динамику дефектов внутри структуры (рис. 3), т.е. процессы диффузии, генерации и рекомбинации можно описать с помощью уравнений для частот процессов [21]:

(1)
${{R}_{G}} = {{f}_{0}}\exp \left( {\frac{{\Delta W - {{\beta }_{1}}{{F}_{{{\text{ext}}}}} - {{\beta }_{2}}\nabla T}}{{{{k}_{b}}T}}} \right),$
(2)
${{R}_{R}} = {{f}_{0}}{\text{exp}}\left( {\frac{{\Delta W - {{\beta }_{1}}{{F}_{{{\text{ext}}}}} - {{\beta }_{2}}\nabla T}}{{{{k}_{b}}T}}} \right),$
(3)
${{R}_{D}} = {{f}_{0}}{\text{exp}}\left( {\frac{{\Delta W - {{\beta }_{1}}{{F}_{{{\text{ext}}}}} - {{\beta }_{2}}\nabla T}}{{{{k}_{b}}T}}} \right),$
где RG, RR и RD – частоты генерации, рекомбинации и диффузии соответственно; ΔW – энергия активации исследуемого процесса; Fext – напряженность локального электрического поля; ∇T – градиент температуры; β1 и β2 – коэффициенты, описывающие снижение энергии активации из-за электрического поля и градиента температуры, соответственно; kb – постоянная Больцмана. Иногда вместо скоростей рассматривают обратное время релаксации [21].

Коэффициент β1 описывает снижение потенциального барьера активации процесса из-за локального электрического поля. Например, для генерации дефекта в качестве коэффициента β1 используют коэффициент поляризации связи [31] или коэффициент симметрии (обычно близок к 0.5), умноженный на постоянную решетки [11]. В некоторых случаях, понижение потенциального барьера рекомбинации из-за локального электрического поля может не учитываться вовсе (β1 = 0) [31]. При рассмотрении диффузии ионов кислорода β1 иногда содержит в себе длину прыжка (расстояние от предыдущего до текущего положения) и зарядовое число иона [12], а при анализе движения вакансий иногда выбирают β1 равным коэффициенту симметрии [29]. Множитель f0 обычно принимают равным характеристической частоте колебаний (1012–1013 Гц) [12].

Коэффициент β2 описывает снижение потенциального барьера из-за градиента температуры. Как правило β2 равен произведению постоянной решетки и постоянной Больцмана и используется при рассмотрении скорости диффузии дефектов за счет термодиффузии ионов вдоль градиента температуры [11, 13]. Для процессов генерации и рекомбинации обычно β2 = 0.

Формирование филамента

Формирование проводящего канала в слое диэлектрика происходит под действием внешнего электрического поля [5]. Распределение электрического потенциала рассчитывается с помощью уравнения Пуассона:

(4)
$\Delta \varphi = - \frac{\rho }{\varepsilon },$
где φ – электрический потенциал, ρ – плотность заряда, ε – диэлектрическая постоянная среды. Для ускорения вычислений иногда предварительно с помощью уравнения Лапласа получают начальное приближение [12, 13]. При этом на границах электродов используются граничные условия Дирихле, а на других границах – периодические граничные условия.

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

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

Существует несколько подходов к моделированию тока, протекающего в слое диэлектрика. Как правило, при моделировании процессов в системе находящейся в высокоомном состоянии (формовка, SET) и процессов происходящих при низком напряжении считают, что ток, проходящий, через диэлектрик описывается моделью тунелирования через ловушки (ТАТ – trap-assisted tunneling) (рис. 4). В TAT модели предполагается, что туннельный механизм движения электронов от одного электрода до другого происходит через разрешенные уровни ловушек в диэлектрике, например, кислородные вакансии. Применение этой модели в слоях HfOx обусловлено экспериментальными результатами говорящими о том, что ток утечки таких структур подчиняется закону Пула-Френкеля [36]. Применение модели ТАТ для систем с резистивным переключением подробно описано в [13]. По мере формирования филамента появляется возможность использования и других моделей тока.

Рис. 4.

Схема TAT процесса и схема изменения заселенности ловушек электронами высокоомного состояния для резистивного переключения на границах HfO2/Hf (1) [11] и TiN/HfO2 (2) [13].

При моделировании процессов в системе находящейся в низкоомном состоянии (RESET) используются дрейфово-диффузионная модель [11, 30, 33], QPC модель [25, 26] и другие, например, сочетание моделей QPC и модели омической проводимости [12]. Под QPC моделью имеется в виду модель квантового точечного контакта, в которой филамент имеет малую площадь сечения и учитывается квантование проводимости [37]. Таким образом, в низкоомном состоянии выбор модели тока зависит от площади сечения филамента.

Условие перехода между моделями для высокоомной и низкоомной систем также может задаваться по-разному. Критерием перехода может быть достижение критической концентрации вакансий в объеме диэлектрика, например, более чем 150 вакансий на 125 нм3 [11]. Другой возможный критерий – достижение определенной площади поперечного сечения филамента, например, для поперечного сечения проводящего канала менее 0.5 нм возможно применение QPC модели [12].

Модель распространения тепла

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

(5)
$P\left( {x,y,z} \right) = {{k}_{{th}}}\Delta T,$
где kth – теплопроводность исследуемой пленки; P(x, y, z) – мощность, локально рассеянная внутри решетки. Для упрощения вычислений обычно используется первое приближение уравнения теплопроводности в предположении, что температура вдоль филамента однородна и растет из-за разогрева электрическим током, как:
(6)
$T = 297\,\,{\text{K}} + \left| {\left( {{{V}^{L}} - {{V}^{R}}} \right)I} \right|{{R}_{{th}}},$
где V L и V R – значения напряжения, приложенного к электродам, I – ток, проходящий в проводящем канале и Rth – эквивалентное тепловое сопротивление. В соответствии с данной моделью было показано, что во время процесса переключения температура может возрастать на 170 K [13].

Возможно использование и второго приближения. Его имеет смысл применять в областях с наибольшей напряженностью электрическое поле [25]:

(7)
$T = {{T}_{{{\text{ambient}}}}} + \frac{{\alpha VI}}{{{{n}_{c}}}}{{R}_{{th}}} - {{f}_{{{\text{dissipation}}}}}{{\left( {T - {{T}_{{{\text{ambient}}}}}} \right)}^{2}},$
где Tambient – температура окружающей среды, nc – концентрация дефектов в проводящем канале, α – коэффициент понижения барьера, fdissipation – коэффициент рассеяния тепла.

Более точное распределение температуры можно получить, если при расчете распространения тепла в материале учесть его теплоемкость [12]:

(8)
${{c}_{p}}\left( r \right)\rho \left( r \right)\frac{{\delta T\left( {r,t} \right)}}{{\delta t}} - \nabla \left[ {\lambda \left( r \right)\nabla T\left( {r,t} \right)} \right] = \frac{{{{j}^{2}}\left( {r,t} \right)}}{{\sigma \left( r \right)}},$
где cp(r) – локальная теплоемкость, ρ(r) – плотность вещества, λ(r) – теплопроводность материала, j(r, t) локальная плотность тока и σ(r) – локальная проводимость внутри структуры.

Метод конечных элементов для моделирования резистивного переключения

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

(9)
$\frac{{\partial {{n}_{D}}}}{{\partial t}} = \nabla \left( {D\nabla {{n}_{D}} - \nu {{n}_{D}}} \right) + G,$
где nD – концентрация вакансий, t – время, D – коэффициент диффузии вакансий, ν – скорость дрифта вакансий, G – скорость генерации вакансий при отрицательном потенциале для структуры [14].

С помощью метода конечных элементов можно оценить влияние температуры окружающей среды на характеристики переключения структуры Pt/Ta2O5/TaOx/W. Оказывается, что при увеличении температуры с 25 до 85°С наблюдается небольшое уменьшение VRESET на величину порядка 0.1 В [14]. Также было показано, что положение разрыва филамента во время перехода в высокоомное состояние может не зависеть от формы образовавшегося проводящего канала – разрыв происходит вблизи верхнего платинового электрода. Помимо этого в случае прямоугольной формы канала при использовании различных материалов верхнего электрода (Pt, Ru, W) наблюдается изменение VRESET на значение порядка 0.1 В, при этом наибольшее значение напряжение переключения имеет структура с верхним электродом из W, а наименьшее структура с электродом из Pt, авторы считают, что это связанно с наибольшей теплопроводностью Pt и наименьшей у W.

Метод конечных элементов позволяет рассмотреть процесс формовки резистивного переключения [38], а также состояние сразу после формовки в слое TaOx с верхним омическим электродом и нижним электродом Шоттки. В таком случае форма образовавшегося проводящего канала зависит от энтальпии миграции вакансий кислорода. При повышении энтальпии уменьшается скорость диффузии вакансий. В результате большее количество вакансии остается у анода, из-за чего филамент имеет коническую форму, а не цилиндрическую.

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

Для задания стохастической (произвольной) компоненты в образовании филамента используется метод Монте-Карло. Чтобы определить вероятность событий генерации, рекомбинации и диффузии дефектов находят их частоты возникновения (1)–(3). Затем определяют вероятности возникновения событий [39]:

(10)
$p = 1 - {\text{exp}}\left( {Rdt} \right),$
где p – вероятность события, R – частота возникновения события, dt – рассматриваемый отрезок времени. На основании данных вероятностей определяются события, которые произошли на этом шагу, частоты которых используются для определения шага по времени. Соответствие между внутренним временем симуляции и временем, за которое процесс происходит в реальности, можно получить при использовании формулы [40]:
(11)
$dt = - \frac{{{\text{ln}}\left( {RND} \right)}}{{\Sigma R}},$
где dt – шаг по времени для следующего шага, ΣR – сумма частот всех произошедших событий, RND – произвольное число от 0 до 1.

Рассмотрим другие примеры использования метода Монте-Карло. Для униполярного переключения в структуре Pt/HfO2/Pt можно представить симметричный филамент в форме песочных часов [28], тогда можно выделить два механизма разложения канала при переключении в высокоомное состояние. Первый связан с электрическим полем в структуре – в этом режиме сопротивление постепенно увеличивается с ростом напряжения. Второй механизм определяется диссипацией тепла в структуре. При этом сопротивление может увеличиваться скачкообразно. Механизм переключения из высокоомного в низкоомное состояния в униполярной структуре Pt/HfO2/Pt можно проанализировать, задав разрыв в филаменте цилиндрической формы [27].

Частота отказов механизма считывания состояния зависит от формы проводящего канала [26], форма которого определяется скоростью нарастания напряжения и ограничением тока (Icomp) в течение процесса переключения. Например, для OxRAM структуры оксидным слоем которой является ультратонкий слой HfOx наибольшее количество отказов происходит для скоростей роста напряжения в пределах 1–10 мВ/с и ограничением тока 1–10 мкА, а наименьшее их количество при скоростях роста напряжения 103–104 мВ/с и ограничением тока 50–100 мкА.

С помощью метода Монте-Карло также можно оценить диаметр проводящего канала и концентрацию вакансий в нем на основе экспериментальных данных о сохранении состояний RLRS и RHRS. Опираясь на оценки изменения диаметра филамента концентрации вакансий в нем, возможно предсказать дальнейшее изменение состояний RLRS и RHRS [41]. Например, для определения параметров филамента данные о RLRS и RHRS, полученные во время измерений, сравниваются с теми же данными, полученными в результате моделирования с помощью критерия согласия Колмогорова. Например, для структуры Ir/Ta2O5 – δ/TaOx/TaN, с размером верхнего электрода 0.5 × 0.5 мкм2 и толщиной диэлектрика порядка 30 нм, диаметр филамента изменяется от 15 до 40 нм. Концентрация вакансий в высокоомном состоянии практически одинакова для разных приборов, а в низкоомном состоянии имеет значительный разброс.

Ниже приведены (табл. 1) значения для напряжения переключения для нескольких циклов различных структур, смоделированные в работах [1114, 30, 32, 34].

Таблица 1.  

Напряжения переключения мемристорных структур

Структура VSET, V VRESET, V Работа
HfO2(5 нм)/Hf 1.6–1.7 –1.7…–2.5 [11]
HfOx(6 нм)/TiOy 1.1 –0.3 [30]
Ti/HfO2(5 нм) 0.36–0.5 –0.58…–0.74 [34]
Ti/HfO2(9 нм) 0.6 –0.5 [12]
TiN/HfOx(10 нм) 2.6–2.7 –1.3 [13]
TiN/HfOx(10 нм) 1.2–1.8 –0.6…–1.4 [32]
Ta2O5(20 нм)/TaOx(40 нм) 1 –0.9…–1.0 [14]

3. КОМПАКТНЫЕ МОДЕЛИ

Компактные модели опираются на готовые выражения для расчета тока и сопротивления, зависящие от электрофизических характеристик структуры. Это позволяет применять их для расчета электрических цепей. Для калибровки таких моделей могут быть использованы данные о ВАХ структуры, полученные из эксперимента, либо характеристики проводящего канала (емкость, сопротивление, количество филаметов), полученные из моделирования с помощью методов конечных элементов и Монте-Карло. Большое количество моделей тока для мемристорных структур представлено в монографии [15].

Также ReRAM могут представляться в виде эквивалентных цепей [42]. Долговечность (время сохранения состояния в ячейке) ReRAM структуры Ta2O5/TaOx [43] моделируют с помощью стохастических дифференциальных уравнений и сети случайных сопротивлений. Такой метод является менее точным, но позволяет предсказать поведение структуры от цикла к циклу намного быстрее, чем метод Монте-Карло.

ЗАКЛЮЧЕНИЕ

Моделирование резистивного переключения в слоях оксидов переходных металлов важно для понимания работы и создания OxRAM и нейроморфных систем. Опираясь на результаты моделирования, можно выбрать наиболее подходящие материалы для мемристоров. Моделирование резистивного переключения обычно разделяют на три основных этапа: ab-initio моделирование материалов, из которых состоит прибор, моделирование с помощью методов конечных элементов и Монте-Карло кинетики процесса резистивного переключения, и моделирование цепей, включающих в себя мемристоры, с помощью компактных моделей.

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

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

Список литературы

  1. Sungjun Kim et al. Neuronal dynamics in HfOx/AlOy-based homeothermic synaptic memristors with low-power and homogeneous resistive switching // Nanoscale. 2019. V. 11. P. 237–245.

  2. Matveyev Yu., Egorov K., Markeev A., Zenkevich A. Resistive switching and synaptic properties of fully atomic layer deposition grown TiN/HfO2/TiN devices // J. Appl. Phys. 2017. V. 117. P. 044901.

  3. Ielmini D. Brain-inspired computing with resistive switching memory (RRAM): Devices, synapses and neural networks // Microelectron. Eng. 2018. V. 190. P. 44–53.

  4. Hu S.G. et al. Associative memory realized by a reconfigurable memristive Hopfield neural network // Nature Communications. 2015. V. 6. P. 7522.

  5. Ielmini D., Waser R. Wiley-VCH Book, ISBN: 978-3-527-33417-9, Wiley-VCH, Weinheim, Germany, 2016.

  6. Yoshio Nishi, Blanka Magyari-Kope. Woodhead Publishing, ISBN: 9780081025840 9780081025857, 2019.

  7. Lanza M. et al. Recommended Methods to Study Resistive Switching Devices // Adv. Electron. Mater. 2019. V. 5. P. 1800143.

  8. O’Hara A., Bersuker G., Demkov A.A. Assessing hafnium on hafnia as an oxygen getter // J. Appl. Phys. 2014. V. 115. P. 183703.

  9. Hao Jiang, Derek A. Stewart. Enhanced oxygen vacancy diffusion in Ta2O5 resistive memory devices due to infinitely adaptive crystal structure // J. Appl. Phys. 2016. V. 119. P. 134502.

  10. Hur J.H. Theoretical studies on oxygen vacancy migration energy barrier in the orthorhombic λ phase Ta2O5 // Comput. Mater. Sci. 2019. V. 169. P. 109148.

  11. Abbaspour E. et al. KMC Simulation of the Electroforming, Set and Reset Processes in Redox-Based Resistive Switching Devices // IEEE Transactions on Nanotechnology. 2018. V. 17. P. 1181–1188.

  12. Dirkmann S., Kaiser J., Wenger C., Mussenbrock T. Filament Growth and Resistive Switching in Hafnium Oxide Memristive Devices // ACS Appl. Mater. Interfaces. 2018. V. 10. P. 14857–14868.

  13. Ximeng Guan, Shimeng Yu, H.-S. Philip Wong. On the Switching Parameter Variation of Metal-Oxide RRAM–Part I: Physical Modeling and Simulation Methodology // IEEE Transactions on Electron Devices. 2012. V. 59. P. 1172–1182.

  14. Kim Sungho et al. Physical electro-thermal model of resistive switching in bi-layered resistance-change memory // Scientific Reports. 2013. V. 3. P. 1680.

  15. Mladenov V. Advanced Memristor Modeling, Memristor Circuits and Networks Memristor Modeling, Memristor Devices, Circuits and Networks; MDPI: Basel, Switzerland, 2019.

  16. McKenna K., Shluger A. The interaction of oxygen vacancies with grain boundaries in monoclinic HfO2 // Appl. Phys. Lett. 2009. V. 95. P. 222111.

  17. Liang Zhao, Sergiu Clima, Blanka Magyari-Köpe, Malgorzata Jurczak, Yoshio Nishi. Ab initio modeling of oxygen-vacancy formation in doped-HfOx RRAM: Effects of oxide phases, stoichiometry, and dopant concentrations // Appl. Phys. Lett. 2015. V. 107. P. 013504.

  18. Brivio S., Frascaroli J., Covi E., Spiga S. Stimulated Ionic Telegraph Noise in Filamentary Memristive Devices // Scientific Reports. 2019. V. 9. P. 6310.

  19. Clima S. et al. First-principles thermodynamics and defect kinetics guidelines for engineering a tailored RRAM device // J. Appl. Phys. 2016. V. 119. P. 225107.

  20. Dan Duncan, Blanka Magyari-Köpe, Yoshio Nishi. Filament-Induced Anisotropic Oxygen Vacancy Diffusion and Charge Trapping Effects in Hafnium Oxide RRAM // IEEE Electron Device Letters. 2016. V. 37. P. 400–403.

  21. Wedig A. et al. Nanoscale cation motion in TaOx, HfOx and TiOx memristive systems // Nature Nanotech. 2016. V. 11. P. 67–74.

  22. Bo Xiao, Xuefang Yu, Satoshi Watanabe. A Comparative Study on the Diffusion Behaviors of Metal and Oxygen Ions in Metal-Oxide-Based Resistance Switches via ab Initio Molecular Dynamics Simulations // ACS Appl. Electron. Mater. 2019. V. 1. P. 585–594.

  23. Azzaz M. et al. Endurance/Retention Trade Off in HfOx and TaOx Based RRAM // 2016 IEEE 8th IMW. 2016. Paris, France.

  24. Linggang Zhu, Jian Zhou, Zhonglu Guo, Zhimei Sun. Synergistic Resistive Switching Mechanism of Oxygen Vacancies and Metal Interstitials in Ta2O5 // J. Phys. Chem. C. 2016. V. 120. P. 2456–2463.

  25. Degraeve R. et al. Modeling RRAM set/reset statistics resulting in guidelines for optimized operation // 2013 Symposium on VLSI Technology. 2013. Kyoto, Japan.

  26. Raghavan N., Frey D.D., Bosman M., Pey K.L. Monte Carlo model of reset stochastics and failure rate estimation of read disturb mechanism in HfOx RRAM // 2015 IEEE International Reliability Physics Symposium. 2015. Monterey, CA, USA.

  27. Yu Li et al. Investigation on the Conductive Filament Growth Dynamics in Resistive Switching Memory via a Universal Monte Carlo Simulator // Scientific Reports. 2017. V. 7. P. 11204.

  28. Shibing Long et al. Voltage and Power-Controlled Regimes in the Progressive Unipolar RESET Transition of HfO2-Based RRAM // Scientific Reports. 2013. V. 3. P. 2929.

  29. Abbaspour E., Menzel S., Jungemann C. The role of the interface reactions in the electroforming of redox-based resistive switching devices using KMC simulations // 2015 SISPAD. 2015. Washington, DC, USA.

  30. Padovani A. et al. Microscopic Modeling of HfOx RRAM Operations: From Forming to Switching // IEEE Transactions on Electron Devices. 2015. V. 62. P. 1998–2006.

  31. Butcher B. et al. Modeling the effects of different forming conditions on RRAM conductive filament stability // 2013 5th IEEE International Memory Workshop. 2013. Monterey, CA, USA.

  32. Shimeng Yu, H.-S. Philip Wong. Characterization and Modeling of the Conduction and Switching Mechanisms of HfOx Based RRAM // Symposium P – Emergent Electron Transport Properties at Complex Oxide Interfaces. 2014. V. 1631. mrsf13-1631-p04-03.

  33. Larcher L., Padovani A. Multiscale modeling of oxide RRAM devices for memory applications: from material properties to device performance // J. Comput. Electron. 2017. V. 16. P. 1077–1084.

  34. Butcher B. et al. Connecting the physical and electrical properties of Hafnia-based RRAM // 2013 IEEE International Electron Devices Meeting. 2013. Washington, DC, USA.

  35. Zhao Y.D. et al. Insights into resistive switching characteristics of TaOx-RRAM by Monte-Carlo simulation // 2015 International Symposium on VLSI Technology, Systems and Applications. 2015. Hsinchu, Taiwan.

  36. Doo Seok Jeong, Cheol Seong Hwang. Tunneling-assisted Poole-Frenkel conduction mechanism in HfO2 thin films // J. Appl. Phys. 2005. V. 98. P. 113701.

  37. Xiaojuan Lian et al. Quantum point contact model of filamentary conduction in resistive switching memories // 2012 13th International Conference on Ultimate Integration on Silicon (ULIS).

  38. Marchewka A., Waser R., Menzel S. Physical modeling of the electroforming process in resistive-switching devices // 2017 SISPAD. 2017. Kamakura, Japan.

  39. Aldana S. et al. A 3D kinetic Monte Carlo simulation study of resistive switching processes in Ni/HfO2/Si-n+-based RRAMs // J. Physics D Applied Physics. 2017. V. 50. № 33.

  40. Kristen A. Fichthorn. Theoretical foundations of dynamical Monte Carlo simulations // J. Chem. Phys. 1991. V. 95. P. 1090

  41. Wei Z. et al. A new prediction method for ReRAM data retention statistics based on 3D filament structures // 2015 IEEE International Reliability Physics Symposium. 2015. Monterey, CA, USA.

  42. Miranda E. et al. Equivalent circuit model for the electron transport in 2D resistive switching material systems // 2017 ESSDERC. 2017. Leuven, Belgium

  43. Zhiqiang Wei, Koji Eriguchi. Analytic Modeling for Nanoscale Resistive Filament Variation in ReRAM With Stochastic Differential Equation // IEEE Transactions on Electron Devices. 2017. V. 64. P. 2201–2206.

Дополнительные материалы отсутствуют.