CC BY
30АПВПМ-2019
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ТУРБУЛЕНТНЫХ ТЕЧЕНИЙ В УСТОЙЧИВО СТРАТИФИЦИРОВАННОЙ СРЕДЕ: ЦИРКУЛЯЦИЯ НАД ГОРОДСКИМ ОСТРОВОМ ТЕПЛА
Л, И, Курбацкая
Институт вычислительной математики и математической геофизики СО РАН, 630090, Новосибирск
УДК 532.517.4, 551.511.32
Б01: 10.24411/9999-016А-2019-10041
Аннотация Представлены результаты численного моделирования турбулентной структуры проникающей конвекции над городским островом тепла малого относительного удлинения в устойчиво стратифицированной покоящейся среде. Использованы представления градиентной диффузии для турбулентных потоков импульса и тепла, зависящие от трех параметров — кинетической энергии турбулентности, скорости её спектрального расходования и дисперсии температурных флуктуаций. Эти параметры находятся из замкнутых дифференциальных уравнений баланса в БАИв приближении описания турбулентности. Распределения осредненных полей скорости и температуры, турбулентных характеристик хорошо согласуются с данными измерений в контролируемом лабораторном эксперименте. Показано, что термическая стратификация вызывает различия в коэффициентах турбулентного обмена между активным скаляром теплом и пассивным скаляром массой. Различие означает, что предположение об идентичности коэффициентов турбулентной диффузии тепла и массы, используемое в обычных моделях турбулентности, дает существенные ошибки в оценивании тепло и массопе-реноса в термически устойчиво стратифицированной среде.
Ключевые слова: турбулентность, планетарный пограничный слой, городской остров тепла, крупномасштабная циркуляция, численное моделирование.
Введение
Циркуляция воздуха над городским островом тепла генерируется за счет энергии антропогенных источников в пределах городской черты. Наиболее интенсивной такая циркуляция наблюдается в ночное время при ясном небе и слабом окружающем ветре. Фундаментальными характеристиками, отражающими структуру ночного городского острова тепла, являются распределения полей скорости и температуры, а также интенсивностей турбулентности. Вертикальный турбулентный тепловой факел от локализованного источника (острова тепла) и, связанная с ним циркуляция, развиваются вследствие различия температуры между островом тепла и его окружением, которое имеет более низкую температуру. Интерес также представляет исследование механизма тепло и массопереноса в течениях окружающей среды. Течения в атмосферном пограничном слое часто стратифицированы по плотности (температуре), и плавучесть оказывает сильное воздействие на перенос скаляра (тепло, масса). Для описания проникающей конвекции, индуцируемой тепловым источником (островом тепла) и воспроизведения полей скорости, температуры и структуры турбулентности используется ИА.ХЯ приближение. Основные уравнения термогидродинамики берутся в гидростатическом приближении [1] в безразмерном виде для осредненных значений скорости и температуры в приближении Буссинеска. Для численного решения основных уравнений необходимо моделировать турбулентные потоки импульса и тепла. Для этих целей используется трехпараметрическая Е — е— < в2 > модель
Работа была выполнена в рамках госзадания Института вычислительной математики и математической геофизики СО РАН № 0315-2019-0004 и частичной поддержки Российского фонда фундаментальных исследований (код проекта № 17-01-00137) и Российского фонда фундаментальных исследований и Правительства Новосибирской области (код проекта 18-48-540005).
ISBN 978-5-901548-42-4
турбулентности. В этой модели кинетическая энергия турбулентности (КЭТ) Е =< щщ >, скорость ее диссипации е и дисперсия флуктуации температуры в2 находятся из дифференциальных уравнений переноса, а вектор турбулентного потока импульса и тепла определяются с помощью явной алгебраической модели, полученной из дифференциальных уравнений переноса в приближении локально равновесной турбулентности. Проникающая турбулентная конвекция, индуцируемая постоянным потоком тепла Но с поверхности пластины диаметра Б моделирует в лабораторном эксперименте прототип городского острова тепла малого относительного удлинения (гг/Е) << 1) [4]. Движение среды рассматривается как обладающее осевой симметрией. Подробности численного моделирования проникающей турбулентной конвекции и связанной с ней циркуляцией могут быть найдены в [2] и [3], и здесь они не приводятся. Ниже будут представлены только результаты численного моделирования, описывающие основные характеристики проникающей турбулентной конвекции и связанной с ней циркуляцией. Основное внимание в статье было уделено описанию турбулентной диффузии пассивного скаляра (примеси).
1 Модель турбулентной диффузии примеси над нагретым источником
Для описания рассеяния пассивной примеси базовая мезомасштабная модель турбулентного атмосферного пограничного слоя должна быть дополнена уравнениями для осредненной концентрации С, вектора турбулентного потока примеси < щс > и корреляции между флуктуациями концентрации и температуры < св >. Поскольку результаты моделирования переноса активной примеси (тепла) над нагретым источником (в данном случае диском) [2], [3] удовлетворительно согласуются с данными измерений [4], для турбулентного потока пассивного скаляра (массы вещества) также может быть сформулирована полностью явная анизотропная алгебраическая модель при упрощении дифференциального уравнения переноса для вектора турбулентного потока скаляра <щ с > и корреляции < св > в приближении локально равновесной турбулентности. Таким образом, определяющие уравнения для поля концентрации пассивного скаляра в цилиндрической системе координат имеют вид [5]:
дС 1 д „л д „л 1 д д
— + -— гСиг + — Сиг = —т-г < игс > — — < щс >, 1)
ел г дг дг г дг дг
д<игс> 1 д д 2 дС дС диг
-—--+ -— г < иг с > иг + 7— < иг с > иг = - <иг > —--< игиг > —--< иг с > ——
ел г дг дг Ог Ог Ог
диг 1 д. Е 2 д<иг с>. 1 д. Е д < иг с>.
- <ихс > —--+ -^-(га^ — <и2 > ---) + -^-(га^— < игиг > ---) +
Ог г Ог £ Ог г Ог £ ог
д Е д < иг с > д Е 2 & < игс > £
+ ^г\гаЛа— < иг иг > ---) + -^-(гаи— <иг > ----а1с — < иг с >, (2)
Ог £ дг дг £ 2 Ог Е
д<и2с > 1 д д 2 дС дС диг +—ТТ— г < и2с > их + — < ихс > их = — < и7 > —--< иги2 > —--< игс >
сЛ г Ог Ог Ог Ог Ог
£ дих 1 д Е 2 д < их с > 1 д Е д < их с >
—«1с <иг с > — < иг с > —--+ - —(таь— <иг > ---) + - —(таь— < иг иг > ---) +
Е Ог г Ог £ Ог г Ог £ Ог
д Е д<иг с> д Е 2 д<иг с> д
+—(гаи— < игиг > ---) + — (гаи— < иг > ----(1 — а.2с) <св > (3)
Ог £ Ог Ог £ Ог в0
д < св > 1 д
т
+ --^[г<св>иА + и: [<св> и*]
д_
дг
«2г
дг
£ 1 д
—«Зст^ <св > +— — Е г Ог
га.2&
Бе*
д < св > дг
Бе*
д <св> ( п дС п дС дв дв*
----<иг в > ——+ < иг в > ——+ < иг с > ——+ <иг с> —
дг дг дг дг дг
(4)
^ Е2 ^дс v/2ДS Гг , 1Л11,ди1 дии — < с >= — с^7 • ^ +(1 — + ) х
^ \ _ 1 / диг ди*. т дС 1 — а2с Е а
х(1 — °2°)в<2(щ + ^Щ + во <>• (5)
В уравнениях (1)-(5): Ur — средняя горизонтальная скорость, Uz — средняя вертикальная скорость (вид уравнений для осредненных скоростей и температуры приведен в [5]), < urc >, < uzc > — вертикальный и горизонтальный потоки массы вещества, ur,uz — горизонтальная и вертикальная турбулентная флуктуация скорости, <Ur >, <UZ > — горизонтальная и вертикальная составляющие КЭТ, < urQ >, < uz0 > — вертикальный и горизонтальный турбулентные потоки тепла, ©о — средняя температура, Dt = 2REZ /'е — коэффициент турбулентной диффузии скаляра, vt = С^Е2/е — турбулентная вязкость, Св = C^/Sct, R = тс/т — параметр отношения характерных масштабов времени скалярного и динамического турбулентных полей. Численные значения констант модели С^ = 0.095, Sct = 0.9, R = 0.6, С10 = 3.28, Cze = 0.5, С1с = 0.4, С в = 0.4, а1с = 4.0, aZc = а^с = 0.4, a1s = aZs = 0.22. Сформулированная диффузионная модель (1)-(5) используется для решения задачи о распространении пассивного скаляра (массы вещества) от непрерывного источника заданной производительностью Q, расположенного на поверхности. На источнике задается постоянный вертикальный поток массы вещества
-Dm(d(:/dz) = нс (6)
где Dm = v/Se — коэффициент молекулярной диффузии, v — коэффициент кинематической вязкости, Se — молекулярное число Шмидта.
2 Коэффициенты турбулентного обмена активного и пассивного скаляра
В обычных градиентно-диффузионных моделях (типа Буссинеска) вертикальные потоки тепла и массы в термически стратифицированном течении имеют вид
дТ
< w0 >= -Kh —,
д
д С
< wc >= -Kc —,
д
(7)
(8)
Где т — турбулентная флуктуация скорости в вертикальном направлении, К и Кс — коэффициенты турбулентного обмена активного и пассивного скаляра.
Главное допущение в таких моделях состоит в том, что турбулентные потоки и массы пропорциональны локальным градиентам средней температуры и концентрации. При использовании таких моделей обычно полагается, что коэффициенты турбулентного обмена теплом К^ и массой Кс равны друг другу. Однако в [7] на основе модели турбулентности второго уровня замыкания для транспортируемых скаляров тепла и влажности показано, что эти коэффициенты не равны друг другу. Покажем различие в коэффициентах К^ и Кс. Для этого запишем уравнение для турбулентного потока тепла
- < uóe >=
СГ Е! ^ - СЯЕ
£ dz С ig £
1 ,дЫ dU,
{[^ + (1 С!9 )а,]-( ^ + ^)
dxi
„ N 1,dUi dUj. , дТ 1 -С!вЕ g ■>
x(1 - С!в)at-(^ + ^+-zl-JL < oZ >.
2 dxj dxi dxj а1в £ ©о
(9)
Для оценки коэффициентов К^ и Кс сохраним лишь ведущие члены с вертикальными градиентами температуры и концентрации, а также члены плавучести в соответствии с (5) и (9).
EZ
- < щв >= Ст — V2R
(1 + Сзв) <ЩЕ> -Сзв_
дТ д
^ ©I *Z >
(10)
- < ujC >= С^
Е ,— Си — V2R
< wZ >
(1 + а3с)^/3Е- -аз\
дС - 1-а^А V2RE<c0>
dz а 1с ©о £
(Н)
Деление (11) на (10) дает отношение коэффициентов
KС = 1 - (1 - a.zc)/(aAcCD) ■ А-1 ^ < св > /Е(дС/дг)-1 Kh = 1 - (1 - Cze)/(СхвСт) ■ А-1^ < 6Z > /Е(дТ/дz)-1
х
где Кс = Кн = <г/д> ~~ коэффициенты турбулентного обмена для массы и тепла, через Дх и Д2
обозначены выражения в квадратных скобках (11) и(10) соответственно. Из (12) следует, что если оба скаляра (тепло и масса) являются пассивными примесями (члены плавучести пренебрежимо малы), то Кн = Кс. Очевидно, что различие между коэффициентами Кн и Кс может возникать, когда один из скаляров (тепло) действует как активная примесь. Величина отклонения зависит от соотношения корреляций < в2 > и < св >. Конечно, о различии в коэффициентах Кн и Кс можно судить по результатам численного моделирования, которое было проведено с использованием полных (неупрощенных) выражений для турбулентных потоков активного скаляра (9) и пассивного скаляра (5). Результаты обсуждаются ниже.
3 Результаты численного моделирования
Представленные результаты, описывающие развитие теплового факела над источником тепла, его турбулентную структуру, соответствуют квазиустановившемуся состоянию термической циркуляции от поступающего тепла от нагревателя, при котором интенсивность острова тепла АТт = Тт — То и поверхностный поток тепла оставались неизменными с течением времени (Тт и Т0 — температура поверхности в центре острова и вне его, соответственно). На рис. 1а показана теневая фотография теплового факела, имеющего верхнюю границу в форме шляпы, что обусловлено воздействием устойчивой стратификации окружающей среды, которая создает запирающую инверсию в верхней части перемешанного слоя при z/zi = 1 (zí — высота перемешанного слоя, ¿ — текущая высота). В лабораторном эксперименте [6] получены картины визуализаций движения воздуха над островом тепла с нагревателем в виде проволоки круглого сечения, а не круглого диска, как в экспериментах [4]. Движение воздуха было, приближенно, плоскосимметричным. (Реальную ситуацию городского острова тепла можно рассматривать, по-видимому, как комбинацию плоско-и осесимметричного движений). Общим для обоих экспериментов является подавление высоты плу-ма устойчивой стратификацией среды, увеличение бокового движения и турбулентности плума. Движение среды над "несущей" вихревой парой (рис. 16) фиксируется картинами визуализаций Noto [6] и может быть связано с локальными циркуляциями над реальным островом тепла. Тепло при этом не транспортируется от урбанизированной поверхности выше в атмосферу, а аккумулируется в локальных циркуляциях.
Рис. 1: Циркуляция над нагреваемым источником (теневая фотография) (а), картина линий тока циркуляционного течения (б), профили средней температуры в различных сечениях (в).
На рисунке 16 видны два крупномасштабных вихря с вращением левого вихря против часовой стрелки, правого — по часовой стрелке, которые формируют основное восходящее движение в центре, простирающееся до слоя вовлечения {z/zi) = 1) и нисходящее движение на периферии. При этом высота теплового факела подавляется устойчивой стратификацией, а боковое движение и турбулентность факела увеличиваются. На рис. 1в хорошо прослеживается область сильной и положительной плавучести в нижней части теплового факела и область слабой отрицательной плавучести вблизи слоя вовлечения {z/zi = 1). Вычисленные вертикальные распределения температуры в центре острова указывают на хорошее перемешивание в нижней и центральной части факела. Такой характер поведения распределения температуры с высотой относится к тем реальным ночным пограничным слоям, в которых преобладают неустойчивые (конвективные) условия вследствие восходящего потока тепла от урбанизированной поверхности при слабом ветре. На рис. 1 в также
видно, что профили температуры внутри факела имеют характерное "вздутие": температура внутри факела оказывается ниже температуры вне его на той же высоте, фиксируя тем самым область отрицательной плавучести, вследствие возвышения факела в центре, связанное с эффектом перекрещивания. Поведение радиальной скорости по вертикали показывает возрастание скорости натекающего потока в нижней части теплового факела, который затухает в зоне вовлечения, скорость оттекающего потока наверху также возрастает от центра острова тепла к периферии (рис. 2 а).
Структура турбулентности теплового факела представлена на рис. 26 в виде распределения среднеквадратичных флуктуаций горизонтальной и вертикальной турбулентных скоростей в зависимости от безразмерной высоты в центре факела. Как данные лабораторных измерений, так и результаты вычислений показывают, что большие значения аи /Ш^ и /Ш^ внутри перемешанного слоя быстро убывают с высотой выше слоя вовлечения. Область вычисленного максимального значения аи/Ш^ простирается по горизонтали дальше, чем значение /Ш^ вследствие горизонтального дивергентного течения, вызываемого устойчивой стратификацией верхней части факела (Шр — конвективный масштаб скорости). Следует отметить, что простая градиентная модель не только правильно описывает характерные особенности распределения аи /Шр и <7ы/Шв-, но и удовлетворительно отражает их анизотропный характер.
N N
1,61,20,80,40,0-
Z / Z. -А— 0.071 -•— 0.14 -♦— 0.28 -▼— 0.42 -ж— 1.05
-0,4
-0,2
0,0 0,2 U /
0,4
0,6
б
Рис. 2: Профили горизонтальной средней скорости в различных сечениях по радиальной координате (а). Распределение горизонтальной ау/Шр и вертикаль ной /Шр дисперсий турбулентной скорости на оси факела над островом тепла (б). Светлые треугольники и квадраты — эксперимент [4], зачерненные треугольники и квадраты — вычисление, зачерненные окружности — вычисленный профиль щи^ >.
а
Турбулентная структура температурного поля представлена на рис. За распределением среднеквадратичных флуктуаций ат/Т^, где Тр — конвективный масштаб температурного поля. (Выражения для W& и Тр см. в [3], [4]). Профиль имеет характерный вид с убыванием от максимального значения вблизи поверхности до минимального значения при z/z{ ~ 0.85, как в эксперименте, так и в RANS приближении с вычисляемыми из уравнений баланса параметрами Е, s, < 62 >. Второй максимум вблизи слоя вовлечения обусловлен большим градиентом скорости генерации < в2 >, подтверждается вертикальным распределением температуры на рис. 1 в. Различие в профилях дисперсии в верхней части перемешанного слоя связано с известным дефектом использования упрощений градиентного типа при аппроксимации процессов турбулентной диффузии (третьих моментов термогидродинамических полей) в RANS приближении моделирования турбулентности с замыканием на уровне вторых моментов. На рис. 36 приведен вертикальный профиль диссипации кинетической энергии турбулентности, вычисленный по трехпараметрической модели, в сравнении с LES данными Moeng и Wyngaard [9] (взятыми из работы Canuto et.el. [8]). Максимум диссипации в области слоя вовлечения связан, по видимому, с высоким уровнем турбулентности вдоль "крышки" термического плума [4], [6], который в свою очередь может быть обусловлен дроблением вихрей, достигающих слоя вовлечения, и генерированием случайных сдвиговых слоев.
На рис. 4а приведено распределение концентрации, полученное в результате численного решения систе-
а б
Рис. 3: Вертикальный профиль дисперсии турбулентных флуктуацпй температуры над островом тепла на оси факела: светлые квадраты — эксперимент [4], зачерненные квадраты — вычисления (а). Вертикальное распределение диссипации энергии турбулентности: светлые квадраты — эксперимент [9], зачерненные квадраты — вычисления (б).
0,0 Л
0,0 -1-'-L
-1,5 -1,0 -0,5 0,0 0,5 1,0 1,5 г D
б
Рис. 4: Изолинии концентрации пассивного скаляра (массы) в диффузионном факеле над поверхностным источником: (а) — анизотропная диффузионная модель, (б) — модель Буссинеска. Профили средней концентрации пассивного скаляра в различных сечениях по вертикальной координате над непрерывным поверхностным источником (б).
мы уравнений (1)-(5), т.е. эйлеровой диффузионной модели, а на рис. 46 показаны изолинии концентрации, вычисленные по модели вихревой диффузии Буссинеска: — < wc >= vt/Sct • dC/dz), в которой плавучесть не принимается во внимание. Рис. 4а и 46 показывают, что под влиянием проникающей турбулентной конвекции, порождаемой двумя вихрями, примесь переносится от поверхностного источника в перемешанный слой, распространяется далее в слой инверсии, претерпевает гравитационное рассеяние в горизонтальном направлении и проникает за пределы инверсионного слоя в область значений z/z{ > 1. Проникновение примеси за пределы инверсионного слоя получено в лабораторных экспериментах по распространению плавучих факелов пассивной примеси от непрерывного и мгновенного источников в конвективном пограничном слое [10].
Однако имеются различия в распределении средней концентрации. Модель Буссинеска без учета эффектов плавучести дает более медленную диффузию примеси в вертикальном (а, следовательно и горизонтальном) направлении по сравнению с моделью (5), учитывающей эффекты плавучести. При использовании модели (5) уровень концентрации оказывается выше как вблизи поверхности, так и по всей высоте перемешанного слоя. В области слоя инверсии (при г/г^ « 1) значение концентрации в центре острова тепла, вычисленное с помощью модели (5), примерно в два раза выше, чем рассчитанные по модели Буссинеска. На верхнем рис. 5а приведены результаты вычислений с использованием модели (1)^(5), на нижнем рис. 5 а — модель Буссинеска.
На рис. 5б приведено поведение коэффициентов турбулентного обмена активного и пассивного скаляра. Можно видеть, что поток тепла меняет знак по высоте и достигает отрицательной максимальной величины в инверсионном слое (при г/г^ = 1). Поток массы по вертикали остается положительным.
Г d <w9> , <wc>10"2
а б
Рис. 5: Профили средней концентрации пассивного скаляра (массы вещества) в различных сечениях по вертикальной координате над непрерывным поверхностным источником (а). Вертикальные профили турбулентных потоков тепла < ыв > (кривая 1) и массы < > (кривая 2) в центре острова: анизотропная диффузионная модель (б).
Заключение
Представленные результаты численного моделирования структуры циркуляции над островом тепла позволяют сделать следующие выводы. Описание турбулентности теплового факела в ИЛ.ХЯ приближении с включением в качестве искомых моментов гидротермодинамических полей второго порядка позволяет воспроизвести структурные особенности проникающей конвекции над источником тепла, включая такие тонкие эффекты, как перекрещивание вертикальных профилей температуры теплового факела с образованием области отрицательной плавучести, свидетельствующей о развитии куполообразной формы у верхней части факела. Вертикальная дисперсия турбулентной скорости и дисперсия турбулентных флуктуаций температуры удовлетворительно согласуется с опытными данными.
Результаты численного моделирования распространения пассивной примеси от непрерывного поверхностного источника показывают, что явные анизотропные выражения для турбулентных потоков тепла и массы с точным учетом эффектов плавучести дают более эффективную диффузию примеси от поверхностного источника в перемешанный слой и выше. Вычисленное отношение коэффициентов турбулентного обмена для пассивной примеси и активной примеси (тепло) в термически устойчивом стратифицированном течении над поверхностными источниками тепла и массы свидетельствуют об их различии.
Список литературы
fl] R. A. Pielke Mesoscale meteorological modeling // Academic Press, New York, Boston 2002. p. 676.
[2] Курбацкий А.Ф., Курбацкая Л.II. Трехпараметрическая модель турбулентности для атмосферного пограничного слоя над урбанизированной поверхностью // Изв. РАН. Физика атмосферы и океана. 2006. Т. 42, № 4, С. 476-494.
[3] Курбацкий А.Ф., Курбацкая Л.И. Турбулентная циркуляция над поверхностным источником тепла в устойчиво стратифицированной окружающей среде // Теплофизика и аэромеханика. 2016. Т.23, № 5, С. 703-719.
[4] Lu J., Arya P. et al. A Laboratory Study of the Urban Heat Island in a Calm and Stably Stratified Environment. Part I and II // J. Appl. Meteor. 1997. Vol.36, No.10, 1377-1402.
[5] A. F. Kurbatskii and L. I. Kurbatskaya Urban aeromechanics: Turbulent circulation and contaminant dispersion above city in stably stratified environment AIP Conference Proceedings 2027, 030024 2018. https://doi.Org/10.1063/l.5065118
[6] Noto K. Dependence of heat island phenomena on stablestratification and heat quantity in a calm environment //Atmospheric Environment. 1996. V.30.№3.P.475-485
[7] Warhaft Z. Heat and Moisture Flux in hte Stratified Boundary Layer // Q. J. Royal Meteor. Soc. 1976. V. 102. P. 703.
[8] Canuto, V.M., Minotti, F., Ronchi, C., Ypma, R.M.,1994: Second-Order Closure PBL Model with New Third-Order Moments: Comparison with LES Data. J.Atmos.Sci. 1605-1618.
[9] Moeng C.-H., Wyngaard J.C. Spectral analysis of large-eddy simulations of the convective boundary layer // J. Atmos. Sci. 1988. V.45. P.3574-3587.
[10] Snyder W.H., Lawson R.E.Je., Shipman M.S., Lu J. Fluid Modeling of Atmospheric in Convective Boundary-Layer // Bound. Layer. Meteor. 2002. V. 102. 3. P.155.
Курбацкая Людмила Ивановна — ст. науч. сотр. Института вычислительной математики и математической геофизики СО РАН;
e-mail: [email protected] Дата поступления — 27 мая 2019 г.
