УДК 004.942
Е. В. Старовойтова, А. Д. Галеев, С. И. Поникаров ЧИСЛЕННЫЙ АНАЛИЗ ВЛИЯНИЯ РАЗМЕРА КАПЕЛЬ НА РЕЗУЛЬТАТЫ РАСЧЕТА ЗОН
ПОРАЖЕНИЯ АММИАКОМ
Ключевые слова: численное моделирование, зоны токсического поражения, аммиак.
Проведен численный анализ влияния размера аэрозольный капель на размеры зон, ограниченных летальной и пороговой токсодозой, при аварийном выбросе сжиженного аммиака. Вычислительная процедура реализована c использованием программы FLUENT.
Keywords: numerical simulation, impact toxic zones, ammonia.
A numerical analysis was carried out of the droplet sizes influence on the characteristics of the toxic impact zones limited lethal and threshold toxic doses, in case of instantaneous release of liquefied ammonia. The computational procedure was implemented by the means of FLUENT program.
Введение
Моделирование и оценка последствий гипотетических сценариев аварий является ценным инструментом для разработки и реализации соответствующих защитных мероприятий. В связи с этим необходимо использовать надежные и проверенные математические модели, описывающие выбросы токсичных и горючих веществ в атмосферу и их распространение в приземном слое атмосферы. Полную математическую модель для прогнозирования последствий аварий можно грубо разделить на три части: моделирование источника, модель рассеивания и моделирование воздействия [1]. Множество работ было сфокусировано на моделировании рассеивания газа, тогда как моделированию источника уделялось недостаточно внимания.
При аварийных выбросах перегретых жидкостей (например, сжиженный аммиак в емкости под давлением) поступление вредных веществ в окружающую среду происходит в результате мгновенного вскипания, испарения аэрозольных капель и парообразования из пролива. На сегодняшний день для описания источника жидкости с поверхности пролива применяют два подхода: упрощенные алгебраические формулы и метод вычислительной гидродинамики (CFD).
Применение метода вычислительной гидродинамики (CFD) и пристеночных функций для описания источника [2,3] позволяет выполнить совместно расчеты по дисперсии газа и парообразованию путем численного решения дифференциальных уравнений переноса массы, импульса, энергии, примесей и турбулентных характеристик.
В данной работе предлагается модель для определения зоны токсического поражения при залповом выбросе сжиженного аммиака на основе программного комплекса FLUENT. Проведен анализ влияния размера аэрозольных капель на результаты расчета зон поражения аммиаком на основе разработанной модели.
Математическая модель
В работе предполагается, что аммиак хранится в резервуаре при температуре выше нормальной температуры кипения. При разгерметизации емкости со сжиженным газом в результате резкого снижения давления и нарушения
термодинамического баланса происходит мгновенное вскипание определенной части жидкой фазы. Доля мгновенно вскипающей жидкости при адиабатическом расширении находится из выражения [4]:
(
x u = 1 - exp
CP,liq (To Tb )
Л
AHg
э у
где Оруд - удельная теплоемкость жидкости, Дж/(кг-К); То - температура сжиженного газа в емкости, К; Ть - температура кипения жидкости, К; ДИэ - удельная теплота парообразования при температуре кипения и атмосферном давлении, Дж/кг.
В результате мгновенного вскипания расширяющиеся пары диспергируют и увлекают часть жидкости, поэтому образующееся облако содержит смесь пара и жидкости в виде аэрозольных капель. Доля капель принимается равной доле мгновенного вскипания жидкости.
Оставшаяся после мгновенного вскипания жидкая фаза разливается на подстилающей поверхности и переходит в состояние кипения. В режиме кипения скорость испарения зависит от разности температур между жидкостью и окружающей средой. Скорость испарения кипящей жидкости регулируется тепловыми потоками от подстилающего твердого слоя и от атмосферы. При температуре жидкости ниже точки кипения процесс протекает в режиме испарения. В этом случае испарение происходит из-за разницы (парциальных) давлений пара на поверхности пролива и в атмосфере. Скорость испарения некипящей жидкости в основном зависит от скорости, при которой пар может перемещаться под действием ветрового потока над проливом. В модели условие перехода между стадиями кипения и испарения определяется как:
W=
q + q,rd+4s - c + Cr
^H
при(Я, + cgrd+C - q, + Ca > Jg5 • Дф
и\, = ть;
приТ < Ть
где W — интенсивность испарения пролива, кг/(м2-с); Ча — тепловой поток из атмосферы, Вт/м2; ЧдгС — поток тепла от грунта к жидкости, Вт/м2; — поток тепла от солнечной радиации, Вт/м2; ЧР — тепловой поток, излучаемый поверхностью пролива, Вт/м2; Чаг — тепловой поток к проливу вследствие излучения атмосферы, Вт/м2; ид5— массовый поток пара при испарении, кг/(м2-с).
В режиме кипения температура жидкости постоянна, при испарении происходит изменение температуры жидкости в соответствии с уравнением:
5TliC
at
Ca + Cgrd + Cs - Cp + Car - Jg,s • ^Hg CP,liC • miiq
Тепловой поток, подводимый от грунта к жидкой фазе, решался, используя аналитическое уравнение. В этом случае Чдгс рассчитывается как функция координат и времени согласно теории теплопередачи через полубесконечное твердое тело:
Cgrd(x,y,t) =
^ grd (Tgrd Tliq ) V^agrd ((tw (Х,У))
(1)
где Хдгс и адгс - теплопроводность и температуропроводность подложки, Т0дгС -начальная температура подложки, ТуЧ - температура пролива и ^ - момент времени, при котором начинается соприкосновение жидкости с твердой поверхностью.
Тепловой поток от подложки рассчитывался по формуле (1), но при этом значения температуры жидкости находились из численного решения.
Тепловой поток из атмосферы вычислялся с помощью пристеночных функций [5]:
Ca = Т + =
(C TP )CpPC°'25kP'5 .
Т+
y + < У Т
|Pr- y+
|Pft(u + + Pt) у + > УТ
PT = 9.24
Г Pr Л3'4
Prt
-1
[1 + 0.28exp(- 0,007 Pr/Prt)]
где Тр - температура в ячейке, прилегающей к стенке, К; Ср - теплоемкость паровоздушной смеси, Дж/(кг-К); р - плотность паровоздушной смеси, кг/м3; Сц - коэффициент; Т+- безразмерная
температура; кр - турбулентная кинетическая
2/2 +
энергия в ячейке, прилегающей к стенке, м /с ; у -безразмерное расстояние от стенки до первого расчетного узла; уТ+ - безразмерная толщина
теплового подслоя; Рг и Рг - молекулярное и турбулентное числа Прандтля соответственно; и -безразмерная скорость; рТ - функция, учитывающая сопротивление теплового подслоя к теплопередаче.
Массовый поток аммиака с поверхности пролива 1в,8 в результате испарения рассчитывался на основе стандартных функций стенки с учетом поправки на конвективный поток:
Jg,s = K
(Yg,s - Yg,P)pC0'25kp5
Y+
Y + =
|Sc•y+
(У + < Ус)
|SCt(u ++ Pc) (y +> yC)
n3/4
Pc = 9.24
Sc Sct
-1
[1 + 0.28exp(- 0.007Sc/Sct )],
K =
41 - Yg,P V(1 - Yg,s))
g,P_
Y - Y
g,s 1 g,P
где Ygs - локальная массовая доля на границе жидкость-газ. Yg,P - массовая доля в ячейке, прилегающей к стенке. Y - безразмерная массовая доля. yc - безразмерная толщина диффузионного подслоя. Sc и Sc - молекулярное и турбулентное числа Шмидта соответственно. Pc - функция, учитывающая сопротивление диффузионного подслоя к теплопередаче. K - поправка на конвективный поток.
Наиболее подробное описание данной модели можно найти в работах [6-9]. В данной работе тепловыми потоками Cs , Cp, Car пренебрегали.
Результаты расчета
Проведено численное исследование влияния размера аэрозольных капель на размеры зон поражения для случая выброса аммиака в атмосферу при двух скоростях ветра 1 м/с и 5м/с. Ранее в работах [4,6-9] распределение капель по размерам в источнике задавалось с использованием соотношения Розин-Раммлера. В настоящей работе диаметры капель принимались постоянными и равными 10 и 100 мкм. Выбор размера аэрозольных капель был сделан на основе результатов экспериментов, приведенных в [10]. Результаты расчета представлены на рис.1 а,б.
При низких скоростях ветра (1м/с) размеры капель влияют незначительно на формы и протяженность как зоны летального воздействия, так и на зоны, ограниченные пороговой токсодозой. Однако, при высоких скоростях ветрового потока наблюдается заметная разница для зоны летального поражения: длина зоны летального поражения больше, а ширина - меньше при daэp=10мкм, чем при daэp=100мкм. Это связано с более быстрым испарением капель меньшего размера и соответственно с более интенсивным рассеиванием первичного облака. Однако, вследствие сильного гравитационного растекания в начальный период формирования облака ширина пороговой зоны при 10мкм больше, чем при 100мкм.
J
»—.
\ \ !
2
if/
li
-200 -100 0 100 200 300 400 500 600 700 300 х,м
б
Рис. 1 - Изолинии токсодоз летальной (1) и пороговой (2) при 1м/с (а) и 5м/с (б): -----0аэр=10мкм; — - 0аэр=100мкм
поражения токсичным газом проводится, в основном, при низких скоростях ветра, то аэрозольные капли можно принимать одинакового диаметра. При высоких скоростях ветра к выбору размера капель нужно подходить более тщательно.
Литература
1. T. Vik, B.A. Petterssxrn Reif . FFI-rappоrt 2010/00254, ^rwegian Defence Research Establishment (FFI) (2010).
2. В. Маршалл, Основные опасности химических производств. Мир, Москва, 1989, 672с.
3. G. Desоutter, C. Habchi, B. Cuernt, Т.РотзЫ;, Int. J. HeatMass Transf, 52, 6028-6041 (2009).
4. A.D. Galeev, E.V. Starоvоytоva, S.I. Роткатоу, Pmcess Saf. Enviran. Pmt, 91, 191-201 (2013a).
5. FluentInc. Fluent 6.1 // User's Guide.—Lebanon, 2003.
6. Е.В. Старовойтова, А.Д. Галеев, С.И. Поникаров, Вестн. Казан. технол. ун-та, 13, 175-179 (2011).
7. Е.В. Старовойтова, А.Д. Галеев, С.И. Поникаров, Безопасность труда в промышленности, 12, 47-51 (2011).
8. . Е.В. Старовойтова, А.Д. Галеев, С.И. Поникаров, Вестн. Казан. технол. ун-та, 15, 6, 207-209 (2012).
9. Е.В. Старовойтова, А.Д. Галеев, С.И. Поникаров, Вестн. Казан. технол. ун-та, 16, 17, 259-261 (2013).
10. С.Е. Якуш, Дисс. докт. физ.-мат. наук, Ин-т проблем механики РАН, Москва, 2000. 336 с.
Вывод
При низких скоростях ветра размеры капель влияют незначительно на размеры зон токсического поражения. В связи с тем, что оценка последствий
© Е. В. Старовойтова - к.т.н., ассистент каф. машин и аппаратов химических производств КНИТУ, [email protected]; А. Д. Галеев - к.т.н., доцент той же кафедры; С. И. Поникаров - д.т.н., проф., зав. каф. машин и аппаратов химических производств КНИТУ.
© E. V. Starovoytova - Ph. D., assistant of Department of Mechanical Engineering, Kazan National Research Technological University, [email protected]; A. D. Galeev - Ph. D., assistant of Department of Mechanical Engineering, Kazan National Research Technological University; S. I. Ponikarov - doctor of technical Sciences, Professor, head of chair Department of Mechanical Engineering, Kazan National Research Technological University.