CC BY
7
Вычислительные технологии, 2020, том 25, № 2, с. 4-21. © ИВТ СО РАН, 2020 Computational Technologies, 2020, vol. 25, no. 2, pp. 4-21. © ICT SB RAS, 2020
ISSN 1560-7534 eISSN 2313-691X
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ
D01:10.25743/ICT.2020.25.2.002
Моделирование распространения выбросов автомобильного транспорта в уличном каньоне
Е.А. Длнилкин1,2, А. В. СТАРЧЕНКО1,2
1 Национальный исследовательский Томский государственный университет, Россия 2Институт оптики атмосферы им. В.Е. Зуева СО РАН, Томск, Россия Контактный автор: Данилкин Евгений А., e-mail: danilkine@gmail.com
Поступила 24 июня 2019 г., доработана 9 января 2020 г., принята в печать 24 января 2020 г.
Работа посвящена построению и апробации микромасштабной математической модели и численного метода ее решения для исследования турбулентных течений и переноса пассивной газообразной примеси в уличных каньонах. Представлены результаты исследования структуры течения в уличном каньоне в зависимости от соотношения ширины улицы и высоты зданий, а также результаты расчета неизотермического турбулентного течения в условиях влияния естественной конвекции.
Ключевые слова: моделирование переноса примеси, моделирование турбулентности, уличный каньон, течения с плавучестью.
Цитирование: Данилкин Е.А., Старченко А.В. Моделирование распространения выбросов автомобильного транспорта в уличном каньоне. Вычислительные технологии. 2020. 25(2):4-21.
Введение
В настоящее время население планеты в основном проживает на урбанизированных территориях и процесс урбанизации продолжается. Как результат — повышение плотности городского населения и масса экологических проблем — от загрязнения атмосферного воздуха до загрязнения почвы и водоемов отходами жизнедеятельности человека. Развитие промышленности и рост среднего числа автомобилей на тысячу жителей сделали задачу улучшения качества атмосферного воздуха в городах одной из важных проблем экологии. Чем крупнее и густонаселеннее город, тем острее стоит обозначенная проблема.
В роли основных источников загрязнения атмосферного воздуха выступают промышленные предприятия, автомобильный транспорт, котельные и ТЭЦ. Особо в этом ряду выделяется автотранспорт, потому как за последние несколько десятилетий именно он стал основным источником загрязнения. Это наиболее заметно в крупных городах, где на его долю приходится до 70 % вредных выбросов в атмосферу. В выхлопных газах автотранспорта содержится сложная смесь более чем из двухсот компонентов, среди которых основными загрязнителями являются оксид углерода, летучие углеводороды, оксид азота и свинец. Постоянный рост общего количества автомобилей и суммарного объема вредных выбросов только усугубляет сложившуюся ситуацию.
У автотранспорта как источника вредных выбросов есть несколько особенностей, выделяющих его на фоне предприятий и котельных. Это низкое расположение источника выбросов, распределенность по всей территории города и близость к среде обитания человека, что оказывает существенное влияние на характер распространения выбросов
автотранспорта. Так, например, если выбросы высоких труб предприятий сравнительно быстро рассеиваются в атмосфере, то в приземном слое наблюдаются наихудшие условия для рассеивания. Отработанные газы автотранспорта, располагающиеся невысоко от земной поверхности, плохо рассеиваются ветром и скапливаются в жилых зонах, попадая непосредственно в зону дыхания человека.
Сказывается и изменение естественного ландшафта, жилая застройка еще больше снижает скорость приземного ветра и приводит к формированию особых условий циркуляции и теплообмена воздуха в приземных слоях. Выбросы от автотранспорта распространяются вдоль дорог, оказывая вредное воздействие на пешеходов и жителей окрестных зданий. Все это приводит к формированию повышенного общего фона загрязнения на территории города.
Сегодня сложно представить себе отказ от использования автотранспорта или одномоментную замену двигателей внутреннего сгорания более экологичными электродвигателями. Поэтому важны контроль и прогноз качества атмосферного воздуха. Для исследования процессов распространения выбросов автотранспорта, происходящих в городской застройке, необходимо развитие и внедрение микромасштабных математических моделей, описывающих движение воздуха и перенос примеси. Используя такие модели, можно прогнозировать распространение выбросов автотранспорта в зависимости от конфигурации городской застройки и метеоусловий.
Моделирование турбулентных течений в окружающей среде осуществляется в основном с использованием осредненных по Рейнольдсу уравнений Навье — Стокса (Reynolds Averaged Navier — Stokes) и переноса скалярных величин, для которых требуется решить проблему замыкания путем привлечения полуэмпирических моделей различного уровня сложности. Этот подход позволяет получить успешные результаты для целого ряда турбулентных течений и является значительно менее трудоемким с вычислительной точки зрения в сравнении с набирающими популярность LES- и DES-подходами [1].
Решению рассматриваемой задачи посвящено большое количество работ как в России, так и за рубежом. За последние несколько десятилетий исследователи значительно продвинулись в части моделирования переноса вредных выбросов в уличном каньоне. В качестве важных параметров, влияющих на структуру течения и механизм распространения примеси, выделяют метеоусловия (скорость и направление ветра) [1-3], геометрические характеристики (соотношение сторон уличного каньона) [4, 5], температурную неоднородность [6-8] и влияние растительности [9, 10].
Применяются новые подходы. Так, в работе [11] моделирование турбулентности и переноса мелкодисперсной примеси в городских каньонах выполнено с использованием LES-модели. За счет явного разрешения крупных вихрей и воспроизведения нестационарной структуры потока изучены особенности вентиляции улиц при одномоментном выбросе примеси. Важно, что в работе сформулированы требования к сеточному разрешению, необходимому для корректного вычисления концентрации примеси.
Активно исследуются процессы переноса для реальных участков городской застройки. В работе [12] выполнено исследование переноса примеси для участка города Ухань (Китай). Рассмотрено течение в трехмерном протяженном уличном каньоне, построенном из зданий различной этажности и протяженности. Показано, как интенсивность источника, сила и направление ветра влияют на концентрацию углекислого газа на высоте полутора метров (зона дыхания).
Однако далеко не все аспекты задачи о моделировании переноса выбросов автотранспорта хорошо изучены. Поэтому целью данной работы является разработка и ве-
рификация микромасштабной математической модели неизотермического турбулентного течения и переноса примеси, опирающейся на ИЛ^-подход моделирования турбулентности, основное назначение которой — выполнение оперативного прогноза и проведение параметрических расчетов уличных каньонов.
Данная работа является продолжением исследований, представленных в [1], и направлена на разработку математической модели и вычислительных технологий для анализа влияния температурной неоднородности в уличных каньонах на движение воздуха и перенос примеси в них. Полученные в работе результаты показывают необходимость учета этого фактора в случаях, например, заметного нагрева солнечными лучами стен зданий, образующих уличный каньон.
1. Физико-математическая постановка задачи
Рассматривается стационарное неизотермическое турбулентное движение в приземном слое воздуха над неоднородной подстилающей поверхностью с элементами крупномасштабной шероховатости. Элементы шероховатости представляют собой неподвижные, непроницаемые для потока прямоугольные препятствия (здания). Размеры зданий соизмеримы с размерами области исследования (рис. 1). Рассматриваются точечные и линейные источники вредных выбросов постоянной интенсивности. Предполагается, что протяженность уличного каньона существенно больше его ширины и высоты, что позволяет использовать двумерную постановку задачи.
Математическая модель включает в себя осредненные по Рейнольдсу уравнения неразрывности, уравнения Навье — Стокса, переноса примеси [1] и теплообмена:
дщ дх;
= 0,
диг и,
дхп
1 др д
+
р dxi дх дТй,
(ьРиЛ
\dxjj
ди'гu'j + дг (Т - То)
д
дх,
дСй, дх,
д
дх,
(о
дх,
дс >
дх, )
j
дТ 'дх,
дх
То
дс' и' _'
дх'
дТЩ
дхj
+ S, i,j = 1, 2.
Рис. 1. Иллюстрация физической постановки задачи Fig. 1. Sketch of the physical problem formulation
Здесь щ, — осредненные и пульсационные проекции вектора скорости на оси координат; р — давление; р — плотность; V — кинематическая вязкость воздуха; Т — температура; С — осредненное значение концентрации примеси; в — функция, описывающая распределение точечных и линейных источников примеси; — тензор напряжений Рейнольдса; а, И — температуропроводность и коэффициент диффузии; 91(0, д) — компоненты ускорения свободного падения.
Замыкание описанной системы уравнений проводится с использованием градиент-но-диффузионной гипотезы Буссинеска [1]:
/ д-щ + диЛ + 2.
\дх, дХг)
<и'з = + ТйГ ) + 3к5ч,
Т V, = -
ут дТ Ргт дх, ш = -ит_ 3 Бот дх
Для вычисления турбулентной вязкости используется двухпараметрическая к — е-модель турбулентности:
деи,•
дки^ дхз д
дх,• дх,•
_д_ дх,•
"т = С,, -—
Рт
у + —
Ут
У + —
о.
—дщ Р = —Щи> д^
\ дк_\ ) дх,)
\ 0е_ N
) дх,)
+ Р + С — £,
+ - (Се1Р + Се2С — С£3е)
к
С = 5.
¿2
Т0 Ргу дх
где ит — турбулентная вязкость, к — кинетическая энергия турбулентности, е — диссипация турбулентной кинетической энергии. Коэффициенты модели турбулентности ак = 1.0, = 1.22, Се1 = 1.44, Се2 = 1апЬ(|ад2|/|мх|), Се3 = 1.92, = 0.09, Бет = 0.5, Ргт = 0.9.
Для задания значений скорости, турбулентных параметров, трения и теплового потока вблизи твердой поверхности используется метод пристеночных функций [13]. В модели не учитывается осаждение примеси на стенах, крышах и подстилающей поверхности, поэтому на границах производные от концентрации по нормали к поверхности задаются равными нулю. Краевые условия на выходе потока из расчетной области и на открытой верхней границе — это равенство нулю производных по нормали. При задании граничных условий на входе используются равномерные профили скорости к, е, Т, С.
При расчете течений вокруг зданий использовался метод фиктивных областей, суть которого заключается в том, что значения векторных и скалярных величин в области преграды равны нулю и на границах фиктивных конечных объемов нет потоков диффузии [1].
2. Аппроксимация дифференциальной задачи и численный метод решения
2.1. Построение конечно-объемной сетки
Дискретизация дифференциальной задачи осуществляется методом конечного объема с использованием разнесенной сетки [1]. При построении конечно-объемной сетки используется алгоритм, в рамках которого вся область исследования делится на подобласти и для каждой области задается свой признак: либо это несжимаемая среда 1, либо непроницаемое препятствие 9. Координаты границ подобластей подбираются так, чтобы разграничить (выделить) в отдельные подобласти непроницаемые участки внутри области исследования. На рис. 2 дан пример двумерной области исследования, представляющий собой модель уличного каньона, где нижний правый 9 и нижний левый 9 квадраты представляют собой здания. Соответственно, в подобластях с признаком 1 будут моделироваться движение воздуха и перенос примеси.
Построение расчетной сетки осуществляется автоматически по данным из файла gгid.еas, который имеет следующую структуру. В первой строке указано общее число подобластей NDX в направлении оси ОХ (рис. 2). В следующей строке указаны ХБХ+1 действительных чисел, являющихся координатами начала подобластей в направлении оси ОХ. Следующие две строки задают число подобластей в направлении оси ОУ и их координаты. Таким образом, область исследования покрывается набором смежных прямоугольников (подобластей). Далее в файле gгid.еas располагается матрица размера NDXxNDY, в которой заданы признаки подобластей: 1 — несжимаемая среда; 9 — непроницаемое препятствие. И в заключение следуют два массива целых чисел длиной NDX и NDY, в которых содержится число узлов в каждой подобласти в направлении осей ОХ, ОУ соответственно. На основании считанных данных строится единая глобальная сетка (сеточные линии), допускающая сгущения в выделенных подобластях, представляющих наибольший интерес для исследования или требующих высокого разрешения с точки зрения моделируемых процессов.
пох[0]=20 пох[1]=20 пох[2]=20
grid.cas
NDX: 3
ХО : 0 15 30 45 NDY: 2
YO : 0 15 45 NTP: 111 9 19 NXO: 20 20 20 NYO: 20 40
ох[0]=0 m ох[1]=15 m ох[2]=30 m ох[3]=45 m
Рис. 2. Разделение области исследования на подобласти и пример файла grid.cas для двумерной модели уличного каньона
Fig. 2. Subdomains of the computational domain and an example of grid.cas file for the two-dimensional model of street canyon
Для построенной конечно-объемной сетки использованы следующие обозначения: х сг, у с] — координаты граней конечных объемов или сеточные линии; Ах сг = х сг+ — х сг, Ау С] = у С]+1 — у сз — расстояния между гранями конечных объемов; хг = (х сг+1 +х сг)/2, Уз = (Усз+1 + Усз)/2 — координаты центров конечных объемов; ¿хг = хг+1 — хг, Ауз = Уз+1 — Уз — расстояния между центрами конечных объемов.
2.2. Аппроксимация и численный метод решения
После того как конечно-объемная сетка построена, необходимо каждое дифференциальное уравнение проинтегрировать по каждому конечному объему Аг ] = [хсг,хсг+1] х [уС], усз+1]. Кратко рассмотрим, как осуществляется дискретизация уравнения переноса для некоторой скалярной величины:
(d—\ + —(d—^
\ дх) ду \ ду)
дС dUC dVC д fdC\ д . . „ + "ТТ" + "ТТ" = ТГ ' D— ) + — ( D— ) + S.
д х д х
Нестационарный член (здесь и далее верхний индекс указывает на шаг по времени, задача решается методом установления):
(f)
(дс\ , , Cf1 - с*
Г I — I axay & I — I ахCidyCj & -—
Aij ij
) dxdy & ( ) dxfydyCj &-—-dxc^dycj.
ij
Конвективное слагаемое
L
Aij
Диффузионное слагаемое д
(^P) dxdy & - UtCkw) dycj.
Aij
дх (d© ^ & (d.- Dw^.
\J kKj \ \J >Xj 1 \ KAjiKji KAjiKji_1 t
Здесь Uk, С*, D. — значения функций, определенные не в центре конечного объема, а на его границе в точке с координатами (хci+1, yj) соответственно; UW, CW, Dw — значения, определенные в точке с координатами (хCi, yj). Проинтегрировав аналогичным образом конвективное и диффузионное слагаемые по другой координате и подставив результат в уравнение переноса, запишем
Cjk+1 _ Ck
%J —— dx Cidy Cj =
(C* _ C* C* _ C* \
D. ^хг ij - Dw ijdXi P'7 1 dycj - (UkCik - UWCW) dyCj+
(C k _ C k C* _ C*
Dn -^ - Ds Cij--^ ) dx Ci - (Vjk Ck - VskCk) dxCi + Si3 dxCi dyc3
dyj dyj-1 1
В полученном выражении значение неизвестной функции на гранях конечных объемов необходимо заменить значениями функции в центре конечного объема. С этой
целью, как правило, используются противопотоковые схемы. Для проведения параметрических расчетов и моделирования течений в областях со сложной геометрией нежелательно использование противопотоковой схемы первого порядка, так как она обладает высокой схемной вязкостью и имеет низкий порядок точности. Поэтому в данной работе применена противопотоковая схема МЬИ Ван Лира [14]. Для ее использования необходимы значения функции в двух предшествующих узлах сетки и в одном последующем по направлению потока. Тогда значение функции на грани определяется по формуле
сп =
Ck + с + dXi cij + сЕ 2 '
X i
rik _ с—_i
ci+1j сЕ 0 '
если Uk > 0, если Uk < 0,
где С+ и функция min mod () определены следующим образом:
С+ = minmod
(Ое + Ь±) 2
2 minmod [ ае, Ье ]
minmod [ае, Ь± ]
<i ье
> | ее
О, если aebf < 0,
е' если | е
е'
если е
е | and aeb±e > 0,
е | and аеЬ±е > 0,
с k с k ci+ij cij
если Ukk > 0,
Ck _ ск
Cij Ci— 1j
_ ) dxi—i с k с k
Ci+2j Ci+1j t тк * c\
—-—' если Ukk < 0.
axi+i
Если сеточная функция на рассматриваемом участке не является монотонной, то схема МЬИ ведет себя аналогично обычной противопотоковой схеме первого порядка, используя подветренное значение. Там, где характер изменения функции монотонный, схема МЬИ обеспечивает второй порядок аппроксимации. В общем случае она остается схемой первого порядка, но является монотонной.
Воспользовавшись схемой МЬИ, приведя подобные и выполнив преобразования с уравнением неразрывности [1], можно получить готовую вычислительную формулу вида
С^ = ар^С^ + ае íjСí+lj + aWijC,í-lj + + аз + bíj,
где коэффициенты вычисляются следующим образом:
^^ ij
D,
dy Cj
& J I ______i ТТП
е , + max (—ЦТ dy Cj, 0)
ap,
i] >
dycj
awij = ( Dw-j—--+ max (UW dycj, 0)
^/XX Ci
anij = I Dn—--+ max (—VT dxc,,, 0)
J x dy, n
ap,
е
(dx Ci \ / 0 Ds---+ max (VI1 dxCi, 0m ар.,,
dyj-i J / J
apij = 1 — (а e ij + awij + anij + as ij),
о d/Xx Ci d'ij сj
aVij A I , bij
о d>xxCi dycj bi3
= At , Ыз = a^,
bij = Sij dx CidyCj —
/ ь> \ _i_ dXi / J \ dXi
— max (Ue dycj, 0) S+—--max (—Ue dycj, 0) S-——+
+ max (U*dyCj, 0) S++ ^T1 + max (—U^dyCj, 0) Sw^i- 1
( к \ —I— 3 ( к \ 3
— max {Vn dx Ci, 0) S+ —--max ( — Vn dx Ci, 0) S— ——+
+ max (Vsk dx a, 0) S+ + max (—Vsk dx a, 0) S—
Полученная явная разностная схема имеет первый порядок аппроксимации по времени и по пространству, является монотонной и условно устойчивой. При расчетах шаг по времени выбирался постоянным и определялся из ограничения на параболическое и гиперболическое числа Куранта:
(h2 h At = 0.25 min —, —
\au' 2Ui
—)
2 Uin J
где h = minij (dxi, dyj) — минимальный шаг по пространству, величина a = max^ (fj/у) оценивалась эмпирически.
Аналогичным образом осуществляется интегрирование уравнений для компонент скорости. Для расчета поля течения использовалась процедура SIMPLE [15], которая на каждом шаге по времени включает в себя последовательность операций (задача решалась методом установления):
— задание начальных приближений для всех зависимых переменных;
— вычисление приближенных значений скорости по готовым явным вычислительным формулам;
— решение СЛАУ для поправки давления методом верхней релаксации;
— коррекция давления и компонент скорости;
— решение дискретного аналога для турбулентных характеристик и концентрации;
— принятие полученных значений зависимых переменных в качестве исходных и переход на следующий шаг по времени с повторением операций, начиная со второй, если не выполнено условие установления, иначе выход из цикла [1].
При решении уравнения Пуассона для поправки давления на выходе из расчетной области задаются граничные условия первого рода (равенство нулю). Это сделано из следующих соображений. Поправка давления определяется через производную скорости, а на выходе из расчетной области производная скорости по нормали к границе равна нулю. При численном решении начальное поле давления принимается равным нулю. Таким образом, давление и поправка давления на выходе из расчетной области всегда равны нулю.
3. Верификация построенной математической модели
и проведение параметрических расчетов течения в уличном каньоне
3.1. Верификация численного метода
Тестирование вычислительного алгоритма проведено на задаче моделирования ламинарного течения несжимаемой среды в начальном участке плоского канала при низких числах Рейнольдса. Течение происходит слева направо, длина Ь и ширина 2Н канала составляют 0.5 и 0.025 м соответственно. Геометрия области исследования представлена на рис. 3.
При моделировании входная скорость задает направление движения и определяет величину перепада давления. В качестве начальных условий использовались следующие: продольная скорость 0.16 м/с, а вертикальная скорость и давление равны нулю. В качестве граничных условий использовались следующие: на входной границе горизонтальная скорость = 0.16 м/с, вертикальная равна нулю; на выходной границе производные скорости по нормали равны нулю; на неподвижных стенках заданы условия прилипания.
Расчет проводился на равномерной сетке 100 х 50 узлов по осям х и у соответственно. Шаг по времени выбран Д£ = 0.000019 с, общее число шагов по времени 300. Из теоретических оценок следует, что при х > 0.18Яе2^ профиль скорости приобретает параболический вид, соответствующий полностью развитому ламинарному течению. В рассматриваемом примере Яе = и\п • 2Н/у = 200, V = 2 • 10-5 м2/с. Аналитически профиль скорости в начальном участке плоского канала можно определить с помощью уравнений, полученных С.М. Таргом приближенным методом [16]:
и(х, у) Uin
=3 i1 - (ПЧ - 2 £ 4 (1 - ехр fя,
2 V ) ^ cos(7„) J PV Ui„ h J'
где и^п — скорость на входе в канал, 2Н — ширина канала, и 7п — последовательные корни уравнения tg(х) = х.
динамическая вязкость,
Рис. 3. Иллюстрация постановки задачи ламинарного течения в начальном участке плоского канала
Fig. 3. Sketch of the laminar flow formulation in the initial section for the flat channel problem
Рис. 4. Распределение продольной компоненты скорости вдоль оси х в семи различных сечениях канала
Fig. 4. Longitudinal velocity values along the x axis at 7 cross-sections of the channel
На рис. 4 представлено распределение продольной компоненты скорости вдоль оси х. Видно, что наблюдается хороший уровень соответствия данных численного моделирования и аналитического решения для рассматриваемой задачи.
3.2. Моделирование течения в уличном каньоне
Предложенная математическая модель была применена для исследования аэродинамики потока в уличном каньоне и переноса примеси, поступающей от точечного источника. Уличный каньон является базовым элементом застройки современного города, поэтому именно он стал объектом экологических исследований, посвященных изучению распространения и рассеивания вредных выбросов в городских кварталах. Исследования показали, что на структуру течения в каньоне в первую очередь влияют геометрия каньона и температурная неоднородность внутри него. Поэтому задачей исследования являлось изучение структуры течения и характера распространения вредных выбросов в зависимости от соотношения высоты и ширины уличного каньона и в условиях заметного влияния естественной конвекции. Также проведена верификация разработанной модели переноса путем сравнения результатов моделирования с данными натурных экспериментов.
Геометрические характеристики уличного каньона (высота Н и ширина Ш) менялись в диапазоне от 5 до 40 м. Расчетная область, представленная на рис. 5, имела следующие размеры: 40 + Ш м в продольном и 50 + Н м в вертикальном направлении. Расчеты
У А
Рис. 5. Моделирование течения в уличном каньоне / Fig. 5. Flow modelling in street canyon
проводились на декартовой сетке 102 х 92. Шаг по времени выбран ДЪ = 0.0025 с, общее число шагов по времени 10 000. Источник поступления примеси постоянной интенсивности располагался в центре уличного каньона вблизи поверхности. Значение вязкости среды принято и = 15 х 10-6 м2/с, что соответствует вязкости воздуха при 20 °С.
Граничные условия задавались следующим образом: на входной границе горизонтальная скорость игп = 5 м/с; на выходной границе равенство нулю нормальных производных скорости; на верхней стенке условия скольжения. Начальные условия: продольная скорость равна скорости на входе (5 м/с), а вертикальная скорость и давление равны нулю.
В соответствии с технологией сгущения конечно-объемной сетки, описанной в под-разд. 2.1, в приграничных узлах размер ячейки в направлении по нормали к твердой границе задан 0.1 м. Это сделано для попадания безразмерной величины у+ в рекомендуемый диапазон от 30 до 300 [13, 17]. Предварительно оценить у+ по входным данным задачи можно, используя следующие формулы [1]:
,+ _ и * ки
У+ = , и
и Щур/у о)
где ур — расстояние от стенки до центра первой ячейки, и * — динамическая скорость, и — вязкость, и — скорость в приграничной ячейке, к = 0.41 — постоянная фон Кармана, уо = 0.000394 — коэффициент шероховатости поверхности.
В рассматриваемой задаче скорость течения у стен и дна уличного каньона значительно меньше скорости основного потока и составляет около 1м/с при скорости основного потока 5 м/с. Поэтому для построенной сетки и заданного коэффициента шероховатости можно использовать приближенное значение динамической скорости и * « 0.0847. Тогда оценочное значение для у+ равно:
+ ур U* 0.2 • 0.0847
у+ = ^-=-— w 282.
у и 15 • 10-6
Выполненная оценка хорошо согласуется с результатами выдачи при работе программы. В расчетах получено, что значения у+ укладываются в 500.
Проведена также серия расчетов на сгущающихся сетках, по результатам которых и выбраны параметры сетки для получения сеточно-независимого решения. В качестве примера, подтверждающего сеточную сходимость, на рис. 6 представлено распределение горизонтальной скорости вдоль прямой х = 30 м для используемой сетки 102 х 92 и более подробной сетки размера 202 х 182.
Анализ результатов расчетов для варианта H = 24м, W = 20 м показал, что основной поток воздуха над крышами зданий приводит в циркуляционное движение воздух, находящийся внутри уличного каньона. Направление ветра над зданиями определяет направление возникающего движения в каньоне и его интенсивность. Возникшее закрученное течение переносит примесь от источника на дне каньона к его подветренной стороне и далее по направлению потока. На верхней границе каньона часть примеси уносится основным потоком, а часть возвращается и циркулирует в нем.
На рис. 7 представлены результаты сравнения рассчитанных нормированных значений концентрации примеси вблизи наветренной и подветренной вертикальных стенок каньона с измерениями Хойдыша и Дэбберти, выполненными в ветровом туннеле [18]. Из рисунка видно, что результаты расчетов удовлетворительно согласуются с экспериментальными данными.
Проведено исследование структуры течения и характера распространения примеси в уличном каньоне в зависимости от соотношения его ширины и высоты. Установлено, что модель воспроизводит основные режимы течения в каньоне, показанные в работе [4], в которой моделирование турбулентности осуществлялось методом крупных вихрей (LES). Так, например, в случае уменьшения высоты уличного каньона H или увеличения ширины W до соотношения H/W = 0.5 центр основного вихря смещается к наветренному зданию и вихрь растягивается во всю длину каньона. Дальнейшее
—80.00 —
70.Ю -:-
60.00
s
га х о л
га (о
и Л
СО
10.01
2.00
Рис. 6. Иллюстрация сеточной сходимости на примере распределения горизонтальной скорости вдоль вертикальной сеточной линии х = 30 м
Fig. 6. Mesh convergence of the computations illustrated by the values of the horizontal velocity along the vertical line x = 30 m
у/Н
0.2
0.4
0.8
0.6
1.0
0
С/С,
0.2 0.4 0.6 0.8 1.0
max
Рис. 7. Рассчитанные и измеренные [18] вертикальные профили нормированной концентрации примеси у наветренной (слева) и подветренной (справа) сторон каньона
Fig. 7. Computed and measured [18] vertical profiles of the normalized concentration of the pollutant near the windward (on the left) and leeward (on the right) sides of the canyon
увеличение расстояния между зданиями до соотношения H/W = 0.125 приводит к образованию двух рециркуляционных зон: большой вихрь у подветренной стороны и малый у наветренной.
Рельеф современного города образован сооружениями и поверхностями из бетона и асфальта. Эти материалы имеют низкую теплоемкость, что вызывает значительные суточные изменения температуры поверхностей. Поэтому важно исследовать влияние температурной неоднородности на рассеивание загрязняющих веществ.
Выполнено исследование влияния нагрева стен зданий и подстилающей поверхности на аэродинамику и перенос примеси в каньоне. Рассматривались три варианта:
а) температура подветренной стенки выше температуры окружающего воздуха на 10 град.;
б) температура наветренной стенки была выше температуры окружающего воздуха на 10 град.;
в) температура нижней стенки была выше температуры окружающего воздуха
Во всех рассмотренных вариантах температура окружающего воздуха 300 К. Максимальное влияние температурной неоднородности будет проявляться в безветренную погоду, поэтому скорость набегающего потока выбрана небольшой — 1.18 м/с. В случае, когда нагрета подветренная сторона (вариант а), восходящий поток вблизи нагретой подветренной стены здания усиливает основной вихрь и увеличивает среднюю скорость кругового движения, вследствие чего примесь интенсивнее выносится за границы уличного каньона (рис. 8, вверху).
В случае нагретой наветренной стороны уличного каньона здесь образуется восходящий поток воздуха. При этом движение образовавшегося восходящего потока направлено в сторону, противоположную движению основного вихря. Это приводит к поджатию основного вихря и образованию нового в зоне между основным вихрем и наветренной стеной здания. На рис. 8 (внизу) представлены векторное поле и распределение температуры, соответствующие рассматриваемому случаю (б), при этом рассеивание загрязняющих веществ ухудшается и примесь накапливается в каньоне.
на 10 град.
Рис. 8. Векторные поля скорости на фоне изолиний концентрации примеси и изолинии температуры для случаев (а) — вверху и (б) — внизу
Fig. 8. Vector velocity fields versus the background of isolines of the pollutant concentration and temperature isolines for cases (а) above and (б) below
В случае, когда нагрето дно уличного каньона (вариант в), образующийся восходящий поток воздуха усиливает основной вихрь и увеличивает среднюю скорость кругового движения, что, как и в случае (а), способствует более интенсивному выветриванию примеси из уличного каньона. Картины течения для случаев (а) и (в) визуально практически идентичны, поэтому приведена иллюстрация векторного поля скорости только для случая (а).
Анализ результатов показывает, что в случае нагрева дна уличного каньона или подветренной стороны вредные выбросы интенсивнее выносятся из каньона, а максимальные концентрации уменьшаются на 10-15% в сравнении с равномерно прогретым уличным каньоном. В то же время в случае нагрева наветренной стороны максимальные концентрации возрастают в 3-3.5 раза.
Заключение
Представлена микромасштабная модель турбулентного неизотермического движения воздуха и переноса примеси в уличных каньонах, учитывающая влияние силы плавучести на характер течения, его турбулентную структуру. Для численного решения уравнений математической модели подробно описаны используемые вычислительные технологии, опирающиеся на метод конечного объема, метод установления, схемы по-
вышенного порядка точности для аппроксимации конвективных членов и итерационные методы для решения сеточных уравнений для поправки давления.
Верификация разработанной модели выполнена на известных экспериментальных данных для переноса примеси при стационарном изотермическом турбулентном пространственном течении в каньоне [18] и на приближенном решении Тарга [16] о задаче ламинарного изотермического движения воздуха в плоском канале.
Для неизотермического турбулентного течения воздуха в каньоне показано, что в случае нагрева его наветренной стороны характер течения существенно изменяется по сравнению с изотермическим потоком [1], что в итоге приводит к гораздо более значительному увеличению значений концентрации примеси в каньоне.
Благодарности. Работа выполнена при финансовой поддержке гранта Президента РФ № МК-1723.2017.5 и Российского фонда фундаментальных исследований и Администрации Томской области № 16-41-700178 ра.
Список литературы
[1] Старченко А.В., Нутерман Р.Б., Данилкин Е.А. Численное моделирование турбулентных течений и переноса примеси в уличных каньонах. Томск: Изд. дом Томского гос. ун-та; 2015: 252.
[2] Balogun A.A., Tomlin A.S., Wood C.R., Barlow J.F., Belcher S.E., Smalley R.J., Lingard J.J.N., Arnold S.J., Dobre A., Robins A.G. In-street wind direction variability in the vicinity of a busy intersection in Central London. Boundary-Layer Meteorology. 2010; 136(3):489-513.
[3] Kim J.J., Baik J.J. A numerical study of the effects of ambient wind direction on flow and dispersion in urban street canyons using the RNG-Turbulence model. Atmospheric Environment. 2004; 38(19):3039-3048.
[4] Danilkin E.A., Starchenko A.V. Large eddy simulation of turbulent flow and of pollutant transport in a street canyon. Proc. SPIE 9680, 21th Intern. Symp. on Atmospheric and Ocean Optics: Atmospheric Physics. 2015:968062-1-968062-6. D0I:10.1117/12.2205490
[5] Chan A.T., So E.S., Samad S.C. Strategic guidelines for street canyon geometry to achieve sustainable street air quality. Atmospheric Environment. 2001; 35(24):5681-5691.
[6] Park S.B., Baik J.J., Raasch S., Letzel M.O. A large-eddy simulation study of thermal effects on turbulent flow and dispersion in and above a street canyon. Journal of Applied Meteorology and Climatology. 2012; 51(5):829-841.
[7] Wang P., Zhao D., Wang W., Mu H., Cai G., Liao C. Thermal effect on pollutant dispersion in an Urban Street Canyon. International Journal of Environmental Research. 2011; 5(3):813-820.
[8] Kim J.J., Pardyjak E., Kim D.Y., Han K.S., Kwon B.H. Effects of building-roof cooling on flow and air temperature in urban street canyons. Asia-Pacific Journal of Atmospheric Sciences. 2014; 50(3):365-375.
[9] Nuterman R.B., Baklanov A.A., Starchenko A.V. Modeling of aerodynamics and pollution dispersion from traffic in the urban sublayer. Mathematical Models and Computer Simulations. 2010; 2(6):738-852.
[10] Fei X., Xiaofeng L. The impact of roadside trees on traffic released PM10 in urban streetcanyon: Aerodynamic and deposition effects. Sustainable Cities and Society. 2017; (30):195-204.
[11] Глазунов А.В. Численное моделирование турбулентности и переноса мелкодисперсной примеси в городских каньонах. Вычислительные методы и программирование. 2018; 19(1):17-37. https://doi.org/10.26089/NumMet.v19r103
[12] Ming T., Fang W., Peng C., Cai C., Richter R., Ahmadi M., Wen Y. Impacts of traffic tidal flow on pollutant dispersion in a non-uniform urban street canyon. Atmosphere. 2018; 9(3):82.
[13] Launder B.E., Spalding D.B. The numerical computation of turbulent flows. Computational Methods in Applied Mechanics and Engineering. 1974; 3(2):269-289.
[14] Van Leer B. Towards the ultimate conervative difference scheme: II. Monotonicity and conservation combined in a second order scheme. Journal of Computational Physics. 1974; (14):361-370. http://dx.doi.org/10.1016/0021-9991(74)90019-9
[15] Патанкар С. Численные методы решения задач теплообмена и динамики жидкости: Пер. с англ. М.: Энергоатомиздат; 1984: 149.
[16] Жукаускас А.А. Конвективный перенос в теплоприемниках. М.: Наука; 1982: 472.
[17] Лойцянский Л.Г. Механика жидкости и газа. М.: Дрофа; 2003: 840.
[18] Hoydysh W.G., Dabberdt W.F. Kinematics and dispersion characteristics of flows in asymmetric street canyons. Atmospheric Environment. 1988; (22):2677-2689.
Вычислительные технологии, 2020, том 25, № 2, с. 4-21. © ИВТ СО РАН, 2020 ISSN 1560-7534
Computational Technologies, 2020, vol. 25, no. 2, pp. 4-21. © ICT SB RAS, 2020 elSSN 2313-691X
MATHEMATICAL MODELLING
DOI:10.25743/ICT.2020.25.2.002
Modelling the transfer of road transport emissions in a street canyon
Danilkin Evgeniy A.1,2, Starchenko Alexander V.1,2 1 Tomsk State University, 634050, Tomsk, Russia
2V.E. Zuev Institute of Atmospheric Optics SB RAS, 634055, Tomsk, Russia Corresponding author: Danilkin Evgeniy A., e-mail: danilkine@gmail.com
Received June 24, 2019, revised January 9, 2020, accepted January 24, 2020
Abstract
The study is focused on developing and testing a microscale mathematical model for the analysis of turbulent flow and passive gaseous admixture transfer in street canyons. Mathematical model is based on Reynolds-averaged Navier —Stokes and continuity equations. Boussinesq approximation and two-parameter k-epsilon turbulence model is used to close the equations.
Numerical solution of the system of differential equations is obtained with the finite volume method on a staggered mesh. Convective terms of the Navier —Stokes equations are approximated with MLU numerical scheme. SIMPLE computational algorithm is used to couple velocity and pressure fields.
Laminar flow on the inlet section of the 2D channel was modelled to test the computational algorithm. Turbulent flow and emission transport in a wind tunnel was modelled to verify the mathematical model.
20
E.А. ^aHH^KHH, A.B. CTap^eHKO
Series of computations of the flow influenced by natural convection in an ideal model of street canyon were performed using the presented mathematical model. Computations were performed for a 24 m height and 20 m wide street canyon. The constant emission source was placed in the center of the canyon near the floor. Analysis of the results has showed that in cases of heating the upwind side or the bottom of the canyon emissions are transported out of the canyon more intensively and maximal concentrations decrease by 10-15 % from the isothermal case. In case of heating the downwind side the structure of the flow changes significantly and maximal concentrations increase by 3-3.5 times.
The structure of the flow in the street canyon was investigated depending on the ratio of the street width to the height of the buildings. Both width of the street and height of the building varied from 5 to 40 m. The results show that increase in the height of the canyon decrease ventilation of the street canyon and increases local maximal concentrations of adverse emissions.
Keywords: admixture transfer, turbulent flow modelling, street canyon, buoyancy flows.
Citation: Danilkin E.A., Starchenko A.V. Modelling the transfer of road transport emissions in a street canyon. Computational Technologies. 2020; 25(2): 4-21. (In Russ.)
Acknowledgements. This work was carried out with financial support of the grant of the RF president No. MK-1723.2017.5 and Russian Fund for Basic Research and Tomsk Region according to the research project No. 16-41-700178 ra.
References
1. Starchenko A.V., Nuterman R.B., Danilkin E.A. Numerical study of turbulent flows and pollution transport in street canyons. Tomsk: Izdatel'skiy Dom Tomskogo Gosudarstvennogo Universiteta; 2015: 252. (In Russ.)
2. Balogun A.A., Tomlin A.S., Wood C.R., Barlow J.F., Belcher S.E., Smalley R.J., Lingard J.J.N., Arnold S.J., Dobre A., Robins A.G. In-street wind direction variability in the vicinity of a busy intersection in Central London. Boundary-Layer Meteorology. 2010; 136(3):489-513.
3. Kim J.J., Baik J.J. A numerical study of the effects of ambient wind direction on flow and dispersion in urban street canyons using the RNG-Turbulence model. Atmospheric Environment. 2004; 38(19): 3039-3048.
4. Danilkin E.A., Starchenko A.V. Large eddy simulation of turbulent flow and of pollutant transport in a street canyon. Proc. SPIE 9680, 21th Intern. Symp. on Atmospheric and Ocean Optics: Atmospheric Physics. 2015:968062-1-968062-6. DOI:10.1117/12.2205490
5. Chan A.T., So E.S., Samad S.C. Strategic guidelines for street canyon geometry to achieve sustainable street air quality. Atmospheric Environment. 2001; 35(24):5681-5691.
6. Park S.B., Baik J.J., Raasch S., Letzel M.O. A large-eddy simulation study of thermal effects on turbulent flow and dispersion in and above a street canyon. Journal of Applied Meteorology and Climatology. 2012; 51(5):829-841.
7. Wang P., Zhao D., Wang W., Mu H., Cai G., Liao C. Thermal Effect on Pollutant Dispersion in an Urban Street Canyon. International Journal of Environmental Research. 2011; 5(3):813-820.
8. Kim J.J., Pardyjak E., Kim D.Y., Han K.S., Kwon B.H. Effects of building-roof cooling on flow and air temperature in Urban Street Canyons. Asia-Pacific Journal of Atmospheric Sciences. 2014; 50(3):365-375.
9. Nuterman R.B., Baklanov A.A., Starchenko A.V. Modeling of aerodynamics and pollution dispersion from traffic in the urban sublayer. Mathematical Models and Computer Simulations. 2010; 2(6): 738-852.
10. Fei X., Xiaofeng L. The impact of roadside trees on traffic released PM10 in urban streetcanyon: Aerodynamic and deposition effects. Sustainable Cities and Society. 2017; (30):195-204.
11. Glazunov A.V. Numerical simulation of turbulence and transport of fine particulate impurities in street canyons. Numerical Methods and Programming. 2018; 19(1):17-37. (In Russ.) https://doi. org/10.26089/NumMet.v19r103
12. Ming T., Fang W., Peng C., Cai C., Richter R., Ahmadi M., Wen Y. Impacts of traffic tidal flow on pollutant dispersion in a non-uniform urban street canyon. Atmosphere. 2018; 9(3):82.
13. Launder B.E., Spalding D.B. The numerical computation of turbulent flows. Computational Methods in Applied Mechanics and Engineering. 1974; 3(2):269-289.
14. Van Leer B. Towards the ultimate conervative difference scheme: II. Monotonicity and conservation combined in a second order scheme. Journal of Computational Physics. 1974; (14):361-370. http: //dx.doi.org/10.1016/0021-9991(74)90019-9
15. Patankar S. Numerical heat transfer and fluid flow. New York: Hemisphere Publ. Corporation; 1980: 214.
16. Zhukauskas A.A. Convective transfer in heat sinks. Moscow: Nauka; 1982: 472. (In Russ.)
17. Loytsyansky L.G. Fluid and gas mechanics. Moscow: Drofa; 2003: 840. (In Russ.)
18. Hoydysh W.G., Dabberdt W.F. Kinematics and dispersion characteristics of flows in asymmetric street canyons. Atmospheric Environment. 1988; (22):2677-2689.
