Full Text

1. Анализ современного состояния описания тонких плёнок применительно к индукционной плавке 1.1. Тонкие плёнки при индукционной плавке В настоящее время технология индукционной плавки широко используется в производственных процессах машиностроения, например в точном литье деталей газотурбинных двигателей. Индукционная плавка заключается в нагреве металла вихревыми токами, которые возбуждаются в нём переменным магнитным полем (далее ПМП) индуктора (водоохлаждаемая медная катушка с переменным током). Кроме того, магнитное поле, взаимодействуя с вихревыми токами, создает объемную силу Лоренца, которая перемешивает расплав. Схема теплового узла индукционной печи показана на рис. 1, геометрические и физические параметры материалов, взятые на основе работ [1-3] приведены в табл. 1. Следует отметить, что в настоящей работе свойства материалов считаются постоянными, так как зависимость физических свойств от температуры следует учитывать при её значительных перепадах внутри объекта или если температура сильно изменяется со временем; в нашем случае перепады температуры невелики, а процессы стационарны, т.е. можно использовать уравнения с постоянными свойствами, взятыми при соответствующей температуре (1500 °С). Технологический процесс происходит в техническом вакууме, однако на поверхности расплава при температурах более 1500 °С может образоваться тонкая (до 0,5 мкм) оксидная пленка, которая может привести к негативным металлургическим эффектам: попадание плёнки в отливку приводит к образованию дефекта, который в металлургии называется «пленой». Кроме того, имея большую, чем чистый металл, излучательную способность, пленка может исказить показания яркостного пирометра. Наиболее склонны к образованию оксидных плёнок сплавы на основе никеля с добавлением алюминия или хрома, а доля дефектов, связанных с оксидными плёнками, для некоторых изделий достигает 10 %. Для удаления оксидной пленки с поверхности расплава применяют «разгон плены» - технологическую процедуру, заключающуюся в интенсификации течений металла, которые разрывают пленку и выталкивает её на стенку тигля. Кусочки плёнки зацепляются за стенку тигля и остаются на ней, затем расплав металла переливают через них форму. Интенсификация потока достигается увеличением тока индуктора от 250 до 300 А. Одновременно с интенсификацией перемешивания расплава увеличивается и тепловыделение, в ходе разгона расплав дополнительно нагревается до 80…100 °С, такой регулярный перегрев снижает срок службы тигля. Рис. 1. Принципиальная схема индукционной печи (а), формализованная схема с указанием границ, соответствующих граничным условиям для расплава (б) и плёнки (в), относительное положение расплава и индуктора (г): при «b - b» совпадают донные (bottoms) части расплава и идеализированного индуктора, при «c - c» - центры (centers) симметрии. Fig. 1. Schematic diagram of induction furnace (a), formalized diagram showing boundaries corresponding to boundary conditions for melt (б) and film (в), relative position of inductor melt (г): "b - b" the bottom parts of the melt and the idealized inductor coincide and "c - c" the centers of symmetry coincide Таблица 1 / Table 1 Геометрические и физические параметры при индукционной плавке Geometric and physical parameters during induction melting Параметр Символ Значение Технологический процесс Внутренний радиус тигля (масштаб длины) Rc 0,1 м Высота цилиндра, занятого металлом Z 0,2 м Высота индуктора Zind 0,35 м Ток индуктора Ic 250 A Число витков в индукторе N 10 Температура расплава Tmelt 1773 K Температура вакуумной камеры Tcham 323 K Свойства расплава Температуропроводность a 1,92×10-6 м2/с Коэффициент температурного расширения b 38,5×10-6 1/K Удельная электропроводность z 0,85 См/м Теплопроводность l 10 Вт/(м×K) Кинематическая вязкость n 5,28×10-7 м2/с Плотность r 7900 кг/м3 Температурный коэффициент поверхностного натяжения aТ 1,5×10-4 Н/(м×K) Излучательные способности расплава / плёнки m 0,1 / 0,2 Некоторые характеристические величины Коэффициент излучения C 3,37×10-10 Вт/(м2×K4) Характерное магнитное поле Hc 2000 A/м Характерная плотность потока излучения qc 3328 Вт/м2 Исходя из вышесказанного, исследование условий формирования и стабильности оксидной плёнки на поверхности расплава при индукционной плавке актуально. Закрытая вакуумная камера, высокая температура, дороговизна материалов и неконтролируемое вероятностное появления плёнок диктуют метод исследования - математическое моделирование. В рассматриваемом случае всю систему можно описывать в осесимметричной постановке, это допустимое приближение, поскольку геометрия тигля - это цилиндр, соосный с катушкой индуктора, магнитное поле которого имеет с хорошей точностью осевую симметрию (отличие наблюдается только вблизи токоподводящих проводов). Этот подход является весьма общим [4-6], однако существуют работы в которых процессы рассматриваются в полностью 3D-постановке [7-9]. В работе [10] показано, что поле скоростей за несколько секунд (в указанной работе - за 6 с) подстраивается под распределение силы Лоренца, определяемое переменным магнитным полем. Поскольку в процессе плавки расплав металла выдерживается до 5 мин для полного плавления всей шихты, распределение скорости у поверхности является установившимся, следовательно, задачу можно решать в стационарной постановке: , т.е. течения, тепловые потоки и состояние упругого напряжения не меняются со временем; процессы, связанные с пластической деформацией, текучестью плёнки и её разрушением не моделируются. Поверхность расплава можно считать плоской и недеформируемой в поперечном сечении, поскольку магнитное поле недостаточно сильно для образования выпуклого мениска, как в технологии плавки с холодным тиглем [11]. Кроме того, на данном этапе не рассматриваются возможные эффекты генерации волн плёнкой, рассмотренные, например, в работах [12-14]. При этом соответственно не учитывается и влияние волн на механическое равновесие плёнки. 1.2. Подходы к описанию тонких плёнок Существуют различные подходы к описанию тонких плёнок, при выборе между ними авторы руководствовались определёнными соображениями. Описание состояния пленки на микроуровне методом молекулярной динамики, как это сделано в [15; 16], неприменимо из-за большой вычислительной ёмкости для рассматриваемого объёма вещества. Модель оксидной плёнки как молекулярного монослоя, упругие свойства которого определяются энергией Гиббса [17], эффективна для медленного потока, где скорости движения жидкости и переноса молекул за счет поверхностной диффузии сопоставимы [18]. Однако наблюдаемая толщина пленки 0,5 мкм не позволяет рассматривать её как монослой, а скорость переноса вещества за счет вынужденной конвекции, значительно превышающая поверхностную диффузию, делает этот подход неприменимым. Существуют работы, в которых тонкая пленка рассматривается как жидкость Максвелла [19]. Учет взаимодействия структурных наночастиц позволяет моделировать упругие и пластические свойства пленки: математически это выражается дополнительными членами в тензоре напряжений и позволяет исследовать устойчивость поверхности в электрическом [20] и магнитном [21] полях. В приближении твердой плёнки из математического описания исключается уравнение движения материи в самой плёнке, однако её механическое равновесие по-прежнему удобно вести в терминах тензора напряжений, как в [22] для магнитной мембраны в переменном магнитном поле. Одним из эффективных и хорошо проработанных методов расчёта напряжённо-деформированного состояния композиционных материалов, содержащих в своей структуре отличные по механическим свойствам тонкие слои или волокна, является градиентная теория упругости, построенная на вариационных методах расчёта. Эта теория обобщает классический закон Гука, включает гипотезу Дюамеля - Неймана о термических напряжениях [23] и уточняет классические представления о напряжённо-деформированных состояниях сложных по структуре материалов, учитывая масштабные эффекты. Так, в работе [24] показаны деформации и напряжения, возникающие на границе контакта двух сред с различающимися на порядок коэффициентами Пуассона. Интересной особенностью указанной теории является возможность нахождения масштабных коэффициентов из моделирования поведения микроструктуры вышеупомянутыми методами молекулярной динамики [25; 26]. Однако при всей привлекательности градиентной теории упругости существуют некоторые сложности её применения к поставленной задаче - неизвестность масштабных параметров: расплавленный металл имеет структуру жидкости, формирующаяся на поверхности движущегося металла плёнка состоит из различных оксидов никеля, хрома и алюминия и также имеет аморфную структуру, поскольку формированию кристаллической структуры окисла на поверхности препятствует перемешивание металла. Определение упругих свойств такой крайне нерегулярной структуры дискретно-атомистическими методами, как в [26], может быть предметом отдельного исследования. Исходя из вышесказанного, макроскопический подход в терминах механических напряжений к описанию оксидной пленки на всей поверхности расплава в условиях интенсивной конвекции представляется наиболее перспективным. 1.3. Теплоперенос в тонкой плёнке В силу малой толщины плёнки пренебрежём её способностью аккумулировать тепло. В этом случае равны плотности потока тепловой энергии: приносимая от расплава к плёнке конвекцией qconv, диффузионная внутри плёнки qdif и уносимая с поверхности электромагнитным излучением qrad - эти плотности потока определяются законами Ньютона - Рихмана, Фурье и Стефана - Больцмана соответственно (1) где aconv - конвективный коэффициент теплоотдачи (Вт/(м2×К)), Tmelt, T1, T2, Tcham - температуры расплава, Рис. 2. Схема, поясняющая теплопередачу через плёнку на поверхности расплава Fig. 2. Scheme explaining the heat transfer through the film on the melt surface нижней и верхней поверхностей плёнки и вакуумной камеры соответственно (К), lfilm - коэффициент теплопроводности плёнки (Вт/(м×К)), Crad - коэффициент излучения (Вт/(м2×К4)); схема, поясняющая систему уравнений (1), показана на рис. 2. Величины Tmelt, Tcham и Crad приведены в табл. 2, температуры T1, T2 подлежат определению. Коэффициент теплоотдачи может быть определён из числа Нуссельта , (2) здесь - среднее значение числа Нуссельта, определённое для пластины длиной l, l - коэффициент теплопроводности расплава (см. табл. 1). В соответствии с работой [27] , (3) где - число Прандтля, характеризующее отношение вязкого и температурного пограничных слоёв, а - температуропроводность (см. табл. 1), - число Пекле, характеризующее относительный вклад конвективного переноса тепла в движущемся расплаве по сравнению с диффузионным. Подставляя в уравнение (3) число Прандтля Pr = 0,282, определённое по свойствам сплава, и число Пекле Pe = 5263, соответствующее средней скорости движения расплава и = 0,1 м/с и размеру пластины равному радиусу тигля l = 0,1 м, получим среднее число Нуссельта . Коэффициент теплоотдачи, полученный из формулы (2), равен aconv = 4682. Решая систему уравнений (1) с lfilm = 6 Вт/(м×К) (это среднее значение между коэффициентами теплопроводности оксидов NiO и Cr2O3 при температуре 1700 К [28]) и толщиной плёнки h = 50 мкм, получим искомые температуры на поверхностях плёнки Т1 = 1699,4 К и Т2 = 1699,38 К. Таким образом, в квазистационарном состоянии температурное поле внутри плёнки можно считать однородным по толщине (перепад Т1 - Т2 = 0,02 К), а температура плёнки в некоторой точке соответствует температуре расплава под плёнкой; погрешность такой аппроксимации составляет не более 1 К или 1/1700 » 0,06 %. 1.4. Действие электромагнитного поля на плёнку Мгновенная объёмная сила, действующая со стороны ПМП на плёнку, выражается из законов Кулона и Ампера (4) где fEM - объёмная сила (Н/м3), rel - объёмная плотность зарядов в плёнке (Кл/м3), Е - напряжённость электрического поля (Н/Кл), J - плотность тока в плёнке (А/м2), В - индукция магнитного поля (Тл). Плотности зарядов и токов можно определить по теореме Гаусса для электрического поля и теореме о циркуляции магнитного поля: (5) здесь D - электрическое смещение (Кл/м2), Н - напряжённость магнитного поля (А/м), Ñ - дифференциальный оператор набла. Для расчёта вклада отдельно магнитного [22] или электрического поля (без учёта электрострикции) [29] подходят соотношения (6) которые могут быть объединены и с учётом вклада электромагнитной индукции дадут (7) где - единичный тензор, Ä обозначает тензорное произведение. Уравнение (7) связывает объёмную силу в плёнке с электромагнитным полем, пронизывающем её, электрические и магнитные свойства плёнки (электрическая и магнитная проницаемости) включены неявно через связи D (E) и B (H). При температурах выше 1500 °C оксидная плёнка парамагнитна (m » 1) и является диэлектриком, поэтому в ней нет токов проводимости (J в (4)). Пленка тонкая и находится на электропроводной поверхности расплава, поэтому она не способна накопить значительный объемный электрический заряд, следовательно, членом можно пренебречь. Проанализируем чисто электромагнитный вклад в объёмную силу (7). Осесимметричное магнитное поле индуктора Н = (Hr, 0, Hz) возбуждает в пространстве вихревое электрическое поле, которое имеет только азимутальную компоненту D = (0, Dq, 0). Векторное произведение в этом случае имеет как радиальную, так и осевую компоненты, однако в соответствие с уравнением Максвелла и гармонически изменяющиеся магнитное поле индуктора и вихревое электрическое поле D колеблются со сдвигом фаз , так что среднее значение за период равно нулю. Таким образом, в применяемом приближении (осесимметричность, температура выше точки Кюри, тонкая диэлектрическая плёнка, гармонические поля) электромагнитное поле не оказывает влияния на напряжения в плёнке. 1.5. Упругие свойства оксидной плёнки Поскольку к образованию оксидных плёнок склонны сплавы, содержащие хром и алюминий, логично предположить, что плёнка образована именно оксидами этих элементов. Однако с выбором свойств плёнки существуют сложности: во-первых, свойства оксидных плёнок, образующихся на поверхности расплава во время индукционной плавки при температуре выше 1500 °С в техническом вакууме, недоступны для непосредственного измерения, во-вторых, после открытия вакуумной камеры и остывания тигля, к которому прицепилась плёнка после процедуры «разгона плены», её химический состав и механические свойства в твёрдом окисленном остывшем состоянии будут иные, нежели у только что образовавшейся. Тем не менее обобщим известные в литературе данные. Для описания упругих свойств плёнок удобно пользоваться эластичностью (Н/м), которая равна силе, необходимой для удлинения плёнки на единицу длины (м): . (8) При всём многообразии: от бета-лактоглобулина (молочный белок) до шлаковых пен и ферромагнитных плёнок эластичность плёнок лежит в некоторых пределах. Так, в работах [12; 13] эластичность плёнки принята в интервале 20 £ £ 40 мН/м, в работе [30] приводятся данные 60 £ £ 100 мН/м и указывается, что разрушение плёнки наступает при уровне деформации 1,5 £ xкр £ 2,5 %. Оценка эластичности по методикам [31], применяемым для шлаковой пены и размерам пузырьков 10…20 мм, даёт значения 29 £ £ 230 мН/м. Коэффициент линейного температурного расширения (далее КЛТР) плёнок, хотя и определяется по косвенным данным, имеет следующие общие черты [32-36]. Во-первых, величина КЛТР в зависимости от химического состава плёнки составляет 7,5×10-6 1/К для Al2O3 [36], 5,4×10-6 1/К для Cr2O3 [36], 10,35×10-6 1/К для Cr-Al2O3 [34]. Во-вторых, и в экспериментальных [32; 33], и в теоретических [35] работах указывается, что с уменьшением толщины плёнки КЛТР увеличивается. Коэффициент Пуассона на основании [32; 37] можно принять g = 0,35. Подводя итог вышесказанному, можно утверждать, что тонкая плёнка является достаточно сложным объектом, и выбор её описания диктуется условиями решаемой задачи. В настоящий момент авторам не известны работы других научных коллективов, в которых рассматривались бы взаимное влияние конвекции металлического расплава в переменном магнитном поле и напряжённо-деформированное состояние тонкой оксидной плёнки на его поверхности. Представленная статья является эволюционным развитием работы [38], в которой рассматривалось напряжённо-деформированное состояние оксидной плёнки, полностью покрывающей поверхность расплава, без учёта термических напряжений для единственной конфигурации расплав-индуктор. Целями настоящей работы являются: 1) формулировка математической модели, описывающей напряженно-деформированное состояние тонкой пленки конечного размера на поверхности расплава в ПМП с учётом термического расширения; 2) оценка точности вычислительного алгоритма; 3) моделирование взаимодействия расплава с упругой пленкой на поверхности, перемешиваемого в ПМП в системах различной геометрией (см. рис. 1, г) при варьировании управляющих параметров; 4) описание влияния ПМП на напряженно-деформированное состояние пленки конечного размера. 2. Математическая модель взаимного влияния движения расплава и напряжённо-деформированного состояния плёнки на его поверхности 2.1. Уравнение упругих деформаций в плёнке на поверхности расплава При формулировке математической модели приняты следующие допущения: 1. Вся задача рассматривается в осесимметричной постановке. 2. Плёнка представляет собой диск равной толщины радиусом R, расположенный в центре поверхности расплава, движение плёнки по поверхности не рассматривается. 3. Анализируются только стационарные состояния. 4. Поверхность расплава предполагается плоской и недеформируемой. 5. Оксидная пленка, покрывающая поверхность расплава, представляет собой непрерывную однородную упругую среду. 6. Сдвиговые деформации отсутствуют: в силу малости толщины плёнки пренебрегаем возможным сдвигом в плоскости r-z, равенство нулю сдвиговых деформаций в плоскостях r-q и q-z следует из осесимметричности. В осесимметричном приближении уравнение баланса сил в плёнке имеет вид [39] (9) где s и t - нормальное и сдвиговое механические напряжения (Па), r, q, z - цилиндрические координаты (м), fV - объемная радиальная сила (Н/м3), вид которой будет обсуждаться ниже. Как уже упоминалось ранее, в силу малости толщины плёнки trz = 0 также не учитываем вес плёнки. В результате во втором уравнении системы (9) отсутствует объемная сила, поэтому система (9) может быть сведена к следующему виду (10) Для учёта температурных воздействий использована классическая гипотеза Дюамеля - Неймана: напряжения в среде связаны законом Гука только с упругой деформацией, которая в осесимметричном приближении имеет следующее представление (11) здесь er и eq - радиальная и азимутальная компоненты нормальной деформации, k - коэффициент температурного расширения (1/К), DТ - перепад температур (К), ¡ - модуль Юнга (Па). Используя определение деформаций (12) где x - смещение (м), уравнения для механических напряжений приводятся к виду (13) Используя выражения (12) и (13), уравнение (10) можно привести к виду (14) В настоящей работе плёнка представлена как неподвижный диск радиусом R, поэтому граничные условия ставятся для центра и края плёнки. Смещения в неподвижном центре отсутствуют: (15) Край плёнки находится в равновесии, если на нём отсутствует радиальное напряжение (16) В терминах смещения, с учётом соотношений (12), (13), условие (16) принимает вид (17) Следует отметить, что граничное условие (15) приводит к следующим равенствам для напряжений и деформаций на оси (18) Как было показано в п. 1.4, переменное электромагнитное поле не оказывает силового действия на плёнку, таким образом, объёмная сила в уравнении (14) может быть вызвана вязким трением со стороны движущегося под плёнкой расплава и термическим расширением. Расплав металла является ньютоновской жидкостью, для которой закон вязкого трения создаёт поверхностную силу (касательное механическое напряжение) (Н/м2) (19) в результате чего плёнка тормозит расплав, здесь r и n - плотность (кг/м3) и кинематическая вязкость (м2/с) расплава, и - радиальная компонента скорости (м/с). Градиент скорости вычисляется вблизи покрытой плёнкой поверхности. Движение расплава, в свою очередь, создает в плёнке касательные напряжения , которые распространяются по всему её объему, причём в соответствие с третьим законом Ньютона (20) В силу малой толщины объёмная сила fV (Н/м3) равномерно распределена по толщине плёнки и быть выражена через поверхностную следующим образом (21) где h - толщина плёнки. Таким образом, объёмная сила вязкой природы, действующая со стороны движущего расплава на плёнку, имеет вид (22) Влияние термических напряжений будет рассмотрено в следующем пункте на основе результатов численного моделирования. 2.2. Движение металла в переменном электромагнитном поле Опишем кратко характер течений металла. В расплаве, помещённом в ПМП, возбуждаются индукционные вихревые токи, что приводит к двум механизмам движения. Первый - свободная тепловая конвекция, которая имеет место при неоднородном нагреве расплава индукционными токами (джоулевы источники теплоты). Второй - вынужденная конвекция, возбуждаемая объёмной силой Лоренца J´B, которая возникает при взаимодействии индукционных токов с магнитным полем. Система уравнений, определяющих течения в расплаве, включает уравнения электродинамики сплошной среды, неразрывности, переноса импульса и тепловой энергии, неразрывности: (23) где m0 - магнитная постоянная (Гн/м), m - магнитная проницаемость, t - время (с), z - удельная электропроводность (См/м), р - давление (Па), g - ускорение свободного падения (м/с2), b - коэффициент объёмного расширения (1/К), Т и TL - температуры (К) локальная и ликвидус сплава соответственно, l - коэффициент теплопроводности (Вт/(м×К)), Ñ2 - оператор Лапласа. Используя масштабы, приведённые работе [40], уравнения электромагнетизма из системы (23) приводятся к безразмерным уравнениям диффузии магнитного поля (24) где Н0 - амплитуда магнитного поля индуктора, Н1 и Н2 - амплитуды гармоник представлении поля индукционных токов , w - циклическая частота тока индуктора (рад/с), (25) параметр диффузии магнитного поля, Rc - характерный размер, в качестве такового взят внутренний радиус тигля. Параметр DH при фиксированных свойствах металла (электропроводность) и геометрии установки (внутренний радиус тигля) пропорционален частоте w поля индуктора, которая является единственным доступным для управления во время плавки параметром. Как это видно из уравнений (24), параметр DH определяет структуру индукционного магнитного поля (Н1 и Н2), а значит и распределение джоулевых источников (26) влияющих на тепловую конвекцию, и распределение объёмной силы Лоренца (27) возбуждающей вынужденную конвекцию, таким образом распределения джоулевых источников теплоты и силы Лоренца для конкретной установки определяются частотой. Уравнения движения и неразрывности можно решать в терминах функции тока y и вихря скорости , в осесимметричной постановке они имеют вид (28) где - число Грасгофа, характеризующее тепловую конвекцию, qc - характерный тепловой поток, - модифицированное число Гартмана, определяющее интенсивность течений, Н - характерное магнитное поле. Следует отметить, что вихрь скорости j и ротор силы Лоренца в осесимметричном случае имеют только азимутальную компоненту. Граничные условия для замыкания системы уравнений (28) включают условие свободного течения на оси симметрии, прилипания на твёрдых границах (плёнка, стенка и дно тигля) и термокапиллярного вклада в конвекцию на свободной от плёнки поверхности расплава: (29) здесь - число Марангони, характеризующее интенсивность термокапиллярной конвекции, aТ - температурный коэффициент (1/К) из разложения коэффициента поверхностного натяжения , индекс L относится к температуре ликвидус. Уравнение переноса тепловой энергии из (23) преобразуется к виду (30) где - параметр тепловыделения, показывающий отношение мощностей тепловыделения и теплового излучения с поверхности, - радиационное число Био, равное отношению характерного перепада температуры внутри металла к разности температур между поверхностью металла и окружающей средой, С - коэффициент излучения, учитывающий излучательные способности материалов и геометрию области радиационного теплообмена. Примеры рассчитанных распределений скорости и температуры, реализующихся при разных значениях параметра DH, приведены на рис. 3, а соответствующие им распределения скорости вблизи поверхности и температуры - на рис. 4. Из рис. 3 видно, что при изменении параметра DH течения изменяют свою структуру от двухвихревого (а) (наблюдается при DH < 18) к четырёхвихревому (в) (DH > 52), что объясняется распределениями силы Лоренца [41]. Следует отметить, что конвекция в условиях индукционной плавки полностью вынужденная [42]. В интервале 18 £ DH £ 52 в зависимости от взаимного расположения расплава и индуктора могут реализовываться режимы, при которых течение в центральной части поверхности направлено к оси симметрии, как на рис. 3 (б) и рис. 4 (а), кривая DH = 40. Такие режимы наиболее «опасны» с точки зрения формирования и устойчивости оксидных плёнок на поверхности расплава. а б в Рис. 3. Распределения скоростей и температур в расплаве при различных значениях параметра диффузии магнитного поля (геометрия с - с) Fig. 3. Velocity and temperature distributions in the melt at different values of the magnetic field diffusion parameter (geometry c - c) а б Рис. 4. Профили скорости вблизи поверхности (а) и температуры на поверхности (б), рассчитанные при различных значениях параметра диффузии магнитного поля. Стрелки на профиле скорости показывают направление течения расплава Fig. 4. Velocity profiles near the surface (a) and temperature at the surface (б), calculated at different values of the magnetic field diffusion parameter. The arrows on the velocity profile show the direction of melt flow Распределение температуры неоднородно (см. рис 3) и соответствует структуре течения, что характерно для интенсивной вынужденной конвекции, при этом перепады температуры невелики по абсолютной величине (не более 1 К) вследствие интенсивного перемешивания. Из рис. 4 (б) видно, что для режимов, где скорость направлена в одну сторону, поверхность расплава с плёнкой можно считать изотермической. Для режимов со сменой направления (18 £ DH £ 52) в центральной части наблюдаются градиенты температуры до 10 К/м. Объёмные силы, создаваемые в плёнке движущимся у поверхности расплавом, имеют вязкую природу. Вынужденные течения, разгоняемые ПМП, определяются его параметрами: в первую очередь, параметром диффузии магнитного поля DH, отвечающим за структуры индукционных токов и силы Лоренца, во вторую - модифицированным числом Гартмана Hm, которое влияет на интенсивность течения. Таким образом, состояние плёнки косвенно определяется параметрами магнитного поля DH и Hm. Для выяснения вклада термического расширения оценим и сравним слагаемые в правой части уравнения (14) на основании рассчитанных распределений скорости и температуры, представленных на рис. 4 для DH = 40. В слагаемое с объёмной силой , как было показано ранее, даёт вклад только вязкое трение, которое определяется соотношением (22), а с учётом определения эластичности (8) приводится к виду . Параметры расплава в соответствие с работами [43; 44] представлены в табл. 1. Значение градиента скорости можно оценить из рис. 4 (а): максимальная скорость в узлах, отстоящих на 2 мм от поверхности, составляет 0,12 м/с, таким образом (1/с). При варьировании в интервале от 20 до 230 мН/м получим интервал значений вязкого слагаемого 0,9511 1/м. Взяв максимальное значение КЛТР для плёнки k = 10,35×10-6 1/К [34] и градиент температуры из рис. 4 (б) около оси симметрии К, имеем вклад температурного расширения 1/м. Таким образом, в случае индукционной плавки, когда развивается вынужденное интенсивное течение, а градиенты температур невелики, влиянием теплового расширения на деформацию плёнки можно пренебречь. Кроме того, при свободном крае плёнки тепловое расширение само по себе не может привести к разрушению плёнки. При таком же градиенте температуры и свойствах материалов слагаемые сравниваются при скорости 1,5×10-3 мм/с, для эластомерных плёнок с КЛТР - при 0,15 мм/с, т.е. тепловое расширение может играть некоторую роль только при слабом ползущем течении. 2.3. Обобщённое уравнение равновесия и интегральная деформация Итак, объёмные силы, создающие напряжения в тонкой плёнке на поверхности, определяются в основном вязкими силами, термическое расширение может проявить себя только в случае неразвитой конвекции и больших температурных градиентов. Хотя связанное с термическим расширением слагаемое при индукционной плавке пренебрежимо мало, сохраним его для сохранения общности модели. В этом случае уравнение баланса (14) с учётом объёмной силы (22) и эластичности (8) принимает вид (31) а с применением масштабов длины Rc, скорости n/Rc и температуры qcRc/l преобразуется к безразмерному виду (32) где K - безразмерный параметр температурного расширения, L - вязкоэластичный параметр, характеризующий отношение напряжения вязкого трения в жидкости к упругому напряжению, возникающему в плёнке, определяемые как . (33) Уравнение (32) замыкается на оси граничным условием (15), вид которого остаётся прежним, для края плёнки уравнение (17) преобразуется к безразмерному виду (34) в котором x и R берутся безразмерными с масштабом Rc, а температура - с масштабом qcRc/l. Изменение размера плёнки приводит к изменению граничного условия для скорости у поверхности (29): Гfilm « Гmelt, размер плёнки можно рассчитать как (35) где R0 - начальный размер плёнки. Рассматривая только стационарные состояния (установившееся движение и размер плёнки), можно рассчитывать движение при неподвижной плёнке, полагая её радиус R финальным состоянием, а начальный радиус R0 восстанавливать из уравнения (35). При этом изменение размера плёнки, интегральная деформация рассчитываются из соотношения (36) 3. Численная реализация, верификация и оценка погрешности Математическая модель реализована в оригинальном авторском коде на языке FORTRAN-95 с применением математической библиотеки Intel MKL. Уравнения диффузии магнитного поля, переноса энергии и импульса, механического баланса дискретизированы методом конечных разностей схемами второго прядка точности [45] по координате и первого по времени. Задача о диффузии ПМП и вынужденного движения расплава в нём верифицированы в работах [46-48]. Численный алгоритм расчета перемещений и деформаций основан на решении конечно-разностного аналога уравнений (31)-(33) с граничными условиями с аппроксимацией всех производных формулами второго порядка точности. Для проверки данного этапа была рассмотрена модельная задача (37) Функция в правой части (37) является по смыслу объемной силой, аналогичной таковой в уравнении (32), и выбрана так, чтобы примерно соответствовать по величине и профилю силе, действующей в реальных потоках расплава. Задача (37) имеет следующее точное решение: (38) Сравнение аналитически и численно рассчитанных распределений смещений и деформаций представлено на рис. 5 (a). Зависимости относительных погрешностей расчета . (39) от количества узлов сетки - на рис. 5 (б). В соответствие с определением (12) погрешность расчёта азимутальной деформации Еq совпадает с таковой для смещения Еx. На рис. 5 представлено решение модельной задачи (37), полученное при R = 0,5 и g = 0,25; видно, что численное решение и для смещения x, следовательно, и азимутальной деформации, и радиальной деформации хорошо совпадает с точным, погрешности расчёта Еx и Еq экспоненциально убывают при сгущении сетки (увеличении числа разбиений N). Однако погрешность Еr на рис. 5 (б) на порядок выше Еx и Eq и убывает немонотонно. Это можно объяснить тем, что смещение x имеет максимум в точке r* » 0,457, в этой точке деформация er, определяемая уравнением (12), равна нулю. Когда знаменатель в формуле (39) мал, и ошибка Еr становится большой, периодическое возрастание погрешности обусловлено совпадением узла сетки и r*, которое происходит, если величина r* кратна шагу сетки. В работе [40] было показано, что оптимальное число разбиения N = 50, так как оно обеспечивает как близость решения к пределу сходимости, так и малую скорости, т.е. устойчивость структуры потока. Таким образом, при таких параметрах сетки погрешности расчета можно оценить по рис. 5 (б): Еx, Еr и Eq < 10-3, т.е. не превышают 0,1 % от абсолютной величины. 4. Результаты моделирования деформации плёнки при индукционной плавке Вычислительные эксперименты проведены с параметрами, приведёнными в табл. 2, в них устанавливалось влияние размеров плёнки и параметров магнитного поля, т.е. фактически течений, реализующихся при указанных условиях, на упругонапряжённое состояние плёнки. Алгоритм экспериментов был следующий. Задавались параметр DH и финальный радиус плёнки R, который (радиус) соответствовал квазистационарному состоянию - установившимся течениям расплава, потокам тепла и напряжённо-деформированному состоянию плёнки. Радиус R определял твёрдую границу (условие полного прилипания) для расплава и условие радиационного теплообмена (излучательную способность и соответствующее число Био), с этим радиусом решалась задача тепломассопереноса в расплаве, помещённом в ПМП, где вычислялись установившиеся распределения скоростей под плёнкой и распределение температуры на поверхности, т.е. в плёнке (см. п. 6). По известным профилям скорости и температуры определялись объёмная сила, действующая на плёнку, и решалась задача упругости (32) с граничными условиями (15) и (17), определялись смещения x, локальные e и интегральная DR деформации. Следует отметить, что смещение x, радиусы R и R0, а также интегральная деформация DR представлены в безразмерном виде (масштаб Rc), локальные деформации e безразмерные по определению (12), рассчитанные при решении безразмерного уравнения (32), не зависят от масштаба и не требуют пересчёта. На рис. 6 показаны профили скорости течения расплава при DH = 40 непосредственно под плёнкой для двух её финальных радиусов, и соответствующие им смещения и радиальная и азимутальная локальные деформации. Распределения скорости и температуры а б Рис. 5. Профили смещения и деформации (а) и зависимости погрешностей расчёта смещения Еx, радиальной Еr и азимутальной Еq компонент деформации от числа разбиений N (б). Штриховая линия - линейная аппроксимация погрешности Еr Fig. 5. Profiles of displacement and deformation (a) and the dependencies of the calculation errors of the displacement Еx, radial Еr and azimuthal Еq deformation components on the grid number N (б). The dashed line is a linear approximation of the error Еr Таблица 2 / Table 2 Управляющие параметры, использованные при расчётах Governing parameters used in the calculations Параметр Символ Значение Числа Био для расплава/плёнки Bi 1,24×10-7 / 2,48×10-7 Параметр диффузии магнитного поля DH 1…150 Число Грасгофа Gr 4,13×107 Модифицированное число Гартмана Hm 4581 Параметр температурного расширения K 4,31×10-4 Вязкоэластичный параметр L 2,57×10-7 Число Марангони Ma 2,31×105 Число Прандтля Pr 0,269 Параметр тепловыделения Q 0,014 Коэффициент Пуассона g 0,35 а б Рис. 6. Профили скорости u, смещения x и локальных деформаций e при DH = 40 и различных финальных радиусах плёнки R: 0,5 (а) и 0,8 (б) Fig. 6. Velocity profiles u, displacement x, and local deformations e at DH = 40 and different final film radii R: 0.5 (a), and 0.8 (б) Рис. 7. Зависимости начального радиуса R0 и интегральной деформации DR от финального радиуса R, рассчитанные при DH = 40 Fig. 7. Dependences of the initial radius R0 and the integral deformation DR on the final radius R calculated at DH = 40 в объёме расплава показаны на рис. 3 (б), а их характерные профили на поверхности - на рис. 4. Из рис. 6 видно, что когда течение расплава под плёнкой направлено к центру (отрицательные скорости - на рис. 6, а, в ней возникают деформации сжатия: x и e < 0, если течение от центра (рис. 6, б), она растягивается в центральной части x и e > 0, а у края наблюдается значительная деформация сжатия er < 0. Связь финального радиуса R в напряжённо-деформированном состоянии с начальным радиусом ненагруженной плёнки R0 показана на рис. 7, там же приведена зависимость интегральной деформации DR. Из рисунка видно, что существует такой финальный радиус плёнки R, при котором её размер не изменяется, а радиус равен таковому в ненагруженном состоянии R = R0. Такое возможно, только когда течение расплава на поверхности направлено в разные стороны, это реализуется в интервале значений параметра диффузии магнитного поля 18 £ DH £ 52 и некоторых геометриях «индуктор - расплав» [47; 48]. Результаты вычислительных экспериментов, в которых рассчитывались локальная радиальная деформация er в центре плёнки и интегральная деформация DR при варьировании финального радиуса в интервале 0,1 £ R £ 1 (при R = 1 край плёнки считается свободным) и параметра диффузии магнитного поля в интервале 1 £ DH £ 150 для двух геометрий «индуктор - расплав» (см. рис. 1, г), представлены на рис. 8. Рис. 8. Карты локальной радиальной деформации er в центре плёнки и интегральной деформации DR на плоскости «финальный радиус плёнки R - параметр диффузии магнитного поля DH», рассчитанные для геометрии “b - b” (а, б) и “c - c” (в, г) Fig. 8. Maps of the local radial deformation er in the center of the film and the integral deformation DR on the plane "final film radius R - magnetic field diffusion parameter DH" calculated for the geometry "b - b" (а, б) and "c - c" (в, г) Из рис. 8 (а) видно, что существует область параметров R - DH, ограниченная изолинией er = 0, в которой плёнка оказывается ненагруженной в центре, другая выделенная изолиния er = 0,025 показывает предел 2,5 %, при котором плёнка рвётся [30]. Наличие области нулевых и отрицательных локальных деформаций объясняется течением расплава, направленным к центру, при этом плёнка сама тормозит течение, снижая его скорость, и тем самым уменьшает деформирующие силы вязкой природы, что делает её ещё более устойчивой. Для геометрии “b - b” в интервале 17 £ DH £ 74 нулевая деформация плёнки определяется её радиусом, а при 45 £ DH £ 60 плёнка не нагружена в центре при любых R. По изолинии er = 0,025 можно сказать следующее: при DH £ 27 радиус R*, при котором плёнка рвётся в центре, зависит от DH, при 28 £ DH £ 70 плёнка не рвётся в центре при любых R, при 71 £ DH £ 82 R* также зависит от DH, а при DH > 83 R* » 0,4 и не зависит от DH. Из рис. 8 (б) видно, что интегральная деформация DR = 0 возможна при 15 £ DH £ 75, при этом радиус не меняющей размер плёнки зависит от DH. Попарное сравнение карт er и DR на рис. 8, построенных для различных геометрий, позволяет выявить влияние взаимного положения индуктора и расплава на состояние плёнки. Видно, что симметричное расположение “c - c” более выгодно с точки зрения противодействия образованию поверхностных оксидных плёнок, чем геометрия с общим дном “b - b”: области параметров, в которых возможны плёнки с недеформированным центром er = 0, в которых не достигается деформация, необходимая для разрыва er = 0,025, и с нулевой интегральной деформацией DR = 0, сократилась как по оси DH, так и по R. Интересно отметить, что для плёнок радиусом R £ 0,4 деформация er в центре не превышает порога разрыва 2,5 % при любых DH, т.е. они могут быть стабильны при любых DH. Кроме того, для геометрии “c - c” в интервале 18 £ DH £ 52 движение расплава на части поверхности направлено к центру (см. на рис. 3 (б)), в этой части плёнка будет удерживаться в центре поверхности течением, не испытывая сколько-нибудь значительных нагрузок. За пределами указанного интервала течения всегда направлены от центра, и плёнка может существовать с точки зрения уровня er, однако она будет находиться в состоянии неустойчивого равновесия, и отклонения от осевой симметрии реальных течений или плёнки может привести к её удалению с поверхности расплава. Заключение Разработана математическая модель, описывающая напряжённо-деформированное состояние тонкой плёнки, находящейся на поверхности расплава металла во время индукционной плавки. Показано, что основной вклад в механическое нагружение даёт вязкое трение, которое возникает на границе плёнки и расплава, перемешиваемого переменным электромагнитным полем. Таким образом, параметры магнитного поля, такие как параметр диффузии магнитного поля DH, пропорциональный частоте и отвечающий за распределения в расплаве источников теплоты и силы Лоренца, т.е. за структуру течения, и число Гартмана Hm, пропорциональное напряжённости поля и отвечающее за интенсивность вынужденной конвекции, через движение расплава влияют на состояние плёнки. Приведены данные о верификации предложенной математической модели как для гидродинамической её части, так и для упругой. Показано, что вычислительная погрешность деформаций не превышает 0,1 %. На основании обзора литературы установлены возможные упругие параметры плёнки и проведены вычислительные эксперименты, в которых установлено влияние геометрии «индуктор - расплав», размера плёнки и параметра диффузии магнитного поля на деформацию плёнки. Показано, что при геометрии “b - b” уровень напряжений в плёнке ниже, чем при “c - c”, что связано с характером течений: при “b - b” движение расплава менее интенсивное и в более широких интервалах DH направлено к центру. Выявлено, что течения расплава во всём исследованном диапазоне DH слабо деформируют плёнки радиусами R £ 0,4, уровень деформации растяжения в центре плёнки для таких плёнок ниже порога разрыва 2,5 %; существуют диапазоны значений DH, при которых плёнка не разрывается при любых радиусах в случае “b - b” и достигает R = 0,85 при “c - c”.Обнаружены интервалы параметра DH с нулевой интегральной деформацией плёнки, т.е. находясь в напряжённо-деформированном состоянии, плёнка не меняет своего радиуса. Таким образом, геометрия «индуктор - расплав» “c - c” менее склонна к формированию и стабилизации оксидных плёнок на поверхности расплава во время индукционной плавки; для интервала 18 £ DH £ 52 для разрыва и удаления плёнки R > 0,4 можно рекомендовать уменьшения DH (понижение частоты магнитного поля), для R £ 0,4 в указанном диапазоне удаление плёнки возможно при нарушении симметрии течения (наклон тигля с расплавом или включение дополнительных катушек, создающих магнитное поле). Таким образом, представленная модель позволяет анализировать напряженно-деформированное состояние пленки на поверхности металлического расплава при индукционной плавке. Это полезно при разработке эффективных технологических режимов и устройств для предотвращения образования плёнок при плавке и повышения качества продукции.

About the authors

I. L Nikulin

Perm National Research Polytechnic University

V. A Demin

Perm National Research Polytechnic University; Perm State University

References

Statistics

Views

Abstract - 147

PDF (Russian) - 142