CC BY
62
МОДЕЛИРОВАНИЕ ПОЖАРОВ
УДК 614.83
МОДЕЛИРОВАНИЕ ФОРМИРОВАНИЯ ВЗРЫВООПАСНОГО ОБЛАКА ПРИ ИСПАРЕНИИ С ПОВЕРХНОСТИ АВАРИЙНОГО ПРОЛИВА НЕФТИ
Представлена модель, описывающая процесс испарения с поверхности аварийного пролива нефти. Модель также позволяет определять изменение во времени и пространстве массы взрывоопасных паров в облаке, ограниченном концентрационными пределами распространения пламени. Ключевые слова: испарение нефти, численное моделирование, взрывоопасное облако.
А. Д. Галеев
канд. техн. наук, доцент, Казанский государственный технологический университет, г. Казань, Республика Татарстан
С. И.Поникаров
д-р техн. наук, профессор, Казанский государственный технологический университет, г. Казань, Республика Татарстан
Обеспечение безопасности взрывопожаро-опасных промышленных объектов предполагает использование процедуры анализа и оценки риска аварий, основу которой составляет определение последствий и вероятности нежелательных событий. Один из важных этапов оценки риска аварий на таких объектах — определение количества опасного вещества, способного участвовать в образовании взрывоопасной смеси. Для решения данной задачи необходимо знание количественных характеристик процесса испарения со свободной поверхности. Многие горючие жидкости, обращающиеся на производстве, представляют собой сложные смеси, состоящие из большого количества компонентов. Процесс испарения с поверхности аварийных розливов таких жидкостей носит нестационарный характер, обусловленный изменением со временем состава жидкой фазы. К тому же значительное влияние на интенсивность испарения будет оказывать динамика потока над поверхностью пролива из-за образования паровоздушной смеси с плотностью, отличающейся от плотности воздуха, и взаимного влияния процессов испарения и дисперсии примеси в атмосфере. Имеющиеся эмпирические корреляции и аналитические формулы не учитывают перечисленных факторов, что обуславливает необходимость совершенствования математической модели испарения.
В данной статье представлен комплекс моделей, позволяющих рассчитывать массовый поток паров при испарении со свободной поверхности пролива нефти и определять массу горючего вещества в об-
лаке с суммарной концентрацией горючих компонентов в диапазоне между нижним (НКПР) и верхним (ВКПР) концентрационными пределами распространения пламени. Модель базируется на численном решении трехмерных нестационарных уравнений гидродинамики и тепломассообмена.
Математическая модель испарения нефти
При разработке модели испарения нефти были сделаны следующие допущения:
• смесь компонентов идеальна;
• взаимное влияние процессов молекулярной диффузии компонентов (приближение независимой диффузии) не учитывается.
Допущение о независимой диффузии является приемлемым для паровоздушных смесей, мольные доли компонентов в которых малы по сравнению с концентрацией воздуха [1].
Нефть представлена дискретной смесью, состоящей из фиксированного числа условных компонентов (псевдокомпонентов), каждый из которых соответствует узкой фракции на кривой истинных температур кипения (ИТК). Разбиение смеси на псевдокомпоненты и расчет для них свойств проводились с помощью программы моделирования химико-технологических систем СИетСАБ.
Концентрация пара г-го компонента на поверхности жидкости определялась исходя из гипотезы о термодинамическом равновесии между жидкостью и ее паром у поверхности раздела. Согласно закону
Рауля мольная доля пара г'-го компонента на меж фазной границе
РщТ)
Yi— Xj
(1)
где Yi,
мольная доля пара i-п примеси на грани-
це раздела фаз;
Хг — мольная доля компонента в жидкой фазе; Рпг (Т1) — давление насыщенных паров г-го компонента жидкости при температуре Т1, Па;
температура жидкоп фазы, К;
Р0 — давление в окружающей среде, Па.
Зависимость давления насыщенных паров компонентов жидкости от температуры аппроксимировалась уравнением согласно данным, полученным с помощью программы СИешСЛБ:
Ры Т) = ЛТ4 + ВТ1 + СТI + ВТ, + Е, (2)
где Л, В, С, В, — коэффициенты.
Массовый поток г-го компонента с поверхности пролива определялся на основе стандартных функций стенки [2] с учетом поправки на стефановский поток (конвективный поток парогазовой смеси, индуцируемый диффузионным потоком компонентов при непроницаемой поверхности раздела фаз и направленный от поверхности жидкости в парогазовую среду):
1 (Сг,™ " Сг,Р )Ри
Ji
C+ —
Cн
при y+ < уС ;
Sciy +
Sct (и + + Pc) при y+ > уС ;
y+ — Ри*Ур .
— - ln( Ey +) -AS; к
AS — - ln(1 + Czz 0l); -
(3)
(4)
(5)
p — плотность паровоздушноп смеси, кг/м3;
и* — скорость трения; и* = (xw/p)0'5;
xw — напряжение трения на стенке;
Sci и Sct — молекулярное и турбулентное числа
Шмидта соответственно;
yP — расстояние по нормали от поверхности испарения (стенки) до соседнего узла расчетноп сетки;
ц — коэффициент молекулярноп динамическоп вязкости;
к — константа Кармана, равная 0,41; E — константа в логарифмическом законе стенки для скорости, равная 9,1; yC — безразмерное расстояние от стенки, определяемое в точке пересечения линепного и логарифмического закона изменения концентрации у стенки;
Cz — коэффициент;
z0l — высота шероховатости твердоп поверхности, прилегающеп к проливу. Вывод уравнения (8) для определения поправки на стефановскип поток представлен в работе [3].
Определенная таким образом интенсивность испарения Jiw использовалась в качестве граничного условия в области источника в задаче эволюции взрывоопасного газа, включающеп численное решение полноп системы трехмерных нестационарных уравненип Репнольдса, замыкаемых уравнением состояния идеального газа и стандартноп &-е-моделью турбулентности. Для дискретизации дифференциальных уравненип применялся метод контрольного объема, реализованньш в пакете FLUENT.
Плотность газовоздушноп смеси p (кг/м3) определяется по уравнению состояния идеального газа с учетом влияния массовых концентрацип компонентов:
Р —
P0
RoTX CJMг
(9)
* pz oiи *
Pc — 9,24
Sci
Sct
0l
3/4
-1
Ks — 1 -X Y,
(6)
[1 + 0,28e -0,007ScilSc ]; (7)
(8)
где Jiw — массовый поток г-го компонента с поверхности аварийного пролива, кг/(м2-с); К8 — коэффициент, учитывающий влияние на интенсивность испарения стефановского потока; С1к — массовая концентрация г'-й примеси на границе раздела фаз, кг/кг; С,Р — концентрация г-й примеси в пристеночном узле расчетной сетки, кг/кг;
где Я0 — универсальная газовая постоянная, Дж/(моль-К);
Т — температура паровоздушной смеси, К; С1 — массовая доля г-го компонента в смеси; М1 — молекулярная масса '-го компонента, кг/моль.
Граничные профили скорости и турбулентных характеристик на входе в расчетную область (рис. 1) определялись из численного решения задачи движения воздушного потока над ровной поверхностью. В силу симметрии задачи рассматривалась только половина области.
Необходимым условием использования пристеночных функций является достаточное удаление расчетного узла Р от w (рис. 2), для того чтобы первый
0
Рис. 1. Геометрия расчетной области
к о Р
............
Рис. 2. Расположение пристеночного узла расчетной сетки
попал в область турбулентного ядра пограничного слоя. Учитывая данное условие, вертикальный размер пристенных ячеек задавали равным 0,2 м.
Изменение температуры жидкости рассчитывалось из уравнения теплового баланса:
йТ , йх
Ча + Чп
Ч и
Е тис1
(10)
рН
где х — время, с;
Ча — тепловой поток из атмосферы, Вт/м2; Чп — поток тепла от грунта к жидкости, Вт/м2; Чисп — теплота, отводимая от жидкости при испарении, Вт/м2;
ти — масса г-го компонента жидкости, отнесенная к единице площади поверхности пролива, кг/м2;
Срц — удельная теплоемкость г-го компонента жидкости, Дж/(кгК).
Тепловой поток, подводимый от грунта к жидкой фазе (чп = Х(дТ/ду)у = 0), определялся из численного решения трехмерного нестационарного уравнения теплопроводности для твердого подстилающего слоя:
С Р дТп .
Сп Р п ИТ =Х'
(д 2Т дх 2
д 2Тп ду2
д 2ТТ д2 2
л
(11)
где Тп (х, у, 2) — распределение температур в слое грунта;
Сп, рп, Хп — теплоемкость, плотность и коэффициент теплопроводности грунта соответственно. При решении уравнения (11) для твердой поверхности, прилегающей к проливу, задавалось условие Тп = Т.
Тепловой поток из атмосферы ча вычислялся с помощью пристеночных функций, описание которых приведено в работе[3].
Теплота, отводимая от жидкости при испарении:
Ч исп -ЕАНг ■ ¿г, *, (12)
где АИг — теплота испарения г-го компонента жидкой фазы, Дж/кг.
Изменение массы г-го компонента в жидкой фазе описывается уравнением
т , *. (13)
ах
В рассматриваемой задаче удельная теплоемкость, коэффициенты молекулярной вязкости и теплопроводности паровоздушной смеси определялись через массовые доли и удельные теплоемкости, коэффициенты молекулярной вязкости и теплопроводности чистых компонентов:
-Е С
(14)
где ф и фг — физические свойства соответственно смеси и ее г-го компонента. Коэффициенты молекулярной диффузии компонентов в воздухе вычислялись по формуле [4]:
-
10 -7 Т1,75
Р( V
1/3
1/3 \ 2 v а )
— + —, (15)
Мг Ма '
где Р — давление, атм;
V — молекулярный объем газа, определяемый как сумма атомных объемов элементов, входящих в состав газа;
М — молекулярная масса газа, кг/кмоль. Индексы "г" и "а" в формуле (15) относятся соответственно к компоненту пара и воздуху. Количество атомов углерода и водорода в суммарной формуле псевдокомпонента определялось по значению средней молекулярной массы.
Пример практического применения модели
Для иллюстрации практического использования разработанной методики был проведен численный анализ процесса формирования взрывоопасного облака при испарении с поверхности аварийного пролива нефти.
Кривая ИТК нефти принималась по лабораторным данным нефтеперерабатывающего предприятия.
Рассматривались два сценария аварийной ситуации:
• сценарий № 1: гильотинный разрыв емкости объемом 100 м3 с выбросом всего содержимого, разлив нефти на подстилающую поверхность, испарение нефти со свободной поверхности, образование взрывоопасного облака;
• сценарий № 2: частичное разрушение емкости с выбросом 10 % содержимого емкости, разлив нефти на подстилающую поверхность, испарение нефти со свободной поверхности, образование взрывоопасного облака. Для расчета принимался следующий состав нефти (% масс.): пропан — 0,34; н-бутан — 0,82; изобу-тан — 0,23; н-пентан — 0,95; изопентан — 0,66; ТК54 — 0,7; ТК70 — 3,23; ТК90 — 1,6; ТК110 — 1,5 5; ТК130 —3,0; ТК150 —4,49; ТК170—2,38; ТК210 — 5,93; ТК280 — 11,3; ТК360 — 10,11; ТК450 — 18,11; ТК594 — 34,6.
Примечание. Аббревиатура ТК... обозначает среднюю температуру кипения узкой фракции.
Температура окружающей среды задавалась равной 38 °С. Состояние атмосферы — изотермия. Скорость ветра на высоте 10 м — 1 м/с.
Принималось, что площадка, на которой размещена емкость, ограничена обвалованием и имеет площадь 900 м2. Исходя из расчета, что 1 л жидкости разливается на площади 0,15 м2 [5], получено, что при полном и частичном разрушении емкости нефть растекается по всей поверхности площадки. Длительность испарения жидкости принималась равной 3600 с [5].
Другие параметры, остававшиеся постоянными, имели следующие значения: высота шероховатости твердой поверхности г0=г01= 0,01м, Зп =1,28 Вт/(м-К), Сп = 1130 Дж/(кг-К), рп = 2300 кг/м3.
Начальная температура жидкости принималась равной температуре окружающей среды.
Проводилось сравнение результатов расчетов скорости испарения, полученных с использованием разработанной методики и полуэмпирического уравнения [6]:
йх
кырп1 (Т,)
(16)
где к — коэффициент массоотдачи, м/с;
М1 — молекулярная масса г-го компонента, кг/моль.
Коэффициент массоотдачи к рассчитывается по формуле [7]:
к = 17,231 и 0,78/
80 -0,115с Г°,67,
(17)
где и — скорость ветра на высоте 10 м, м/с;
0 — диаметр пролива, м.
На рис. 3 представлены значения удельной массы испарившейся жидкости тисп в зависимости от времени испарения х, рассчитанные с использованием уравнения (16) и разработанной модели.
Из данных рис. 3 видно, что при одинаковой площади пролива масса паров, поступивших в окружающее пространство, в значительной степени зависит от толщины слоя пролитой жидкости. Так, при ре-
Рис. 3. Зависимость удельной массы тисп испарившейся жидкости от времени х при полном разрушении (а) и частичной разгерметизации (б) емкости: 1 — численный метод; 2 — численный метод, без учета влияния силы тяжести; 3 —расчет по уравнению (16)
ализации сценария № 1 количество испарившейся жидкости за рассматриваемый промежуток времени в 3,2 раза превышает значение, полученное для сценария № 2. Удельная масса испарившейся жидкости, определенная по уравнению (16), при реализации сценария № 1 примерно в 1,5 раза, а сценария № 2 в 1,2 раза превышает соответствующее значение, рассчитанное по численной модели. Уравнение (16) не учитывает изменения гидродинамических характеристик потока над поверхностью пролива в результате образования паровоздушной смеси с плотностью, отличающейся от плотности воздуха. Компоненты, образующиеся при испарении, имеют плотность, превышающую плотность воздуха, и, соответственно, при смешении с воздухом образуют тяжелую паровоздушную смесь. Отрицательный вертикальный градиент плотности (устойчивая стратификация) вызывает подавление турбулентности в облаке. Уравнение (16) не учитывает также падения движущей силы процесса диффузии вследствие образования над поверхностью пролива высоких концентраций паров.
Рис. 4 демонстрирует влияние тяжелого газа на турбулентный обмен. В начальный период испарения происходит резкое подавление турбулентности вследствие интенсивного испарения и образования над проливом высоких концентраций тяжелых примесей. Далее происходит постепенный рост коэффициента турбулентной вязкости, обусловленный уменьшением содержания легколетучих компонентов в жидкой фазе и, соответственно, снижением интенсивности испарения. Из рис. 4 видно, что при развитии сценария № 2 эффект устойчивой стратификации выражен слабее, так как концентрации примесей над проливом ниже, чем в случае реализации сценария № 1, вследствие более резкого уменьшения содержания легколетучих компонентов
Ц, Кг/(М'С)
0'00090 400 800 1200 1600 2000 2400 2800 3200 т, с
Рис. 4. Динамика турбулентной вязкости над поверхностью пролива во времени х: 1 — сценарий № 1; 2 — сценарий № 2
в жидкости и, соответственно, снижения интенсивности испарения.
Для сравнения было выполнено численное моделирование процесса испарения без учета влияния силы гравитации. Из рис. 3 видно, что без учета влияния гравитационных эффектов значение массы испарившейся жидкости будет завышено на 32 % для сценария № 1 и на 15 % для сценария № 2.
К изменению интенсивности испарения приводит также снижение температуры за счет испарения компонентов. При реализации сценария № 1 жидкость за рассматриваемый промежуток времени охладилась на 4,1 °С, сценария № 2 — на 5,8 °С. Снижение температуры жидкости в результате испарения компенсируется в большей степени подводом тепла от твердой поверхности.
Горючие пары, смешиваясь с воздухом, образуют взрывоопасные смеси. Взрывоопасной называется смесь, для которой значения концентрации горючего газа лежат в диапазоне, ограниченном НКПР и ВКПР. Обозначим массу газа, способного участвовать в горении, как тв, а отношение тв к общему количеству испарившейся жидкости как 2. Масса тв определяется путем интегрирования:
тв - |Ц С(х, у, г, £) йхйуйг, (18)
2НКПР < У <2ВКПР
где С (х, у, 2, £) — суммарная концентрация компонентов паров бензина, кг/м3; 2НКПР и 2ВКПР — поверхности, ограничивающие область между НКПР и ВКПР. Значения НКПР и ВКПР принимались равными соответственно 1 и 6 % об.
Нестационарное воздействие источника является причиной изменения тв с течением времени (рис. 5). Максимальное значение тв достигается по истечении примерно 8 мин при развитии сценария № 1 и 3 мин — сценария № 2.
Масса горючих паров, способных участвовать в горении, снижается со временем. При частичной раз-
0 600 1200 1800 2400 3000 т, с
Рис. 5. Зависимость тв от времени х: 1 — сценарий № 1; 2 — сценарий № 2
0,4 - V \
0,3- \ \/
0,2 - \ X.
0,1 - \2
0 600 1200 1800 2400 3000 т, С
Рис. 6. Зависимость относительного количества газа, способного участвовать в горении, 2 от времени х: 1 — сценарий № 1; 2 — сценарий № 2
герметизации емкости (сценарий № 2) по истечении 25 мин взрывоопасные концентрации не образуются.
На рис. 6 представлена зависимость относительного количества газа, способного участвовать в горении, 2 от продолжительности испарения. Из данных рис. 6 видно, что в случае реализации сценария № 2 значения 2 ниже, чем при реализации сценария № 1, в силу того, что при меньшей интенсивности испарения выше доля паров, разбавленных до безопасных значений концентраций.
Отношение максимального значения тв к суммарному количеству паров, образующихся за время 3600 с, составило 0,1 для сценария № 1, 0,11 — для сценария № 2. Данные значения согласуются с рекомендуемыми нормативными методиками при оценке последствий аварийных взрывов на открытых пространствах (2 = 0,1).
Результаты численного моделирования показали, что динамика формирования взрывоопасного облака при испарении с поверхности аварийного пролива нефти существенно зависит от толщины слоя разлитой жидкости. Корректное определение массы реагирующей смеси важно при оценке потенциала взрывоопасности технологических блоков и определении ожидаемых последствий аварий в рамках процедуры анализа риска.
СПИСОК ЛИТЕРАТУРЫ
1. Франк-Каменецкий, Д. А. Диффузия и теплопередача в химической кинетике. — М. : Наука, 1987. — 502 с.
2. Fluent Inc. Fluent 6.1. User's Guide, Lebanon, 2003.
3. Галеев, А. Д. Прогнозирование зон токсической опасности и пожаровзрывоопасности при авариях на объектах хранения нефтепродуктов / А. Д. Галеев, С. И. Поникаров // Безопасность жизнедеятельности. — 2009. — № 5. — С. 29-34.
4. Рид, Р. Свойства газов и жидкостей : справочное пособие / Р. Рид, Дж. Праусниц, Т. Шервуд. — Л. : Химия, 1982. — 592 с.
5. НПБ 105-03. Определение категорий помещений, зданий и наружных установок по взрыво-пожарной и пожарной опасности : утв. Приказом МЧС РФ от 18 июня 2003 г. № 314 : ввод. в действие 1 августа 2003 г. — М.: ГУГПС, 2003.
6. Reed, М. A coastal zone oil spill model: development and sensitivity studies / M. Reed, E. Gund-lach, T. A. Kana // Oil and Chemical Pollution. — 1989. — № 5. — P. 441-449.
7. Kawamura, P. I. The evaporation of volatile liquids / P. I. Kawamura, D. Mackay// Journal of Hazardous Materials. — 1987. — Vol. 15. — P. 343-364.
Материал поступил в редакцию 16 декабря 2009 г.
© Галеев А. Д., Поникаров С. И., 2010 г.
(e-mail: galeev_ainur@mail.ru).
О DQ
О
_о
СО
«
Представляем новую книгу
П0ЖНАУКА
»
СВОДЫ ПРАВИЛ. Системы противопожарной защиты. - 2009. - 618 с.
С мая 2009 г. введен в действие Федеральный закон № 123 «Технический регламент о требованиях пожарной безопасности» (полный текст закона опубликован в журнале «Пожаровзрывобезопасность». — 2009.—Т.18,№1).
С вступлением в силу указанного закона теряют свое значение многочисленные Нормы пожарной безопасности (НПБ), Строительные нормы и правила (СНиП), регламентировавшие требования пожарной безопасности к зданиям и сооружениям. В качестве нормативных документов добровольного применения введены Своды правил (СП) и Государственные стандарты.
Настоящий сборник включает Своды правил, которые рекомендуются для применения проектными, строительными и эксплуатирующими строительные объекты организациями при решении вопросов обеспечения пожарной безопасности.
121352, г. Москва, ул. Давыдковская, д. 12, стр. 7; тел./факс: (495) 228-09-03; e-mail: firepress@gmail.com
