CC BY
32
МОДЕЛИ И МЕТОДЫ ДЛЯ ОЦЕНОК ЭКОЛОГИЧЕСКИХ ВЗАИМОДЕЙСТВИЙ В СИСТЕМЕ ГОРОД - РЕГИОН*
В. В. Пененко, Е.А. Цветова Институт вычислительной математики и математической геофизики СО РАН, Новосибирск, Россия e-mail: penenko@sscc.ru, E.Tsvetova@ommgp.sscc.ru
The ecological interconnections, risk assessment and prospects in the city-region system for Siberian conditions are considered. The global atmospheric dynamics is taken into account. To describe these interconnections, the models of transport and transformation of pollutants are combined with the hydrodynamics models at different time and space scales. Two kinds of objects are considered: a region as a part of the global climatic system and a big city as a part of the region. Connections of the source-receptor type do always exist between such objects. The methods of decomposition and complexing provide a convenient and consistent tool for a description of the processes at different scales. Results of the scenarios on forward and inverse problems are presented in which some of the big industrial cities of the Siberian Federal District are taken as the aggregated sources and receptors in these scenarios.
Введение
Динамика климатоэкологической системы является результатом взаимодействия широкого спектра процессов естественного и антропогенного происхождения с различными пространственно-временными масштабами. Получить их адекватное описание с помощью одной модели проблематично, поэтому приходится создавать многофункциональные комплексы. Естественно, возникают вопросы согласования элементов комплекса и эффективности его реализации. На разных этапах технологии моделирования требуются и методы декомпозиции для детализации описания процессов, и методы комплексирования с целью объединения различных элементов в единую систему.
Успешно справиться с этой проблемой помогают вариационные принципы. При специальной организации определяющих функционалов методы декомпозиции, расщепления и комплексирования выступают как инструменты алгоритмической реализации этих принципов. Методы декомпозиции и расщепления дают возможность построить простые и эффективные численные схемы для сложных многомерных моделей [1-4]. В свою очередь,
* Работа выполняется по Программам фундаментальных исследований РАН (№ 13), ОМН РАН (грант № 1.3.2), при финансовой поддержке Российского фонда фундаментальных исследований (грант №04-0564562) и Интеграционными грантами СО РАН (гранты №03-130, №03-137, №03-138).
© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2005.
вариационные принципы обеспечивают согласование различных элементов моделей и алгоритмов на всех этапах технологии моделирования, организацию взаимодействия между различными моделями и данными наблюдений, построение прямых и обратных связей в соответствии с заданными критериями и ограничениями. Ключевую роль вариационные принципы играют также при организации методов обратного моделирования с использованием прямых и сопряженных задач [3] и в формировании различных вариантов распараллеливания вычислительных процессов.
Аэрозоли и газовые примеси, являющиеся загрязнителями атмосферы, оказывают влияние на ее радиационный баланс и качество воздуха. Наиболее значимые источники загрязнений в регионе — это крупные города. В то же время города являются рецепторами, т. е. получателями загрязнений как от своих источников, так и за счет дальнего переноса. Процессы в окружающей среде тесно связаны с процессами в социальной сфере. Это значит, что необходимы адекватные средства моделирования для оценок экологических перспектив и рисков с целью планирования природоохранной и хозяйственной деятельности с позиций безопасности и устойчивого развития. И прежде всего необходимо понимание основных физических и химических процессов, управляющих полями концентраций примесей. Исследования по проблемам окружающей среды активно ведутся во всем мире. Обзоры различных подходов к решению природоохранных задач можно найти, например, в серии регулярных выпусков трудов "Air pollution modelling and its application" (N.Y.: Kluver Acad. / Plenum Publ.) и в списках литературы, имеющихся в этих трудах.
В настоящей работе рассматриваются лишь некоторые аспекты развиваемой нами технологии моделирования на базе вариационных формулировок, относящиеся к построению и комплексированию дискретных моделей для прямых и сопряженных задач динамики атмосферы и охраны окружающей среды.
1. Математические модели
Для исследования масштабов долгопериодных экологических взаимодействий в регионе за счет импорта-экспорта и трансформации загрязнений между различными источниками и рецепторами требуются модели нескольких типов: мезомасштабные, с детальным воспроизведением специфики городской атмосферы; региональные (для описания прямых и обратных связей в регионе) и глобальные (для обеспечения региональных моделей информацией о фоновых процессах). Все эти модели организованы иерархически по масштабам и определены во вложенных областях. Они основаны на одних и тех же физических предпосылках, но имеют различные пространственно-временные масштабы и, как следствие, различные параметризации. Возникает проблема их взаимной согласованности. Как было упомянуто выше, вариационный принцип позволяет обеспечить это согласование на алгоритмическом уровне. Для унификации построения численных схем и алгоритмов модели представляются единообразно, в виде интегральных тождеств с использованием функционалов энергетического типа [4, 5].
Запишем вариационную формулировку комплекса региональных моделей гидротермодинамики атмосферы совместно с моделями переноса и трансформации влаги и примесей с помощью интегрального тождества
I(<р, Y, = а^к(Я(ЩкФ
к
*
- + б»Ж* + 7г(и*ёгас1Ф - иёг&аФ*) + (ф* - Ф^тг-)
¿=1 V о /
+ -(а* -аТ*) + т 1-ст
дп дп
(и* - иТ*)^ + (V* - ьТ*)^-дх ду
+ аххх*> dDdt+
+
дЬ ]
0
dSdt + пФ*undОdt = 0;
Пг
(Лу*,у*) = { I |0.5(у*М(пу) — уИ(пу*)) + п^г^ у grad у* |> dDdt+ КРг
ду
Пг
0.5срип — оп) жу*<Ш1<И
г = 1, 9 + п.
Здесь М(пу) — антисимметричный оператор переноса
кшр/т) + -^-(пиу/т) х у
+
дпуо
(1)
(2)
X = па + divs(пu)dа — о divs(пu)dа , ц = пЯТу/(п + рт/о).
Вектор Ф = = 1,5 + гг} = {(//?,/3 = 1,5, С«, ск = 1,п, п > 1} представляет совокуп-
ность оптически и химически активных субстанций в газовом и аэрозольном состоянии; и = (и, V, о) — вектор скорости в системе координат х, у, о; о = dа/dt, о = (р—рт)/(рв —рт), 7Г = р3 — рт, <у> = (и, V, а, Т, Ф, Ф, х, 7г) = {Уг, г = 1,п+ 12} € <3(А;) — вектор-функции состояния, & ср* £ — сопряженные функции такой же структуры; V = {У^, г = 1, А} £ R(Dt) — вектор параметров модели и источников; Т — температура; Ть — виртуальная температура; Ф — геопотенциал; ^ — характеристики влажности атмосферы (водяной пар, облачная и дождевая вода, ледяные и снежные кристаллы); Са — концентрации загрязняющих примесей; п — число различных субстанций: р, рт, рв — давление, давление на верхней границе воздушной массы и на поверхности Земли; ^ = (р,^х, , ) — коэффициенты турбулентного обмена для субстанции у в направлении координат х, у, о соответственно; (R('Ф))k — операторы трансформации субстанций; т — масштабный коэффициент карты; о — координата, следящая за рельефом поверхности Земли; Qi — источники; е = {е*} — функции, учитывающие неопределенности моделей; а^ — весовые коэффициенты скалярного произведения; О — граница области изменения пространственных переменных D; S — горизонтальная проекция области D; ип — нормальная составляющая вектора скорости к границе. Формы (Лу *, у*) соответствуют операторам переноса и турбулентного обмена с учетом краевых условий, функции х вводятся так, чтобы учесть уравнение для расчета вертикальной составляющей вектора скорости как дополнительное ограничение на класс функций состояния и чтобы упростить алгоритмы построения сопряженных задач. Граничные и начальные условия для функций состояния учитываются в тождестве (1), (2) через интегралы по границам О^ и St области Dt и при Ь = 0. Условия для сопряженных функций являются следствием этого тождества и конкретной
1
1
а
постановки задачи. Свойство энергетического баланса хорошо видно на примере формы (1) при подстановке = <рг,Т* = 1}.
Для решения практических задач экологии и климата требуются оценки обобщенных характеристик изучаемых процессов, определенных на множествах функций состояния и параметров. Их можно всегда записать в виде функционалов
Фк(<р) = ! У)Х*(х, к=1,К, К > 1, (3)
где (ф, V) — функции, подлежащие оценке; Хк(х, ^ > 0 — весовые функции; Хк(х, t)dDdt порождаемые ими меры Радона или Дирака в области Dt. Здесь мы главным образом рассматриваем задачу об оценке рисков для областей-рецепторов. В этом случае функции Хк € Я* ) используются для описания конфигурации областей-рецепторов, а функции ^к € Q(Dt) — для оценок качества воздуха в атмосфере рецепторов.
Под риском мы подразумеваем вероятность поступления в область рецептора загрязнений от действующих и потенциально возможных источников эмиссии и образование опасных высокотоксичных продуктов их трансформации. Количественные оценки рисков получаются с помощью соотношений теории чувствительности для функционалов (3).
Для комплексирования моделей разных масштабов введем целевой функционал, связывающий функции состояния, рассчитанные с помощью модели (1), (2), с полями, рассчитанными по данным наблюдений или с помощью "фоновой" модели с более грубым разрешением:
Фо(ф) = I(ф - фь)ТWо(ф - фь)Хк(х, t)dDdt + ^(еТW1 ei)xldDdt. (4)
В В
Здесь весовые функции Хг(х, t) > 0 и меры выбираем так, чтобы описать векторы фь и е в сеточной области Dh € Dí; Wг, г = 0, 1, — диагональные матрицы с положительными элементами, определяемые структурой "внутренних" скалярных произведений в составе тождества (1), (2). Таким образом, тождество (1), (2) объединяет различные по физическому содержанию модели, а функционалы (3), (4) предназначены для совместного анализа результатов этих моделей и для связей их с глобальными процессами. Изучать функционалы и их вариации в пространстве параметров и источников будем с помощью вариационного принципа в предположении, что модели описывают связи между функциями состояния и параметрами. При таком предположении для построения необходимых алгоритмических конструкций введем семейство расширенных функционалов
Фк(ф) = Фк(ф) + [/(ф, У, <р*)]в„ к = 0,К. (5)
2. Дискретные модели и декомпозиция
Для дискретизации интегрального тождества и функционалов (5) воспользуемся методом слабой аппроксимации с дробными шагами по времени [4]. Для удобства представления формул перепишем тождество (1), выделяя интегралы по времени:
г
/(ф, У, Ф*) = / {0.5 [(В%, ф% - (Вд-£, ф)в] + а(ф, У, Ф*)} Я + о
+0.5(Вф, ф*)г>|о = 0, (6)
где а(ф, У, ф*) = /в (С(ф, У), ф*) йВ — функционал модели по пространственным переменным; С(ф, У) — совокупность "пространственных" операторов модели с учетом функций источников и неопределенностей; В — блочно-диагональная матрица. Следуя идее метода расщепления [1, 3], используем представление
г
С(ф, У) = £ Ск(ф, У), г > 1; (7)
к= 1
а(ф, У, ф*) = ^ ак(ф, У, ф*), ак(ф, У, ф*) = [(Ск(ф, У), ф*)йВ. (8)
к=1 в
Необходимое требование к выбору представления (7), (8) состоит в том, чтобы антисимметричные пары слагаемых в тождестве (1) не разрывались, т.е. чтобы сохранялось свойство энергетического баланса. Теперь можно непосредственно воспользоваться техникой вариационного исчисления для построения дискретных аналогов модели [4]. Сначала строятся сумматорные аналоги функционалов (1)-(6), для чего вводятся сеточные области В^, Б^, и определяются дискретные аналоги функциональных пространств Qh(Вí),Q*h(Вt). Затем интегралы аппроксимируются кубатурными формулами, а подынтегральные выражения — дискретными представлениями.
Для приближения адвективно-диффузионных операторов переноса в составе форм (2) строятся монотонные схемы, обладающие свойством транспортивности по отношению к потоку воздушных масс. Для этих целей вводятся конструкции типа соотношений баланса с использованием аналитических решений локальных сопряженных задач по координатным направлениям [6].
Кубатуры для аппроксимации функционалов (3), (4) в составе (5) выбираются в соответствии с аппроксимациями интегрального тождества. При необходимости, когда детально рассматривается атмосфера городов-рецеторов, производится декомпозиция Dí на подобласти мезомасштабов. В каждой из таких подобластей вводятся мезомасштабные версии моделей в форме тождества (1), (2), их аппроксимации строятся формально независимо, но так, чтобы в совокупности получалась единая конструкция в области Вг. При таком подходе все согласования между частями модели на подобластях происходят автоматически. В силу свойств аддитивности функционалов (1)-(6) при выборе представлений (7), (8) и их дальнейших преобразований естественно использовать принцип распараллеливания. В этом случае можно организовывать параллельные ветви алгоритмов, начиная с уровня обобщенного описания моделей.
Основная идея в построении схем расщепления с помощью метода слабой аппроксимации по времени выражается следующим соотношением:
г
/' Д 4- '
к=1 2 к=1
а(ф, У, ф*)^ = / V а,(ф, У, ф*)^ « ^ > а^1 + ^ \, (9)
имеющим точность 0(Д£3); 7 = (2к — 1)Д^/4г. Окончательно, аппроксимируя производные по времени и выполняя дискретизацию по пространственным переменным, получим
сумматорный аналог тождества (6)
ь
*\Ь\.1
+ Е К3 +7 + а^1 дЛ + 0.5(£ф, ф*)ь|. = о, (10)
к=1
аппроксимирующий тождество с точностью 0(А^2). Индекс к обозначает дискретный аналог соответствующих объектов.
Численые схемы для модели получаются из условий стационарности сумматорных функционалов /ь(ф, У, ф*) и Фь(ф) к вариациям функций ф* и ф в узлах сетки Для прямых задач эти условия выражаются с помощью операций
д
{/"(¥>, у, Ф*) + = 0, Ф* е д^д71), (11)
дф
а для сопряженных задач
д_ дф
— (Лф,у,ф*) + ф£(ф)} = 0, <ре<У(1%), к = 0, К. (12)
В результате этих преобразований получаются двуциклические схемы расщепления [1, 4] для решения основных и сопряженных задач. Эти схемы взаимно согласованы и энергетически сбалансированы.
Решения задач (11), (12) {<р, к = 0, К} являются входными данными для реализации алгоритмов теории чувствительности [4]
д_
6Фк(<р) = ^гас1тФ^),£У) = — Г(<р,<р*к, Уь + £<ПГь)|*=о, к = 0,К, (13)
где £ — вещественный параметр; ^гас1тФ&(ф) = дФк((р/дУг, г = 1, А} — набор функций чувствительности функционалов Фк (ф) к вариациям параметров в окрестности Уь основного состояния параметров, который задавался для решения задач (11), (12).
3. Комплексирование моделей
Для целей комплексации решений задачи (11) с решением фь фоновой задачи, т.е. задачи для моделей более крупного масштаба, аппроксимации строим так, чтобы функции неопределенности модели включать только на последнем этапе схемы расщепления на дробных шагах \Ъ^-1/2г]. Система дискретных уравнений этого этапа имеет вид
*—^-+ = + 4+1, (14)
где Агр = д[а!Ъ(ф, У, ф*)]/др* — дискретный аналог оператора схемы расщепления на этапе с номером г; /^ — функции источников, а — функции, учитывающие неопределенности исходной модели и ее схемы расщепления на шаге . Процедура комплексиро-вания сводится к аддитивной процедуре усвоения данных [7]. Расчетная схема получается из условия минимума функционала (4) относительно функций при ограничении,
выраженном системой уравнений (14). Такая модификация осуществляется только на отдельных этапах расщепления, не затрагивая общую структуру системы в целом.
Фоновый поток представляющий собой пространственно-временную конструкцию, играет роль направляющего фазового пространства, на фоне которого реализуются региональные процессы, описываемые моделями (1), (2) и (11), (12). Его учет в режиме комплексирования позволяет на базе моделей более мелких по иерархии масштабов рассчитывать долговременные сценарии развития ситуаций. Расчеты по моделям на ограниченной территории с формально задаваемыми краевыми условиями без учета динамики фонового движения не имеют практического смысла.
4. Пример сценария моделирования
Развиваемая система моделирования предназначена для решения в рамках сценарного подхода задач, связанных с оценками экологической перспективы и безопасности индустриального региона при различных вариантах антропогенных нагрузок. Так, решаются прямые задачи распространения загрязнений от городов как агрегированных источников и обратные задачи по оценке рисков загрязнения атмосферы в городах-рецепторах. Приведем пример сценария моделирования на базе описанной системы по изучению масштабов взаимодействия в Сибирском федеральном округе (СФО). Конкретно, приведены результаты решения обратной задачи по оценке рисков получения загрязнений несколькими крупными городами. Названия этих городов-рецепторов приведены на рисунках. Для оценок используется функционал вида (3), в котором оцениваемая функция представляет собой тождественный оператор, а носители весовых функций указывают на местоположение этих городов в регионе. В расчетах использовались две модели — глобальная, в полусферном варианте, и региональная. Общая структура комплекса и схема взаимодействий между моделями строятся по аналогии с системой моделирования, предназначенной для исследования Байкальского региона, описанной в [8]. Восстановление общей циркуляции атмосферы в Северном полушарии осуществляется с помощью полусферной модели по данным реанализа [9]. Способ построения сценариев изложен в [10]. Циркуляция в регионе рассчитывалась по региональной модели с учетом данных полусферной модели в режиме комплексирования — усвоения с функционалом качества типа (4). На фоне полученной циркуляции атмосферы в регионе решались прямые и сопряженные задачи на базе модели переноса примесей. Для оценок рисков использованы функции чувствительности целевого функционала (3) к вариациям источников примесей, распределенных в регионе.
Выбранный регион представляет собой часть СФО и расположен между 80-95° в.д. и 40-50° с.ш. На рисунках ось ординат показывает дополнение до широты. Разрешение по горизонтали 15' х 15' в сферических координатах, по вертикали 19 уровней в гибридных (о — р) координатах до уровня 10 мб. Перепад высот рельефа в регионе составляет около 3 км. Шаги по времени во всех моделях адаптивные и выбираются в зависимости от интенсивности процессов в соответствии с критериями аппроксимации. Длительность сценария — 30 суток, 1-30 июня 2002 года. Летний сценарий интересен тем, что июнь является типичным для атмосферной циркуляции в период действия летнего Саяно-Алтайского циклогенеза. На рисунках показаны результаты сценария обратного моделирования, соответствующие 12 и 8 июня 2002 года. Для иллюстрации выбраны два фрагмента, представляющие поля циркуляции атмосферы и функции риска по отношению к источникам, расположенным на поверхности Земли. Функция риска рассчитывается через решения
.У 12 6.2002 11: 30
1с ть
У а В .2002 18: 0
ютй
Двумерный разрез функции риска на уровне приземного слоя, относящийся к моментам времени 11:30 12.06.2002 (а) и 18.00 08.06.2002 (б)
сопряженных задач. Двумерные разрезы полей соответвуют уровню приземного слоя. Атмосферная циркуляция в регионе имеет высокую степень изменчивости, и, соответственно, функции риска следуют движениям атмосферы. Для иллюстраций выбраны два момента времени из месячного интервала, когда риск получения загрязнений для рецепторов наибольший, если в качестве меры риска рассматривать площади, занимаемые носителем функции риска заданного уровня значимости.
Заключение
Представленная система моделирования позволяет проводить сценарные оценки экологических перспектив для индустриальных регионов и городов при одновременном воздействии естественных и антропогенных факторов. Наряду с традиционными прямыми
задачами по изучению процессов распространения и трансформации примесей в системе активно используются методы обратного моделирования для оценок экологических взаимодействий городов и промышленных объектов в регионе с позиций коллективной безопасности и устойчивого развития. Вариационный принцип в сочетании с методами декомпозиции и комплексирования предоставляет удобный формализм для организации исследований и построения многофункционального комплекса моделей с заданными свойствами и алгоритмов их реализации, согласованных на всех этапах технологии моделирования. Объединение моделей с помощью оптимизационых процедур последовательного усвоения данных дает возможность рассчитывать долговременные сценарии развития климатоэкологических ситуаций в регионе с учетом динамики глобальных процессов и трансграничных переносов загрязняющих веществ.
Список литературы
[1] Марчук Г.И. Методы вычислительной математики. М: Наука, 1977.
[2] Марчук Г.И. Сопряженные уравнения и анализ сложных систем. М.: Наука, 1992. 335 с.
[3] ЯнЕнко Н.Н. Метод дробных шагов решения многомерных задач математической физики. Новосибирск: Наука, 1967.
[4] Пененко В.В. Методы численного моделирования атмосферных процессов. Л.: Гидроме-теоиздат, 1981.
[5] Пененко В.В., Алоян А.Е. Модели и методы для задач охраны окружающей среды. Новосибирск: Наука, 1985. 254 с.
[6] Пененко В.В. Численные схемы для адвективно-диффузионных уравнений с использованием локальных сопряженных задач. Новосибирск, 1993 (Препр. РАН . Сиб. отд-ние. ВЦ СО РАН. № 984). 49 с.
[7] Цветова Е.А., Пененко В.В. Быстрое усвоение данных в атмосферных и океанических исследованиях // Вычисл. технологии. 2002. Т. 7. Спец.выпуск. С. 141-146.
[8] Пененко В.В., Цветова Е.А. Математические модели для исследования взаимодействий в системе озеро Байкал — атмосфера региона // ПМТФ. 1999. Т. 40, № 2. C. 137-147.
[9] Kalnay E., Kanamitsü M., Kistler R. et al. The NCEP/NCAR 40-year reanalysis project // Bull. Amer. Meteor. Soc. 1996. Vol. 77. P. 437-471.
[10] Пененко В.В., Цветова Е.А. Формирование сценариев атмосферной циркуляции для климатоэкологических исследований // Оптика атмосферы и океана. 2001. Т. 14, № 6-7. C. 491-495 .
Поступила в редакцию 2 июня 2005 г.
