CC BY
73
Раздел II. Математическое моделирование в аэро- и гидрофизике
УДК 519.86
В.К. Гадельшин, Д.С. Любомищенко, А.И. Сухинов
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ПОЛЯ ВЕТРОВЫХ ТЕЧЕНИЙ И РАСПРОСТРАНЕНИЯ ЗАГРЯЗНЯЮЩИХ ПРИМЕСЕЙ В УСЛОВИЯХ ГОРОДСКОГО РЕЛЬЕФА МЕСТНОСТИ С УЧЕТОМ k - £ -МОДЕЛИ
ТУРБУЛЕНТНОСТИ
Предложен подход моделирования поля ветровых течений на основе системы уравнений Навье—Стокса с учетом неровностей подстилающей поверхности. В качестве численного метода используется SIMPLE-метод для неразнесенных сеточных значений давления и компонент скорости. Предложено описание учета сжимаемости и турбулентности .
Моделирование поля ветровых течений; моделирование распространения примеси; SIMPLE-метод; неразнесенные сетки; k - £ -модель турбулентности.
A.I. Sukhinov, V.K. Gadelshin, D.S. Lubomishenko
MATHEMATICAL MODELING OF AIR VELOCITY AND POLLUTION DISTRIBUTION IN CITY ATMOSPHERE BOTTOM LAYER CONSIDERING k - £ -TURBULENT MODEL
In the paper there is an approach for wind field modeling in city conditions based on Navie-Stoks equations system. SIMPLE method on collocated grid is used as basic for solving the system. Compressibility and k - £ -turbulence model have been considered.
Wind field modeling; air pollution modeling; SIMPLE-method; collocated grids; k - £ -turbulence model.
1. Структура приземного слоя атмосферы города. Приземны й слой атмосферы города имеет сложную трехмерную структуру (рис. 1). Часто при описании данного слоя используют широкий спектр параметров шероховатости (здания, трава и т.д.). Для корректного описания процессов в приземном слое необходимо выделять подстилающий приземный слой, в котором детально учитывать такие элементы шероховатости, как здания, деревья, кустарники и т.д.
При рассмотрении подстилающего приземного слоя необходимо выделять горизонтальное и вертикальное направления. В горизонтальном направлении за счет изменения коэффициентов шероховатости и термических свойств поверхностей его структура может меняться в большей или меньшей степени. Для описания явлений на этом уровне вводится понятие внутреннего приземного слоя. Во внутреннем приземном слое должны учитываться свойства шероховатости и термические свойства, форма и размер уличных каньонов и т.д.
Вертикальное направление подстилающего приземного слоя описывается с помощью параметра шероховатости zH , который варьируется в диапазоне от 0
до средней высоты элементов шероховатости. В данном слое можно применять модели высокой размерности для описания движения среды, в которых должны учитываться форма зданий, термические свойства поверхностей, условия движения транспортных средств. Такое моделирование может дать существенную информацию о градиентах величин между наветренными и подветренными сторонами улицы. Затем, по соотношению высоты и ширины уличных каньонов можно судить о процессах, протекающих в подстилающем приземном слое.
Подстилающий приземный слой города включается в шероховатый подслой.
*
Высота данного подслоя обозначается г и зависит от высоты и плотности эле.
параметров ветровых течений от свойств шероховатости. Для простоты вводится ,
приземного слоя города:
г* = агн,
где а - параметр, изменяющийся в диапазоне от 2 до 5.
Выше шероховатого подслоя находится внутренний подслой. В случае, если исследуется стационарный поток среды, слабо зависящий от граничных условий и , -
.
слоями, расположенными выше, и является тонким.
Над внутренним подслоем находится городской граничный слой. Данный слой представляет собой зону перемешивания и связывает все предыдущие слои с метеорологическими.
Гарсцкхэв г/г.н:пхігі слан
ППҐЯ
Рис. 1. Схема слоев городской среды
2. Моделирование приземной аэродинамики. Моделирование поля ветровых течений в условиях городской застройки является важным этапом в задаче мониторинга загрязнения окружающей среды. В данном разделе предложен вариант 81МРЬЕ-метода для построения поля течений ветра над неровной подстилаю, , а также предложена неявная процедура Стоуна для решения семиточечных СЛАУ .
Запишем базовую трехмерную систему уравнений для сжимаемого вязкого газа в тензорном виде [1, 2, 3]:
др д(Рип)
дї
- + -
дх
0,
д(Ри) + д{риБЩ)
дї
дх-
др д
дх дхі
і ]
,н
дщ
дх]
р — рЯТ,
д(рТ) + д(ри}Т) =^Э_
дї
дх:
дх:
ке// дТ Срт Ц
+
о
С
(1)
(2)
(3)
(4)
рт
Здесь р - плотность среды, /Чед - эффективная вязкость, р -давление, щ - компоненты вектора скорости в направлении трех координатных осей, Ы/ - компоненты вектора скорости в направлении трех координатных осей, удов-
(1).
(1), (2), (4)
в форме уравнения переноса для переменной Ф:
д(РтФ)_+д(РтЫи Ф)= д
дї
дх:
дх:
е//
дФ
дх:
+ Б— .
(5)
Значения коэффициентов для разных уравнений приведены в табл. 1.
1
Значения коэффициентов для различных типов уравнений
№ уравнения Ф Г 1 е// Бе//
1 1 0 0
2 иі №е// -др/дх, +рВі +
4 Т ке// / /Срт ОТ/ Срт
Применим к уравнению (5) метод дискретизации, основанный на интегрировании по контрольному объему (ЮСУ) [4, 5]. Контрольный объем, по которому
, . 2.
дх.
(6)
Далее, расписывая интегралы по объему и поверхности, получим:
рАй
2Аї
(3Ф "+1 - 4Фп + Ф ”-1) — ) Ф Р+1 - 2—
(7)
где
а'в — 3рАй и О — рАй(4Фп -Фп-1).
2Аї
2Аї
ч
w 1 1 1 1 п р1 -• е
•- - — — -•
5» 1 1 1 1
ЬА
Рис. 2. Контрольный объем
В предыдущих формулах индекс п означает п -й временной слой.
Поверхностный интеграл ¥с = |ры{]ФёБ в трехмерном пространстве пред-
я
ставляется в виде суммы шести интегралов по соответствующим плоскостям кон:
, ,
¥с — ¥с + ¥с + ¥с + ¥с + ¥с + ¥с
'№ е Б п Ь ї
гс — [ риг—йБ ~ т —
е у л е е
(8)
(9)
где
те = |ры • ПйБ ~ (ры )т 1 Я,.
Аналогично расписывается поверхностный интеграл для диффузионных членов:
(10)
где
дх
Б их}
рй — рй + рй + рй + рй + рй + рй
(11)
(12)
(2)
компоненты скорости ы можно представить в виде
(13)
Б
Б
В общем случае интеграл по объему для источникового члена уравнения можно представить в виде
03 = 15фйа- 5фДУ, (14)
а
где
3Ф = 3и + 3рфр. (15)
Для аппроксимации конвективных членов используется схема задержанной коррекции. Суть данной схемы сводится к тому, что в процессе осуществления итераций вначале конвективные члены аппроксимируются по противопоточным схемам, а когда происходит установление итерационного процесса, то аппроксимация осуществляется с помощью более точной центрально-р^ностной схемы. Например, на плоскости е контрольного объема поверхностный интеграл (9) можно записать так:
F; = т ф^+т (ф ( - фГ )"-1. (16)
С учетом вышесказанного дискретную форму для переменной Ф мож но записать в виде
4 Ф р + ¥с = ¥й + О + О. (17)
Или, если речь идет об уравнении движения (2), то в правой части появляется , :
а‘рФр + ¥с = ¥й + 0р + О + О. (18)
, (17)
ар Ф р + X а1 ФI = Ор , (19)
где
аР = а{Р - X Щ , l = E,W, N, 5, T, B, Qp = Q + Qp + QS + Qc, (20)
Qc =
i \UDS / \CDS
(r) ~(r)
Г S . / . Г S„
m-1
aE = min (rne ,0)-e—^~, aW = min (me ,0)--------------------,
xE - xP xE - xP
■ ( ■ n) reSe ■ ( ■ (\) reSe
aN = min (me ,0)-e—^ aS = min (me ,0)---------------------------e—e—
yE - Ур , yE - Ур ,
• ( ■ n) reSe ■ ( ■ n) reSe
aT = min (me ,0)-e—^ aB = min (me ,0)---------------------------e—e—
ZE - ZP ZE - ZP
Уравнение неразрывности можно записать в виде
m* + mtw + m* + m* + m* + m*b =Am*P, (21)
или
(pSu')e - (pSu)w + (pSv')n - (pSv')s + (pSw')t - (pSw')b + &m*P = 0 . (22)
Звездочки означают, что потоки массы считаются с использованием прогнозных компонент вектора скорости, не удовлетворяющих в общем случае уравнению .
Алгоритм SIMPLE осуществляется в несколько этапов:
♦ производится вычисление прогнозных компонент вектора скорости (в общем случае полученные значения не удовлетворяют дискретному аналогу уравнения неразрывности);
♦ по прогнозным значениям компонент вектора скорости вычисляется поправка к давлению;
♦ с использованием поправки к давлению вычисляются компоненты обновленного вектора скорости, удовлетворяющего дискретному аналогу урав-
.
*
Прогнозное значение для ue компоненты вектора скорости можно записать так:
*
и, =■
Qp - Qp - Т
aiui
(23)
a
Остальные компоненты вычисляются аналогично. Далее решается уравнение для поправки к давлению:
a
lp'P + £app; = -ДmP, l = E,W, N,5,T,B,
(24)
где
aE =-
( pS2 ap =- > UW f pS2 ] ap =- > UN ( pS2 , ap =- f ps 2 ^
V ap e ^ aP w ^ aP n ^ ap / s
aTp = ■
r pS_2 Л
v aW Л
ґ pS2 Л
V aP Jb
, aP
ap , l = E,W, N, S,T, B.
После того, как поправка к давлению найдена, по явным схемам корректируем прогнозное значение вектора скорости и находим давление на текущей итерации:
> * Se i ' ' \
ие = ие-[Pe - Рр ),
aP
* . '
Р = Р + Р .
(25)
(26)
Для численного решения уравнения используется 81Р-метод. Сеточная задача представляется в виде разреженной СЛАУ АФ = Ь .
Далее ищется приближенное разложение вида
м = ьи = А + N, (27)
где м - разреженная матрица, приближенно равная А , N - разреженная матрица с малыми элементами, Ь - нижнетреугольная матрица с тремя диагоналями, и - верхнетреугольная матрица с единичными элементами на диагонали и тремя побочными ненулевыми диагоналями.
Основной задачей является построение таких матриц Т и и , чтобы их произведение М аппроксимировало оператор А как можно лучше. В результате перемножения матриц Т и и указанного вида матрица М будет иметь шесть дополнительных диагоналей, вклад которых должен компенсироваться с помощью матрицы N .
,
(М Ф)р = МР Ф р + МЕ Ф Е + Мш Фш + N Ф N + М5 Ф 8 + Мт Фт + Мв Ф в +
+ М NW Ф NW + М' БЕФ БЕ + МТШ ФТШ + М ЕВФ ЕВ + МТБ ФТБ + М т Ф т
.(28)
Матрица N должна содержать шесть дополнительных диагоналей матрицы М . Необходимо задать элементы оставшихся диагоналей так, чтобы NФ « 0 . В результате получаем уравнение
М (ф -Ф* ) + М ( -Ф ) + М (ф -Ф* ) +
NW \ NW ш) БЕ\ БЕ БЕ ! ТШ \ ТШ ТШ )
+М (ф -Ф* \ + М (ф -Ф* \ + М (ф -Ф* ) = 0
ЕВ \ ЕВ ЕВ ) ТБ \ ТБ ТЕ ) Ш \ Ш Ш ) ’
(29)
где ф* - аппроксимации для соответствующих Ф . Например,
= +Ф N -ФР ), а< 1.
Формулы разложения имеют вид
Т1 = Тр =
Т = - А
ТВ ~
'1 , (л1 -NI■N: Т1\-N¡N¡
1 + а[им ‘ ( + и
Т1 =
Б
и1-1** + и1Т-^
і - А:
(і + а(-1 + и1Т-1)
Л1р + а[Ьви1~№т + 1^и1~Н]) + а[Ьви1~№т + Ё$и1~1) +Л, (30)
+а(и1т-т + ЬиТ-1) - )1т-тщ - 1^и1Е~т - 4и-1 и1м = (-(т -а(и1-№т + Ь^и1~т ), иЕ = (- аЕ - а(Е-тт+АиЕ1), и‘т = (-( -а(иТ-т + 4и1т-1 ).
Далее осуществляется последовательное решение исходной задачи обычным методом ьи -р^ложения:
" (31)
Ьи8п+1 =рп, и8п+1 = Г1рп = Г,
(32)
Rl =( -8l = Rl
llsRl -1
- UN 8
LwR
і - N
l pl - NiNj
■IBR
)/)p
l+1 UlE8l+Nj -U‘T8l+NiNj
(33)
(34)
3. Неразнесенные сетки. Записывая уравнения модели на неразнесенных
, , .
*
Прогнозное значение компоненты ие должно быть получено с учетом последую-
щей интерполяции и , для получения информации о компоненте скорости в центре контрольного объема. Также необходимо избежать осцилляций в решении уравнения для корректировки давления, с которыми обычно сталкиваются в ходе непосредственного применения формул для разнесенных сеток к случаю разнесенных. Для того чтобы избежать этих проблем, формула (23) модифицируется:
и =
(и*)е -ДЦ
V ap J
8p
-і \
V
8x
Je
8 p
-і \
V
8x
Je
(36)
означает соответствующую интерполяцию,
где черта сверху ДЧ =( - хр )ДУ .
, ие
сать в виде
ч\ ( х'\
а
и' = -ДЦ
1
\ap j
V
8p
8x
= -S
У e
1
\ap j
(cE- pP ).
(37)
Уравнение для поправки давления получается способом, аналогичным случаю с разнесенными сетками.
4. Учет сжимаемости атмосферы. Задача моделирован ия поля аэродинамических течений в приземном слое атмосферы города является важным этапом в проблеме описания распространения загрязняющих веществ. Сложная форма подстилающей поверхности существенно влияет на трансформацию потоков в около. .
В сжимаемом случае уравнение неразрывности можно записать в виде
Дї
+m* + К + m* + m*+m*+m* = Q.
(38)
Член Q‘m отражает дисбаланс, вызванный массовыми потоками m и плотностью pm-1, не скорректированными процедурой SIMPLE. В случае сжимаемой
атмосферы массовый поток через грани контрольного объема зависит не только от нормальной компоненты скорости, но и от плотности среды. Поэтому, для того чтобы компенсировать дисбаланс Q*m , необходимо ввести корректировки для скорости v'n и плотности среды р .
Для скорректированного массового потока на грани e контрольного объема можно записать выражение вида
т„
СРт-1 +P')e ( + V'n )e .
(39)
А корректировка массового потока запишется в виде
те = (pm-1 Sv'n ) + (* Sp')e + (p'v'nS )e.
(40)
Последний член правой части уравнения (40) можно не учитывать, так как он обладает вторым порядком малости по сравнению с другими членами (40). Выражения для оставшихся двух членов нужно получить. Первый член правой части (40) , . -
ражение для нее в точности совпадает с выражением для несжимаемой среды:
1
Apv
ґ 8р'л Sn
(41)
где п - внешняя нормаль к грани контрольного объема.
Запишем выражение для р , связанное с корректировкой давления р :
/ /-Г /
Р = Срр ,
(42)
где Cp .
Таким образом, второй член правой части выражения (40) запишется в виде
СрШ
V
Ре .
(43)
В итоге, подставляя (41) и (43) в (40), мы можем записать все выражения для корректировки массового потока:
те = (pm-1S АУ )е
АрУ’
S р' Sn
(44)
Для вывода уравнения для корректировки давления необходимо записать корректировки массового потока для всех шести граней контрольного объема, применяя необходимую аппроксимацию для компонент давления, находящихся на гранях контрольного объема.
Записывая уравнение неразрывности для поправки массовых потоков и плотности, получим:
(Рр )АУ
At
+ тЄ + mw + m'n + ms + mt + mb - Qm = 0.
(45)
Подставляя выражение для рр из (42) в (45), можно записать итоговое уравнение для поправки давления, учитывающее сжимаемость среды:
app'+X А/Р/ =-Q*
(46)
Уравнение (46) решается в рамках полунеявной итерационной процедуры SIMPLE так же, как и в случае несжимаемой среды. Данный подход позволяет бо-
лее точно описывать движение воздушных масс в области со сложной подсти-.
5. Модель распространения вредных примесей от автотранспорта. У равнение транспорта вредной примеси представляет собой нестационароное трехмерное уравнение с параметризуемыми коэффициентами турбулентного обмена и постоянной деструкции:
с;+исх+мсу+(* - * )с;+ас = пс)+(Пксу) у+(пс;) + /с, (47)
где с(X, у, г, г) - концентрация вредной примеси, и , V, * - компоненты вектора скорости движения воздушной среды, - скорость оседания загрязняющего ве-
щества под действием силы тяжести, 7 - член, ответственный за деструкцию вещества, /с (X, у, г, г) - функция ИСТОЧНИКОВ загрязняющих веществ, Пн И П - коэффициенты соответственно горизонтального и вертикального турбулентного обмена.
Задача распространения вредной примеси формулируется на базе универсальной формы записи уравнений переноса при Ф = с .
[(ір + Бр] ср = о^с^ + аЕсЕ + я5с5 + а^с^ + овсв + отст + Би, (48)
р образует контрольный объем с центром в точке Р.
Область численного интегрирования ограничена подстилающей поверхностью, входной и выходной границами, боковыми плоскостями, а также верхней границей. Граничные условия для входной и выходной границе определяются в виде
эс, .
0. (49)
На боковых границах:
На верхней границе:
дп
дс
дп
дс . _____і
дп
= 0.
п - .
На подстилающей поверхности в общем случае граничное условие имеет вид
-С
к -С- -%С. = 0, (50)
дп
где С; - концентрация I -го загрязняющего вещества, ¡3- - коэффициент, зависящий от свойств подстилающей поверхности.
Для различных подстилающих поверхностей коэффициент /3- будет отличаться. При моделировании распространения примеси следует иметь в виду, что возможны случаи взаимодействия с почвой и водной поверхностью. В случае взаимодействия с почвой необходимо отметить, что примеси слабо взаимодейст-
0
вуют с ней. Загрязняющее вещество, попав на поверхность почвы, обычно не на, .
В связи с этим граничное условие на подстилающей поверхности принимает вид:
к — = 0. (51)
дп I
Водная поверхность поглощает примесь, поэтому концентрацию примеси у поверхности воды можно принять равной нулю. Граничное условие на поверхности воды имеет вид
С = 0 . (52)
Граничные условия для точечного источника применяются к контрольному объему, в котором этот точечный источник расположен. Необходимо также учитывать, что газовоздушная смесь обладает скоростью —, и в соответствии с этим граничное условие для вертикальной составляющей скорости принимает вид
V (( Уo, ¿0 ) = V*. (53)
Примесь обладает температурой Т , поэтому граничное условие для темпе-
ратуры определяется в виде
Т ( Уo, ¿0 )= Т . (54)
Для определения концентрации примеси в устье источника по величине мощности выброса можно воспользоваться соотношением:
\ М
С (Хр Уo, ¿0 ) = V"’ (55)
где М - мощность выброса в кг / с , а V1 - величина расхода газовоздушной
V 3 / ~
смеси В I / П .
Величины VI и V* связаны очевидной формулой
VI = Л- ,
1 *7
где 5 - площадь устья источника.
6. Моделирование распространения примесей от линейного источника.
Рассмотрим область численного интегрирования задачи распространения примеси. Основным источником загрязнения воздушной среды города является автотранспорт. Вклад транспортных средств в общую картину загрязнения удобно описывать с помощью линейных источников загрязнения. При расчете загрязнения от линейного источника он представляется в виде совокупности точечных источни-.
, :
Лт 5ь4й
N =---------, (56)
х
где X - наименьшее расстояние в метрах от источника до расчетной точки на местности; и - расчетная скорость ветра в м/с.
В соответствии с формулой (56) количество точечных источников увеличивается с увеличением протяженности линейного источника. Для проведения расчета необходимо, чтобы точечные источники находились в центре контрольных .
В качестве граничных условий для линейных источников выбираются граничные условия для совокупности точечных источников.
Тогда функцию подвижного источника примеси можно представить в виде
N
/ = Ё Еп3{7 - ГВ (57)
1=1
где Е = Е (Г, I) Г1 - количество выбросов примеси от источников в узле дискретной сетки Г = ( , , XI ) в момент времени t .
Примесь обладает температурой Т, поэтому необходимо определять граничное условие для температуры.
Для определения концентраций в устьях точечных источников, в случае аэра-
,
вещества М (кг/с) принимается равной суммарной мощности источника:
М
Мы = —. (58)
N
Тем самым аэрационный фонарь эквивалентен при расчете группе одинако-.
В случае проведения расчета распространения концентраций загрязняющих веществ от автомобильных дорог необходимо учитывать условия движения автомобилей в транспортном потоке. В общем случае величина выбросов автомобилей I -го загрязняющего вещества М{ на участке улицы длиной I за единицу времени может быть определена по формуле
М1 = Ми + Д., (59)
где М. - выброс I -го загрязняющего вещества в кг/с при непрерывном движении транспортного потока;
Д - дополнительный выброс I -го загрязняющего вещества в кг/с, связанный с задержкой транспортных средств.
Величина Мц характеризует мощность выброса при непрерывном движении автотранспорта на участке автомагистрали. Она определяется видами транспорт, , , интенсивностью движения. Величина Д связана с условиям и движения по автомагистрали. Она отражает мощность выброса загрязняющих веществ, который
складывается из выбросов при торможении и разгоне транспортных средств,
а также работе двигателя на холостом ходу.
Граничные условия для концентрации С. для совокупности точечных источ-
, ,
С
Сі (*0, Уо, ^о) = N ,
(60)
где С. - концентрация г -го загрязняющего вещества для линейного источника.
7. Модель турбулентности приземного слоя атмосферы для несжимаемого случая. Запишем осредненную по Рейнольдсу базовую трехмерную систему уравнений для несжимаемого вязкого газа в тензорном виде
д(риі) = 0 Эх ’
(61)
д(ри)
дї
д /______ ~гт\ др дТ]
+----(рии; + рии ; I =------+ —-
дхЛ 1 ' дх дх;
• + ■
дх дх;
(62)
(63)
д(рфг) д / “ТГП д
-------+—(щ ф (+р«1(Ф, )=дх-
дї
V
,дФ
дХ;
і У
(64)
Здесь величины с чертой означают осредненные по необходимому масштабу времени переменные исходной системы уравнений, а осредненное по времени
произведение флуктационных компонент скорости рии ■ нуждается в дополни.
Будем исходить из предположения, что явление турбулентности в модели описывается дополнительной турбулентной вязкостью. Тогда выражение рии'■
можно записать в виде
ди, ди;
■ + ■
дх дх
-- Рк. з 1
(65)
В случае, если речь идет об уравнении пе реноса для скалярной величины Ф, соответствующий член ри 'Ф'. представим в виде
ри,Ф1 =Г, £
(66)
(65) -
ставимо в виде
7 1л'' 1/'',
к = — ри( и; = —\иих + и уиУ + ии„ I.
г\ • і 1 г\ у * * уу гг}
(67)
Выражения и Гt - соответствующие коэффициенты турбулентной вязкости. В простейшем случае явление турбулентности можно описывать с помощью двух характеристик: кинетической энергии к и масштабам длиной L :
(68)
42к, си - безразмерная
где д = V 2К , - Оезразмерная константа.
Для достаточно простых течений значение величины д можно получить по известному полю распределения скорости
ди
д = Ь
ду
(69)
Данная формула достаточно точна для течений с мало выраженным трехмерным характером. Определение L для сложных трехмерных течений является сложной задачей, поэтому возникает необходимость к переходу на более сложные .
8. К - £ -модель турбулентности. К - £ -модель является логическим продолжением модели на основе двух характеристик. Здесь используются дифференциальные уравнения для определения кинетической энергии и энергии диссипации.
Уравнение для кинетической энергии имеет вид
д( рк ) +д(ри,к ) = д
дґ
дх;
дх;
і V
дк
дх:
і /
дх;
р
/ / / , / /
и иі иі + Р и :
(70)
/ / диі ди' ди'
-риі и; —L - и—і------------------.
1 дхі дхк дхк
Второй член в правой части уравнения описывает турбулентную диффузию . -диента кинетической энергии:
Р ' ' ' . •
—и и и + р и
2 1 і і И
/ /
дк_ ак дх} ’
(71)
где СГК - турбулентное число Прандтля.
Третий член в правой части уравнения описывает уровень кинетической , .
в виде
Рк =-рщ и1 ^ «м
' диі
дх;
Гди: ди 1
—'- +—1
, дх; дх
V 1 1 J
диі
дх;
(72)
Для определения длины масштаба турбулентности требуется другое уравнение. Наиболее подходящим соотношением между энергией диссипации и кинетической энергией турбулентности является:
є = •
к2 Ь
Уравнение для диссипации имеет вид
3
зрі+зір1^ ср Є — рСє2є-+4-
дґ дх к к дх.
ґ М дє є
ає дх;
є 1 У
(74)
Здесь и = рСм4кЬ = рСм — .
Данная модель включает в себя пять параметров. Для задачи моделирования поля ветровых течений в условиях городского рельефа местности они принимают вид
С = 0,09, Г, = 1,44, С2 = 1,92, ок = 1,0, о; = 1,3.
9. Комплекс программ и результаты моделирования. Для опис анных математических моделей был построен программный комплекс расчета приземной аэродинамики в условиях городской застройки. Структурная схема комплекса . 3.
Рис. З.Структурная схема организации программного комплекса
В рамках данной организации выделяются следующие элементы:
♦ ГИС - геоинформационная система (набор карт города высокого разрешения в электронном виде);
♦ Препроцессор - блок обработки входных данных для использования в
;
♦ Процессор - основная программа, рассчитывающая поле ветровых течений в приземном слое атмосферы;
♦ Постпроцессор - программа обработки выходных данных процессора и визуализации скалярных и векторных полей.
В программном комплексе ГИС играет роль набора необработанных карт, которые описывают геометрию области, тип и свойства объектов области моделирования и т.д.
Основной задачей препроцессора является подготовка входных данных, полученных из ГИС, для удобного использования процессором. После получения данных о моделируемом объекте из ГИС препроцессор формирует карты типов и расположения контрольных объемов для всех физических величин, рассчитываемых про.
для всех компонент поля скорости ветра (пустые, входные, выходные и т.д.).
- , -та. Из этой процедуры последовательно вызываются другие процедуры, которые рассчитывают необходимые физические величины. В простейшем случае вычисляются компоненты скорости и поле давления.
Постпроцессор осуществляет визуализацию результатов моделирования. В нем предусмотрена возможность отображения скалярных полей в виде изолиний
и векторных полей в виде стрелок, в каждой точке сетки, указывающих направления ветрового потока.
На рис. 4 и 5 представлены результаты моделирования при следующих условиях: интенсивность подвижных источников (автотранспортного потока) -1400 единиц/ч, загрязняющая примесь (однокомпонентный газ СО) - скорость движения - 50 км/ч, мощность выбросов - 5,5 мг/(с/км), ветер восточный - 5 м/с, температура - 20° С, давление - 760 мм. р. ст.
Рис. 4.Векторное тле скорости движения среды
Рис. 5.Поле загрязняющих веществ
07125179
В работе предложен алгоритм моделирования поля ветровых течений в условиях городской застройки. Численная реализация основана на универсальном дис.
Навье-Стокса используется SIMPLE-метод. Предложено описание неявного алгоритма Стоуна для решения СЛАУ специального вида, получающихся в ходе аппроксимации сеточных задач. Представлены результаты моделирования поля ветровых течений и распространения загрязняющих веществ от транспортных магистралей города Таганрога.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Берлянд М. Е. Прогноз и регулирование загрязнений атмосферы. - Л.: Гидрометеоиз-дат, 1985. - 271 с.
2. Марчук Г. И. Математическое моделирование в проблеме окружающей среды. - М.: Наука, 1982. - 319 с.
3. Матвеев Л.Т. Физика атмосферы. - СПб.: Гидрометеоиздат, 2000. - 779 с.
4. Патанкар С.В. Численные методы решения задач обмена и динамики жидкости. - М.: Энергоатомиздат, 1984. -152 с.
5. Ferziger J., Peric M. Computational Methods for Fluid Dynamics - 3., rev. ed. - Berlin; Heidelberg; New York; Barcelona; Hong Kong; London; Milan; Paris; Tokyo: Springer, 2002. - 423 p.
Сухинов Александр Иванович
Технологический институт федерального государственного образовательного учреждения высшего профессионального образования «Южный федеральный университет» . .
E-mail: sukhinov@gmail.com.
347928, г. Таганрог, пер. Некрасовский, 44.
Тел.: 88634310599.
Гадельшин Валерий Камельянович
E-mail: gadelshin@mail.ru.
Тел.: 88634601461.
Любомищенко Денис Сергеевич
E-mail: dexusint@gmail.com.
Тел.: 88634601461.
Sukhinov Alexander Ivanovich
Taganrog Institute of Technology - Federal State-Owned Educational Establishment of Higher Vocational Education “Southern Federal University”.
E-mail: sukhinov@gmail.com.
44, Nekrasovskiy, Taganrog, 347928, Russia.
Phone: +78634310599.
Gadelshin Valeriy Kamelianovich
E-mail: gadelshin@mail.ru.
Phone: +78634601461.
Lubomishenko Denis Sergeevich
E-mail: dexusint@gmail.com.
Phone: +78634601461.
