Оценка полей концентрации метана над северным полушарием по данным измерений и модели переноса и диффузии пассивной примеси Текст научной статьи по специальности «Науки о Земле и смежные экологические науки»

ОЦЕНКА ПОЛЕЙ КОНЦЕНТРАЦИИ МЕТАНА НАД СЕВЕРНЫМ ПОЛУШАРИЕМ ПО ДАННЫМ

_ __и ______

ИЗМЕРЕНИЙ И МОДЕЛИ ПЕРЕНОСА И ДИФФУЗИИ ПАССИВНОЙ ПРИМЕСИ*

Н. В. КиллновА, Е. Г. Климова Институт вычислительных технологий СО РАН, Новосибирск, Россия e-mail: [email protected], [email protected]

A passive pollution data assimilation algorithm based on the Kalman filter theory is investigated The concentration fields forecast is given with the help of semi-Lagrangian model of passive pollution transport and advection for the Northern Hemisphere. Kalman filter algorithm requires calculation of very high order matrixes, which currently makes its practical realization impossible. Therefore, an approach, which relays on the simplified model for calculation of the forecast error covariance matrixes, is considered. The algorithm for estimation of the methane emissions based on the assimilation procedure is proposed.

1. Полусферная модель переноса и диффузии пассивной примеси

Предварительная оценка полей концентрации примеси дается с помощью модели переноса и диффузии пассивной примеси:

ds udsvds ds

тгг н--ттт + -ТГ- + {w - Wg) -Trot a cos yd A ady dz

1 d ds\ 1 d Í1 ds\ d Í1 ds'. „

a2 cos2 yd A \ dAj a2 cos ydy \ dA) dz \ dz

s(A, i, z, t) = s(A + 2n, i, z,t), s = sG при y = 0,

ds , , , , ds

v— = a\s — s0) ПРИ £ = b + n, v— = 0 при z = ti,

dz dz

где A — долгота; y — широта; z — высота над уровнем моря; s — концентрация примеси; V = (u,v,w — wg) — вектор скорости ветра с компонентами в направлениях A, y, z соответственно; wg — скорость гравитационного оседания; a — средний радиус Земли; k\,k2 — коэффициенты турбулентного обмена в горизонтальном и k3 — вертикальном направлении; f — функция, описывающая размещение и мощность источников; (A , y) — функция,

* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (гранты № 04-05-64481, № 03-05-65279).

© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2005.

ОЦЕНКА ПОЛЕЙ КОНЦЕНТРАЦИИ МЕТАНА

133

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

Задача рассматривается в области Dt = О х [0, Т], где О = Б х [Ь + К, И]; Б = {0 < Л < 2п;0 < <р < п/2}.

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

2. Усвоение данных наблюдений в задаче переноса и диффузии пассивной примеси

Используемый в задаче оценивания полей концентрации пассивной примеси алгоритм усвоения данных основан на теории оптимальной фильтрации Калмана. Постановка задачи оптимальной фильтрации подробно изложена в [5], описание алгоритма фильтра Кал-мана можно найти в [6-11]. Дискретный алгоритм фильтра Калмана в линейном случае изложен в [12].

Предположим, что модель переноса и диффузии (1) в дискретном по времени виде записывается следующим образом:

хк = Лк-х хк-1, (2)

где Хк — п-вектор прогнозируемых значений концентрации примеси в момент времени Ьк, его размерность равна количеству узлов сетки: п = п\п(рпх, где п\ — количество узлов сетки по долготе, п^ — по широте, пх — по высоте.

Пусть "истинное" значение концентрации в момент времени Ьк описывается векторным уравнением

хк = Лк-1Хк-1 + £к-1, (3)

где £к — случайный вектор "шумов" модели, удовлетворяющий условиям

Е£к = 0, Е (£к )(£,)Т = Як 5кг ■

Вектор данных наблюдений представим в виде

Ук = Мк Хк + %, (4)

где Ук — т-вектор наблюдений в момент времени Ьк; Мк — (пхт)-матрица; £к — случайный т-вектор ошибок наблюдений, удовлетворяющий условиям

Е£к = 0, Е (£кЖ )Т = Як 5кг,, Е (ф(£г)т = 0.

134

п.о. Киланова, &.±. Климова

Оптимальная оценка концентрации x{ по данным наблюдений (4) и модели переноса и диффузии примеси (2) ищется из условия минимума следа матрицы ковариаций ошибок оценивания. Решением этой задачи является алгоритм фильтра Калмана:

xk = Ak-ixa-i,

pf = Ak-iPk-iAl-i + Qk-i, Kk = Pf MT (Mk Pf MT + Rk )-l,

Pk = (I - Kk Mk )PÍ, xak = x{ + Kk (yk - Mk x{).

f f Здесь x{ — прогноз концентрации по модели, Pf, Pk — матрицы ковариаций ошибок

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

Реализация алгоритма фильтра Калмана на ЭВМ в полной постановке невозможна, так как для современных глобальных моделей порядок ковариационных матриц достигает миллиона. Одним из способов решения этой проблемы является использование упрощенных моделей для расчета матриц ковариаций ошибок прогноза. Такой алгоритм называют субоптимальным алгоритмом фильтра Калмана [9-11, 13].

3. Численные эксперименты

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

xX = AX{-i + nk-i + £k-i, По = ñ, (5)

где A — оператор модели; xf — "истинная" концентрация метана в момент времени tk; nk-i — "истинная" эмиссия метана в момент времени tk-i.

Прогноз концентрации по модели x{ (предварительная оценка поля концентрации) задавался следующим образом:

x{ = Axk-i + nL^ nf = 0,

f

где nf-i — эмиссия метана в момент времени tk-i. В моменты наблюдений моделировались данные наблюдений

yO = xX + Со, (6)

где С0 — случайная величина, распределенная по нормальному закону с нулевым математическим ожиданием и дисперсией afi. В этом случае имеем систему уравнений для ошибок прогноза концентрации Axk и эмиссии Ana

AXX = AAxk-i + Ank-i + £k-i,

Дщ = Ank-i, Ano = Axo.

Оценки полей концентрации и эмиссии в момент усвоения рассчитывались по формулам

X{ = X{ + Pf MT(MkPf MT + Rk)-i(yo - MkX{); (7)

ОЦЕНКА ПОЛЕЙ КОНЦЕНТРАЦИИ МЕТАНА

135

т}% = 4 + Ах(Аг])Т(МкР^М/: + Я,)-1 (у°к - М,4>, (8)

При вычислении матрицы ковариаций использовалось предположение, что вероятностное осреднение можно заменить осреднением по времени:

1 "

Р/ = Л.Г/,(Л.Г/, = ^ Ахг(Ахг)Т.

г=1

По вычисленным ошибкам полей концентрации и эмиссии оценивались ковариации ошибок прогноза полей концентрации метана и кросс-ковариации ошибок концентрации и эмиссии, а затем производилась оценка по данным наблюдений поля концентрации по формуле (7) и поля эмиссии по формуле (8).

Для проверки свойств предлагаемого алгоритма были проведены численные эксперименты с моделируемыми данными наблюдений. В численных экспериментах расчеты проводились на сетке 2.5 х 2.5° по горизонтали и 15 уровнях по вертикали, шаг по времени Д£ = 15 мин. Для расчета распределения пассивной примеси были использованы реальные данные о полях скорости ветра, температуры и давления и модельные данные о концентрации метана [14]. Расчеты проводились на двое суток с усвоением каждые 12 часов, при этом данные наблюдений моделировались по формуле (6) так, что в каждый момент усвоения задавалась информация в полосе [12]. На рис. 1 приведено распределение данных наблюдений. При этом каждая вертикальная полоса соответствует области, в которой задаются данные в 12, 24, 36 и 48 часов.

В процессе усвоения производилась оценка по данным наблюдений полей концентрации и эмиссии метана. Результаты численного эксперимента приведены на рис. 2. На рис. 2,а представлена относительная ошибка концентрации, полученной в процессе усвоения данных с помощью субоптимального алгоритма фильтра Калмана. По оси абсцисс отложены шаги по времени. Оценка поля эмиссии производится в процессе усвоения данных по формуле (5). "Истинное" значение эмиссии задается постоянным на каждом шаге по времени в каждой точке сетки на двух нижних уровнях по вертикали, равном 1.063438-10-5 ррш/е. В начальный момент времени прогностическое значение эмиссии задавалось равным нулю.

Рис. 1. Распределение моделируемых данных наблюдений.

130

н.в. киланова, е. г . Климова

Ù -I-1-1-1- о -1-1-1-

i 30 9? u6 19г i 50 99 us 19?

а б

Рис. 2. Относительная ошибка концентрации метана (а) и оценки эмиссии (б).

На рис. 2,б приведена относительная ошибка оценки эмиссии в субоптимальном алгоритме. Видно, что в процессе усвоения данных происходит восстановление поля эмиссии, при этом ошибка оценки эмиссии после трех процедур усвоения составила около 30 %.

Заключение

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

Список литературы

[1] Марчук Г.И., Алоян А.Е. Глобальный перенос примеси в атмосфере // Изв. РАН. Физика атмосферы и океана. 1995. Т. 31, № 5. С. 597-606.

[2] Марчук Г.И. Математическое моделирование в проблеме окружающей среды. М.: Наука, 1982. 319 с.

[3] Толстых М.А. Полулагранжева модель атмосферы с высоким пространственным разрешением для численного прогноза погоды // Метеорология и гидрология. 2001. № 4. С. 5-16.

[4] Bermejo R., Staniforth A. The conversion of semi-Lagrangian advection scheme to quasimonotone scheme // Mon. Wea. Rev. 1992. Vol. 120. P. 2622-2632.

[5] Jazwinski A.H. Stochastic processes and filtering theory. N.Y.: Acad. Press, 1970. 377 p.

[6] Ghil M., Malanotte-Rizzolli P. Data assimilation in meteorology and oceanography // Advances in Geophysics. 1991. Vol. 33.

оценка нолей концентрации метана

137

[7] Menard R., Cohn S.E., Chang L.-P., Lyster P.M. Assimilation of stratospheric chemical tracer observations using a Kalman filter. Pt 1: Formulation // Mon. Wea. Rev. 2000. Vol. 128. P. 2654-2671.

[8] Menard R., Chang L.-P. Assimilation of stratospheric chemical tracer observations using a Kalman filter. Pt 2: %2-validated results and analysis of variance and correlation dynamics // Mon. Wea. Rev. 2000. Vol. 128. P. 2672-2686.

[9] КлимовА Е.Г. Упрощенные модели для расчета ковариационных матриц в алгоритме фильтра Калмана // Метеорология и гидрология. 2000. № 6. C. 18-30.

[10] Климова Е.Г. Модель для расчета ковариаций ошибок прогноза в алгоритме фильтра Калмана, основанная на полных уравнениях // Метеорология и гидрология. 2001. № 11. С. 11-21.

[11] КлимовА Е.Г. Численные эксперименты по усвоению метеорологических данных с помощью субоптимального фильтра Калмана // Метеорология и гидрология. 2003. № 10. С. 54-67.

[12] КлимовА Е.Г., КилАновА Н.В. Усвоение данных наблюдений в задаче переноса и диффузии пассивной примеси // География и природные ресурсы. 2004. Спец.выпуск: Тр. Меж-дунар. конф. ENVIR0MIS-2004. Новосибирск, 2004. С. 175-180.

[13] Todling R., Cohn S. Suboptimal shemes for atmospheric data assimilation based on the Kalman filter // Mon. Wea. Rev. 1996. Vol. 124. P. 2530-2557.

[14] КрупчАтников В.Н., Крылова А.И. Численное моделирование распределения метана по данным наблюдений на поверхности Земли // Оптика атмосферы и океана. 2000. Т. 13, № 6-7. С. 622-626.

Поступила в редакцию 2 июня 2005 г.