Численное моделирование испарения летучей жидкости Текст научной статьи по специальности «Физика»

УДК 614.8

А. А. Салин, А. Д. Галеев, С. И. Поникаров

ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ИСПАРЕНИЯ ЛЕТУЧЕЙ ЖИДКОСТИ

Ключевые слова: аварийный пролив, испарение жидкости, численное моделирование, проверка модели.

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

Keywords: emergency spill, pool evaporation, numerical simulation, model validation.

The paper presents the results of the validation of the developed pool evaporation model using literature and our own experimental data.

Введение

Корректное прогнозирование последствий аварийных событий является ключевым элементом управления промышленной безопасностью на опасных производственных объектах, в части проведения процедуры анализа риска и разработки мероприятий по локализации, ликвидации последствий аварий, вследствие чего становится возможным снижение материального ущерба, а также вреда жизни и здоровью персонала при реализации гипотетических сценариев. Точность прогнозирования определяется достоверностью используемых математических моделей, описывающих высвобождение токсичных химических или легковоспламеняющихся веществ из пролива и их последующее распространение в атмосфере.

При авариях, связанных с разгерметизацией емкостного оборудования со стабильными жидкостями, поступление вредных веществ в окружающую среду обусловлено испарением с поверхности пролива. Существует много работ, касающихся моделирования распространения опасных веществ в окружающей среде, но в то же время мало внимания уделяется моделированию источника выброса.

Описание математической модели испарения

Массовый поток газовой фазы с поверхности пролива вследствие испарения Jg,s определялся с помощью стандартных пристеночных функций [1] с учетом поправки на стефановский поток [2]:

Jg,s = К

(Yg,s - Yg,p)• p• с0-25 .k°-5

Y + =

|Sc•y+

|sct • (u++ PC)

(y +< yC ) . (y +> yC )'

p C0.25 k0.5 y y + = p • C|j • kp • yp ;

J '

u + = - • ln(Ey+)-ЛБ'

(1)

(2)

(3)

(4)

ЛБ = 0 при(К+ < 2.25)'

ЛБ = 1ln(l + Cs • К+)при(К+ > 90)

в иномслучаеЛ Б =

= "lnfК+- 225 к 1 87.75

+Cs • К

• sinjD.4258 (|пК+ - 0.81")}

p К C0.25 k0.5 + p ^s^j kp

(5)

Кэ =

PC = 924

J

f л 34 Sc

Sct

V t.

-1

[l + 0.28 • exp(- 0.007 Sc/Sct ) x ln((1 - Yg,P V(1 - Yg,s))

Yg,s - Yg,P

(6)

(7)

(8)

Предполагается, что переход между полностью турбулентной областью и вязким пристеночным слоем происходит при значении у+ = 11, независимо от концентрации испаряющегося компонента.

Распространение паров в атмосфере определяется решением трехмерных уравнений переноса энергии, массы и импульса осредненных по Рей-нольдсу. Уравнения Навье-Стокса замыкаются к-е моделью турбулентности. Влияние плавучести на генерация кинетической энергии турбулентности (Щ) включены в уравнения переноса в виде:

Gb =-gi

J t dp

(9)

pPrt axj

В случае устойчивой стратификации силы плавучести приводят к подавлению турбулентности (Gb<0). В программном комплексе вычислительной гидродинамики Fluent влияние плавучести на генерацию кинетической энергии турбулентности включается в расчет при наличии гравитационного поля и градиента плотности или температуры. В данной работе влияние сил плавучести на скорость диссипации s турбулентности не учитывается.

В уравнениях (1) - (6) коэффициент CM принимается постоянной величиной, равной Cu=0.09. В модели распространения паров данный коэффициент определяется отношением CM к средней деформации потока, как принято в Realizable k-e модели турбулентности [3].

Изменение температуры испаряющейся жидкости определялось из выражения:

qa + qgrd + qs - qp + qar - Jg,s •ЛН

aTliq

■g

at

CP,liq • mliq

(10)

где да - тепловой поток из атмосферы, Вт/м2; - тепловой поток от грунта к жидкости, Вт/м2; д8 - тепловой поток за счет солнечной радиации, Вт/м2; др - тепловой поток, обусловленный длинноволновым излучением поверхности жидкости, Вт/м2; дЗГ - тепловой поток, поглощаемый проливом вследствие длинноволнового излучения атмосферы, Вт/м2.

Тепловой поток от грунта к жидкой фазе, д^ = ^м(ЭТ6й/ду)у=0 определялся из численного решения трехмерого нестационарного уравнения теплопроводности для подстилающего слоя.

Тепловой поток от атмосферы да определялся с помощью стандартных пристеночных функций:

qa

(Tliq - TP )•

Tp )• C

„0.25 . 0.5 • P • C|j • kr

P

T +=■

PT = 9.24

T

Pr • y+ Prt(u+

+P

У +< УT ) У+ >yT '

Pr

Prt

V t /

f

1 + 0.28exp-0,007-

V

Pr

Pr

(11)

(12)

(13)

Изменение массы жидкой фазы определялось из уравнения:

amliq

at

= -j

g,s-

(14)

Уравнения (10) и (14) были включены в CFD код FLUENT в качестве определяемых пользователем скаляров (UDS) без конвективных и диффузионных членов в уравнениях переноса.

Массовый поток Jg,s использовался в качестве граничного условия при решении задачи распространения в области источника пролива.

Для сравнения проводился расчет изменения массы жидкости по формуле, приведенной в нормативной методике [6]:

amli 6

—IH = -(10-6•Vm• u• p (T,. )), at Hv liq"

(15)

где М - молярная масса опасного вещества, кг/моль; и - скорость ветра, м/с; рн(Т¡¡д) - давление насыщенных паров вещества, при температуре жидкости Тйд. При этом в формуле (10) тепловой поток из воздуха да определялся по формуле [6]: ,0,

• Л0 • ( I „ — I ,. _ )

' (16)

0,035 • u0

Ча =■

'• Ла • (Т0

Tliq)

,0,8 .0,2 • d '

где Т0 - температура воздуха, К; - коэффициент теплопроводности воздуха при температуре Т0; и -скорость воздушного потока над поверхностью испарения, м/с; d - характерный диаметр пролива, м; vа - кинематическая вязкость воздуха при Т0, м2/с.

Проведение эксперимента и сравнение

опытных и расчетных данных

Для подтверждения достоверности модели испарения использовались собственные экспериментальные данные по испарению ацетона. Жидкость наливалась в емкость с размерами 0.6 м х 0.4 м х 0.035 м. Скорость ветра измерялась на высоте 2 метров акустическим анемометром с погрешностью. ± 0.2 м/с. Показания прибора усреднялись каждые две минуты. Для измерения изменения массы емкости с жидкостью использовались весы МеШег-То^о ХР80028 с точностью измерения 0.01 г. Для измерения температуры жидкости использовалась термопара (хромель-копель) с точностью 0.1°С.

3300

3150

4 3000 -

s

I

п :"оо

2550 -

2400

ч X

2i<

30

0 ? 10 15 20 Время, мин

Рис. 1 - Изменении массы ацетона при испарении: штриховая линия - эксперимент, сплошная - расчет по разработанной модели, с маркерами - расчет по нормативной методике [6]

Эксперимент проводился на открытой площадке при температуре 23° С и влажности воздуха 40%. Начальная температура ацетона составляла 23.5°С. Скорость ветра в процессе проведения эксперимента изменялась, и в среднем расчетное значение скорости ветра было принято равным 1.05 м/с. Тепловой поток от солнца определялся по уравнению из работы [4], с учетом широты и долготы местности и облачности его значение принято равным 430 Вт/м2. Величина теплового потока вследствие длинноволнового излучения атмосферы и поверхности пролива определялись по уравнениям из того же источника [4]. Тепловой поток от грунта принимается равным нулю, поскольку чаша с жидкостью была установлена на весы с теплоизоляционной подложкой. Аэродинамическая шероховатость поверхности вокруг пролива принимается равной 0.01 м (^=0.098 м), что соответствует поверхности с невысокой травой [5].

1

На рисунках 1 и 2 представлены результаты моделирования, из которых видно, что расчетные данные по интенсивности испарения, полученные с использованием разработанной модели испарения, занижены на 17,5 % по сравнению с результатами эксперимента. Также стоит обратить внимание, на то что масса испарившегося вещества, определенная по нормативной методике [6], оказывается в три раза меньше значений, полученных экспериментальным путем, как в случае, рассмотренном в работе.

Различие между экспериментальными и прогнозируемыми с помощью разработанной модели зависимостями температуры от времени испарения объясняется расхождением расчетной и экспериментальной интенсивностей испарения [7]. Стоит отметить, что в течение первых десяти минут, экспериментальные и расчетные значения испарившейся массы и температуры жидкости близки, хотя расчетная скорость ветра (1,05 м/с) была примерно в 2,5 раза больше, чем в эксперименте. 300

Ьй 295

и

0

3 :90

в

Л

рц

н

сз

&

<1)

§

Н 280 275

О 5 10 15 20 25 30 Время, мин

Рис. 2 - Изменение температуры ацетона при испарении: штриховая линия - эксперимент, сплошная - расчет по разработанной модели, с маркерами - расчет по нормативной методике [6]

Это можно объяснить, во-первых, погрешностью при измерении скорости ветра, а, во-вторых, тем фактом, что при очень низких значениях скорости ветра и при небольшом размере пролива может преобладать

молекулярный механизм [8] отвода паров от поверхности испарения, при котором влияние скорости ветра на интенсивность испарения ниже, чем при турбулентном режиме. Условно принято считать, что число Шервуда, Sh, (безразмерный параметр, описывающий массоперенос) пропорционален Re-' 5 для ламинарных течений над гладкой поверхностью, и в тоже время Sh~Re-08 для турбулентных течений [9].

Заключение

В результате проверки разработанной модели испарения, интегрированной в программный комплекс вычислительной гидродинамики Fluent, c собственными экспериментальными данными было установлено удовлетворительное совпадение измеренных и прогнозируемых интенсивностей испарения, в отличие от расчетов по нормативной методике, рекомендуемой для проведения процедуры анализа риска.

Литература

1. А.А. Салин, А.Д. Галеев, С.И. Поникаров Основы прогнозирования последствий аварийных проливов бинарных растворов. Изд-во академии наук РТ, Казань, 2015. 134 с.;

2. A.D. Galeev, A.A. Salin, S.I Ponikarov, Journal of Loss Prevention in the Process Industries, 38, 39-49 (2015);

3. T.H. Shih, W.W. Liou, A. Shabbir, Zhu, Comput. Fluids, 24, 227-238 (1995);

4. P.I. Kawamura, D. Mackay, J. Hazard. Mater, 15, 343364 (1987);

5. I. Troen, E. L. Petersen, European Wind Atlas. Riso National Laboratory, Roskilde, Denmark, 1989;

6. Руководство по безопасности "Методика моделирования распространения аварийных выбросов опасных веществ" (утв. Приказом Ростехнадзора от 20.04.2015 № 158)

7. А.Д. Галеев, Г.С. Дьяконов, С.И. Поникаров, А.А. Салин, Вестник Казанского технологического университета, 17, 24, 238-241 (2014);

8. А.А. Салин, А.Д. Галеев, С.И. Поникаров, Вестник Казанского технологического университета, 16, 2, 208-211 (2013);

9. F.P. Incropera, D.P. DeWitt, T.L. Bergman, A.S. Lavine, Fundamentals of Heat and Mass Transfer. Wiley, 2006.

© А. А. Салин - канд. техн. наук, асс. каф. машин и аппаратов химических производств КНИТУ, c888aa@mail.ru; А. Д. Галеев - канд. техн. наук, доц. той же кафедры, galeev_ainur@mail.ru; С. И. Поникаров - д-р техн. наук, проф., зав. каф. машин и аппаратов химических производств КНИТУ, mahp_kstu@mail.ru.

© A. A. Salin - Candidate of Technical Sciences, Assistant, Department of Mechanical Engineering for Chemical Industry, KNRTU, c888aa@mail.ru; A. D. Galeev - Candidate of Technical Sciences, Associate Professor, Department of Mechanical Engineering for Chemical Industry, KNRTU, galeev_ainur@mail.ru; S. I. Ponikarov - Doctor of Engineering, Professor, Chair of the Department of Mechanical Engineering for Chemical Industry, KNRTU, mahp_kstu@mail.ru.