CC BY
309
ПРОЦЕССЫ И АППАРАТЫ ХИМИЧЕСКОЙ ТЕХНОЛОГИИ
УДК 614.83
А. Д. Галеев, Е. В. Старовойтова, С. И. Поникаров
ДИНАМИКА ФОРМИРОВАНИЯ ВЗРЫВООПАСНОГО ОБЛАКА ПРИ АВАРИЙНОМ ВЫБРОСЕ СМЕСИ СЖИЖЕННЫХ УГЛЕВОДОРОДНЫХ ГАЗОВ В АТМОСФЕРУ
Ключевые слова: сжиженный углеводородный газ, интенсивность парообразования, распространение примеси,
численное моделирование.
Проведено численное моделирование аварии с выбросом смеси пропан-бутан. Показано влияние метеорологических условий на изменение во времени массы горючего газа в облаке, ограниченном концентрационными пределами воспламенения.
Keywords: liquefied hydrocarbon gas, evaporation rate, pollutant spread, numerical simulation.
The numerical simulation of accident with propane-butane mixture emission was conducted. The influence of meteorological conditions on changing in time the mass of explosive vapours in cloud which bounded with concentration limits of flame propagation was shown.
Объекты, на которых хранятся, перерабатываются и транспортируются сжиженные углеводородные газы (СУГ), являются источниками повышенной опасности. Опасность таких объектов связана с интенсивным испарением СУГ при его выбросе из оборудования и образованием взрывоопасных облаков, способных перемещаться под действием ветра на значительные расстояния и воспламеняться от различных источников.
Выбор способов и методов повышения безопасной эксплуатации объектов использования СУГ должен проводиться на базе комплексного анализа и прогноза вероятных аварийных ситуаций и сопровождающих их физических процессов. Важной задачей при этом является расчет массы опасного вещества, способного участвовать в горении.
При разгерметизации емкости с перегретой жидкостью в результате резкого снижения давления и нарушения термодинамического баланса происходит мгновенное вскипание определенной части жидкой фазы за счет высвобождения внутренней энергии. При этом устанавливается новое состояние равновесия, а температура оставшейся жидкой фазы понижается до температуры кипения при атмосферном давлении. При мгновенном испарении жидкости формируется первичное облако пара. Расчет доли мгновенно вскипающей жидкости при адиабатическом расширении проводился с использованием уравнения состояния Соаве-Редлиха-Квонга.
В результате мгновенного вскипания расширяющиеся пары диспергируют и увлекают часть жидкости, поэтому образующееся облако содержит смесь пара и жидкости в виде аэрозольных капель. Масса выброшенных капель жидкости принималась равной массе мгновенно образовавшейся паровой фазы [1].
Для описания движения аэрозольных капель используется лагранжев подход, согласно которому вся дисперсная фаза разбивается на большое число групп капель, в пределах каждой из которых параметры всех физических капель считаются одинаковыми, и, следовательно, эволюция каждой группы может быть прослежена путем расчета движения лишь одной представительной капли [2].
Фазовые переходы в объеме аэрозольного облака влияют на его состав, температурный режим и рассеяние в атмосфере.
Обратное влияние дисперсной фазы на несущий поток, обусловленное межфазным обменом теплотой, импульсом и массой, учитывается включением соответствующих источниковых членов в уравнения переноса энергии, импульса, примеси и в уравнение неразрывности для сплошной фазы.
Оставшаяся после мгновенного вскипания часть сжиженного газа, образуя пролив, испаряется вследствие интенсивного теплообмена с подстилающей поверхностью и атмосферой.
Изменение содержания ьго компонента в жидкости описывалось уравнением:
_ дп + да _ з •м,
dx ДИ •^у.з • м., (1)
I
2 2 где дп - тепловой поток от грунта, Вт/м ; да - тепловой поток из атмосферы, Вт/м ; ДИ - теплота
парообразования, Дж/кг; У.з - равновесная мольная доля компонента в паровой фазе при
температуре жидкости Т|; М. -молекулярная масса 1-го компонента, кг/кмоль; Л-шаг по времени.
В модели предполагалось, что при пузырьковом режиме кипения теплообмен между твердым подстилающим слоем и жидкостью лимитируется переносом тепла в твердой фазе [3]. Исходя из этого, тепловой поток, подводимый от твердой поверхности к жидкой фазе, дп = ^(дТп/^Оу=о определялся из численного решения трехмерного нестационарного уравнения теплопроводности для твердого подстилающего слоя:
С Р ЗТ^
Сп Рп ~ ^ п
ЗХ
йх2 Зу2 Зг2
(2)
где Тп(х,у,2) - распределение температур в подстилающем слое; Сп, рп, ^п - теплоемкость, плотность и коэффициент теплопроводности подстилающего слоя соответственно.
При решении уравнения (2) на твердой поверхности, прилегающей к проливу, задавалось условие Тп=Т|.
Вторая составляющая теплового потока да рассчитывалась при решении задачи распространения облака неизотермичного тяжелого газа в приземном слое атмосферы по найденному температурному полю в области источника с помощью пристеночных функций [2].
Изменение температуры жидкой фазы со временем происходит в соответствии с изменением состава жидкости согласно кривой фазового равновесия. Данная кривая для смеси пропан-бутан изображена на рис. 1. Представленная диаграмма получена с помощью программы СИетСАО, где константы парожидкостного равновесия рассчитаны из уравнения состояния Соаве-Редлих-Квонга. На рис. 1 также приведена кривая фазового равновесия (показана в виде полых точек), полученная из закона Рауля. Температура кипения смеси определялась последовательным приближением по уравнению изотермы жидкой фазы [4]:
N р
X Хк_ 1; к_ -р", (3)
М Р0
где к - константа фазового равновесия; X. - мольная доля 1-го компонента; Ру,| - давление насыщенных паров при температуре Т|, Па; Р0 - давление окружающей среды: Р0 = 101325 Па.
Давления насыщенных паров пропана и бутана определялись из уравнения [5]:
( В ^
Р,,. _ ехР А, + Т + С1|п(Т|) + 0;Т|Е. , (4)
V Т| У
где А., В., С., й., Е. - коэффициенты. Значения коэффициентов даны в таблице 1.
Таблица 1 - Значения коэффициентов для расчета давления насыщенных паров [5]
Компонент Д| В| С| й| Е|
Пропан 54,276 -3368 -5,261 8,6-10-6 2
Бутан 67,164 -4386,2 -7,1805 9,7145-10-6 2
Рис. 1 - Кривая фазового равновесия смеси пропан-бутан: сплошная линия - константы равновесия определены из уравнения SRK; полые точки - из закона Рауля
Модель распространения примесей в атмосфере включает численное решение полной системы трехмерных нестационарных уравнений Рейнольдса, переноса массы и энергии, замыкаемых Realizable k-s моделью турбулентности [2].
Для учета дополнительного нагрева воздуха вследствие конденсации водяного пара использовалась функция источникового члена:
в уравнении энергии:
Q _ V(T)' Cv 'Р' rw (5)
Qv _ Ат ,
в уравнении переноса компонента (паров воды):
S ._—. (6,
где «(Т) - функция, определяющая долю сконденсированной влаги:
0, T > Ts
С v — Cv s
Ф(Т) _j vC v,s, 273,15 < T < Ts , (7)
Cv
1, T < 273,15 К
где Ts - температура точки росы при данном парциальном давлении паров воды в паровоздушной смеси; Cv - концентрация паров воды в воздухе, кг/кг; Cv,s - концентрация водяного пара, при которой он достигает насыщения при данной температуре, кг/кг; rw -теплота парообразования воды, Дж/кг.
Результаты расчетов
В целях оценки применимости изложенной модели рассматривалась аварийная ситуация, связанная с полным разрушением автоцистерны, содержащей сжиженный углеводородный газ (СУГ). Масса жидкости в цистерне принималась равной 18000 кг. Состав СУГ: пропан - 50% (мас.), бутан - 50% (мас.). Температура СУГ в цистерне принималась равной 303 К.
Расчетные характеристики источников первичного и вторичного облаков в начальный момент времени представлены в таблице 2.
Первичное облако представлялось в виде цилиндра с радиусом равным высоте.
2
Площадь пролива (м ) жидкости рассчитывалась по формуле [6]:
F* _ fP ' V*, (8)
1 1
где fP - коэффициент разлития, м- (150 м- при проливе на бетонное или асфальтовое покрытие); V* - объем жидкости, м3.
Таблица 2 - Характеристики источников выбросов
Параметр выброса Значение
Масса бутана в первичном облаке, кг 1197
Масса пропана в первичном облаке, кг 4144
Масса аэрозолей в облаке, кг 5341
Масса жидкой фазы, пролитой на твердую поверхность, кг 731S
2 Площадь пролива, м 1S2Q
Диаметр капель, мкм 1QQ
Начальная температура капель и жидкости в проливе, К 242,5
Массовая доля бутана в проливе и каплях Q,6S
Массовая доля пропана в проливе и каплях Q,32
Радиус и высота первичного облака, м 12
Расчеты проводились для значения скорости ветрового потока в невозмущенной атмосфере 1,25 м/с при неустойчивых, нейтральных и устойчивых условиях стратификации атмосферы.
Устойчивая стратификация (инверсия) соответствует повышению потенциальной температуры с высотой, неустойчивая (конвекция) - понижению температуры с высотой, равновесная (изотермия) - градиент температуры по высоте равен адиабатическому (примерно
0,6 К на 100 м высоты). Поскольку характерные значения вертикальных градиентов температур в пограничном слое атмосферы значительно превосходят величину адиабатического градиента, то в приземном слое можно не различать обычную и потенциальную температуры.
Степень устойчивости атмосферы учитывалась путем задания соответствующего распределения температуры в слое грунта.
Процесс испарения капель имеет короткую продолжительность и завершается по истечении примерно 30 с. При конвекции скорость испарения капель ниже, чем при изотермии и инверсии. В последних случаях наиболее сильно сказывается влияние отрицательной плавучести облака, приводящей к увеличению скорости смешения с теплым внешним воздухом и уменьшению времени испарения.
На рис. 2 показаны временные зависимости отношения текущей массы жидкости к первоначальной. Для сравнения приведены данные, полученные из эмпирического выражения [7]:
m _ 20-Г>
ДH
( \------------------------ с -л
2 I т + 5,WReX aт
2Л"Л “Т +--------^----
У ла d
v у
(9)
2
где тисп - удельная масса испарившейся жидкости, кг/м ; То - начальная температура
материала, на поверхность которого разливается СУГ, соответствующая расчетной
температуре, К; а - эффективный коэффициент температуропроводности материала, на
поверхность которого разливается СУГ, м /с; Ре=(и-Ь)/уа - число Рейнольдса (и - скорость воздушного потока, м/с; Ь — характерный размер пролива СУГ, м; va— кинематическая вязкость воздуха при расчетной температуре, м /с); Ха —коэффициент теплопроводности воздуха при расчетной температуре, Вт/(м^К).
Рис. 2 - Изменение относительного количества жидкости со временем: 1— изотермия, 2— конвекция, 3—инверсия, 4—значения, вычисленные согласно [7]
Полное испарение жидкой фазы, согласно формуле (9), происходит в течение 180 с, а исходя из численной модели - приблизительно в течение 600 с. Такое расхождение объясняется тем, что формула (9) не учитывает изменение температуры кипения жидкой фазы.
Горючие пары, смешиваясь с воздухом, образуют взрывоопасные смеси. Обозначим массу газа во взрывоопасных пределах, как тв, а отношение тв к общему количеству горючего газа как 7. Согласно [8], зону опасности определяли по половине нижнего концентрационного предела воспламенения (НКПР) вследствие турбулентного характера пролива и флуктуационного изменения параметра течения. Значение нижнего предела воспламенения НКПР принималось равным 1,8 % (об.). На рис. 3 представлены временные зависимости доли горючего газа, способного к воспламенению.
При инверсии и изотермии величина 7 увеличивается со временем в результате смешения горючего газа с воздухом до концентраций ниже ВКПР (верхнего концентрационного предела воспламенения), достигая максимума, а затем снижается. При конвекции характер изменения доли горючего газа, способного к воспламенению, существенно отличается: величина 7 начинает резко снижаться. В целом полученные результаты согласуются с существующим представлением о влиянии степени устойчивости атмосферы на процесс рассеяния примеси в атмосфере. При неустойчивой стратификации интенсивное перемешивание примеси обуславливает более резкое снижение концентраций горючего газа со временем, и по истечении примерно 900с концентрации в пределах облака снижаются до безопасных значений. При инверсии происходит более медленное снижение 7 со временем, чем при конвекции и изотермии, что обусловлено более неблагоприятными условиями для рассеяния газа. Ухудшение рассеяния при инверсии связано не только с подавлением турбулентности, но и с более низкой скоростью ветра в непосредственной близости от подстилающей поверхности, чем при изотермии и конвекции. Кроме этого, при инверсии доля горючего газа, способного участвовать в горении, превышает значение рекомендуемое для расчета (7=0,1) и по истечении 3600 с составляет 0,55.
z
0.8
0.7
0.6
0.5
0.4
0.3
0.2
0.1
0
0
600
1200 1800 2400 3000
т, с
Рис. 3 - Временные зависимости доли горючего газа в облаке, ограниченного 0,5НКПР: 1— изотермия, 2—конвекция, 3 — инверсия
Представлена математическая модель, описывающая основные физические явления, которые определяют динамику формирования взрывоопасного облака при аварийном выбросе смеси сжиженных углеводородных газов:
• мгновенное вскипание жидкости вследствие падения давления при полном разрушении оборудования;
• испарение аэрозоля в облаке;
• парообразование при кипении пролива с учетом изменения состава жидкой фазы, теплообмена с подстилающей поверхностью и атмосферой;
• рассеяние тяжелого газа в атмосфере.
Достоверное прогнозирование массы горючего газа в облаке, ограниченном концентрационными пределами воспламенения, имеет большое значение при оценке потенциала взрывоопасности технологических блоков и определения ожидаемых последствий аварий в рамках процедуры анализа риска.
1. Маршалл, В. Основные опасности химических производств / В. Маршалл. - М.: Мир, 1989. - 672 с.
2. Fluent Inc. Fluent 6.1. User’s Guide, Lebanon, 2003.
3. Старовойтова Е.В. Математическая модель процесса парообразования при аварийном проливе сжиженного газа / Поникаров С.И., Галеев А. Д., Старовойтова Е.В. // Вестник Казан. технол. ун-та. -№ 1. - 2011. - С. 323-325.
4. Брусиловский А.И. Фазовые превращения при разработке месторождений нефти и газа. - М.: Грааль, 2002. - 575 с.
5. ChemCad 5.0. User’s Guide, Chemstations Inc.
6. Методика определения расчетных величин пожарного риска на производственных объектах (утв. приказом МЧС РФ от 10 июля 2009 г. № 404).
7. ГОСТ Р 12.3.047-98. ССБТ. Пожарная безопасность технологических процессов. Общие
требования. Методы контроля.—Введ.—2000—01—01. —М.: Изд-во стандартов, 1998.
8. РД-03-26-2007. Методические указания по оценке последствий аварийных выбросов опасных веществ (утв. Приказом Федеральной службы по экологическому, технологическому и атомному надзору от 14 декабря 2007 г. №859).
© А. Д. Галеев - канд. техн. наук, доц. каф. машин и аппаратов химических производств КГТУ,
mahp_kstu@mail.ru; Е. В. Старовойтова - асп. той же кафедры; С. И. Поникаров - д-р техн. наук,
проф., зав. каф. машин и аппаратов химических производств КГТУ.
Заключение
Литература
