ГОУ «Санкт-Петербургский государственный политехнический университет»
Факультет технологии и исследования материалов
Кафедра моделирования металлургических процессов
Отчёт
о лабораторной работе №1 «Построение детерминированной программной модели».
Работу выполнила студентка группы 4064/1 Петрова С.С.
Преподаватель Вяххи И.Э.
Санкт-Петербург
2009 г.
Цель работы – знакомство с методикой построения детерминированной модели и разработка численной модели металлургического процесса для системы с сосредоточенными параметрами.
Система – целесообразная совокупность взаимодействующих элементов, которая ориентирована на выполнение той или иной функции. Системы классифицируются по различным параметрам. По пространственной структуре бывают системы с сосредоточенными и распределёнными параметрами.
Система с сосредоточенными параметрами имеет равномерное распределение выходных параметров по объёму и характеризуется их осреднёнными величинами.
Система с распределёнными параметрами имеет более сложную структуру, в ней выходные параметры неравномерно распределены по объёму.
Модель – объект, находящийся по отношению к натурному объекту в отношении подобия (т.е. взаимно-однозначного соответствия); приближённое описание процессов, происходящих в системе, ориентированное на выполнение определённых функций.
Детерминированные модели описывают процессы с известным механизмом с помощью физико-химических уравнений, влияние случайных возмущений не учитывается.
Постановка задачи:
Необходимо разработать математическую модель процесса кристаллизации сплава системы Fe – 0,16 %
C с учётом равновесного выделения твердой фазы по диаграмме состояния. Изменение температуры данного сплава происходит в интервале от tз=1800°С (температуры заливки в форму) до tк=800°С (конечной температуры охлаждения слитка). Удельная теплоёмкость сплава С=444 Дж/(кг·К), его плотность ρ=7000 кг/м3
, коэффициент теплоотдачи α=126,5 Вт/(м2
·К), скрытая теплота фазового превращения L=277 кДж/кг.
Поскольку мы имеем дело с системой с сосредоточенными параметрами, то модель, с помощью которой мы будем описывать процессы, происходящие в системе, так же будет с сосредоточенными параметрами, то есть будет характеризоваться осреднёнными параметрами.
Построение физической модели:
На рисунке 2 показан характер изменения температуры во времени при охлаждении сплава заданного состава. Разобьем температурный интервал охлаждения слитка сплава Fe – 0,16 %
C на участки:
I. Участок жидкого состояния.
При заданной температуре заливки tз=1800°С сплав Fe – 0,16 %
C находится в жидком состоянии; жидкость охлаждается до пересечения с линией А-В при температуре ликвидуса (tliq
) в точке А1
.
II. В точке А1
начинается выпадение из расплава кристаллов феррита,
Ж → Ф
при этом концентрация жидкости меняется по линии А-В, а концентрация феррита – по линии А-С.
III. В точке С1
,
когда сплав достигает температуры перитектики (tp=1499°С), выпадение кристаллов феррита заканчивается. При данной температуре происходит перитектическое превращение:
Ж + Ф → А
С-С1
-В – линия нонвариантного перитектического превращения.
IV. Ниже точки С1
охлаждение идёт в области диаграммы, соответствующей твёрдому состоянию сплава, до заданной температуры охлаждения слитка tк=800°С.
Для дальнейшего построения физической модели примем следующие допущения:
- для простоты модели будем рассматривать систему как объём с сосредоточенными параметрами, т. е. не учитывая перепада температур по сечению слитка и принимая его в качестве материальной точки, имеющие постоянные поверхность теплообмена F, объем V и плотность ρ;
На рисунке 3 изображён слиток заданного сплава:
(м2
), (м3
).
- теплообмен слитка со средой происходит по закону конвективного теплообмена Ньютона с постоянным коэффициентом теплоотдачи α, т.е. пренебрежём лучистым (радиационным) теплообменом;
- основные параметры системы (плотность р, теплоемкость С,
теплота кристаллизации L) являются постоянными, не зависящими от температуры и состава выделяющихся фаз, в том числе не учитываются объёмные изменения и физико-химическое взаимодействие между сплавом и окружающей средой;
- при расчете предполагаем линии диаграммы состояния отрезками прямых (линейные зависимости) и рассчитываем соответствующие концентрации как линейные функции от температуры.
Формулировка математической модели
.
Основным соотношением рассматриваемой математической модели является уравнение баланса энергии: изменение внутренней энергии слитка dQc
равно количеству денного в среду тепла dQB
.
Теплообмен слитка со средой происходит по закону конвективного теплообмена Ньютона:
, (1)
где α – коэффициент теплоотдачи от поверхности слитка площадью F
и объёмом V в окружающую среду;
tcp
- температура окружающей среды;
t - температура слитка;
τ - время.
В общем случае изменение внутренней энергии сплава в зависимости от этапа кристаллизации имеет вид
|
|
|
(2)
где dm – изменение относительного количества твердой фазы ( т
≤ 1). Рассмотрим участки охлаждения и кристаллизации сплава Fe – 0,16 %
C
I. Для первого участка охлаждения жидкой фазы от температуры заливки (tз) до температуры ликвидуса этого сплава (tn
) уравнение баланса энергии имеет вид обыкновенного дифференциального уравнения первого порядка:
, откуда .(3)
Условие окончания первого этапа охлаждения жидкого металла имеет вид:
Значение tл
(С0
) рассчитываем по диаграмме состояния в виде линейного соотношения
, где С0
– состав сплава.
Зависимость температуры ликвидуса от состава сплава при допущении о прямолинейности линий диаграммы состояния находим в виде уравнения прямой
по двум точкам. Как видно из диаграммы состояния (рис. 1), на линии А-В при С=0% tл
=1539°С; при С=0,51% tл
=1499°С. Составим систему уравнений с использованием данных значений:
(4)
Решив эту систему, получили уравнение зависимости температуры ликвидуса от состава сплава:
tл
=1539–78,43С (5)
II. Для второго участка охлаждения в интервале температур от перитектики (tp
) уравнение баланса энергии с учетом теплоты фазового превращения при выделении феррита имеет вид:
a (t
–tep
)Fdτ =
–V
С
pdt + VLpdm,
(6)
где в правой части содержатся дифференциалы двух взаимозависимых переменных – t и m.
Для исключения «лишней» переменной и преобразования этого уравнения к виду
необходимо связать дополнительным соотношением t и m ,
а затем выразить dm через количество твердой фазы – феррита при температуре t в условиях равновесной кристаллизации найти из диаграммы состояния по правилу отрезков:
, (7)
где С1
(t) и С2
(t) – концентрация углерода С в жидком сплаве (линия А-В) и в феррите (линия А-С) при температуре t..
Зависимости С1
(t) и С2
(t) находим, используя допущение о том, что линии диаграммы состояния являются прямыми. Как видно из диаграммы состояния (рис. 1), на линии А-В при t=1539°С С=0%; при t=1499°С С=0,51%. Составим систему уравнений с использованием данных значений: (8)
а = – 0,01275;
Решив эту систему, полу
(t)=19,62 – 0,01275t .
Аналогично находим уравнение для прямой С2
:
На линии А-С при t=1539°С С=0%; при t=1499°С С=0,1%. Составим систему уравнений: (9)
Решив эту систему, получили уравнение зависимости С2
(t)=3,85 – 0,0025t.
После соответствующих преобразований и подстановки в соотношение (7) определим темп выделения твёрдой фазы путем дифференцирования уравнения (7) для последующего использования этого выражения в математической модели (6):
, тогда
Введём обозначение dm/dt = A, тогда из (6) получим следующее дифференциальное уравнение:
a (
t
–tep
)
F
dτ
= –VС
p
dt +
VLp
Adt,
откуда получаем
(10)
Условие окончания второго этапа: (11)
III Для третьего участка кристаллизации (перитектическое превращение при постоянной температуре tр
) уравнение баланса энергии принимает вид:
a(t
–tep
)F dτ
= VLpdm
, откуда (12)
(13)
Условие окончания третьего этапа: (14)
где mф
– количество твёрдой фазы, которая остаётся после завершения перитектической реакции. m(1) легко определяемое по правилу отрезков на диаграмме состояния (ниже линии tр
).
IV. Для четвертого участка охлаждения твердой фазы в интервале температур от tр
до tк
уравнение баланса энергии по виду аналогично первому участку:
, (15)
Условие окончания четвёртого этапа: (16)
После определения всех уравнений модели переходим к выбору метода и составлению алгоритма численного решения задачи.
Выбор метода и разработка алгоритма решения задачи:
Преобразуем дифференциальные уравнения модели в разностную форму методом Эйлера:
(17)
(18)
(19)
(20)=(17)
Программирование задачи:
Оформленная в системе MATLAB программа по описанному выше итерационному алгоритму:
%Programma of coling Fe-0.16C alloy
clear;
c0=0.16;C=444;R=7000;L=277000;a=126.5;F=0.025;V=0.00025;
t(1)=1800;tp=1499;tk=800;tsr=20;dtau=1;m(1)=0.854;n=10000;
tliq=1539-78.43*c0;
%step_1
for i=1:n; if t>tliq;
t(i+1)=t(i)-a*F*dtau*(t(i)-tsr)./(V*C*R);
end;end; s1=length(t);
%step_2
for i=s1:n; if t>tp;
t(i+1)=t(i)+a*F*dtau*(t(i)-tsr)./(-V*C*R-0.0016*V*L*R/(15.77-0.01025*t(i)).^2);
end;end; s2=length(t);
%step_3
m(s2)=0.854; for i=s2:n; if m<1;
m(i+1)=m(i)+a*F*dtau*(t(i)-tsr)./(V*L*R);
t(i+1)=t(i);
end;end; s3=length(t);
%step4
for i=s3:n; if t>tk;
t(i+1)=t(i)-a*F*dtau*(t(i)-tsr)./(V*C*R);
end;end; sk=length(t);
h=0:dtau:dtau*(length(t)-1);
plot(h,t); grid
title('temperature as function of time')
ylabel('temperature,C')
xlabel('time,sec');
Полученная в результате запуска программы математическая модель процесса кристаллизации сплава системы
Fe – 0,16 %
C
Проверка адекватности модели:
Проверка адекватности программной модели состоит в сопоставлении расчетных ур
и экспериментальных уэ
значений, которые выдает преподаватель, и в оценке величины коэффициента парной корреляции rУрУз
по формуле:
и коэффициента регрессии a
=tgφ
где yр
i
и yэ
i
- пары соответственных значений расчетных и экспериментальных данных; и — средние значения ур
и уэ
по совокупности всех п
сопоставляемых величин.
Таблица 1
Анализ согласованности расчётных и экспериментальных данных.
| время τ, сек |
50 |
100 |
150 |
200 |
250 |
300 |
| t расчётная, °С |
1520 |
1510 |
1490 |
1215 |
1000 |
810 |
| t эксперим., °С |
1600 |
1550 |
1450 |
1100 |
850 |
730 |
Рис. 5 Проверка адекватности численной модели с помощью графического образа.
rкр
=0,707 (число степеней свободы n=6, уровень доверительной вероятности р=0,95)
Из графического и численного анализа согласования полученных расчетных ур
и экспериментальных уэ
данных (рис. 3) можно заключить о принципиальной адекватности полученной модели на основе тесной зависимости между расчетными и экспериментальными данными (так как точки на графике располагаются вблизи биссектрисы угла, образованного координатными осями и r> rкр
(0,95;6)
и а≈
1).
Исследование процесса с помощью модели
.
Поскольку результат проверки адекватности модели положителен, проведём серию расчётов по данной модели, варьируя ряд исходных данных.
Изменим коэффициент теплоотдачи α со 126,5 до 1000 Вт/(м2
·К):
Увеличим температуру среды до 200 °С (tsr
начальная = 20°С):
В 2 раза увеличим размеры слитка с 0,05×0,05×0,1 м3
до 0,1×0,1×0,2 м3
Вывод:
Мы ознакомились с методикой построения детерминированной модели и разработали численную модель металлургического процесса для системы с сосредоточенными параметрами. С помощью коэффициента корреляции мы доказали, что построенная модель является адекватной и отражает реальную картину охлаждения и кристаллизации заданного расплава.
При изменении некоторых параметров (коэффициента теплоотдачи α=1000 Вт/(м2
К), температура среды tsr
=200°C и размеры слитка 0,1×0,1×0,2 м3
) налицо изменения самой модели.
Получили, что при увеличении коэффициента теплоотдачи α примерно в 7,9 раз охлаждение расплава протекает быстрее во столько же раз. А также площадка, соответствующая перитектической реакции, стала ровной, что в большей степени отвечает реальному процессу, поскольку при нонвариантном превращении температура сплава остается неизменной.
При увеличении температуры среды III и IV стадии протекают дольше. Несмотря на то, что tsr
является параметром, который присутствует во всех стадиях, её изменение не оказало видимого влияния на замедление скорости протекания I и II стадий.
При увеличении размеров слитка значительно замедлился процесс кристаллизации сплава на всех стадиях. Здесь столь сильное влияние объясняется тем, что увеличив размеры, мы изменили сразу два параметра: поверхность теплообмена F и объем V, которые фигурируют на каждом шаге (на втором шаге V фигурирует дважды ).