УДК 532.534
B.C. СУЗДАЛЬ, канд. техн. наук, Л.В. ДЕРБУНОВИЧд-р техн. наук, проф. каф. АУТС НТУ "ХПИ", Л.И. ГЕРАСИМЧУК, канд. техн. наук, Ю.М. ЕПИФАНОВ, канд. техн. наук, Институт сцинтилляционных материалов НАН Украины (г. Харьков)
МОДЕЛИРОВАНИЕ ПРОЦЕССА ВЫРАЩИВАНИЯ КРУПНОГАБАРИТНЫХ ЩЕЛОЧНО-ГАЛОИДНЫХ КРИСТАЛЛОВ МЕТОДОМ ЧОХРАЛЬСКОГО
Задача управління процесом вирощування монокристалів вирішена на основі використання тривимірних математичних моделей гідродинаміки бінарного розплаву, тепло- і масопереноса в рідкій і твердій фазах. Обґрунтована й отримана стохастична модель процесу кристалізації шляхом представлення рівнянь Навье-Стокса, теплопровідності і дифузії системою звичайних диференціальних рівнянь з випадковими коефіцієнтами і випадковою правою частиною. Для переходу до системи лінійних рівнянь використана підстановка Кармана, лінеаризація і поділ перемінних.
Process control of the single crystal growth was based on the three-dimensional mathematic modeling of the binary melt hydrodynamics and mass transfer in liquid and solid states. Stochastic model of the crystallization was worked out and discussed in terms of the Navier-Stokes equations with random coefficients and random right part of the equation. The Carman substitution, linearization and separation of the variables were used in transition to the system of the linear equations.
Постановка проблемы. Необходимость повышения качества управления технологическими процессами (ТП) в производстве щелочно-галоидных кристаллов (ЩГК) больших размеров связана с тем, что при значительном увеличении объемов ростовой камеры и выращиваемого образца применение известных физико-математических моделей ТП [1 - 5] недостаточно эффективно.
Качество управления производством крупногабаритных ЩГК может быть повышено с учетом присущих ТП кристаллизации как объекту управления (ОУ) свойств стохастичности, нестационарности и нелинейности. Эти свойства, в частности, определяются тем, что для рассматриваемого ОУ характерны: изменяющаяся зона нечувствительности по управляемому
параметру (диаметру кристалла), высокий уровень внутренних и внешних возмущений, продолжительный период реализации (7 - 12 суток). Игнорирование этих факторов в системе управления (СУ) может привести к выходу характеристик ТП из нормативно-технической области даже при незначительных вариациях технологических параметров на стадиях подготовки сырья, роста или охлаждения кристалла. В таких условиях применение алгоритмов, основанных на детерминированных (баланс масс), либо одно-двумерных моделях ОУ, не гарантирует воспроизводимость характеристик готовой продукции.
Для формализации и решения задачи управления ТП выращивания крупногабаритных ЩГК из расплава следует перейти к пространственным моделям физических процессов, с учетом неоднородных характеристик физических процессов, внутренних и внешних возмущений в ОУ, обладающих стохастическими свойствами.
Анализ литературы. Математическое моделирование процесса выращивания монокристаллов методом Чохральского основано на численном решении нестационарных уравнений Навье-Стокса совместно с уравнениями переноса тепла и массы в приближении Буссинеска [2, 3]. Математическая модель, представленная в [2], составлена с учетом осевой симметрии системы кристалл - расплав.
Численное решение уравнений физической системы кристалл - расплав осуществлено в [4], получена математическая модель для кристалла CsI(Tl), описывающая теплоперенос в цилиндрической части кристалла, а также распределение температурного поля в объеме кристалла и в плоскости Х02. Решение исследовано при различных значениях степени черноты поверхности образца, что позволило определить диапазон изменения параметров модели теплопереноса и учесть его при синтезе каналов СУ. Исследования движения расплава в зависимости от скорости вращения тигля и образца, размеров выращиваемого кристалла по отношению к размерам тигля, от соотношения высоты тигля к его ширине, числа Рейнольдса проведены в [5]. Однако, разработанные в [2 - 5] модели недостаточно учитывают особенности процесса выращивания ЩГК больших размеров, в частности, стохастический характер ОУ.
Цель статьи. Анализ существующих моделей кристаллизации из расплава и результатов математического моделирования позволяет сделать вывод об актуальности задачи разработки моделей процесса выращивания кристаллов из расплава в пространстве Я3 с учетом стохастических свойств системы кристалл - расплав, что и определило цель настоящей статьи.
Общий подход к моделированию процесса Чохральского. При
выращивании из расплава методом Чохральского растущий кристалл 2 в форме кругового цилиндра с радиусом Я вытягивают на затравочный кристалл 1 со скоростью V со свободной поверхности расплава 4 из тигля 5 (Н - уровень расплава в тигле, к - высота фронта кристаллизации (ФК)) и одновременно вращают его вокруг оси г с угловой скоростью О (рис. 1).
В процессе выращивания кристалла происходит перераспределение масс жидкой и твердой фазы в рабочем пространстве ростовой установки, изменение характеристик теплового поля, параметров тепломассопереноса, положения ФК и величины градиента температуры в области поверхности раздела кристалл - расплав.
При моделировании ТП выращивания монокристалла учитывают
Рис. 1. Схема ТП выращивания МК
параметры, характеризующие вращение образца и тигля, гравитационную, тепловую и концентрационную конвекции, а также поверхностные механизмы движения - термокапиллярную и капиллярно-концентрационную конвекции [2 - 4]. Для оценки влияния этих параметров на результат процесса кристаллизации, а также на устойчивость ТП анализируют систему уравнений, включающую модель мениска расплава [6]. При совместном решении уравнений непрерывности, диффузии и движения расплава (Навье-Стокса) используют граничное условие на свободной поверхности мениска расплава в виде уравнения капиллярности Лапласа в дифференциальной форме с учетом стабильности угла роста кристалла. Численное решение задач движения границ раздела областей кристалл - расплав и газ - расплав осуществляют с применением методов энтальпии и преобразования координат [7].
Для осесимметричных менисков приближенное решение уравнения Лапласа описывает образующую свободной поверхности мениска расплава и имеет наибольшую точность возле ФК. Погрешность приближенного решения зависит от радиуса вытягиваемого кристалла [5]. Наиболее часто применяют одномерное приближение, например, при анализе динамики ФК [8], при описании поглощения тепла на ФК [9].
Для численного решения уравнений, описывающих стационарные и нестационарные явления переноса в расплаве, применяют метод конечных элементов [6, 10]. Основные приближения связаны с предположениями о стационарности граничных условий, двухмерности ФК и свободной
поверхности расплава, а также об отсутствии течений, вызванных линейным перемещением кристалла.
Положение ФК, а также массообмен между расплавом и твердой фазой существенно зависят от тепловых условий выращивания кристаллов [8, 11]. Математическое моделирование теплопередачи в кристаллизационной установке включает в себя расчеты тепловых процессов в ростовой камере и выявление механизмов конвективного тепло- и массообмена в тигле. Основа анализа тепловых условий процесса кристаллизации - решение уравнения теплопроводности Лапласа для разных стадий выращивания с учетом теплофизических свойств материалов ростовой системы и изменяющихся граничных условий. Существенным фактором, влияющим на результаты выращивания из расплава крупногабаритных кристаллов, является степень стационарности тепловых процессов в объеме кристаллизационной камеры [12], что определяет необходимый переход к постановке и решению задачи пространственного моделирования ОУ.
Стохастическая модель ТП в пространстве Я3. Математическая модель процесса кристаллизации при выращивании МК методом Чохральского, представляющая собой систему нелинейных дифференциальных уравнений Навье-Стокса, переноса тепла и примеси, в приближении Буссинеска приведена в [4].
Для формирования математических моделей физических процессов кристаллизации при выращивании МК из расплава с учетом стохастичности параметров ТП следует перейти от исходной нелинейной системы дифференциальных уравнений в частных производных к системе линейных обыкновенных дифференциальных уравнений (ОДУ). С этой целью использована подстановка Кармана совместно с методами линеаризации и разделения переменных [13].
Входящие в состав модели уравнения Навье-Стокса для стационарного движения расплава в тигле в цилиндрической системе координат (г, ф, X) и уравнение непрерывности (в предположении осевой симметрии распределения компонент вектора скорости и давления в расплаве, когда все производные по координате ф равны нулю) имеют вид [4]:
дм V2 дм
и-------------Ъ м— =
дг г дх
1 др
дw дw
и--------Ъ м-------=
дг дх
— -----ъ V —-—\— ----------------\-— ;
р дх дг г дг дх
1 др д 2w 1 дм д2
ди и дм
— + - + — = 0, (2)
дг г дх
где и, V, м - радиальная, осевая и окружная компоненты вектора скорости в цилиндрической системе координат, р - плотность расплава, р - давление.
Система дифференциальных уравнений (1) - (2) приводится к
безразмерному виду подстановкой Кармана [13]: и=—-г-Н’(є),
^2 г2
v=— -гО(є), м = —2у1 — • V • Н(є), р = С —-— р - 2р • — • V • Р(є), где — - угловая
скорость вращения, 1/с; V - кинематическая вязкость расплава, м2/с; Н, О, Р -
_ , — дє —
безразмерные функции от є = х • Л1— ; — =./— .
V V дх V V
После стандартных преобразований [13] получим систему ОДУ:
й3Н ^ /йН\2 ^2 „ гтй2Н
= С + 1------I — О2 — 2 Н
йє3 І йє ) йє2
= 2 Г йНо — Н — I; (3)
йє V йє йє
йР й2Н „ гг йН
- + 2 Н-
йъ йъ йъ
Уравнение непрерывности (2) при подстановке Кармана удовлетворяется тождественно. Модель теплопроводности в системе кристалл - расплав в случае, когда движением расплава можно пренебречь [4]:
д0
— = %вД6, (4)
д/
где %е - константа, Д - оператор Лапласа, 0 = Т ^ - безразмерная величина
температуры, ДТ=Тт- Тк, ТК и ТТ - максимальная и минимальная величина
температур на поверхности кристалла и стенках тигля, соответственно.
Уравнение (4) преобразуем в систему из трех ОДУ, используя метод Фурье разделения переменных [13], что в случае осевой симметрии задачи дает: 0(г,х,/) = 01(г)02(х)0з(/). Пусть 0(г,х,/) = 0п(г)0ъ(/). Тогда из (4) следует:
= Хе<е,(П%2 + 0,0)-^ + е,<,,%^. (5)
ш дг г дг дх
Умножая полученное равенство на------------1-----, получим
х)0Ъ(/)
й03(ґ)
Л Х0 ґд 2012(г, х) 1 д012(г, х) д2012(г, х)^
( 12 +1). (6)
03(/) 012(г, х) д г г д г д х
Равенство двух функций от различных аргументов возможно тогда и только тогда, когда обе они равны одной и той же постоянной. Обозначая эту
2
постоянную Vе , переходим к двум уравнениям:
.. ,д 012(Гх) . 1 12(Гх) . д 012(г, х)ч 2д , , „• (7)
Х0(—--2—+---------2-----+—7-2—)^e0l2(г,х) = 0 ; (7)
д г2 г д г д х2
й03(О „2,
йґ
^203(0 = 0. (8)
Зависимость (7) является уравнением в частных производных, для разделения переменные также используем метод Фурье. Пусть
0п(г,х)=01(г)02(х). (9)
Подставляя в уравнение (7) выражение (9), приходим к равенству
Х0(02(х) й ^12(г) +02(х)1 Л01(г) + 01 (г) й 022(х)) - Ve0l(г)02(х) = 0, (10)
йг г йг йх
разделив которое почленно на произведение 0;(г)02(х) и перенеся вправо член, зависящий от х , находим:
й201(г) Л01(г) й202(х)
+ь-V 2 ==-,л~Лх- . (11)
01(г) г 01(г) 0 2 (х)
Выражение (11) определяет равенство двух функций от различных аргументов. Приравнивая обе части этого равенства постоянной ^, получаем два ОДУ и, добавляя к ним (8), окончательно - систему ОДУ:
й201 (г) 1 й01 (г) о . , ч Л
%0 ----7"2-----+Х0----- -----(^0 +V 0 )01(г) = 0;
йг2 г йг
Х0 + 402( х) = 0; (12)
йх
^ 203(0 = 0,
йґ
где У0, Х0 - константы, которые находят из граничных условий.
Для стационарного режима остаются только первые два уравнения системы (12):
й201(г) 1 й01(г) 2
Ь Д +Ьг-^ — ^^(г) = 0;
ь-
г йг
(13)
йх
+ Г2 02 (х) = 0.
Уравнение диффузии имеет вид (С-концентрация) [5]:
дС дС 1 дС дС _ 1 д
—+ V —+ V.------------------+ Ух —=--------(-
дґ дг ф г дф дх ЯеРГд дг
д2С 1 дС
2' + _^“ + -2 г дг
1 д 2С д2С
22 г 2 дф 2
Если движением расплава можно пренебречь, то
дС кг
1Г = Х ‘ ДС,
-) . (14)
(15)
где Хс - константа.
Используя метод Фурье разделения переменных, дифференциальное уравнение в частных производных (15) представим системой трех ОДУ, и в случае осевой симметрии задачи после преобразований, аналогичных рассмотренным выше, система ОДУ для стационарного режима (Хс -константа) получим:
1 СП — Х5С1(г) = 0;
X,
йг2
й 2С2( хК,2
г йг
(16)
+ Х2сС2( х) = 0.
При линеаризации нелинейной системы Навье-Стокса (3) представим уравнения, соответственно, третьего, второго и первого порядков, в виде нелинейной системы шести уравнений первого порядка. Введя обозначения:
сіє ёв
є - безразмерная величина,
Су2А ^ СуЧ г>
= У2Й, ~~г~ = Уза , О = У^ —7~ = У2я, Р = У^ и Учитывaя, что
—
є = х„1— .
преобразуем систему (3) и представим результат в векторной форме:
+
йУ = / (У, є)
йє
где У = У(є) - вектор с компонентами (У1Н, У2Н, У3Й, У^, У2g, Уз*, У1р), /(у, є) -векторная нелинейная функция с компонентами /1н, /2н, /н, /lg, ^ /з*, /1р).
Если у0 (є) - известное решение системы
ЛУп - _
-7° = / (У0, є), є0 < є < є1, йє
а у(є) представимо в виде
У(є) = Уо(є) + ~(є), У(єо) = Уо(єо) + ~(є0) , є0 < є <є1,
где у (є), у (є) - малые возмущения, то систему (17) можно линеаризовать, подставляя в нее у(є) и применяя разложение в ряд Тейлора с точностью до линейных членов включительно
й~1
1Н
йє
ЛУ2Н
йє
й~3Н
= У2н ;
= Узн;
йє
Лу1*
= 2 У0зн -у1Н + 2yoeнУeн 2yoeнyзн 2Усі^ Уlg ;
(18)
йє
йє й~1 р йє
= —2У02*у1Н + 2У01*У2Н ;
= 2 У02Н у1Н + 2 У01Н У 2Н + узн.
Результат математического моделирования физических процессов кристаллизации МК из расплава в стационарном случае позволяет перейти к стохастическим моделям: Используя представление процессов
гидродинамики, теплопроводности, диффузии в виде систем ОДУ со случайными коэффициентами и такой же правой частью, построим общую стохастическую модель физических процессов кристаллизации в стационарном случае:
йу ~ _
-------А(е; ю) у = а (ю);
йе
,4^0^) ^1 й0Дг) 20 , , , ч
Хе(ю)—+Х0И—-—^001(г) = ^1(ю); йг2 г йг
й 2 0 2(г) 2
Х0 (ю)--+ Х20 2(2) =л 2(ю); (19)
йг
й2 С- (г) ^ 1 йС- (г) 2 ^ , ч к , ч
х с (ю)-2—+х с (ю)— -----^сС1(г) = ^1(ю);
йг 2 г йг
й 2 С 2( г) 2
Хс (ю)--2--+ ХсС2 (2) = £, 2 (ю),
йг 2
где А(е; ю) = [./(у0(е); е;ю)] - случайная матрица Якоби размерности 6x6, а(ю) - случайный вектор-столбец размерности 1x6, у(е) = у0 (е) + у(е). Х0(ю), %(ю), ^2(ю) - случайные величины; юеП, (О, Е, Р) - вероятностное пространство; О - пространство элементарных событий, Е - ст-алгебра событий, Р - вероятностная мера на ст-алгебре событий %с(ю), ^(ю), §2(ю) -случайные величины.
Практическое использование стохастической модели (19) физических процессов кристаллизации из расплава в стационарном случае позволит более точно определить влияние различных факторов на результат кристаллизации и оценить вклад погрешности каждого технологического параметра в ошибку моделирования.
Выводы. Для повышения качества крупногабаритных ЩГК необходимо дальнейшее повышение точности управления процессами их получения. Это возможно путем разработки моделей реальных ОУ, обладающих свойствами нелинейности, нестационарности, несимметричности физических полей, с учетом высокого уровня производственных шумов.
В методе Чохральского форма и размеры растущего кристалла определяются капиллярными силами, формирующими мениск в зоне межфазной границы, и условиями тепло- и массообмена в системе кристалл -расплав. Исходная для моделирования система включает в себя уравнения движения участвующих в процессе сред, их непрерывности и является математической моделью неоднородного расплава в приближении Буссинеска.
С целью синтеза системы управления ТП выращивания МК на основе стохастических моделей этого ОУ осуществлен переход от исходной нелинейной системы дифференциальных уравнений в частных производных к системе линейных ОДУ с использованием подстановки Кармана. В результате получена общая математическая модель физических процессов
кристаллизации при выращивании МК для стационарного случая с учетом стохастических свойств ОУ.
Стохастический подход позволяет учитывать естественную неопределенность исходных данных о процессе кристаллизации, вследствие чего применение полученной модели эффективно для синтеза СУ промышленных установок типа "РОСТ" с заданными геометрическими размерами тигля, в частности, для решения задачи определения законов программно-логического управления параметрами ТП.
Список литературы: 1. Авдонин Н.А. Математическое описание процессов кристаллизации. -Рига: Знание, 1980. - 178 с. 2. Тевяшев А.Д., Суздаль В.С., Бородавко ЮМ., Пелипец А.А. Математические модели физических процессов при выращивании монокристаллов методом Чохральского // Радиоэлектроника и информатика. - 2001. - № 4 (17). - С. 33-43. 3. Crochet M.J., Wouters P.J., Geyling F.T., Jordan A.S. Finite-element simulation of Czochralski bulk flow // J. Cryst. Growth. - 1983. - 65, № 1-3. - Р. 153-165 4. Тевяшев А.Д., Суздаль В.С., Бородавко ЮМ., Пелипец А.А. Математическое моделирование и численный анализ физических процессов при выращивании крупногабаритных щелочногалоидных монокристаллов // Новые технологии. -2004. - № 1-2 (45). - С. 135-142. 5. Тевяшев А.Д., Суздаль В.С., Бородавко Ю.М., Пелипец А.А. Численное моделирование процессов тепло- и массообмена при выращивании монокристаллов методом Чохральского // Радиоэлектроника и информатика. - 2002. - № 1 (18). - С. 12-16. 6. Татарченко В.А. Устойчивый рост кристаллов. - М: Наука, 1988. - 310 с. 7. Crowley A.B. Mathematical modeling of heat flow in Czochralski crystal pulling // IMA J. Appl. Math. - 1983. -Vol. 30. - № 2. - Р. 173-189. 8. Простомолотов А.И., Сидельников СА., Черкасов А.В., Чернышенко О.В. Математическое моделирование тепловых процессов при выращивании монокристаллов: математические модели и их программная реализация. - М.: Ин-т пробл. мех. АН СССР, препр, 1990. - 44 с. 9. Kuroda Ekyo, Kozuka Hirotsugu. Influence of growth conditions on melt interface temperature oscillations in silicon Czochralski growth // J. Cryst. Growth. - 1983. - Vol. 63. -№ 2. - Р. 276-284. 10. Virzi A. Finite elements analysis of the thermal history for Czochralski growth of large diameter silicon single crystals // J. Cryst. Growth. - 1989. - Vol. 97. - № 1. - P. 152-161. 11. Бакирова О.И., Пелевин О.В., Соколов А.М. Численное исследование тепло- и массопереноса в системе расплав - кристалл при получении полупроводниковых материалов // Изв. АН СССР. Сер. физ. - 1983. - Т. 47. - № 2. - С. 334-337. 12. Virzi A. Finite elemente analysis of the thermal history for Czochralski growth of large diameter silicon single crystals // J. Cryst. Growth. - 1989. - Vol. 97. - № 1. - Р. 152-161. 13. Тихонов А.Н., Самарский А.А. Уравнения математической физики. - М.: Наука, 1966. - 724 с.
Поступила в редакцию 11.10.2005