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

114

Вестник СамГУ — Естественнонаучная серия. 2006. №6/1(46).

УДК 001.891.57

О МАТЕМАТИЧЕСКОЙ МОДЕЛИ МИГРАЦИИ И ОСАЖДЕНИЯ ПОЛИДИСПЕРСНОГО АЭРОЗОЛЬНОГО ОБРАЗОВАНИЯ1

© 2006 М.В.Меньшов2

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

Введение

Математические модели, описывающие осаждение аэрозольных образований при обработке сельскохозяйственных полей, разнообразны. К наиболее значимым отличиям стоит отнести:

• полноту описания протекающих при этом физических процессов

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

Практика применения сельскохозяйственной авиации в настоящее время эффективна для полетов на предельно малых высотах над достаточно однородной поверхностью, при малых скоростях ветра и, как правило, при наличии нейтральной стратификации пограничного слоя атмосферы (ПСА). При этом осуществляется распыление жидких рецептур с диаметром капель 80-600 мкм с нормой плотности выпадений 300 г/га и более.

Имеющие место несистематические случаи применения метода уль-трамалообъемного (УМО) распыления рецептур при нормах плотности

1 Представлена доктором технических наук, профессором И.С. Загузовым.

2Меньшов Максим Владимирович (ssaa@samara.ru), кафедра высшей математики Самарской государственной сельскохозяйственной академии, 44642, Россия, Самарская обл., пос. Усть-Кинельский, ул.Учебная 2.

10-300 г/га с диаметром капель до 30 мкм математически не проработаны по условиям эффективного ввода аэрозольных образований. Это зачастую становится причиной негативных последствий в виде непрогнозируемых выносов частиц аэрозоля далеко за пределы обрабатываемых территорий, неравномерности обработки и т.д.

В настоящей работе приводятся результаты верификации модифицированной математической модели [1] движения и осаждения полидисперсного аэрозольного облака, образованного при распылении самолетом сельскохозяйственной авиации жидких рецептур, путем поточечного сравнения с результатами полевых экспериментов, проведенных в 1986-1991 годах в условиях равнинной и среднепересеченной местности.

Известны работы К.П. Куценогого [2], В.Ф.Дунского и других ученых [3] по подобной верификации моделей осаждения аэрозольных образований, базируюшихся на математическом моделировании турбулентной струи Гаусса [4]. В этих работах однозначно выявлена необходимость более тщательного описания состояния приземного слоя атмосферы как определяющего фактора особенностей осаждения частиц аэрозоля на подстилающую поверхность местности. Опыт проведения полевых экспериментов подсказывает, что положительные результаты верификации модели в [2] получены только для условий достаточной устойчивости динамических процессов в приземном слое атмосферы.

Следует особо отметить, что измерение метеопараметров при проведении полевых экспериментов, используемых для верификации, как в работах [2, 3], так и в настоящей работе, проводилось дискретно по времени и по положению на экспериментальном поле. Это неизбежно приходит в несоответствие с продолжительностью миграции аэрозольного образования (30 мин и более) и большими линейными размерами экспериментального поля (более 4 км), ввиду возможной изменчивости метеопараметров, регистрируемых в конкретных точках, за непродолжительное время.

1. Описание математической модели

Вводится декартова прямоугольная система координат (х, у, г), в которой ось г направлена вверх, ось х — на восток, а ось у — на север. Область решения определяется в виде параллелепипеда с неровной нижней границей, отражающей неоднородность рельефа подстилающей поверхности

0 ^ х ^ Ьх, 0 ^ у ^ Ьу, 6 ^ г ^ Н, (1)

где Ьх, Ьу — размеры области по горизонтали, функция г = 6(х, у) задает рельеф местности, Н — положение верхней границы рассматриваемой области.

Динамические и термические параметры воздушного турбулентного по-

тока представлены автомодельными зависимостями вида [5]:

ит(г) = -ш, т = Т0 + —М';), К&) = км,фв(§), (2) к к

где ит = Ум2 + V2 — модуль скорости, д = —, Ьт — масштаб турбулентности,

кии „ к(Гк - То)

и* =--динамическая скорость, У* =--масштаб температу-

/и(Як) !т (Як)

ры, щ, Гй—скорость и температура на уровне z = 1г, ^ = —, фи, /„, /т —

Рш

безразмерные универсальные функции переменной к — постоянная Кармана.

Согласно теории подобия в слое постоянных потоков выполняются соотношения

ди Зу

К7— = саити, К7— = саиту, ж = О, (3)

дг дг

где с^ — коэффициент трения потока о неровности подстилающей поверхности.

Трехмерная структура метеополей аппроксимируется в виде представления

Ф(X, у, г, г) = Фх (х, у, г)Фг(г, г). (4)

Это действительно оправдано, если учесть факт подобия процессов по вертикале и гипотезу о соотношении горизонтальных и вертикальных масштабов.

При этом вектор-функция Ф2 считается безразмерной величиной со структурой, определенной соотношениями (3). Толщина устойчивого слоя к принимается искомой переменной к = к(х, у, г). Уравнение для определения компонентов Фх получаются путем подстановки (3) в исходные уравнения гидротермодинамики и их осреднения по вертикали от г = б(х, у) до г = к:

дки дкии дкиу —I" «1—---1- ах—-— =

дг дх ду

д(к + 6) д ди д ди = -1АТк~——- + -т-НКх— + —ЪКу— + ^ - сл\ф, дх дх дх ду ду

дку дкиу дкуу

-дГ+ал—+ал—= (5)

д(к + 6) д Зу д Зу = -ХАТк \ + —кКх— + — ЪКу— + ^ - са\ф, ду дх дх ду ду

дН дик дук ^ дг дх ду

где ДТ > 0 — вертикальный перепад температур между устойчивым приземным слоем и вышележащей свободной атмосферой, X — параметр плавучести , |и| = V«2 + V2 , ¥Ру — внешние силы. Выражения для внешних сил представляются в виде

Рх = Рх, — Рх2, Ру = Ру, — Ру2, (б)

где Ех1, Еу1 — динамическое воздействие на поток со стороны верхних слоев атмосферы, ¥Х2, Еу2 —проекции силы сопротивления движению воздуха в проницаемых массивах растительности. Выражения для ¥х2, Еу2 имеют вид:

Гх2 = цНfCfa\u\u, Гу2 = цНfCfa\u\v, (7)

где п —доля поверхности, покрытой деревьями, Cf — коэффициент сопротивления, Нf — средняя высота деревьев, а — плотность растительности в лесном массиве.

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

1) при задании произвольных как угодно гладких начальных полей в задаче обтекания рельефа на горизонтальных масштабах порядка нескольких километров в расчетной области генерируется система фронтальных волн, распространяющихся как вдоль, так и против течения и заметно трансформирующих начальное распределение скорости. Следовательно, исходные поля ио, Но не должны задаваться произвольно, а иметь согласованный характер. Невыполнение этого требования приводит к неадекватной картине течения при проведении рабочих расчетов;

2) формулировка краевых условий должна обеспечивать свободный выход волн из области через границы (в том числе и на границах втока). Это означает, что задание крупномасштабной скорости как внешнего параметра путем определения ио на входной границе неприемлемо.

В качестве альтернативы жесткого условия рассматриваются условия вида

ди ди дН дН

ТГс«лГа л«5Г° при (8)

ди ди дН дН

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

Боковые граничные условия формулируются в виде:

и = иг, V = vг, при (х,у) е Г+,

ди ди дv д^ ^ / ^ (10)

= = = Т = при 0>.у)ег дх ду дх ду

где Г+ — граничные участки втекания потока в область, Г- —участки, где поток направлен из области.

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

К = а,АхАу + 1)1, (11)

ди дv дv ди

где = —--—, = ——I- — , а^ ~ 1 —множитель пропорциональности,

дх ду дх ду

Ах, Ду — размеры элементарного сеточного бокса.

Для описания миграции облака примеси рассматривается уравнение полуэмпирической теории переноса и турбулентной диффузии для концентрации С(х, у, г, г)

дС диС дvC д^ — wc )С — +-+-+ —-— =

дг д х ду дг ( )

д дС д дС д дС (12)

= ас—Кх— + ас—Ку— + ас—Кг— + Яс - у/С, д х д х ду ду дг дг

где и, v, w — компоненты вектора скорости, wc —скорость гравитационного оседания примеси, Кх, Ку, Кг — коэффициенты турбулентного обмена в направлениях х, у, г соответственно, ас = 1/Бт, 8т —число Шмидта, Яс — интенсивность эмиссии вещества (г-м—3-с-1), уf — скорость осаждения примеси на кронах деревьев и/или других проницаемых препятствий, подлежащих явному описанию в принятом сеточном разрешении.

Базовый набор краевых условий для уравнения (12) формулируется в следующем виде. На боковых границах области задаются условия

дС

— = 0 при х = О, х = Ьх, * „ „ <13)

—— = 0 при у = 0, х = Ьу. ду

Концентрация на подстилающей поверхности определяется с помощью соотношения

дС

К7 = — = 0 при г = 6 + гм, (14)

дг

где ги — параметр шероховатости поверхности. На верхней границе полагается

дС

К7— при г = Н. (15)

дг

Для получения плоской модели переноса и диффузии примеси, согласованной с базовыми динамическими уравнениями (5), проводится вертикальное осреднение (12) в предположении, что вертикальная структура концентрации близка к гауссовому факелу. Исходя из представления (4), для С2 задается следующее выражение [6]:

С2(г) = N

(г — Нс)2\ / (г + Нс)2

ехр---— + ехр

2о!т / \ 2о2

(16)

где Нс — высота источника эмиссии над землей, ог —стандартное отклонение (дисперсия), Ыс — нормировочная константа, рассчитываемая из соотношения

оо

J С2Шг = 1. 6

Величина ог рассчитывается в соответствии с теоретическими асимптотиками для классов устойчивости Е, F по Пэскуилу. На расстояниях, не превышающих 10 км от источника, для класса Е используется формула

0,03г

1 + 0,0003г'

(17)

где г — длина радиус-вектора, соединяющего источник с расчетной точкой на плоскости (х, у).

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

дкС дкиС дкуС / д , дС д , дС \

-+-+-= ас\—НК— + —НК—\ + Rc-w„c. (18)

дг дх ду \дх дх ду ду)

Краевые условия для этого уравнения задаются в виде (13).

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

Пространственная аппроксимация дифференциальных операторов основана на монотонных схемах и схемах с невозрастанием полной вариации. Для уравнения (16) она определяется соотношением второго порядка точности

^ = (<?<•+1 - <?0 - - <?г-1)> (19)

где г — номер сеточного узла, Сг+1/2, А^щ — неосциллирующие интерполян-ты, определенные на 5-точечном шаблоне с узлами г, г± 1, г±2 и зависящие от значений ц в них.

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

2. Краткое описание условий проведения полевых экспериментов

Экспериментальные поля размерами 6 X 8 км были оборудованы средствами технического контроля плотностей концентрации распыляемой жидкой рецептуры по неравномерной сетке. Эксперименты проводились в летнее время при температуре воздуха 11—20°С на высоте 0,5 м, скорости сносящего ветра 0,5-5,6 м/с на высоте 1 м и 4,2-6,8 м/с в среднем по слою.

Эксперименты проводились в условиях блокирующей инверсии и при нейтральной стратификации ПСА.

Пролет самолета сельскохозяйственной авиации выполнялся, как правило, перпендикулярно вектору ветра, измеренного шар- пилотным способом. Высота полета варьировалась от 30 до 100 м.

Средства технического контроля располагались на поверхности земли или примятой травяной растительности, как правило, параллельно земной поверхности. Ввиду полидисперсного характера оседающего аэрозольного образования, они аккумулировали частицы аэрозоля от 10 мкм и более, фиксируя фактическую плотность осевшего аэрозоля за все время осаждения. С целью контроля материального баланса массы распыленной и осевшей на экспериментальное поле рецептуры обработке подвергались все средства контроля, зафиксировавшие плотность осаждения аэрозоля с отличием до 0,1% от ожидавшегося расчетного равномерного значения. Итоговая суммарная ошибка материального баланса не превышала в опытах 5%.

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

Верификация математической модели, изложенной в разделе 1, осложняется необходимостью учета динамики осаждения частиц любого диаметра более 10 мкм. С целью решения этой проблемы фактическая концентрация на техническом средстве контроля представлялась в виде интегральной характеристики по всему спектру частиц аэрозоля. Источник аэрозоля представлялся в виде суперпозиции частиц диаметром до 30 мкм и далее до 160 мкм с интервалом в 20 мкм с заданием плотности распределения частиц аэрозоля в каждом интервале.

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

Из двенадцати полевых экспериментов на равнинной местности приняты зачетными все двенадцать. Поточечное сравнение результатов математического моделирования в десяти экспериментах показало, что в 80% расчетные значения практически совпадают с экспериментальными, а в оставшихся 20%—отличаются на 10-18% от фактически наблюдавшихся.

В двух полевых экспериментах имело место совпадение (отличие не более 5%) в 58% точек размещения средств технического контроля. Отличие до 20% наблюдавшихся в эксперименте и расчетных результатов математического моделирования имело место в 31% точек размещения средств технического контроля. Превышение более чем на 50% наблюдавшихся концентраций от расчетных имело место в 11% точек контроля. Более тщательный анализ в ходе математического модельного варьирования метеоусловий при проведении эксперимента показывает, что в указанных двух полевых

опытах могло присутствовать случайное изменение скорости и направления ветра (при этом имеет место совпадение результатов полевого мониторинга и результатов математического моделирования на уровне 85-88%), которое, видимо, не было зафиксировано средствами измерения профиля и скорости ветра в ходе эксперимента.

Из шести полевых экспериментов на пересеченной местности с общим числом равномерно расположенных четырехсот точек контроля зачетными приняты пять. Из-за отсутствия технических средств непрерывного мониторинга метеоусловий над экспериментальным полем размером 2 X 2 км измерение скорости и направления ветра производилось в одной точке шаро-пилотным методом. Очевидно, что распространение измеренного профиля ветра на все расчетное поле имело достаточную условность, особенно для случаев блокирующей инверсии ПСА.

Поточечное сравнение результатов математического моделирования с наблюдавшимися в эксперименте в указанных пяти опытах показало отклонение в пределах 5% в 64% контролируемых точек, отклонение не более чем на 20% в 28% контролируемых точек, а в 10% контролируемых точек — отклонение до 50%.

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

Вывод

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

Литература

[1] Шлычков, В.А. Численная модель пограничного слоя атмосферы с детализацией конвективных процессов на основе вихреразрешающего подхода / В.А. Шлычков // Аэрозоли Сибири. - Новосибирск: изд-во СО РАН, 2005.

[2] Куценогий, К.П. Рассеяние аэрозолей в приземном слое атмосферы / К.П. Куценогий // Прикладная механика и техническая физика. -1970. - №4.

[3] Оседание грубодисперсного аэрозоля на подстилающую поверхность земли / В.Ф.Дунский [и др.] // Труды ГГО. - 1966. - Вып. 185.

[4] Бызова, Н.Л. Экспериментальные исследования атмосферной диффузии и расчеты рассеяния примеси / Н.Л. Бызова, В.Н. Иванов, Е.К. Гаргер // - СПб.: Гидрометеоиздат, 1991. - С. 8-23.

[5] Монин, А.С. Статистическая гидромеханика / А.С. Монин, А.М.Яглом // - СПб.: Гидрометеоиздат, 1992. - Т. 1. - 694 с.

[6] Ханна, С.Р. Применение исследований в области турбулентности для моделирования загрязнений воздуха / С.Р. Ханна // Атмосферная турбулентность и моделирование распространения примесей. - Л.: Гидро-метеоиздат, 1985. - С. 281-314.

Поступила в редакцию 31/ VIII/2006;

в окончательном варианте — 31/^11/2006.

A MATHEMATICAL MODEL OF POLYDISPERSION AEROSOLIC COMPOUND MIGRATION AND SEDIMENTATION3

© 2006 M.V.Menshov4

The paper is aimed at verification of mathematical model of polydispersion aerosolic cloud migration and sedimentation. The cloud is formed from the liquid compounds spraying by means of the aero planes used for fields spraying treatment from the average height. The paper describes the environment and the methodology of field experiments. It also presents the results of a comparative analysis, which are obtained by measurements of land control and data computation according to the proposed model.

Paper received 31/ VIII/2006. Paper accepted 31/VIII/2006.

3Communicated by Dr. Sci. (Tech.), Prof. I.S. Zaguzov.

4Menshov Maxim Vladimirovich (ssaa@samara.ru), Samara State Agriculture Academy, 446409, Kinel, Russia.