CC BY
15
1.2. УСТОЙЧИВОСТЬ АТОМНОЙ СТРУКТУРЫ НАНОСИСТЕМЫ ПРИ МОЛЕКУЛЯРНОМ ПЕРЕКЛЮЧЕНИИ
Попов Александр Михайлович, профессор, д-р физ.-мат. наук, профессор, факультет вычислительной математики и кибернетики, Московский Государственный Университет им. М. В. Ломоносова, e-mail: [email protected]
Шумкин Георгий Николаевич, докторант, канд. физ.-мат. наук, факультет вычислительной математики и кибернетики, Московский Государственный Университет им. М. В. Ломоносова, e-mail: [email protected]
Аннотация: В статье разрабатывается математическая модель фазового перехода из аморфного слабопроводящего состояния в кристаллического проводящее состояние. Исследование проводится методом квантовой молекулярной динамики. Проводится анализ структуры полупроводника, и оцениваются изменения проводимости и удельной теплоемкости при фазовом переходе. Найденный процесс моделируется с помощью модели сплошной среды с джоулевым источником нагрева. Изучается устойчивость решения нелинейной задачи теплопроводности в условиях, когда вольт-амперные характеристики находятся в области отрицательной дифференциальной проводимости. Обсуждаются особенности работы устройств памяти на фазовых переходах, связанные с неустойчивостью молекулярной системы.
Ключевые слова: фазовый переход, память на фазовых переходах, квантовая молекулярная динамика, нелинейная задача теплопроводности.
NANOSYSTEM ATOMIC STRUCTURE STABILITY DURING THE MOLECULAR SWITCHING
Popov Alexander Mikhailovich, professor, Doctor of physical and mathematical sciences, Faculty of Computational Mathematics and Cybernetics, Moscow State University named after M. V. Lomonosov, e-mail: [email protected]
Shumkin Georgiy Nikolaevich, doctoral student, PhD, Candidate of Physical and Mathematical Sciences, Faculty of Computational Mathematics and Cybernetics, Moscow State University named after M. V. Lomonosov, e-mail: [email protected]
Abstract: In the article the mathematical model of the phase transition from high-resistance amorphous state to low-resistance crystalline state is developed. The study is made, employing quantum molecular dynamics. The semiconductor structure is analyzed and changes of conductivity and heat capacity during the phase transition are estimated. Found process is simulated by the model of a continuous medium with Joule heating source. Stability of the nonlinear heat conduction problem solutions is studied, when the current-voltage characteristics are in the region of negative differential conductivity. The phase-change memory performance, associated with the molecular system instability, is discussed.
Index terms: phase transtion, phase-change memory, quantum molecular dynamics, nonlinear heat conductivity problem.
1. Введение.
При создании новых технологий в наноэлектронике используют эффект переключения фазового состояния молекулярной системы - фазового перехода между плохо проводящим и хорошо проводящим состоянием системы. Фазовые переходы в молекулярных системах используются для создания сверхмалых устройств наноразмеров в электронике [1], [2].
К таким системам относится, например, сверхмалая память (основная часть системы - проводящая точка диаметром 20 нм на тонкой пленке аморфного углерода) на фазовых переходах [3]. Теоретическое изучение фазовых переходов и переключений в таких молекулярных системах показывает возможность развития неустойчивостей тока с различными исходами [4], [5], [6].
Понимание работы таких систем существенно зависит от возможностей их моделирования. Система состоит из проводящей наноточки, которая является частью электрической цепи. Такой прибор является сложной нелинейной системой. Моделирование должно содержать несколько этапов описания на разных масштабах. В начале надо выяснить механизм происхождения фазового перехода в проводящей точке. Од-
ной из главных причин фазового перехода считается изменение молекулярной структуры под действием температуры. Такое моделирование для случая аморфного углерода проведено в [7] методом квантовой молекулярной динамики [8].
В случае аморфного углерода было показано, что под действием температурных эффектов молекулярная система переходит в графитовую структуру, обладающую высокой проводимостью. Начало теории фазовых переходов второго рода положено в [4]. Была предложена теория теплоемкости при фазовом переходе второго рода. Фазовые переходы второго рода называют также точками Кюри или точками. Переход может быть связан с изменением упорядоченности кристалла. Скачок теплоемкости при фазовом переходе моделируется в данной работе при определении точки фазового перехода алмаз - графит. В этом случае моделирование требует вычислений на суперкомпьютере. Многомасштабная модель переключения предложена в [20] и предполагает решение уравнения теплопроводности с джоулевым источником, меняющим форму в зависимости от развития процесса фазового перехода. Получающиеся электротермические уравнения содержат Вольт-Амперную характеристику (ВАХ), имеющую участок с
отрицательной дифференциальной проводимостью (S форма ВАХ). Этот режим может быть связан с развитием неустойчивости тока в проводящей точке.
Следующим шагом в процессе моделирования является определение реакции электронной составляющей системы на структурные изменения. В настоящей работе предлагается соответствующая модель временной квантовой редуцированной молекулярной динамики Эренфеста. Здесь определяется реакция электронной системы на внешнее электрическое поле.
Главной характеристикой системы является полученная Вольт-Амперная характеристика (ВАХ) прибора. Эта характеристика измеряется в эксперименте. В переключательных экспериментах на ВАХ имеется участок с отрицательной дифференциальной проводимостью. Этот режим связан с развитием неустойчивости тока в проводящей точке [6]. Моделирование должно учитывать развитие молекулярной структуры точки и пространственную форму электрической проводимости. Причины и результат неустойчивости могут быть различными. В настоящей работе изучается перегревная неустойчивость и электротепловая модель, параметры которой берутся из молекулярной модели - геометрическая форма проводящей точки, полученный в процессе молекулярной динамики структурный параметр, форма теплоемкости и положение ее особенности на температурной оси.
Процессы в наносистеме с отрицательной дифференциальной проводимостью могут приводить к новым эффектам -движению и колебаниям шнуров тока или доменов. От этих процессов может зависеть эффективность процесса молекулярного переключения.
Проведенные эксперименты по изучению полупроводников с отрицательной дифференциальной проводимостью [5] показывают, что все наблюдавшиеся S-образные ВАХ и токовые шнуры связаны либо с пробоем, либо с инжекцией из контактов. Отсутствует удовлетворительная теория пробоя с возникновением S-образной ВАХ и особенно теория шнурования тока при пробое и инжекции из контактов. Почти не исследовано влияние магнитного поля на явления неустойчивости в полупроводниках с отрицательной дифференциальной проводимостью. Следует ожидать, что при определенных условиях в полупроводнике могут возникать магнитные страты, аналогичные магнитным стратам в слабо ионизованной газовой плазме.
До настоящего времени не удалось наблюдать бегущие слоевые волны, домены слабого поля, а также (в двухдолинных полупроводниках) широкие домены с плоской вершиной.
Исследования, представленные в настоящей работе, направлены на выяснение механизма неустойчивости переключения сопротивления в аморфном углероде на основе вычислений из первых принципов с помощью квантовомеханической молекулярной динамики Кара-Парринелло, представленной в численном коде CPMD [18], [19].
Вычисления проводились на суперкомпьютере IBM Blue Gene/P.
2. Постановка задачи.
При построении моделей молекулярной динамики из «первых принципов» исходят из уравнения Шредингера и кулонов-ского взаимодействия заряженных частиц. Мы используем математическую модель молекулярной динамики Кара-Парринелло [8] на основе приближения Борна-Оппенгеймера. В приближении Борна-Оппенгеймера считается, что движение
массивных ядер может быть описано в классическом приближении. Силы, действующие на ядра со стороны электронов, рассчитываются на основе решения уравнения Шредингера для электронов в основном состоянии. Мы будем использовать квантовомеханическую формулировку, основанную на уравнениях функционала плотности («DFT» [21]). Уравнения модели Кар-Парринелло записываются в виде
M,R,(t) = -щEKS[{ipt UR,}],
= -^*EKS[m,{R,}] + Zj^j ipj (i)
где М1 и RI есть масса и координата ядра I, ß — фиктивная масса электрона, — орбитали Кона-Шэма и — энергия Кона-Шэма. Последний член в уравнении (1) есть голономные связи, которые учитываются требованием ортогональности орбита-лей. Это приводит к необходимости вычислять силы, действу-
д и
ющие на ионы — — EKS на каждом шаге моделирования, орбитали ——;EKS и множители Лагранжа Aij.
Энергия Кона-Шэма определяется выражением:
осс
EKS[{^i},{R,}]= ^ fif tf (r)[-2 V2Tpi(r)]dr + f Uext (r)p(r)dr i
+ ^Я^^ drdr' + E*c[p] + E,[{R,}1 (2)
где Uext есть ионный потенциал и
P(r)= Zifi №(r)\2 (3)
есть плотность заряда электронов, ft — число заселенности, Ехс
— обменно-корреляционная энергия и Е1 — электростатическое ион-ионное взаимодействие.
Используется разложение волновых функции по плоским волнам. Рассматривая представление первой зоны Бриллюэна, орбитали Кона-Шэма записываются в виде:
Vi (г) = -=LeCt (G)elGr, (4)
где П есть объем ячейки, G — вектор обратной решетки и ct (G)
— Фурье-коэффициенты орбитали i. В CPMD используется быстрое преобразование Фурье для преобразования величин между двумя пространствами.
3. Кластеризация существующих областей графита sp2 внутри алмазной матрицы sp3 как механизм переключения сопротивления.
Начальная конфигурация атомов углерода в супер-ячейке состоит из области графита sp2 и области алмаза sp3. Атомная конфигурация в супер-ячейке состоит из 192 атомов:
44 атома с 3 связями, обозначаемые как Nsp2 и представляющие графитовые области;
78 атомов с 4 связями, обозначаемые как Nsp3 и представляющие алмазные области;
70 атомов с двумя связями, обозначаемые как Nsp и представляющие линейные области.
Был проведен ряд расчетов молекулярной динамики при разных температурах. Во всех расчетах при поднятии температуры наблюдался эффект замещения алмазных областей графитовыми. Сравним установившиеся атомные структуры после процесса молекулярной динамики. Две конфигурации соответствующие расчетам молекулярной динамики приведены на рис. 1(а) и на рис. 1(б). Видно, что при высокой температуре атомная конфигурация представлена в основном графитовыми областями.
Рисунок 1. Геометрия молекулярной системы после процесса молекулярной динамики, проведенной при (а) Т=2000К и (б) Т=3000К.
Найдены два различных типа конечной атомной структуры. Первый тип структуры получен при температурах ниже Т=2500К. Он имеет сложную структуру, включающую смесь Бр2 и Бр3 подобных атомов. Во всех случаях повышения температуры число Бр2 подобных атомов возрастает, а число Бр3 подобных атомов убывает по сравнению с начальной конфигурацией. Особый тип конфигураций возникает при температурах больше Т>2500К. Постоянный рост числа Бр2 атомов приводит к формированию графитоподобной послойной структуры. В таких структурах число атомов Бр3 очень мало. Можно утверждать, что структуры второго типа обладают малым сопротивлением.
Можно сделать главный вывод о том, что изменения в кластерной структуре существуют ниже точки плавления аморфного углерода. Это изменение в кластерной структуре можно рассматривать фазовый переход в наноструктуре материала, приводящий к пороговой смене сопротивления.
Сценарий с нагревом приводит к порогу по температуре, выше которого возникают высокопроводящие области. Температура является основной причиной неустойчивости кластерной структуры и появления высокопроводящих кластеров графита внутри алмазной матрицы.
Джоулев нагрев в эксперименте приводит к возрастанию температуры до пороговых значений. Этот сценарий соответствует нижней части вольт-амперной S-образной характеристики, получаемой в эксперименте [3].
При больших температурах (больше 2500К) число атомов составляющих графитовые области возрастает, что приводит к падению сопротивления в конечном состоянии. Главный вывод заключается в том, что замещение графитовыми областями алмазных приводит к высокопроводящим атомным конфигурациям.
Важно получить противоположный сценарий перехода от жидкого состояния к твердому аморфному состоянию.
Переход от состояний с высокой проводимостью к состояниям с низкой проводимостью проводился мгновенным охлаждением пленки аморфного углерода, полученной при высокой температуре. Этот сценарий предложен в работе [13] для аморфного алмаза. Вычисления начинались с конфигурации, полученной при Т=3700К, и далее проводилась тер-мализация при температуре Т=5000К. Уже в этом процессе число атомов Бр2 графитовых областей сокращается, а число атомов Бр3 алмазных областей возрастает. Затем атомная структура (в состоянии жидкости) использовалась для проведения перехода жидкость-аморфное состояние. Жидкость мгновенно охлаждалась до температуры Т=300К и атомная конфигурация релаксировала к твердому аморфному состоя-
нию (с большим сопротивлением) в процессе молекулярной динамики. Этот основной эффект наблюдался как дальнейший рост числа атомов Бр3 алмазных областей, и в установившемся стационарном состоянии число таких атомов составляло 10 % от общего числа.
4. Структурные изменения в аморфном углероде. Структурная переменная S(t). Эффект стратификации области проводящей точки.
Мы изучали динамику атомов, вычисляя среднеквадратичное отклонение (MSD) для каждой атомной конфигурации следующим образом:
«а
<г2(г) >а= -1 V <\г1 (о-г;(0)|2 >,
где величина в скобках <> вычисляется как статистическое среднее по атомам типа . Вычисления проводились при различных температурах. Основываясь на кривых, мы вычисляли коэффициенты диффузии, используя соотношение Эйнштейна для процесса самодиффузии:
< \г£(0 - Г;(0)\2 >= 6£> где D есть коэффициент самодиффузии.
Некоторые параметры следующие: энергия связи между слоями есть 15.1 Дж/моль, расстояние между слоями равно 1=0.3354 нм, энергия связи в слое 167.6 Дж/моль (1118 С). Расстояние между атомами равно 0.1418 нм.
Введем координационную переменную, чья эволюция будет характеризовать фазовый переход. Для этого вычислим плотность атомов с использованием техники гистограмм. Пример плотности атомов атомной конфигурации при Т>2500К, аппроксимируемой гауссианами, показан на рис. 2. За координационную переменную выберем полуширину гауссианов, которые аппроксимируют плотность атомов.
6 7 г (А)
Рисунок 2. Гистограммы плотности ионов, аппроксимируемые гауссианами.
Эволюция координационной переменной показана на рис. 3. Видно, что координационная переменная от одних значений переходит к другим значениям, т.е. она характеризует фазовый переход. Осцилляции координационной переменной характеризуют процессы диффузии, вызванные тепловым движением атомов.
Time (ps)
Рисунок 3. Эволюция координационной переменной.
Рисунок 4. Диффузия атомов в различных направлениях.
Диффузия атомов в различных направлениях показана на рис. 4. Хорошо виден анизотропный характер процесса диффузии. Это обусловлено тем, что расстояние между графеновыми листами I больше, чем расстояние между атомами. Тепловое движение атомов диктуется сильными связями внутри слоя и большим свободным пространством между слоями.
Удельная теплоемкость в зависимости от температуры показана на рис. 5. Расчет удельной теплоемкости проводился по формуле:
< и2 > -< и >2 С' (Т)=-^-'
где и есть потенциальная энергия системы, скобки <> обозначают усреднение по времени за определенный промежуток времени. В формуле также используется среднее мгновенных значений температуры за определенный промежуток времени. В расчете температура растет ступеньками (используется термостат Нозе-Хувера), проводится термализация на каждом уровне, или ступеньке, температуры. Согласно теории Ландау, теплоемкость должна резко измениться при подходе к точкам фазового перехода при фазовом переходе 2-го рода. Удельная теплоемкость тесно связана с числом степеней свободы системы. Переход из аморфного состояния к графитовым слоям означает, что меняется число степеней свободы, что отражено двумя пиками теплоемкости при Т = 2500К и Т = 3600К. Первый пик соотносится с ситуацией, когда затравочные графитовые области захватывают алмазные области, но остаются атомы между графитовыми слоями. Второй пик соответствует ситуации, когда атомы, связывающие 2 графитовых слоя, разрывают связи с одним из слоев и устанавливают сильные связи с другим графитовым слоем.
Рисунок 5. Удельная теплоемкость в зависимости от температуры. 5. Электронная структура атомных конфигураций. Модель уровня II. Следующие вычисления были выполнены для проведения анализа электронной структуры эволюции в процессе структурного фазового перехода. Вычисленные возбужденные состояния могут дать эффекты фазового перехода, обусловленные электрическим полем. Это может быть сделано в рамках многомасштабной модели на основе молекулярной динамики Эренфеста [10]. Рассмотрим конструкцию редуцированной модели молекулярной динамики Эренфеста.
Из модели уровня I мы можем вычислить электронные орбитали (не только для основного состояния, но и для возбужденных состояний) и построим редуцированную модель молекулярной динамики Эренфеста на основе временной версии DFT [10]. Целью является получение возмущения электронной плотности из-за внешнего электрического потенциала и термически обусловленного переключения между HOMO орбита-лью и возбужденными орбиталями.
Разложим электронные волновые функции по полному набору орбиталей: V(r,R, t) = lmcm (Wm(r,R). (5)
Здесь ipm(r, R) удовлетворяет стационарному уравнению: Щт (r, R) = Em^m(r, R) (6)
Молекулярная динамика Эренфеста имеет дело с временными зависимостями коэффициентов cm(t) во время эволюции системы. Соответствующие уравнения движения запишутся с учетом введения координационной переменной S следующим образом: описывающие фиктивное движение координационной переменной S
d2S Vi г<-лидЕ™
dT2 = -L1 c-(t)2HsT-
m
(t)ck (t)(Em-Ek )dmk (S) (7) и уравнения движения электронов
= cm(t)Em(S) -iZkck (t)f dmk + V& ck(t), (8)
где г есть положения электронов и ядер. Уравнения включают следующие интегралы перекрытия электронных орбиталей:
атФ к (S) = iipm^dr. (9)
Допустим, что Гамильтониан имеет вид: H = Н0 + Vext (10)
где Н0 есть полный молекулярный стационарный Гамильтониан, и Vext есть зависящий от времени внешний потенциал.
Положим, что мы знаем положения ядер из расчетов атомной структуры при различных температурах (модель уровня I). Также рассчитаны орбитали основного и возбужденных состо-
яний. Целью является вычисление интегралов перекрытия и использование системы уравнений (8) для получения эволюции электронной системы. Интегралы перекрытия определяются следующим образом:
vmk = i^Prn (х,у,z)Vext(z)xpk(x,у,z)dxdydz (11)
Структурный параметр S(t) рассчитывается на основе статистики из CPMD вычислений на основе положений атомов при определенной температуре. В течение эволюции структура локализуется в окрестности графитовых слоев, приводя к фи-ламентации тока через образец и его локализации. Этот структурный параметр в поведении плотности тока учитывается в модели теплопроводности III.
Молекулярные орбитали для основного состояния (HOMO) и низшего незаполненного возбужденного состояния (LUMO) системы, полученные при T=600K, показаны на рис. 6: (а) HOMO, (б) LUMO.
Рисунок 6. Молекулярные орбитали системы (а) HOMO, (б) LUMO.
Заметим, что орбиталь HOMO распределена в основном по алмазным областям аморфного углерода, тогда как орбиталь LUMO распределена по графитовым затравочным областям. Это означает, что при возбуждении системы перекрытие ор-биталей должно привести к распространению затравочных графитовых областей по всей системе за счет замещения алмазных областей.
Модель теплопроводности III. Макроскопическая электротепловая транспортная модель.
Структурный фазовый переход в молекулярной структуре обусловлен температурным эффектом. Однако для модели молекулярной динамики температура является заданной термодинамической постоянной, как по времени, так и по пространству. Этот факт в макроскопическом смысле противоречит эксперименту, где фазовый переход осуществляется в проводящей точке, диаметром порядка 20 нм. Таким образом, учет границ проводящей точки и пространственная зависимость температуры могут оказывать существенное влияние на результаты моделирования. Мы предполагаем, что температура для расчетов молекулярной динамики должна диктоваться макроскопическим уровнем модели с учетом граничных условий на температуру. Для моделирования на макроскопическом уровне мы используем уравнение теплопроводности с источником Джоулева тепла, генерируемого плотностью тока j и напряженностью поля Е. Требуется согласование распределения температуры в модели теплопроводности и передача в модель молекулярной динамики. В то же время структурные параметры, определяющие электропроводность, должны быть переданы в макроскопическую модель теплопроводности. Эта процедура повторяется через определенные промежутки времени.
= 7±(К±7±Т) + + <12>
Закон Ома записывается в виде j = а(Е, Т, A(t)(E), (13)
где тензор проводимости (продольная и поперечная составляющие) а(Е, Т, A(t)) учитывает зависимость от температуры и структурного параметра, полученных при моделировании на микроскопическом уровне.
Для графитовых подобластей принимаются следующие транспортные коэффициенты: плотность = 2 • 103, теп-
m3
лоемкость С\—-—] = 700 (для алмаза С\—^~] = 500), тепло-
1кдк' ^ кд^К ''
проводность параллельно графитовым слоям W
к»
m-К.
= 278.11 - 2435
и перпендикулярно слоям К± = -йц. Базовые значения элек-
1
тропроводности принимаются равными: а»\-1 = 2
11 П-т
а± =8 - 102
106 и
Характерные значения поля и тока в эксперименте равны:
Ч-п J
= 108 , л4] =
50-10 А
= -- 1011
4
2.010-8т '-"-т2* 400 10-1ау
Времена релаксации имеют следующие характерные значения
Ц = 1Л2104 [Б Т± = 6.0[Б
Время задается в единицах t° = 0.28ps
Ко/(СоРо)
= 2.8 - 10-13s =
ч° = И°Е°)
L2° 5
к°т°
В безразмерных переменных уравнение теплопроводности принимает вид:
S=+^-Е)
(14)
где
Q = q°(J-E)= a(x, y, z, T)E2
и в приведенных расчетах полагается
-(г-Гр)2
Q = Чо T(x,y,z,t)e а22
Расчеты показывают наличие эффекта локализации тепла в определенной конечной области. В то же время зависимость электропроводности от температуры показывает существование границы по параметрам мощности и структурного параметра, при которых наступает температурная неустойчивость и неограниченный рост температуры в центре. Назовем этот эффект температурным пробоем или температурным срывом. Механизм неустойчивости определяется тем, что вкладываемая источником мощность становится больше, чем диссипация тепла за счет теплопроводности. При линейной зависимости проводимости от температуры порог срыва определяется как граница температурной неустойчивости. При более сложных зависимостях может наблюдаться нелинейный эффект локализации температуры, когда локализованное решение может существовать неограниченно долго.
Нормализация. Пороговые параметры для теплового пробоя.
Выберем Джоулев источник в виде
-( г-Гр)2
е ао
f (х, у, z, Т) = q0- Т (х, у, z, t)---(г-гр)2—
С е а2 dr
или
-(г-Гр)2 е ао
Q (х, у, z, Т) = Е0 а(х, у, z, Т) = Е$ ■ Т(х, у, z, t)--<r^—
S1 е ао dr
В таком виде источник тепла нормализован на интеграл от экспоненты, чтобы исследовать влияние формы локализованного источника.
Найденный порог после нормализации теплового пробоя равен q0 = 4.0
Уравнение включает электропроводность, пропорциональную температуре T(x,y,z,t), структурный параметр и величину приложенного напряжения. Структурный параметр идет из CPMD вычислений. Изменения в структуре происходят в атомной конфигурации и зависят от температуры. В принципе, член с источником включает вольт-амперную характеристику.
Уменьшение полуширины структурного параметра в течении подпорогового режима.
a(t) = afin + (aini - afin) ■ e-yt
Приближение к срыву по температуре хорошо видно по температурным профилям в радиальном направлении по времени, когда q0 = 50, показанным на рис. 7. В радиальном направлении диссипация тепла значительно уступает процессу росту температуры за счет локализации источника.
Рисунок 7. Температурные профили в радиальном направлении по времени, когда д 0 = 5 0. и а0 = 0.1 .
T(x,y,z,t)
0-0 O.î OA 0.6 0.8 J IЛ
Рисунок 8. Температурные профили в направлении z по времени, когда q0 = 5 0.
Температурные профили в направлении z по времени, когда q0 = 50, показанные на рис. 8, говорят о том, что процесс диссипации в направлении z превосходит процесс роста температуры в отличие от случая с радиальным направлением.
По теории Валеса [22] электрическая проводимость а в направлении, параллельном основной плоскости в кристалле графита, дается следующим выражением
16пе2 т(Т) ° = "Т2-ЬвТ log (2)
где т(Т) - время релаксации электронов проводимости, и L -константа решетки в направлении, перпендикулярном базовой плоскости. Из формулы Валеса хорошо виден линейный вклад температуры в температурную зависимость проводимости, что важно для модели теплопроводности. Оценка температурной зависимости т(Т) подлежит теоретической или вычислительной оценке.
6. Критический радиус проводящей точки. Тепловая неустойчивость.
Уравнение теплопроводности включает источник:
дТ(М, t)
—^^- = К ■ АТ(М, t)+ji ■ Т(М, t) Решение уравнения имеет следующую форму:
Т(М, t)= t ■ Фп(М)
п=1
где Лп и Фп(М) - собственные значения и собственные функции задачи.
Если один из показателей экспоненты будет положителен, то будет иметь место рост решения в экспоненциальной форме. Критерий роста можно выразить в форме: Р-КХ1 >0
Когда р фиксирован, тогда критический размер определяется первым собственным значением задачи Л1 < р
Критический диаметр проводящей точки для цилиндрической области будет иметь форму:
, =
сг = S ,
где у.0 = 2.4048.
Для выбранных параметров dcr — 100нм. Критерий может быть записан с использованием порогового значения напряжения. Если а = а ■Т, то
Sz W
Vth = ~т Vо I— d0
где Sz — толщина пленки, d0 — радиус точки, полученный на основе вычислении структурного параметра S(t) в молекулярной динамике.
T(x,y,z,t)
LL
Винтовая структура графитовых поверхностей в атомистическом моделировании.
По данным [3] экспериментальная Вольт-Амперная характеристика была дополнена до S-формы с использованием полинома третьей степени, что показано на рис. 9.
< Ё
S
0.0 0.5 1.0 ' 1
Voltage, V
Рисунок 9. Вольт-Амперная характеристика S-формы.
Одним из возможных объяснений появления S-образной характеристики является появление винтовых графитовых поверхностей, что показано на рис. 10. Атомы трех графитовых поверхностей, изображенных справа на рис. 10, пронумерованы числами от 1 до 3 так что, единички принадлежат только первой поверхности, двойки только второй, тройки только третьей поверхности. Соответственно, рисунок слева содержит ту же нумерацию атомов, только для начальной конфигурации. Хорошо видно, что графитовые поверхности из «2» из «3»(зеленая и синяя) искривлены винтовым образом в отличие от относительно плоской поверхности из «1»(красной).
Рисунок 10. Винтовая структура графитовых поверхностей.
Хорошо видны винтовые графитовые поверхности в той же атомной конфигурации, которая была показана предыдущем рисунке, на рис.11. Вместе с изображенными атомами на этом рисунке представлена поверхность уровня электронной плотности. Электронные плотности графитовых нижней и средней плоскостей должны отталкиваться, что препятствуют полному соприкосновению этих графитовых плоскостей.
Рисунок 11. Винтовая структура электронных поверхностей. Внутри поверхностей показана структура атомов в виде графитовых колец.
7. Заключение. Предложена многомасштабная модель наносистемы позволяющая изучать процессы переключения на разных уровнях детализации.
Главное в модели - это согласование нелинейных параметров, участвующих в процессе переключения. Сама нелинейная молекулярная система (проводящая тока) характеризуется координационными переменными - конфигурацией атомов, полученной при эволюции и являющейся главным параметром при фазовом переходе. Проведено квантовое молекулярно-динамическое моделирование из первых принципов кластерной структуры аморфного углерода, направленное на объяснение механизма смены сопротивления в проводящей наноточке в исследованиях памяти, основанной на фазовых переходах.
Квантовое моделирование процессов переключения проведено в работе с использованием суперкомпьютеров. На уровне квантовой молекулярной динамики найдены координационные параметры, характеризующие транспортные свойства, которые характеризуют изменение фазового состояния. Расчеты проведены при эволюции температуры в процессе расчета. Этот подход использован для нахождения скачка теплоемкости при фазовом переходе алмаз - графит. Тщательное рассмотрение поведения теплоемкости позволяет точнее определить температуру фазового перехода и изменение коэффициентов теплопроводности.
Использование электротепловой макроскопической модели позволило получить Вольт-Амперную характеристику проводящей наноточки. Джоулев нагрев приводит к локализации проводимости около графитовых плоскостей. Показано, что локализация проводимости и тока приводит к возрастанию температуры вплоть до теплового срыва. Найдены условия теплового срыва и начала шнурования тока.
Показано возникновение винтовой неустойчивости плотности тока и графитовых плоскостей за порогом условий неустойчивости.
Для понимания возможных эффектов изменения структуры, индуцированных электрическим полем, были получены и исследованы возбужденные орбитали проводящей графитовой точки внутри алмазной матрицы.
Список литературы:
1. Handbook of Theoretical and Computational Nanotechnology, Third Edition by Michael Rieth and Wolfram Schommers, Karlsruhe, Germany, 10-Volume Set, 200б, S000 p.
2. Nano and Molecular Electronics, Handbook, Edited by Sergey Edward Lyshevski, .- CRC Press, Taylor and Francis Group, (2007)
3. Abu Sebastian, Andrew Pauza, Christophe Rossel, Robert M Shelby, Arantxa Fraile Rodriguez, Haralampos Pozidis and Evangelos Eleftheriou, .New Journal of Physics (2011)
4. Л.Д. Ландау и Е.М. Лифшиц, Статистическая физика, том V, Изд. Наука, Москва, 56s c., 19б4.
5. А.Ф. Волков, UI.M. Коган, Физические явления в полупроводниках с отрицательной дифференциальной проводимостью, УФН, том 96, вып.4,1968.
6. Б.Б. Кадомцев, Коллективные явления в плазме, Изд. Наука, Москва, 304 c., 1988.
7. G.N. Shumkin, F.Zipoli, A.M. Popov, A.Curioni, Multiscale quantum simulation of resistance switching in amorphous carbon, Procedia Computer Science, 9, 2012, pp. 641-650.
S. R. Car and M. Parrinello, Unified approach for molecular dynamics and density-functional theory. Phys.Rev.Lett., l985,pp. 2471-2474.
9. В.И. Иванов-Омский, А.Б. Лодыгин, С.Г. Ястребов, Сканирующая туннельная спектроскопия аморфного углерода: модель туннелирова-ния. Письма в ЖТФ, том 25, вып. 24 , с. 66-71, 1999.
10. G.N. Shumkin, A.M. Popov, A.Curioni, T. Laino, A multiscale modelling of naphthalocyanine-based molecular switch, Procedia Computer Science 1 (1) pp. 185-192, 2010.
11. Standley B., Bao W., Zhang H., et al. Graphene-based atomic-scale switches, Nano Lett. , pp. 3345-3349, 2008.
12. Meijer G., Who wins the non-volatile memory race? Science , pp.1625-1626, 2008.
13. Marks N.A., McKenzie D.R., Pailthorpe B.A., Bernasconi M., Parrinello M., Ab initio simulations of tetrahedral amorphous carbon. Phys. Rev. B , N.14, pp.9703-9714, 1996.
14. Takai K, Oga M, Sato H, et. al. Structure and electronic properties of a nongraphitic disordered carbon system and its heat-treatment effects. Phys. Rev. B , 214202, 2003.
15. Yu He, Jinyu Zhang, Ximeng Guan, et. al. Molecular Dynamics Study of the Switching Mechanism of Carbon-Based Resistive Memory. IEEE
TRANSACTIONS ON ELECTRON DEVICES, Vol. 57, N. 12, pp.3434-3441, 2010.
16. Ivanov-Omskii V.l., Lodygin A.B., Yastrebov S.G., Chelnokov V.E., Journ. of Chemical Vapour Deposition. V.5, p. 198, 1997.
17. Averin D.V., Likharev K.K., Single-electron Tunnelling and Mesoscopic Devices. Ed. by H. Koch, H. Lubbig, Springer-Verlag, Berlin, p.3, 1992.
18. W.Andreoni and A.Curioni New Advances in Chemistry and Materials Science with CPMD and Parallel Computing Parallel Computing, , pp. 819842, 2000.
19. The CPMD consortium: http://www.cpmd.org.2010.
20. Popov, A. M., Shumkin, G. N., Nikishin, N. G., Multiscale simulation of thermal disruption in resistance switching process in amorphous carbon, Journal of Phys.: Conf. Series, Volume 640(1), 012027, 2015
21. W. Kohn, Density Functional and Density Matrix Method Scaling Linearly with the Number of Atoms, Phys.Rev.Lett., ,p. 3168, 1996.
22. P. R. Wallace, The Band Theory of Graphite. Phys. Rev. , p. 622, 1947.
