ФИЗИКА ТВЕРДОГО ТЕЛА
УДК 548.4 + 519.24
МЕХАНИЗМ РОСТА ГРАНИ КРИСТАЛЛА С ЗОНАРНОЙ МИКРОСТРУКТУРОЙ
© 2008 г. С.С. Гуськов, М.А. Фаддеев, Е.В. Чупрунов
Нижегородский госуниверситет им. Н.И. Лобачевского [email protected]
Поступрла виедакцрю 27.05.2008
Предложен механизм образования зонарной микропериодической структуры монокристаллов с неоднородным распределением примеси. Для анализа условий возникновения зонарности использован алгоритм стохастического моделирования роста двухкомпонентного кристалла, позволяющий исследовать закономерности как нормального, так и послойного роста грани. В модели предусмотрена реализация слоисто-спирального дислокационного роста и послойного роста за счёт образования двумерных зародышей. Полученные в результате моделирования зависимости нормальной скорости роста грани от пересыщения аппроксимируются функциями, известными из теории и подтверждёнными экспериментами. Показано, что в зависимости от численных значений параметров модели возможно формирование как равномерного, так и квазипериодического распределения примеси, количественные характеристики которого согласуются с зонарной микроструктурой реальных кристаллов.
Ключевые слова: рост монокристаллов, двумерные дефекты, зонарная микроструктура, примесь, математическое моделирование.
Введение
Для многих природных и синтетических монокристаллов характерны дефекты, известные как зонарная (в минералогической литературе используется термин «зональная») микроструктура. Дефекты этого типа проявляются в виде чередующихся слоёв неодинакового качества, параллельных фронту кристаллизации. Характерный период зонарной структуры составляет единицы и десятки микрометров. Согласно современным представлениям [1, 2], различие физических характеристик отдельных слоёв обусловлено, в первую очередь, неравномерным распределением примеси, причем достаточно малой концентрации (~0.01%). Причины и механизмы неравномерного распределения примеси изучены слабо, общий метод описания зонарной структуры не разработан. Считается, что зонарная микроструктура образуется либо из-за квазипериодических изменений внешних условий роста кристалла (температуры, пересыщения) за счёт конвекционных потоков [1, 3], либо как результат внутренних процессов роста кристалла, т.е. путём самоорганизации [2, 4].
Зонарная микроструктура содержит ценную информацию о процессах роста кристалла и
изменениях в среде кристаллизации. Исследование зонарности позволяет оценивать качество ряда кристаллов, важных в лазерном приборостроении [5], и в надлежащем направлении совершенствовать технологию их выращивания. Наличие микрослоистой структуры может приводить к возникновению новых свойств кристаллов [6].
В данной работе предлагается механизм образования зонарной микроструктуры за счёт накопления примеси в среде перед фронтом кристаллизации и последующего встраивания примеси при превышении некоторого критического уровня её концентрации. Проведён анализ условий, приводящих к неравномерному распределению примеси.
Математическая модель роста кристалла с примесью
Разнообразие процессов, происходящих при росте кристаллов, обусловливает сложность уравнений, описывающих кинетику роста для разных значений параметров физической системы. Как следствие, широко используется исследование роста кристалла путем прямого статистического моделирования методами Монте-
Карло [7]. Это позволяет учитывать вероятностный характер присоединения и отрыва атомов в разных положениях на растущей грани кристалла.
Рассмотрим кристаллизацию как переход ростовых единиц (РЕ) - атомов, молекул и т.п. -из маточной среды в кристалл. В данной работе мы ограничимся моделированием роста грани (001) кристалла с простой кубической решёткой (кристалл Косселя). Пусть имеется два типа РЕ - основное вещество и примесь. РЕ обоих типов могут встраиваться только в узлы решётки, причем используется приближение ближайших соседей, т.е. рассматривается шесть соседей каждой РЕ.
Пусть узел на грани кристалла имеет координаты (х, у, z). Процесс роста кристалла представим как последовательность элементарных событий четырех типов: 1) присоединение РЕ
1-го типа в узел (х, у, z+1); 2) присоединение РЕ
2-го типа в тот же узел; 3) отрыв от поверхности кристалла РЕ, находящейся в узле (х, у, z); 4) пустое событие. Таким образом, множество элементарных событий включает 4 следующих элемента:
{Ж} = {Г1+,Г2+, Ж- Щ}. (1)
Вероятности элементарных событий определяются физико-химическими процессами на границе фаз. В [7] показано, что вероятность Щ12 перехода системы из состояния 1 в состояние 2 может выражаться через потенциальные энергии системы Е1 и Е2 в состояниях 1 и 2 соответственно следующим образом:
Ж12-----------------------------------1-, (2)
12 1 + ехр[( Е2 - Е1)/ кТ ]
где к - постоянная Больцмана; Т - абсолютная температура.
Пусть потенциальная энергия РЕ в кристалле Е5, а в маточной среде - Еь. Обозначив
АЕ = Еь - Е8, (3)
получим из (2) следующие выражения для вероятностей присоединения и отрыва:
1
Ж = С ------
;+ ; 1 + ехр(-АЕ / кТ)’
Ж. = Ж- =--------1------г,
1 + ехр(АЕ / кТ )
С + С2 = 1. (5)
Влиянием концентрации на вероятность отрыва Ж. пренебрегаем.
В нашей математической модели роста кристалла энергия данной РЕ в среде полагается постоянной, а энергия данной РЕ в кристалле определяется количеством и типом ближайших соседей:
2
(6)
ЕЬ = E0, Е5] £ 1ке Д :
(4)
где индекс ] принимает значения 1 и 2 для РЕ основного вещества и примеси соответственно.
В выражениях (4) учитывается, что вероятность присоединения при росте двухкомпонентного кристалла зависит от концентрации С данного типа РЕ в слое среды у поверхности кристалла. Концентрации РЕ двух типов связаны соотношением:
к=1
где 1к - число соседей типа к; е;к > 0 - энергия связи]-й и к-й РЕ.
Поведение реальной физической системы «кристалл - среда» зависит от отклонения этой системы от равновесия. Отклонение от равновесия принято характеризовать разностью химических потенциалов среды и кристалла Дц [1]. Химический потенциал данной фазы равен работе, которую надо совершить, чтобы изменить число частиц в этой фазе на единицу. Химический потенциал цг, где индекс г принимает значения 5 для кристалла и Ь для маточной среды, можно записать в виде
ц, = Е, - ТБ, + РО,. (7)
Здесь Б, - энтропия, Р - давление, Ц - удельный объем фазы в расчете на одну частицу. Количественные оценки показывают, что в (7) первое (энергетическое) слагаемое заметно больше остальных [1]. Таким образом, в некотором приближении можно считать, что ц, = Ег. Как известно, химический потенциал кристалла равен химическому потенциалу частицы в изломе (в положении полукристалла) цу2 [1]. Число ближайших соседей частицы в положении полукристалла равно половине числа таких же соседей в объеме кристалла. Как следствие, при условии С2 << С (т.е. в случае малой концентрации примеси) энергия частицы в изломе для рассматриваемой грани (001) примитивной кубической решетки Е1/2 = - 3ец. Тогда разность химических потенциалов может быть записана в следующем виде:
ДЦ = ЦЬ - = ЦЬ - Ц1/2 = ЕЬ - Е1/2 = Еь + 3е1Ь (8)
Следовательно,
Еь = Дц - 3е„. (9)
Путем несложных преобразований выражений (4) с учетом (3), (6) и (9) получим
С,.
Ж =-
Ж}+
1 + f (Ац)ехр[(-£ 1к е }к)/кТ ]
к=1
Ж =
1 + (I(Ац)) 'ехр[(£ 1ке]к)/кТ]
к=1
(10)
1
где введено обозначение
I (Ац) = ехр[(3е11 - Дц)/кТ ]. (11)
Связь энергии РЕ в среде и Дц (9) обеспечивает соответствие разработанной модели с принципом детального равновесия: при нулевом отклонении системы «кристалл - среда» от равновесия кристалл расти не должен. Это эквивалентно выполнению условия Ж1+(Дц = 0) = = Ж_(Дц = 0). Легко видеть, что выражения (10) удовлетворяют данному требованию. При возрастании Дц значение Еь увеличивается, что соответствует увеличению скорости роста кристалла, а при Дц < 0 кристалл постепенно разрушается (испаряется, растворяется, расплавляется).
Отклонение физической системы «кристалл - среда» от равновесия определяется в случае роста из паровой фазы давлением пара, в случае роста из раствора - концентрацией кристаллизующегося вещества в растворе, при росте из расплава - температурой расплава [1, 5]. Например, для пара - идеального газа и идеального (разбавленного) раствора выполняется соотношение
Дц = кТ !п(о +1).
Рассмотрим так называемый приповерхностный слой (ПС) - слой маточной среды у поверхности кристалла. Перед началом расчетов задается средняя концентрация примеси в среде С20. После каждого МКШ производится пересчет концентрации основного вещества и примеси в ПС.
Пусть ПС содержит фиксированное число п0 РЕ, причем из них п типа у. Тогда относительная концентрация примеси определяется соотношением:
Г = Пк
Г2 = .
П0
(15)
Пусть за время Дt из ПС к кристаллу присоединилось ко РЕ, причем среди них к РЕ типа у. Присоединившиеся РЕ покинули ПС, на смену им из среды кристаллизации в ПС поступило к0Г20 РЕ примеси и к0(1-Г20) РЕ основного вещества. Таким образом, новая концентрация в ПС определяется через старую следующим рекурсивным соотношением:
Г
(12)
В зависимости от типа среды параметр о представляет собой относительное пересыщение пара (Р - Р0)/Р0 или относительное пересыщение раствора (С - С0)/С0 , где Р0 и С0 - равновесные значения давления пара и концентрации раствора соответственно. Для переохлажденной жидкости (расплава)
Дц = ДЯо, (13)
где ДН - удельная теплота плавления в расчете на одну частицу, параметр о = (Т0 - Т)/Т0 - относительное переохлаждение расплава, Т0 -температура плавления.
Вероятность пустого события Ж0, в соответствие с (1), равна
Ж0 = 1 - (Ж1+ + Ж2+ + Ж_). (14)
Базовый алгоритм моделирования заключается в повторении следующих операций: 1) выбор узла на поверхности кристалла; 2) расчет нормированных вероятностей возможных событий (10) и (14); 3) выбор элементарного события в зависимости от значения генерируемого псевдослучайного числа; 4) выполнение выбранного события.
При использовании численного стохастического моделирования время t процессов принято измерять в шагах алгоритма Монте-Карло (МКШ). Один МКШ реализуется при переборе всех узлов решетки на поверхности рассматриваемого участка грани. Расчет проводится для заранее заданного числа МКШ.
(і+1) = П2 к2 + к0Г20
1 П0
кГ - к
= Г (і) + к0Г20 к2
Г(0) = Г Г2 = Г20 5
где і - порядковый номер МКШ.
(16)
Реализация различных механизмов роста грани кристалла
Под механизмами роста кристаллов понимают способы присоединения частиц к кристаллу. Различают два основных механизма [1]: 1) механизм нормального роста; 2) механизм послойного роста (дислокационный или на двумерных зародышах). Нормальный рост грани происходит на несингулярных поверхностях, при этом частицы присоединяются к поверхности практически повсеместно. Данный механизм обычно проявляется на стадии регенерации граней или в случае равновесного состояния шероховатой поверхности грани растущего кристалла. Послойный рост представляет собой образование гладких граней путем последовательного отложения слоев вещества, при этом ростовые ступени распространяются от центров роста, в качестве которых могут выступать выходы на поверхность грани винтовых дислокаций или двумерные зародыши.
В рамках нашей модели в зависимости от значений безразмерных параметров ец/кТ и Дц/кТ могут осуществляться ярко выраженные
п
0
Н, отн. ед.
1
1
! 1
/Т
1 1
1
* / /
1 У
1 I'
1
1 '
0.05
0.10
0.15
Рис. 1. Зависимости нормальной скорости роста грани от пересыщения для разных вариантов модели: 1 - дислокационный рост, 2 - нормальный рост шероховатой грани, 3 - рост путём образования двумерных зародышей. Точки - результаты моделирования, линии - графики аппроксимирующих функций. Видно, что в случае дислокационного роста зависимость R(ст) при малых значениях ст аппроксимируется квадратичной функцией, а при больших - линейной. Переходная область между этими функциями показана штрихпунктирной линией
0
послойный и нормальный механизмы роста грани. При малых значениях параметра еп^Т наблюдается нормальный рост грани кристалла, при больших - послойный рост. Характерное граничное значение параметра ец^Т лежит в диапазоне от 1 до 5 в зависимости от величины отклонения от равновесия Дц/^Т. Переход от послойного к нормальному росту происходит при уменьшении значения е\\/¥Т или увеличении Дц/^Т. Путем выбора граничных условий в области расчета возможна реализация послойного роста путем образования двумерных зародышей, путем движения ступени или за счет выхода на поверхность грани одной или нескольких винтовых дислокаций. В последнем случае наблюдается слоисто-спиральный рост.
Значения энергий связи РЕ можно оценить по теплоте испарения (плавления) кристалла [8, 9]. Для молекулярных (ван-дер-ваальсовых) кристаллов эта энергия связи сравнительно мала, порядка 1 ккал/моль, ионные и ковалентные кристаллы обладают наибольшей энергией связи - около 200 ккал/моль, в металлах энергия связи порядка 50 ккал/моль. Следовательно, с учетом различия температур фазовых переходов, значение безразмерного параметра еп^Т ~ 10 для роста ионных кристаллов из раствора и е\\/¥Т ~ 1 для кристаллизации металлов из расплава. Таким образом, результаты моделирования согласуются с известными экспериментальными данными: рост металлических и полупроводниковых кристаллов из расплава,
характеризующийся малыми значениями энергии ец и высокими температурами, происходит в основном по нормальному механизму, а ионные кристаллы из раствора при низких температурах растут слоями [1].
Полученные с использованием для Дц выражения (12) зависимости нормальной скорости роста грани R от пересыщения о могут быть аппроксимированы при нормальном механизме роста линейной функцией при любых значениях параметра о; при дислокационном механизме роста - параболической функцией при низких значениях о и линейной при высоких о; при образовании двумерных зародышей - экспоненциальной функцией (рис. 1). При образовании двумерных зародышей наблюдается характерное критическое значение пересыщения о0, такое, что при о < о0 зародыши не образуются и скорость роста равна нулю. Такие функциональные зависимости хорошо согласуются с теоретическими представлениями о связи скорости роста с отклонением от равновесия [1, 5, 10] и результатами работ [11, 12] по моделированию роста кристаллов методами Монте-Карло.
Распределение примеси в растущем кристалле
Систематическое моделирование роста двухкомпонентного кристалла показало, что возможна реализация четырех вариантов распределения примеси. Вариант 1 - равномерное
0.05
0.04
0.03
0.02
0.01
(
0.25
0.20
0.15
0.10
0.05
С
0
500
1000
1500
а
0
500
1000
1500
Рис. 2. Распределение примеси: а - равномерное; б - квазипериодическое. Слева - распределение примеси в кристалле, а0 - параметр кубической элементарной ячейки. Справа - зависимость концентрации примеси в ПС от времени, которое измеряется в шагах алгоритма Монте-Карло (МКШ)
распределение примеси в кристалле без накопления и без дефицита примеси в среде у фронта кристаллизации (рис. 2а). Варианты 2 и 3 - это равномерное распределение примеси в кристалле с накоплением и с дефицитом примеси соответственно в ПС среды. Принципиальное отличие варианта 4 заключается в квазипериодиче-ском распределении примеси в кристалле с соответствующими осцилляциями концентрации примеси в среде (рис. 2б).
Реализация того или иного варианта распределения примеси определяется значениями энергии связи «примесь - основное вещество» є12 и энергии связи «примесь - примесь» є22. Для наглядности результаты моделирования представлены в виде энергетических диаграмм (рис. 3). Каждой точке диаграммы соответствует определенный вариант распределения примеси. В ходе расчетов выяснилось, что четырем рассмотренным выше вариантам распределения примеси соответствуют четыре односвязные области на энергетических диаграммах (рис. 3). Размер и конфигурация данных областей довольно чувствительны к изменению параметра Єц^Т, отклонения от равновесия о и средней концентрации примеси в маточной среде С20.
На рис. 3 видно, что характер распределения примеси существенно зависит от соотношения безразмерных величин е^ец и е22/ец при постоянных значениях параметров ецЛТ, о и С20. При е12 > ец, независимо от значения е22, присоединение к кристаллу частиц примеси энергетически более выгодно, чем присоединение частиц основного вещества, поэтому имеет место дефицит примеси в ПС и равномерное распределение примеси в кристалле (вариант 3). Если значения е12 и ец близки, т.е. е^/ец ~ 1, то частицы примеси и основного вещества можно считать равноценными. Естественно, в этом случае примесь равномерно распределена в кристалле, а её концентрация в ПС близка к средней концентрации примеси в маточной среде С20 (вариант 1). При е12 < ец примесь накапливается в ПС и на поведение ростовой системы, в отличие от случаев, рассмотренных выше, начинает оказывать влияние значение е22. Когда е22 существенно (на порядок) меньше, чем ец, концентрация примеси в ПС возрастает до такого значения, что поток примеси из маточной среды в ПС становится равным потоку примеси из ПС в кристалл (вариант 2). Когда е22 значительно (на порядок) больше, чем е1Ь то первые
Єц/АТ = 2 Єц/АтТ = 4
Рис. 3. Области с разным характером распределения примеси при о = 0.2, С20 = 0.01 и различных значениях є11/АТ. На осях отложены значения отношений энергий є12/є11 и є22/є11. Цифровые обозначения областей соответствуют номерам вариантов распределения примеси
же присоединившиеся к кристаллу частицы примеси способствуют переходу примеси из ПС в кристалл и уменьшению С2 до С20 (вариант 1). Самое интересное поведение ростовой системы наблюдается, когда е22 и еп - величины одного порядка (вариант 4). В этом случае примесь в среде у поверхности растущей грани накапливается и при превышении некоторого критического уровня концентрации С0 присоединяется к кристаллу, в результате чего образуется неравномерное (слоистое) распределение примеси в кристалле. Величина С0 зависит от соотношения энергий взаимодействия РЕ разных типов.
В результате моделирования выявлены следующие закономерности поведения параметров слоистого распределения примеси. Во-первых, при увеличении параметра о уменьшается пространственный период слоистой микроструктуры и толщина легированных примесью слоев. Это связано с тем, что при повышении о частицы примеси активнее адсорбируются на поверхности кристалла, тем самым способствуя дальнейшему переходу примеси из ПС в кристалл. Во-вторых, снижение концентрации примеси в среде С20 ведет к уменьшению содержания примеси в слоях и к увеличению их пространственного периода, т.к. примесь медленнее накапливается в ПС. Как следствие, пока концентрация примеси увеличивается до критического значения С0, успевает вырасти более толстый чистый слой.
В зависимости от значений энергии связи е12 и е22 возможно получение грубой зонарности с
большим периодом или тонкой зонарности с меньшим периодом и менее выраженными максимумами концентрации примеси. Таким образом, при изменении параметров физической системы возможно возникновение квазиперио-дического распределения примеси с разными количественными характеристиками. Например, при С20 = 0.02% период слоев с повышенным содержанием примеси равен 104а0 (где а0 -параметр кубической элементарной ячейки), что составляет 10 мкм, если положить а0 = 10 А. Средний период зонарной структуры кристаллов КОР, выращенных в присутствии 0.03% А1(К03)3-9Н20, по результатам измерений методами оптической микроскопии в сочетании с травлением, составляет около 8 мкм [13]. Таким образом, результаты моделирования хорошо согласуются с экспериментальными данными. Естественно, в случае модели косселевского кристалла мы можем говорить лишь о совпадениях по порядку величины, т.к. количественные детали заметно зависят от кристаллографической симметрии и дальности взаимодействий между атомами.
Заключение
Разработанная модель роста двухкомпонентного кристалла позволяет наблюдать и исследовать закономерности как нормального, так и послойного роста грани кристалла. Поведение модели при изменении параметров совпадает с наблюдаемым экспериментально поведением
реальных кристаллов и хорошо укладывается в рамки современных теоретических представлений. Предложен механизм, объясняющий формирование квазипериодической зонарной микроструктуры кристаллов исходя из предположения о возможности накопления примеси перед фронтом кристаллизации. При этом образуется кристалл, в котором высокосовершенные слои основного вещества разделяются более тонкими прослойками с повышенным содержанием примеси. Показана возможность применения модели для исследования условий возникновения зонарности и расчета характеристик зонарного распределения примеси.
Список литературы
1. Современная кристаллография. Т.З. Образование кристаллов / Под ред. А.А. Чернова. М.: Наука, 1980. 408 с.
2. Smolsky I.L., Voloshin A.E., Zaitseva N.P. et al. // Philosophical Transactions of the Royal Society of London. 1999. V. 357. P. 2631-2649.
3. Holten T., Jamtveit B., Meakin P. et al. // American Mineralogist. 1997. V. 82. P. 596-606.
4. Merino E., Wang Y. // Non-Equilibrium Processes and Dissipative Structures in Geoscience, Yearbook Self-Organization. 2001. V. 11. P. 13-45.
5. Портнов В.Н., Чупрунов Е.В. Возникновение и рост кристаллов. М.: Физматлит, 2006. 328 с.
6. Зайцева Е.В., Портнов В.Н., Фаддеев М.А., Чупрунов Е.В. // Кристаллография. 1997. Т. 42. № 6. С. 969-971.
7. Методы Монте-Карло в статистической физике / Под ред. К. Биндера. М.: Мир, 1982. 400 с.
8. Киттель Ч. Элементарная физика твердого тела. М.: Наука, 1965. 368 с.
9. Лейбфрид Г. Микроскопическая теория механических и тепловых свойств кристаллов. М.: Физматлит, 1963. 312 с.
10. Элементарные процессы роста кристаллов / Под ред. Г.Г. Леммлейна и А.А. Чернова. М.: Издательство иностранной литературы, 1959. 300 с.
11. Александров Л.Н. Моделирование роста и легирования полупроводниковых пленок методом Монте-Карло. Новосибирск: Наука, 1991. 125 с.
12. Gilmer G.H. // Journal of Crystal Growth. 1977. V. 42. P. 3-10.
13. Ким Е.Л. Дис... канд. физ.-мат. наук. Н. Новгород: ННГУ, 2003. 133 с.
THE MECHANISM OF GROWTH OF THE CRYSTAL FACE WITH AN OSCILLATORY ZONING PATTERN
S. S. Guskov, M.A. Faddeev, E. V. Chuprunov
A mechanism for forming monocrystal oscillatory zoning patterns has been proposed. A stochastic modeling algorithm for two-component crystal growth has been used to analyze the conditions of oscillatory zoning origin. The algorithm allows studying the regularities of normal and level-by-level growth of crystal face. The model provides the realization of layer-spiral dislocation and level-by-level growth mechanisms due to two-dimensional nuclei formation. The model dependences of face normal growth rate on supersaturation are approximated by theoretical functions confirmed by experiments. The model simulation has shown a possibility to form both uniform and quasi-periodic impurity distributions whose parameters agree with oscillatory zoning patterns of real chips.