Курсовая работа - Математическое моделирование в системе конечноэлементностных расчетов Femlab

Подождите немного. Документ загружается.

1.2. Практическая часть

Для моделирования выберем теплообменник типа «труба в трубе». Смоделируем

его в системе FEMLAB и получим кривые переходных процессов. В качестве

реального прототипа выберем процесс нагрева сырья колонны К-101 блока ВТ

установки У-731. В качестве теплоносителя выберем перегретый водяной пар.

При моделировании приняты следующие допущения:

1) теплопроводность стенок труб велика по сравнению с теплопроводностью

среды внутри теплообменника

2) пар достаточно перегрет для того, чтобы вклад процесса кондесации греющего

пара в теплообмен был мал по сравнению с процессом охлаждения

3) конвективные потоки, вызванные изменением плотности жидокости при

нагреве, слабы по сравнению с ее движением вдоль теплообменника

Создадим модель в системе FEMLAB, для чего выберем из предложенного списка

моделей heat transfer, convection and conduction, transient analysis. Укажем, что

пространство двумерное (2D).

Примем следующие геометрические размеры теплообменника: длина 12 м, радиус

внутренней трубы 0,5 м, внешней трубы 1,5 м. После построения получаем

следующее:

Рис. 1. Модель теплообменника типа «труба в трубе»

Арабскими цифрами помечены границы, римскими – области с одинаковыми

параметрами среды, или домены.

Переходим к заданию параметров областей (physics, subdomain settings). Во

внешней рубашке (области I и III) идет перегретый пар, начальная температура

которого 250

С. Смоделируем процесс, когда в разогретый теплообменник подают на

нагрев холодное сырье с температурой 20

С. Физические параметры сырья

(стабильного гидрогенизата фракция НК – 350

С) предположим близкими к

параметрам смеси предельных и ароматических углеводородов С

– С

. Учитывая, что

время задаем в минутах для удобства, вводим следующие параметры:

Таблица 1.

Параметр Размерность Область I Область II Область III

Время-масштабный

коэффициент 

безразмерный 1 1 1

Теплопроводность k

(изотропная)

Кминм

Дж



1,75 7 1,75

Плотность 

кг

8 800 8

Теплоемкость C

Ккг

Дж



2000 1880 2000

Внутреннее тепловыделение

Дж

0 0 0

Скорость в горизонтальном

направлении u

мин

0,27 0,42 0,27

Скорость в вертикальном

направлении v

мин

0 0 0

Начальная температура T(0) К 523 523 523

Скорости потоков рассчитываются исходя из расхода пара 50 м

/час, а сырья 10

/час.

Далее зададим граничные условия (physics, boundary settings). Учтем потери тепла

в окружающую среду через стенки внешней трубы (границы 4 и 5). Принимая

температуру окружающей среды 20

С, потери с достаточной точностю могут быть

описаны формулой:

минм

Дж

TTq

tSQ

Км

Вт











60)293(75,014,3))273(058,03,9(

058,03,9



Таблица 2.

Граница Вид условия T

, K q

минм

Дж



1, 3 Температура 523 -

2 Температура 293 -

4, 5 Тепловой поток -

)293())273(058,03,9(3,141  TT

6, 7, 8 Конвективный теплообмен - -

9, 10 Расчетная граница - -

Зададим параметры решения задачи (solve, solver parameters). Выберем решение,

зависимое от времени (time dependent) и время до 50 мин с шагом 0,1 мин.

Относительную и абсолютную точность зададим 0,0001. Все остальные параметры

оставим без изменения.

Инициализируем сетку (mesh, initialize mesh). Предложенная сетка весьма грубая,

поэтому сделаем ее более мелкой (mesh, refine mesh). Чем мельче элементы сетки, тем

точнее решение, но тем больше ресурсов оно требует. Я остановился на сетке из 3456

областей. Статистику разбиения можно увидеть в меню mesh, mesh statistics. Внешний

вид сетки показан на рисунке.

Рис. 2. Разбиение модели на расчетные области

Запускаем решение (solve, solve problem). Через некоторое время мы можем

увидеть результат решения – распределение температур в теплообменнике.

Рис. 3. Распределение температур внутри теплообменника

Как видно, FEMLAB предоставляет возможность наглядной демострации

результатов расчета. Однако, необходимо представить результаты моделирования в

виде кривой разгона. Для этого воспользуемся сопряжением со средой MATLAB 6.5.

Запустив FEMLAB в связке c MATLAB, экспортируем структуру с параметрами

модели в MATLAB (file, export, FEM structure).

Выходом системы будем считать среднее значение температуры сырья в

выходном сечении трубы. Для анализа результатов моделирования составим

небольшую программу на встроенном языке MATLAB.

filename: get_h.m

function T = get_h (fem, t, step_t, x, y, r, step_r)

T = zeros (length (t(1):step_t:t(2)), 1);

for i = 0:step_r:r

if (i == 0) | (i == r)

k = 0.5;

else

k = 1;

end

T = T + postinterp (fem, 'T', [x y+i]', 'T', t(1):step_t:t(2))*i*k*step_r;

end

T = 2 * T / r^2;

plot (t(1):step_t:t(2), T);

Запустим ее командой

T = get_h (fem, [0 50], 0.1, 6, 0, 0.25, 0.01);

Получим массив T, содержащий в себе переходной процесс в системе, заданный с

шагом 0,1 мин. Построим переходной процесс.

Рис. 4. Переходной процесс в системе при подаче сырья

Как видно, на графике есть всплески температуры, не могущие иметь места в

реальности и обусловленные ошибками моделирования. Для их сглаживания

воспользуемся фильтром скользящего среднего.

filename: mean_filter.m

function Signal_filtered = mean_filter (Signal_noised, m)

Signal_filtered = [];

for i = 1:(length (Signal_noised) - m)

s = [];

for j = 0:m - 1

s = [s Signal_noised(i+j)];

end

Signal_filtered = [Signal_filtered mean(s)];

end

for i = (length(Signal_noised) - m + 1):length (Signal_noised)

s = [];

for j = 0:(length(Signal_noised) - i)

s = [s Signal_noised(i+j)];

end

Signal_filtered = [Signal_filtered mean(s)];

end

После фильтрации получаем следующий переходной процесс.

Рис. 5. Переходной процесс в системе при подаче сырья (фильтрованный)

Как видно, процесс в теплообменнике при подаче сырья без большой

погрешности можно считать процессом идеального вытеснения, поскольку чистое

запаздывание на порядок превышает длительность переходного процесса без

запаздывания.

Смоделируем переходной процесс, происходящий в системе при подаче большего

количества греющего пара. Для этого изменим время решения на 250 мин. с шагом 0,1

мин., а также скорость движения в доменах I и III (расход пара). Зададим повышение

расхода пара с 50 до 250 м

/час на 150-й минуте моделирования выражением:

u = 0.27 + 1.08 * flsmhs(t – 150, 1)

После завершения моделирования и экспорта в MATLAB, выполним команду

T = get_h (fem, [150 250], 0.1, 6, 0, 0.25, 0.01);

Получаем следующую кривую.

Рис. 6. Переходной процесс при подаче дополнительного греющего пара

Как видно, переходной процесс при воздействии на этот канал в корне отличается

от предыдущего, что неудивительно. Здесь, напротив, транпортное запаздывание мало,

а переходной процесс более длительный. Передаточная функция может быть

аппроксимирована, например, по методу Симою. После аппроксимации получаем:

17,262,31504,2062

0071,0

)(





ppp

Построим исходную и аппроксимированную кривые разгона. На рисунке

исходная кривая показана красным.

Рис. 7. Исходный и аппроксимированный переходной процесс при подаче

дополнительного греющего пара

2. Моделирование химической реакции

2.1. Теоретическая часть

Рассмотрим моделирование химического реактора, в котором происходит процесс

частичного окисления о-ксилола до фталевого ангидрида. Процесс происходит в

многотрубном реакторе при температуре порядка 400

С на катализаторе из оксида

ванадия и сульфата калия, наеасенных на основу из диоксида кремния. Время

протекания реакции колеблется от 0,5 до 5 с. Наиболее важным параметром в процессе

является температура в реакторе, поскольку реакции, проходящие в реакторе, сильно

экзотермичны, и при нарушении охлаждения реактора он теряет устойчивость.

Распределение температуры внутри реактора влияет на его КПД и выход фталевого

ангидрида, который, естественно, необходимо повысить. Температуру можно

контролировать путем варьирования диаметра трубок реактора, скорости потока в

трубках, температуры стенок трубки (степени охлаждения).

Поскольку трубки имеют круглое сечение и идентичны друг другу, мы имеем

право воспользоваться осесимметричной моделью лишь одной трубки.

Химизм процесса, протекающего в реакторе, достаточно сложен, и не может быть

здесь детально рассмотрен, но упрощенно его можно представить как совокупность 3

реакций:

1) окисление о-ксилола до фталевого ангидрида

2) горение о-ксилола с образованием смеси оксидов углерода

3) горение фталевого ангидрида с образованием смеси оксидов углерода

Как видно, целевая реакция – первая. Поскольку процесс ведется при

значительном избытке кислорода, можно считать, что реакции первого порядка.

Обозначив константы скорости реакций k

, k

и k

соотвественно, запишем уравнения:

 

BAAC

BAAB

AAA

xkxkyyr

xkkyyr

2300

2100

3100

)1(

)]1)((







, где

– молярная доля кислорода

– молярная доля о-ксилола

– общая скорость расхода о-ксилола

– скорость превращения о-ксилола во фталевый ангидрид

– скорость образования оксидов углерода

Коэффициенты скорости реакций в соответствии с законом Аррениуса

определяются из выражения:

eAk









, где

– температура стенки реактора

TTT 



- перегрев среды относительно стенки

Рис. 8. Модель трубки реактора

Определение уравнения в частных производных

Уравнения, описывающие данный процесс, выглядят следующим образом.















iibPgSeff

Atot

iSieff

rHTCuT

xuxD





)(

, где

– поверхностная скорость



– плотность газа



– плотность катализатора

tot

– общая концетрация



eff

– эффективная теплопроводность среды

Учтя, что x

= x

+ x

, решаем уравнения только для x

и x

. После подстановок и

упрощений получаем:

)()(

)(

CBbPgSeff

Atot

CSCeff

Atot

BSBeff

rHrHTCuT

xuxD

















Уравнения скорости реакций:

 

BCBAC

BCBAB

xkxxkyyr

2300

2100

)1(





Граничные условия

eff







































Полагаем, что на выходе реактора определяющее влияние имеют конвенкция и

теплопередача.

Задание модели

Ввиду того, что длина трубки реактора значительно превосходит ее диаметр,

масштабируем уравнения для большей наглядности моделирования. В данном случае,

длина превышает диаметр более чем в 200 раз – трубки имеют в длину 3 метра при

диаметре 1,27 см. Масштабирование также избавит нас от проблем при формировании

сетки.

scale

r 

После масштабирования необходимо внести изменения в коэффициенты

уравнений. Наиболее сильно это отразится на коэффициентах диффузии и

теплопередачи, которые становятся анизотропными.





























scale

