Оценка влияния лесозащитных полос на перенос диоксида серы (SO2) в приземном слое атмосферы с помощью гидродинамической модели
М. С. Рыжова1, Ю.В. Мухартова1^, Н. Е. Шапкина1, А. В. Ольчев26
1 Московский государственный университет имени М. В. Ломоносова, физический факультет, кафедра математики. Россия, 119991, Москва, Ленинские горы, д. 1, стр. 2.
2 Институт проблем экологии и эволюции имени А.Н. Северцова РАН.
Россия, Москва, Ленинский просп., д. 33.
E-mail: a muhartova@yandex.ru, 6 aoltche@yandex.ru Статья поступила 16.06.2015, подписана в печать 04.08.2015.
Численная двумерная гидродинамическая модель применена для оценки влияния лесозащитных полос различных размеров на турбулентный перенос диоксида серы (SO2) в приземном слое атмосферы. Результаты расчетов показали, что защитная лесополоса приводит к существенному ослаблению горизонтального переноса SO2 за счет уменьшения скорости воздушного потока и поглощения примеси кронами деревьев. Коэффициент ослабления потока примеси растет с увеличением размеров защитной лесополосы и уменьшением высоты источника загрязнения.
Ключевые слова: двумерная гидродинамическая модель, турбулентный перенос, атмосферное загрязнение, сухое осаждение SO2, полуторное замыкание системы усредненных уравнений гидродинамики, лесозащитные полосы.
УДК: 519.63. PACS: 02.30.Jr, 92.60.Fm.
Введение
Задаче изучения процесса переноса атмосферных примесей в приземном слое атмосферы в настоящее время посвящено огромное число экспериментальных и теоретических исследований [1-4]. Это обусловлено прежде всего усилившейся в последние годы антропогенной нагрузкой на окружающую среду и необходимостью разработки и проведения комплексных мероприятий по снижению негативного влияния антропогенного загрязнения воздуха на природные экосистемы и улучшению качества и условий жизни населения.
Априори хорошо известно, что растительность, и прежде всего леса, активно поглощает аэрозоли и газовые примеси из атмосферного воздуха, снижая уровень атмосферного загрязнения. Количественные оценки подобного влияния пока очень немногочисленны, что связано как с отсутствием необходимого оборудования для проведения мониторинговых измерений концентрации газовых примесей в воздухе, так и с отсутствием репрезентативных подходов для адекватного количественного (модельного) описания процесса переноса газовых примесей от источников загрязнения в приземном слое воздуха в условиях пространственно-неоднородного антропогенного ландшафта.
В рамках настоящей работы на основе двумерной гидродинамической модели были проведены расчеты переноса воздушным потоком пассивной загрязняющей газовой примеси (диоксида серы Э02) от двух видов антропогенных источников загрязнения (автомобильной дороги и заводской трубы) и получены
оценки влияния на горизонтальные потоки защитной лесополосы различных размеров. Диоксид серы — достаточно вредное вещество, особенно опасное для людей с болезнями органов дыхания, например, астмой. Его выбрасывают в атмосферу топливосжигаю-щие предприятия и заводы, работающие на мазуте и угле (например, ТЭЦ или ТЭС). Также выделение Э02 характерно для двигателей внутреннего сгорания автомобилей.
1. Двумерная модель турбулентного движения воздуха в приземном слое атмосферы при наличии растительности
Для решения задачи переноса пассивных газовых примесей в приземном слое воздуха в случае нейтральной стратификации в работе была использована двумерная модель [16], основанная на усредненных уравнениях гидродинамики. Существенным преимуществом данной модели является ее способность описывать процессы переноса над горизонтально неоднородной поверхностью, и в том числе на границе растительных сообществ, рассматриваемых как сплошная проницаемая среда, оказывающая сопротивление воздушному потоку.
Будем рассматривать задачу в некоторой прямоугольной области с размерами х е [-Ь, Ь], г е [Н0, Н]. Система уравнений для усредненной по пространству и времени скорости ветра V(х, г, 0 = {и(х, г, 0, №(х, г, Щ, где и и № — горизонтальная и вертикальная компоненты скорости ветра соответственно и 5Р — усредненного отклонения давления от гидростатического, имеет
31 ВМУ. Физика. Астрономия. № 6
вид
dU + и —U + W —U _
dt дх dz
__ д_ ( p + л +o_
dx \ p0 ) dx
>«dUi)+
+ Fu
+
(1)
d fК dU\ d (RdW
dz \ dz J dz у dx dW UdW W dW _ t x z
__d(P + E) + d (KdW
dz \ p0 J dx \ dx
+ — (2KW) + — (KU
dz у dz J dx у dz dU dW _ о x z
где p0 — плотность воздуха, E — кинетическая энергия турбулентного движения, К — коэффициент турбулентной диффузии, величины Fu, Fw описывают изменение скорости воздушного потока при взаимодействии с растительностью:
F _ {Fu, Fw } _ -cd LAD \V \V.
Здесь LAD — суммарная односторонняя площадь поверхности листьев деревьев в единице объема, cd — коэффициент аэродинамического сопротивления элементов растительности (в нашем исследовании он был принят равным 0.2).
Для замыкания системы уравнений (1) будем использовать следующее выражение для коэффициента турбулентной диффузии: К _ CpE2е-1, где Cр _ 0.09 — безразмерный коэффициент пропорциональности [7], е — скорость диссипации турбулентной кинетической энергии E. Система уравнений для E и ф _ еЕимеет вид [7-9, 16]:
увеличение диссипации турбулентной кинетической энергии за счет взаимодействия с растительностью [7]:
Дф _ 12С1/2(Сф2 - CvX)cd LAD \V ф.
Боковые границы x _ ±L расчетной области являются свободными. Обычно на свободных границах используют так называемые условия сноса [17], которые заключаются в том, что производная искомой функции по нормали к свободной границе полагается равной нулю. В нашем случае эти условия будут иметь вид
dU dW dE дф
+ FW, x x_±L x x_±L x x_±L dx
x_±L
dA + uddE + w—e _
t x z
x
d i К dE\ d f \оф dx ) dz \
d (к dE\
^^ + ^ ^^ + Pe - е,
оф dz )
(2)
+
dx ^P
_ 0. (3)
x_±L
dW dE дф
z z=H z=H z „ dz z=H
_ 0. (4)
z_H
На верхней границе г = Н также будем использовать условия сноса для всех функций, кроме избыточного давления:
ди_
дг
При условии, что верхняя граница расчетной области находится достаточно далеко от расположенных на поверхности земли элементов шероховатости, можно сделать предположение, что избыточное давление 5Р при г = Н равно нулю.
На нижней границе области для горизонтальной компоненты скорости ветра воспользуемся условием [18, 19]:
k—U
z
KClJ4 Е1/2
г_
u
(5)
z_h0
г=Н0 1п((г " а)/г0)
где к = 0.4 — эмпирическая постоянная фон Кармана, го — параметр шероховатости, й — высота слоя вытеснения. Для остальных функций используем следующие условия:
^ + и ^ + Г ^ = д1 дх дг
= д ( К дф \ д ( К дф \
дх \аф дх ) дг \аф дг )
+ ф (Сф1РЕ - Сф2£} + Дф,
где а'Е = 2 и аф = 2 — числа Прандтля для Е и ф соответственно [10], РЕ — сдвиговая генерация турбулентной кинетической энергии. Выражение для РЕ в двумерном случае записывается следующим образом:
П п„(( ди )2 (дж )2\ / ди дж )2
рЕ=2К(Ы )+КЫ+ Ж)
Множители Сф1 = 0.52 и Сф2 = 0.8 в правой части уравнения для функции ф являются модельными константами [10], а слагаемое Дф характеризует
ф\ z_h0 _
W\z_h0 _ 0,
C..E
дЕ
z
_ 0,
z_h0
_ ^^
К
z_h0
dz ^P
(6)
_ 0.
z_h0
Будем решать систему уравнений (1)-(2) с граничными условиями (3)-(6) как задачу на установление, используя в качестве начального условия для и логарифмическое распределение скорости ветра по высоте [5, 6], справедливое в случае однородной подстилающей поверхности. Также предположим, что Ж и 5Р равны нулю в начальный момент времени, а при I = 0 зависимость функций Е и К от г описывается в соответствии с полуэмпирическими формулами, справедливыми для однородной подстилающей поверхности [6].
2. Перенос SO2 воздушным потоком и его поглощение растительностью
После того, как в ходе решения начально-краевой задачи (1), (2), (3)-(6) получено установившееся поле скорости ветра V = {и, Ж} и распределение
коэффициента турбулентной диффузии К, решается начально-краевая задача для усредненной концентрации С переносимого воздушным потоком вещества. Уравнение переноса для С имеет вид
—+U —+Г — - — (— + — (— dt dx dz dx \ Sc dx) dz \ Sc dz J
+Fc ,
где Бе — турбулентное число Шмидта [20, 21], которое в нашей работе принято равным 0.75, а функция РС описывает источники и стоки переносимого вещества:
РС = ^оигее + ^¡пк.
В качестве начального условия в задаче для С возьмем некоторое фоновое значение концентрации Б02 в атмосфере: С|^=о = Со = 6 мкг • м _3.
В качестве граничных условий для С используем условия сноса на свободных границах и условие равенства нулю потока С на нижней границе области:
дС дх
x-± L
f dz
0.
z—h0, z-H
Лабораторные и полевые исследования показывают, что скорость поглощения Б02 растениями (сухое осаждение) определяется интенсивностью турбулентной и молекулярной диффузии в слое воздуха внутри растительности, а также устьичной проводимостью, характеризующей пропускную способность устьиц, расположенных на поверхности листьев растений, для атмосферных газов [22-24]. В случае наличия влаги на поверхности листьев (например, роса или задержанные листвой капли дождя) атмосферный Б02 может вступать в химическую реакцию с водой с образованием сернистой кислоты Н2Б0з (влажное осаждение). Также в виде Н2Б0з он может попадать на растения и с атмосферными осадками.
В рамках нашего исследования будем рассматривать только сухое осаждение, скорость которого будет определяться исключительно устьичной проводимостью, и в наиболее простом случае может быть представлена как Уа = С, где Q — поток поглощаемого Б02, С — его концентрация. Опубликованные в литературе значения скорости сухого осаждения Уа для Б02 в зависимости от типа и видового состава растительности изменяются в пределах от 0.1 до 2 см • с _1 [22-24].
Таким образом, растительность является стоком Б02 из атмосферы и скорость поглощения Б02 растительностью может быть представлена как
-5;пк(х, 2, ^ = ЬЛО(^(х, 2, ^ = Уа •ЬДБ^К (х, 2, Г).
В рамках нашего исследования в качестве источника Б02 будем рассматривать автомобильную дорогу шириной /5оигее, расположенную перпендикулярно плоскости 0x2, и трубу небольшого промышленного предприятия высотой ^5оигее = 20 м.
При рассмотрении переноса Б02 воздушным потоком от автомобильной дороги в двумерной задаче
источником загрязнения служит поперечное сечение дороги, которое моделируется в виде прямоугольника шириной /5оигее и высотой 25оигее, расположенного у нижней границы расчетной области. При моделировании процесса переноса Б02 от трубы предприятия также будем использовать эквивалентный источник, расположенный на высоте 20 м над поверхностью земли.
3. Результаты модельных расчетов
При расчетах переноса Б02 от автомобильной дороги было сделано предположение, что она имеет ширину /5оигее = 10 м, и мощность источника Б02 составляет 329.76 мкг• с_1. На первом шаге были проведены расчеты для случая открытой ровной поверхности земли без каких-либо препятствий (рис. 1).
Z, м
50 45 40 35 30 25 20 15 10 5 0
-300 -200 -100
0
х, м
100
200
300
Рис. 1. Распределение концентрации SO2 [мкг• м-3] вблизи автотрассы в направлении движения воздушного потока при условии отсутствия вокруг автомобильной дороги каких-либо препятствий. Предполагается, что ветер дует слева направо. Скорость ветра на высоте 20 м принята равной 3.5 м • с-1
На втором шаге в качестве возможного препятствия воздушному потоку были рассмотрены защитные лесополосы различной ширины. Было сделано предположение, что исследуемая лесополоса состоит преимущественно из мелколиственных пород деревьев (береза, осина, ольха) высотой около 20 м. Предполагается, что край лесополосы расположен в 20 м от автомобильной дороги. Индекс листовой поверхности растительности в лесополосе LAI — 5. Средняя скорость сухого осаждения SO2 в кронах деревьев Vd — 1.0 см • с-1, что соответствует малооблачным погодным условиям с температурой воздуха в диапазоне 20-25 ° C при достаточных условиях почвенного увлажнения. При динамической скорости u* — 0.4 м • c-1 и выбранных параметрах лесополосы средняя скорость ветра непосредственно над кронами деревьев составляет 2 м • c-1.
Рис. 2. Распределение концентрации SO 2 [мкг-м 3] вблизи автотрассы в направлении движения воздушного потока при наличии защитной лесополосы шириной 150 м (а) и 330 м (б). Ветер дует слева направо
Результаты расчетов распределения концентрации БОг в приземном слое воздуха при наличии двух лесополос шириной 150 и 330 м (рис. 2) показывают, что наличие лесополосы существенно меняет пространственную картину распределения концентрации БОг (при сравнении с открытой безлесной поверхностью). Это происходит как за счет прямого поглощения БОг растительностью (сухого осаждения), так и благодаря изменению скорости и направления воздушного потока при его взаимодействии с элементами растительности в лесополосе.
На рис. 3 приведена зависимость плотности горизонтального потока БОг (Щс) от высоты при разной ширине лесополосы на расстоянии 510 м от источника загрязнения. При наличии установившегося воздушного потока с усредненной скоростью и
выражение для Wx может быть записано как
Рис. 3. Зависимость горизонтального потока БО? от ширины лесополосы на расстоянии 510 м от источника загрязнения (автотрасса)
г)Г
^ = 1С + = и (С -Со)-Кс где =
В табл. 1 приведены средние значения горизонтальных потоков на расстоянии 510 м от источника:
о
Таблица 1
Зависимость интегрального потока 80г от ширины защитной лесополосы для случая, когда в качестве источника загрязнения была рассмотрена автомобильная дорога
Ширина лесополосы, м Усредненный поток, мкг-м_2-с_1 Ослабление потока, %
0 10.86 —
30 6.85 37
70 4.07 63
150 2.96 73
230 2.65 76
330 2.12 80
380 1.83 83
Расчеты показывают, что наличие лесополосы даже минимальной ширины позволяет снизить значение горизонтального потока БОг на расстоянии 510 м от источника почти на 40% с 11 до 7 мкг-м_2-с_1. Для лесополосы шириной 150 м и более снижение потока БОг вблизи поверхности достигает 73% и более (рис. 3 и табл. 1).
Для оценки влияния высоты источника загрязнения на скорость переноса загрязняющей примеси в работе в качестве источника БОг была рассмотрена труба небольшого промышленного предприятия.
При расчетах было сделано предположение, что источник загрязнения расположен на высоте 20 м над земной поверхностью и имеет одинаковую с автомобильной дорогой мощность.
Результаты расчетов показывают, что расположение источника БОг над поверхностью земли приводит к существенному увеличению площади загрязнения и скорости распространения загрязняющей примеси. Это связано, прежде всего, с более высокими скоростями ветра вблизи источника загрязнения и огибанием воздушным потоком с примесью препятствия по его верхней границе без существенного ослабления (рис. 4, 5).
х, м
-3
Рис. 4. Распределение концентрации БСЬ [мкг-м вблизи заводской трубы в направлении движения воздушного потока при отсутствии препятствий. Ветер дует слева направо
Анализ зависимости горизонтальных потоков БОг от ширины защитной лесополосы показывает, что в случае приподнятого положения источника
загрязнения величина потока БОг существенно слабее зависит от ширины лесополосы, чем в случае расположения источника загрязнения у поверхности земли (рис. 6, табл. 2). В частности, вблизи поверхности земли коэффициент ослабления усредненного потока БОг для лесополосы шириной 380 м составляет 52%, а для лесополосы шириной 30 м — всего 18% (табл. 2).
Рис. 6. Зависимость горизонтального потока БСЬ от ширины лесополосы на расстоянии 510 м от источника загрязнения (заводская труба)
Заключение
На основе разработанной двумерной гидродинамической модели получены данные о распределении концентрации БОг вокруг источников загрязнения, а также значения горизонтальных потоков загрязняющей примеси при отсутствии и наличии лесозащитных полос. Было показано, что лесополоса существенно снижает горизонтальный поток БОг как за счет того, что, являясь механической преградой,
х, м
х, м
Рис. 5. Распределение концентрации SO? [мкг-м 3] вблизи заводской трубы в направлении движения воздушного потока при наличии лесополосы шириной 150 м (а) и 330 м (б). Ветер дует слева направо
Таблица 2
Зависимость интегрального потока SO2 от ширины защитной лесополосы для случая, когда в качестве источника загрязнения была рассмотрена заводская труба
Ширина лесополосы, м Усредненный — 2 —1 поток, мкг-м -с Ослабление потока, %
0 41.58 —
30 34.11 18
70 30.27 27
150 27.14 35
230 25.48 39
330 21.88 47
380 20.01 52
она перераспределяет воздушные потоки, так и за счет прямого поглощения SO2 кронами деревьев.
В ходе расчетов было установлено, что защитная лесополоса с минимальной шириной 30 м позволяет снизить значение горизонтального потока SO2 на расстоянии 510 м от источника загрязнения, расположенного на поверхности, почти на 40% и в случае его расположения над поверхностью земли на высоте крон деревьев — лишь на 18%.
Результаты расчетов показали, что разработанная математическая модель может служить эффективным инструментом, позволяющим находить пространственное распределение концентрации примеси в приземном слое воздуха без привлечения измерительного оборудования. Она может также быть пригодным инструментом для разработки защитных мероприятий (в том числе проектирования лесозащитных полос) и строительства сооружений, снижающих потоки загрязняющих веществ от источников загрязнения в окружающем воздухе.
Исследование выполнено в рамках проекта Российского научного фонда (грант 14-14-00956).
Список литературы
1. Алоян А.Е., Йорданов Д.Л., Пененко В.В. // Метеорология и гидрология. 1981. № 8. С. 32.
2. Austin J., Brimblecombe P., Sturges W. // Elsevier Science. 2002. P. 569.
3. Kampa M., Castanas E. // Environ. Pollut. 2008. 151. P. 362.
4. Anderson J.O., Thundiyil J.G., Stalbach A. // J. Med. Toxicol. 2012. 8. P. 166.
5. Wyngaard J.C. Turbulence in the Atmosphere. Cambridge, 2010.
6. Garratt J.R. The Atmospheric Boundary Layer. Cambridge, 1992.
7. Sogachev A., Panferov O. // Boundary-Layer Meteorol. 2006. 121, N 2. P. 229.
8. Sogachev A. // Boundary-Layer Meteorol. 2009. 130, N 3. P. 423.
9. Sogachev A., Menzhulin G.V., Heimann M., Lloyd J. // Tellus. 2002. 54B, N 5. P. 784.
10. Wilcox D.C. Turbulence modeling for CFD. DCW Industries, Inc. La Canada, CA, 1998.
11. Белоцерковский О.М. Численное моделирование в механике сплошных сред. М., 1994.
12. Дубов А.С., Быкова Л.П., Марунич С.В. Турбулентность в растительном покрове. Л., 1978.
13. Бояршинов М.Г., Горемыкин В.Д. // Матем. модел. 2004. 16, № 7. С. 31.
14. Ольчев А.В., Радлер К. // Изв. Самарского научн. центра РАН. 2009. 11, № 1(7). С. 1538.
15. Olchev A., Radler K., Sogachev A. et al. // Ecological Modelling. 2009. 220. P. 3046.
16. Мухартова Ю.В., Левашова Н.Т., Ольчев А.В., Шап-кина Н.Е. // Вестн. Моск. ун-та. Физ. Астрон. 2015. № 1. С. 15. (Mukhartova Yu.V., Levashova N.T., Ol-tchev A.V., Shapkina N.E. // Moscow University Phys. Bull. 2015. 70. N 1. P. 14.)
17. Елизарова Т.Г. Математические модели и численные методы в динамике газа и жидкости. М., 2005.
18. Krzikalla F. Numerical Investigation of the Interaction between Wind and Forest under Heterogeneous Conditions: Diploma Thesis. Institute for Hydromechanics, University of Karlsruhe, 2005.
19. Foudhil H., Brunet Y., Caltagirone J.P. // Environ. Fluid Mech. 2005. 5. P. 247.
20. Stull R.B. An Introduction to Boundary Layer Meteorology. Dordrecht, 1988.
21. Flesch T.K., Prueger J.H., Hatfield J.L. // Agricult. and Forest Meteorol. 2002. 111. P. 299.
22. Fowler D., Unsworth M.H. // Quart. J. R. Met. Soc. 1979. 105. P. 767.
23. Gravenhorst G., Oltchev A., Sogachev A. et al. // Meteorologische Zeitschrift. 2005. 14, N 2. P. 117.
24. Erisman J.W., Baldocchi D. // Tellus. 1994. 46B. P. 159.
An impact assessment of forest belts on the SO2 transport within the atmospheric boundary layer using a hydrodynamic model
M.S. Ryzhova1, Yu.V. Mukhartovala, N.E. Shapkina1, A.V. Oltchev2b
1 Department of Mathematics, Faculty of Physics, Lomonosov Moscow State University, Moscow 119991, Russia.
2 Severtsov Institute of Ecology and Evolution, Russian Academy of Sciences, Moscow, Russia. E-mail: a muhartova@yandex.ru, b aoltche@yandex.ru.
A numerical two-dimensional hydrodynamic model was used to describe the influence of forest belts of different sizes on turbulent transport of SO2 within the atmospheric surface layer. The results of the model calculations showed that the presence of a forest belt results in an substantial reduction of the horizontal
SO2 flux due to the decrease of the wind speed and the absorption of SO2 by tree crowns. The extinction coefficient of SO2 flux increases with an increase in the forest belt size and decrease with the pollution source height.
Keywords: two-dimensional hydrodynamic model, turbulent transport, atmospheric pollution, dry deposition of SO2, one-and-a-half closure of averaged hydrodynamic equations, forest belts. PACS: 02.30.Jr, 92.60.Fm. Received 16 June 2015.
English version: Moscow University Physics Bulletin 6(2015). Сведения об авторах
1. Рыжова Мария Сергеевна — аспирант; тел.: (495) 939-10-33, e-mail: maria.msu.ff@gmail.com.
2. Мухартова Юлия Вячеславовна — канд. физ.-мат. наук, ст. науч. сотрудник; тел.: (495) 939-10-33, e-mail: muhartova@yandex.ru.
3. Шапкина Наталья Евгеньевна — канд. физ.-мат. наук, доцент; тел.: (495) 939-10-33, e-mail: neshapkina@mail.ru.
4. Ольчев Александр Валентинович — канд. геогр. наук, ст. науч. сотрудник; e-mail: aoltche@yandex.ru.