ВЕСТНИК
ПРИАЗОВСКОГО ГОСУДАРСТВЕННОГО ТЕХНИЧЕСКОГО УНИВЕРСИТЕТА 2000 г. Вып.№9
УДК 532.5:669.18
Жук В.И.1
КОМПЬЮТЕРНОЕ МОДЕЛИРОВАНИЕ ТЕПЛОПЕРЕНОСА, ДИФФУЗИИ И ДВИЖЕНИЯ ЧАСТИЦ В ЗАТВЕРДЕВАЮЩЕМ РАСПЛАВЕ
На основе применения конечно-разностных схем для уравнений конвективного тетомассопереноса и гидродинамики частиц в кристаллизирующемся расплаве построена компьютерная модель, имеющая удобный интерфейс для пользователя. С помощьц) модели может быть проведен анализ процессов охлаждения жидкого ядра слитка, диффузии радиоактивного источника, перемещения кристаллов твердой фазы, удаления неметаллических включений в завиашости от способа перемешивания расплава.
Память и быстродействие современных вычислительных устройств позволяет на современном этапе развития теории металлургических процессов выйти на уровень непосредственного моделирования многих явлений, происходящих при кристаллизации слитков и отливок, с помощью компьютера ("компьютерного моделирования"). Вследствие больших энергозатрат и технических трудностей, а порой и невозможности прямых измерений на натурных слитках, многие исследователи пошли по пути "вычислительного эксперимента", как показано в работе [1], пытаясь найти адекватное описание полученным закономерностям. Время, затраченное на проведение такого "вычислительного эксперимента", становится сравнимым со временем протекания процесса, что позволяет проверить многие гипотезы без значительных капиталовложений в научные разработки. Однако очень мало исследований ставят своей целью проведение целевых и контрольных экспериментов для проверки правильности выдвигаемых гипотез. Некоторые опытные данные до сих пор широко не используются в расчетах. Например, полученные в работе [2] методом радиоактивных изотопов распределения эффективных коэффициентов диффузии позволили выдвинуть определенную схему естественной конвекции в слитке. Однако их использование затруднительно из-за отсутствия соответствующей методики. В работе [3] показано, что вынужденное перемешивание расплава влияет на коэффициент теплопередачи, и, как следствие, на формирование кристаллической структуры слитка. Однако механизм этого влияния и оценку степени перемешивания трудно представить без соответствующей компьютерной модели. Сопоставление математических моделей с уже накопленными опытными данными и выдвинутыми схемами явно недостаточно именно по причине отсутствия удобных в пользовании программ. Целью настоящей работы является создание компьютерной модели2, позволяющей на основе современного математического аппарата произвести расчет распределения температур, концентраций и частиц в условиях конвективного движения расплава и проверить некоторые существующие предположения о схемах течения в слитках и отливках.
Математическая модель процессов тетомассопереноса.
Дифференциальные уравнения, выражающие законы сохранения энергии и массы и описывающие распределение температур и концентраций растворенной в расплаве примеси, имеют следующий вид:
1 ПГТУ, канд. техн. наук, доцент
2 В разработке приняли участие студенты ПГТУ Белов А., Делиев А., Мищенко С.
~ + (УУ)Т = аАТ, (1)
81
^- + (УУ)С = ОАС, (2)
где / - время, Т - температура расплава, С - концентрация примеси (вещества трассера), V -вектор скорости, а - коэффициент температуропроводности, О - коэффициент диффузии (самодиффузии).
Введем понятие рабочего объема как некоторого п-мерного пространства, в котором исследуется процесс тепломассопереноса. Дальнейшая задача заключается в моделировании процесса конвективной теплопроводности и диффузии с учетом того, что у нас уже имеется
значения скоростей V в каждой точке рабочего объема в произвольный момент времени. Это могут быть данные расчетов на модели, например, в задаче о естественной или вынужденной конвекции [3], либо данные, заложенные исследователем в виде схемы, полученной при обобщении результатов эксперимента [2]. Для решения уравнений (1,2) требуется задать начальные и граничные условия в конкретном случае.
Рассмотрим, например, задачу теплопереноса при охлаждении жидкой сердцевины слитка. Пусть в начальный момент времени температура расплава постоянна. Т |<=0 = Т0, где
Т0- начальная температура. В следующий момент времени температура на границах
становится равной температуре кристаллизации Тк. Верхнюю границу будем полагать теплоизолированной. С помощью ряда преобразований уравнение (1) приводится к безразмерному виду в области с неподвижной границей (так как размеры рабочего объема
могут изменяться). Введем безразмерные координаты области 7] 1 = ——, 7]2 — ——, где х0-
х0 х0
Т-Т
характерный размер области, Ро=М/х0-безразмерное время (число Фурье), (9 =
т0-тк
безразмерная величина температуры. Тогда уравнение (1) преобразуется к виду (3)
80 „ 80 „ 80 Г 820 820
- + V; - + V, -
8Ро 8 т], 8т]2
(3)
дг], дт]2
где V], - безразмерные компоненты скорости вдоль координат XI и х2.
Уравнение (3) - основное уравнение конвективной теплопроводности, приведенное в безразмерном виде. Добавив начальные и граничные условия (4), получаем математическую модель охлаждения жидкой сердцевины слитка:
30
= 0 01 = 0, г/к - координата фронта кристаллизации. (4)
6|Г =7
{Чк
Пг-о
Аналогично формулируется задача конвективной диффузии примеси, например, радиоактивного изотопа. Пусть в расплав, залитый в изложницу конечного объема, в определенный момент времени помещается источник диффузии определенной концентрации. В реальных условиях это может быть радиоактивный источник из тех же атомов (самодиффузия) или атомов другого вещества (диффузия). Источник помещен в ампулу, которая быстро растворяется, и начинается процесс конвективной диффузии источника (трассера), о концентрации которого в данной точке расплава можно судить по интенсивности радиоактивного излучения. Так как источников может быть несколько, в начальный момент
времени концентрация в каждой ^той точке объема равна С
^=С](х,у), где С)-
начальная концентрация в j-той точке.
Приведем уравнение (2) к безразмерному виду. Для этого введем величины:
£ = С /С0 -безразмерная величина концентрации, С0 -масштабная величина (принять равной
максимальной величине концентрации в объеме), Рга В/а - диффузионное число Прандтля. Тогда уравнение (2) преобразуется к виду (5)
д8 „ Э5 „
■+У,
5Fo дт]1 Уравнение
= РГ;
дд28
+
дг)1 дт)
(5)
(5) - основное уравнение конвективной диффузии, приведенное в безразмерном виде. Добавив начальные и граничные условия (6), получаем математическую модель конвективной диффузии источника (трассера):
дБ
дг1 „=.
О, т]к- координата фронта кристаллизации.
(6)
Для получения численного решения уравнений (3-6) воспользуемся методом конечных разностей или методом сеток [4]. Для этого двумерную рабочую область разобьем горизонтальными и вертикальными линиями с шагами Ь] и Ь2 соответственно (рис. 1).
Ц2 Ь1
м 1 Г6*
т
1 к 1
т 1 Л Ь, У
0 1 I
Рис. 1 - Схема дискретизации расчетной области.
Таким образом, рабочий объем представляется набором узловых точек 1*М с координатами 0<К1; 0<т<М. Расчет концентрации ведется только в узловых точках, в связи с этим заменяем производные функции концентрации от времени и координаты в точке в уравнении (3) на соответствующие приращения функции концентрации от времени и координаты. Переписывая уравнение (3) с использованием конечно-разностных формул получим •
в'.
1, т
1~&к -
0 * . - 0 к ,
I + /, т_/ - /, т +у к
2 к, 2
~0К
I, т + /_I, т - 1
2И
0
/ - 1, т
■ & к + 1 + & к I, т
. * 0к I + 1, т + I, т
2-0к+} +&к 1 1,т I, т + I
(7)
где т - выбранный шаг по времени, - температура в (1, ш) точке, к - индекс по времени. С
помощью уравнения (7) можно рассчитать температуру в (1. т) точке на к+1 шаге по времени.
Особое внимание необходимо уделить краевым условиям, т.к. в этой задаче они играют определяющую роль в поведении системы. Раскладывая в ряд Тейлора функцию температуры (концентрации) в окрестности некоторой граничной точки (0,т) и принимая во внимание граничные условия (4), получим:
0
= 0
1т?2=0
+
■К:
0
к+1 - /Л к+' , ^ и 2 1.0 -&1.1
д0 д20
дРо 11=0 дц] 1, = 0 _
.о - 20 к+! + 0 к ""^1.0 ^ ^ 1 + 7,0
К
(8)
(9)
2 г А/
Используя формулы (8-10) совместно с уравнением (7) можно вычислить температуру в любой точке объема в произвольный момент времени. Отдельного рассмотрения требуют так называемые "угловые" точки, которые в расчетах не участвуют, однако, для компьютерной модели их значение также определяется путем разложения в ряд Тейлора по двум переменным.
Для уравнения диффузии (5) с граничными условиями (6) получены аналогичные соотношения с той лишь разницей, что в формулах (7-10) присутствует диффузионный критерий Прандтля Ргд, и в начальный момент времени картина распределения источников диффузии задается произвольным образом. Используя соответствующие выражения, можно вычислить концентрацию примеси в любой точке объема в произвольный момент времени.
Дтя численной реализации представленных математических моделей и алгоритмов созданы программы расчета, имеющие удобный для пользователя интерфейс, в котором данные представлены в виде картин изотерм и изолиний концентраций. Такую программу в дальнейшем будем называть компьютерной моделью соответствующего процесса. Ниже приведены два типа моделей для практически важных случаев исследования кристаллизационных явлений.
Компьютерная модель охлаждения жидкой сердцевины слитка
Программа разработана для OS Windows 95 с помощью алгоритмического языка Borland С++ Builder 4.0. Для работы необходима машина не ниже 486 DX с 8Мб RAM. Внешний вид главного окна программы приведен на рис.2.
Рис.2 - Главное окно программы "Компьютерная модель охлаждения жидкой сердцевины слитка".
Рассмотрим назначение управляющих кнопок и расшифруем приведенное изображение:
• Кнопка "Настройка" - при нажатии этой кнопки открывается окно, в которое
пользователь вводит критерии моделирования, размер рабочей области, шаг сетки, шаг по времени и значение времени окончания процесса моделирования.
• Кнопка "Запуск"- запускает расчет. Кнопка "Стоп"- останавливает расчет.
• Кнопка "Пауза"- приостанавливает расчет для того, чтобы просмотреть температуру в отдельных точках и, в случае необходимости, распечатать экран.
• Кнопка "Сохранить картинку"- позволяет сохранить изображение области в файле
для сопоставления результатов.
Слева видна полоска с цветовой гаммой, позволяющая визуально представить численную картину в масштабе длин волн цвета от красного до фиолетового (в черно-белом изображении контраст может быть не заметен). Справа, в масштабе изложницы, который может изменяться пользователем, рассмотрена рабочая область, в которой с помощью цветовой гаммы представлена картина распределения температур в жидкой сердцевине слитка. Чередование глобального и локального масштаба температуры позволяет судить об автомодельности задачи. С помощью "мыши" можно определить текущую температуру в любой точке области.
Компьютерная модель диффузии радиоактивного источника
Программа разработана для OS Windows 95 с помощью алгоритмического языка Pascal в среде Delphi 3.0. Для работы необходима машина не ниже 486 DX с 8Мб RAM. Внешний вид главного окна программы приведен на рис.3:
Длима изложницы : 400 Ш ирина изложницы: 3Q0
П уск.
П арамеп"ры
И нициализация
Рис.3 - Главное окно программы "Компьютерная модель диффузии радиоактивного источника".
Рассмотрим назначение управляющих кнопок и расшифруем приведенное изображение:
• Кнопка "Инициализация" - при нажатии этой кнопки открывается окно, в которое
пользователь вводит размер рабочей области и шаг сетки. В дополнительных окнах вводятся координаты точек, в которых задается уровень начальной концентрации радиоактивного источника, и координаты точек, за которыми будет производиться наблюдение в процессе счета.
• Кнопка "Параметры" - при нажатии этой кнопки открывается окно, в которое
пользователь вводит параметры моделирования: шаг по времени, значение времени окончания процесса моделирования, диффузионное число Прандтля. Пользователь также может изменять граничные условия и отключать конвекцию. Предусмотрена возможность регистрировать промежуточные результаты вычислений в файлах отчета.
• Кнопка "Пуск" - при нажатии этой кнопки программа строит рабочую область и
запускает расчет.
Слева в окне в цветовом разрешении изображаются линии изоконцентраций, соответствующие распределению радиоактивного изотопа в жидком ядре с течением времени. Справа выводятся графики изменения концентрации с течением времени в определенных точках области.
Математическая модель гидродинамики частиц.
Математическое моделирование гидродинамики частиц в затвердевающем расплаве основывается на алгоритме, предложенном в работе [5]. Сформулируем следующую гидродинамическую задачу: пусть в расплаве появляются частицы (кристаллы,
неметаллические включения, инокуляторы, модификаторы и т. п.) определенного размера, формы и плотности и в некоторый момент времени известны скорость и координаты каждой из них. Дальнейшее поведение любой частицы, т.е. зависимость скорости и координат от времени, а также траекторию частицы, можно предсказать на основе анализа независимого движения частицы в несущем потоке. Полагая частицу сферической, уравнение ее движения с учетом только гидродинамических сил запишем в виде (11):
dW 1 d(W-V) - 9v -Г
Ps~lT+-2Pl—t-= <PS-Pl>'*-^2-<W-V> <И>
Здесь W - скорость движения частицы, V - скорость движения несущей среды, ps -
плотность частицы, р1 - плотность несущей среды, v - вязкость несущей среды, g - ускорение
свободного падения, R - радиус сферы. Вводя компоненты скорости частицы W, и W2, число
Пекле Ре = а / Хд, согласующее безразмерное время Fo с предыдущими задачами,
коэффициент Е — рн /pt — 1, получим уравнения:
dW j dV 9v
(s+1,5)--l-------=--(W,-V.)
dFo 2Pe dFo 2PeR2 1 1
dW j dV 9 (12)
(e + 1,5)--L--i---= ----~r(W.-V.) v '
dFo 2Pe dFo Pe 2PeR2 2 2
Уравнения (12) - основные уравнения движения частицы, приведенные в безразмерном виде. Из уравнений (12) конечно-разностным методом можно найти компоненты скорости частицы в следующий момент времени:
W*1 = Wj +—-------(W,k -V,k+J)
1 ' £ + 1,5 2R Ре 11 1 '
urk+1 rrrk _£__(F—___—_CWk -Vk + 1 ))
2 £ + 1,5 ( Pe 2R2Pe 2 })
(13)
Зная скорость на следующем шаге можно найти новые координаты частицы: £к + 1- £к+уук + 1.х gk + l- gk + wk + l.T (14)
Для тех случаев, когда частица не попадает в узел, применяется интерполирование методом Ньютона для функции 2-х переменных.
Компьютерная модель движения частиц
Программа разработана для OS Windows 95 на языке С++, при помощи пакета Borland С++ Builder 1.0. Для работы необходима машина не ниже 486 DX с 8Мб RAM. На жестком диске необходимо иметь место для хранения результатов расчета из расчета один шаг модели -20 байт для каждой частицы, рассматриваемой в опыте. Для работы программы необходимо иметь файл с уже посчитанными скоростями. Результаты моделирования могут быть сохранены в файле, отдельном для каждой частицы. Внешний вид главного окна программы приведен на рис.4:
Рис.4 - Главное окно программы "Компьютерная модель движения частиц".
На панели присутствуют две кнопки:
• Кнопка "Инициализация" - служит для задания параметров модели. Параметры
могут относиться как ко всей модели в целом, так и к конкретной частице.
• Кнопка "Моделирование" - запускает процесс моделирования, при этом в левой
части экрана строится траектория движения частицы, а внизу показывается процент завершения расчетов.
По результатам расчета на экране строится траектория движения частицы или нескольких частиц, которые могут появляться в расплаве. В принципе в программе можно задать изменение размеров частицы в ходе движения и таким образом полностью имитировать поведение кристаллов и неметаллических включений при кристаллизации слитка.
Выводы
Разработана компьютерная модель, основанная на применении конечно-разностных схем для уравнений конвективного тепломассопереноса и гидродинамики частиц и имеющая удобный интерфейс для пользователя.
С помощью модели может быть проведен анализ процессов охлаждения жидкого ядра слитка, диффузии радиоактивного источника, перемещения кристаллов твердой фазы, удаления неметаллических включений в зависимости от способа перемешивания расплава.
Представленная компьютерная модель может оказаться полезной для металлургов-исследователей и научных работников для анализа конвективных течений и согласования экспериментальных и теоретических результатов.
Перечень ссылок
1. Жук В.И. Проблемы математического моделирования кристаллизации слитков в изложницах.//Вестник Приазов. гос. тех. ун-та: Сб. науч. тр. - Мариуполь, 1996. - Вып.2. -С.48-51.
2. Скребцов А.М. Конвекция и кристаллизация металлического расплава в слитках и отливках,- М. Металлургия, 1993 - 144 с.
3. Ефимов В.А., Элъдарханов А С. Современные технологии разливки и кристаллизации сплавов. - М.: Машиностроение, 1998.-360с.: илл.
4. Самарский A.A. Введение в теорию разностных схем.- М.: Наука, 1971. - 552 с.
5. Жук В. И. Математическое моделирование движения кристаллов в затвердевающем расплаве.//Вестник Приазов. гос. тех. ун-та: Сб. науч. тр. - Мариуполь, 1998. - Вып.6. -С.358-360.
Жук Виктор Иванович. Канд. техн. наук, доцент кафедры физики ПГТУ, окончил Донецкий государственный университет в 1972 году. Основные направления научных исследований - математическое и физическое моделирование теплообменных процессов при кристаллизации металлов.