УДК 519.6
DOI: 10.33764/2618-981 Х-2019-6-1-94-101
МЕТОД УСВОЕНИЯ ДАННЫХ ДЛЯ ЗАДАЧИ РАСПРОСТРАНЕНИЯ ПАССИВНОЙ ПРИМЕСИ В АТМОСФЕРЕ, ОСНОВАННЫЙ НА ДИНАМИКО-СТОХАСТИЧЕСКОМ ПОДХОДЕ
Марина Владимировна Платонова
Новосибирский национальный исследовательский 630090, Россия, г. Новосибирск, ул. Пирогова, 2, e-mail: [email protected]
Екатерина Георгиевна Климова
Институт вычислительных технологий СО РАН, 630090, Россия, г. Новосибирск, пр. Академика Лаврентьева, 6, доктор физико-математических наук, доцент, старший научный сотрудник, тел. (383)332-42-57, e-mail: [email protected]
В данной работе рассмотрена методика усвоения данных для задачи распространения концентрации пассивной примеси в атмосфере. Описаны классические подходы к решению подобных задач, выделены особенности применения алгоритмов, их минусы и плюсы.
Рассмотрены два алгоритма: ансамблевого фильтра Калмана и ансамблевого сглаживания Калмана. Рассмотрены различные способы улучшения сходимости этих алгоритмов, такие как локализация и увеличивающий множитель.
Ключевые слова: усвоение данных, ансамблевый фильтр Калмана, ансамблевое сглаживание Калмана, увеличивающий множитель, прогнозирование, перенос и диффузия, концентрация пассивной примеси, локализация, оценка, модель.
DATA ASSIMILATION METHOD FOR THE TASK OF PASSIVE IMPURITY PROPAGATION IN THE ATMOSPHERE BASED ON A DYNAMIC-STOCHASTIC APPROACH
Marina V. Platonova
Novosibirsk National Research State University, 2, Pirogova St., Novosibirsk, 630073, Russia, Graduate, phone: (996)382-07-29, e-mail:[email protected]
Ekaterina G. Klimova
Institute of Computing Technology SB RAS, 6, Prospect Аkademik Lavrentiev St., Novosibirsk, 630090, Russia, D. Sc., Senior Researcher, phone: (383)332-42-57, e-mail: [email protected]
In this paper, we consider the method of data assimilation for the problem the propagation of the concentration a passive impurity in the atmosphere. Classical approaches to solving such problems are described, features of the application of algorithms, their minuses and pros. Two algorithms are considered: the ensemble Kalman filter and the ensemble Kalmans moother. Various ways to improve the convergence of these algorithms, such as localization and inflation factor, are considered.
Key words: data assimilation, ensemble Kalman filter, ensemble Kalman smoother, inflation factor, forecasting, transfer and diffusion, concentration of passive admixture, localization, assessment, model.
государственный университет, магистрант, тел. (996)382-07-29,
Введение
Под усвоением данных принято понимать совместное использование наблюдений и математической модели для оценки состояния окружающей среды. В работе рассмотрена методика усвоения данных для задачи распространения пассивной примеси в атмосфере. Описаны классические подходы к решению подобных задач, выделены особенности применения алгоритмов, их минусы и плюсы. Рассмотрены два алгоритма: ансамблевый фильтр Калмана и ансамблевое сглаживание Калмана. Получены результаты сравнительного анализа этих алгоритмов. Рассмотрено влияние различных значений входящих параметров алгоритмов на получаемые результаты.
1. Ансамблевый фильтр Калмана, ансамблевое сглаживание Калмана Будем рассматривать линейную модель системы:
где уп - вектор наблюдений;
гп - вектор ошибок наблюдений, являющийся гауссовской случайной переменной.
Будем считать хгп «истинным» состоянием системы [1]. Задается ансамбль значений моделируемой переменной в начальный момент времени. Тогда для каждого элемента ансамбля производится шаг прогноза и шаг анализа. Сам ансамблевый фильтр Калмана (EnFK) состоит из двух шагов [2]: шаг прогноза: х^) = Апхлп (г);
шаг анализа: хА(') = х^) + Кп(уП') - Их^}), где номер элемента ансамбля I = 1.. .к; к - количество ансамблей.
Матрица Кп - матрица усиления, она задается формулой:
Хп+1 = Апх1
.а
п '
где Ап - оператор модели;
х{+х - модельный прогноз состояния системы; хпа - найденная оценка, п - номер шага по времени. Наблюдения представим в следующем виде:
Уп = хп +£ п ,
кп = РпРИт (ИРпрИт + Я)-1,
1 ^
матрица ковариаций ошибок про-
гноза;
М'
—р ¿=1
хп+1 = —--среднее значение прогноза вектора состояния;
к
Т
Я = е-8 - матрица ковариаций ошибок наблюдений.
При этом среднее по ансамблю значение вектора состояния системы является оптимальной оценкой в задачи фильтрации. Она достигается при минимуме следа матрицы ковариаций ошибок [2]:
РА = V (ХА°) - ХА )(ХА°) - ХА V п ~ к -1 ^ п+1 п+1 п+1 п+1 '
Задача ансамблевого сглаживания заключается в получении оптимальной оценки в некотором временном интервале, при наличии данных в этом интервале. Причем состояние на шаге п оценивается по данным на временных шагах с номерами больше п. Как показано в работе [3], при наличии гауссовского распределения у ошибок наблюдения и отсутствии корреляции между ошибками наблюдения и прогноза, алгоритм позволяет получать оценку состояния системы по мере поступления данных. Обозначим через АЕпКР (х, ti) ансамбль оценок
состояния системы, полученный после применения алгоритма фильтрации для шага времени ti; а через АЕпКЯ (х, ti) - ансамбль оценок, получаемый в алгоритме сглаживания для шага времени ti. После применения алгоритма сглаживания полученный результат можно представить в виде:
АЕпК8 (^ Ь ) = П ХА
]=к
к+1
V А
'ЕпКР'
где X, - матрицы коэффициентов, вычисляемые для каждого шага по времени: X, = I + Ц; О] = У - Нх; С] = + & -1)* Я,,
У] - вектор наблюдений в момент времени];
Я, - матрица отклонений от среднего для ансамбля;
Я, - матрица ковариаций ошибок наблюдений [3, 4].
]
В данной работе алгоритм усвоения данных рассматривался для задачи с неизвестным параметром. Для решения данной задачи вводим расширенный
Л
, который содержит информацию о концен-
V 8 )
трации пассивной примеси и эмиссии. Будем использовать расширенный ансамблевый фильтр Калмана [5]:
вектор состояния системы: х
шаг прогноза: х^п+1 = АхЩ,
шаг анализа: ха = х1 + Кх (У - Их1); ga = gf + К& (У - Их1);
Кх = Рхх(Рхх + Я)-1;кg = (Рхх + Я)-1;Рх - кроссковариационные матрицы.
Маленький размер выборки является причиной низкой точности оценки, могут существовать ложные ковариации или появиться расходимость алгоритма. Важными элементами корректировки данных алгоритмов являются увеличивающий множитель и локализация. Выборки небольших размеров являются причиной появления ложных ковариаций. Локализация состоит в поэлементном умножении матрицы на функцию от расстояния:
о о -'-У Р/2 ¿2 Р . = Р .ое 1 1 '
где Р у - элемент матрицы ковариации.
Здесь масштаб локализации, который подбирается исходя из размеров сетки по пространству.
Другая проблема возникает, когда теоретическая ошибка оценки (след матрицы ковариаций) и реальная не совпадают. В таком случае прибегают к помощи увеличивающего множителя (тйа1;юп£ас1;ог [6]). Ниже представлен классический вариант увеличивающего множителя - мультипликативный:
Р = Р-а2,
где а - постоянный коэффициент, который берется близким к единице (часто используют значения 1,1; 1,04). В таком случае на каждом шаге по времени матрица ковариаций ошибок оценки умножается на заранее подобранную константу.
Другой подход к заданию коэффициента увеличивающего множителя -адаптивный [7]. В этом случае на каждом шаге по времени вычисляется коэффициент увеличивающего множителя А по формуле[8]:
т
. d d - №асе(Я) 7 — .
Д =-—; d = у. - Их. ,
Шсе{ ИРИ) у у
где /гасе(...) след матрицы, d - вектор невязок.
Основной задачей увеличивающего множителя является предотвращение расходимости алгоритма. Величину увеличивающего множителя принято ограничивать [9]: 0,9 <а< 1,5. Для того чтобы сгладить значения по времени предлагается использовать формулу: Д1 = кДм + (1 - к)Д^, где к = 0,05.
В данной работе рассматривается способ вычисления увеличивающего множителя, состоящий в осреднении вектора невязок по времени. Осреднение производилось на некотором временно интервале, который выбирался, исходя из параметров эксперимента.
2. Численные эксперименты
С использованием рассматриваемых методов были проведены численные эксперименты с моделью переноса-диффузии пассивной примеси. Использовался метод расщепления задачи по физическим процессам. Решение уравнения переноса находилось с помощью полулагранжева метода. Уравнение диффузии решалось с помощью циклической прогонки. Во всех экспериментах область определения: 0 < х < 1, 0 < t < 1. В качестве граничного условия было взято условие периодичности q (0, t) = q(1, t). Эмиссия была задана
&(t) =
0,1 при 0,2 < х < 0,5; 0 иначе.
Для оценки качества работы алгоритмов были рассмотрены следующие параметры [10]: след матрицы ковариаций ошибок оценки и среднеквадратиче-
t а 2
ская ошибка оценки, которая считалась следующим образом: тш8п = хп - хп , в качестве нормы берется вторая норма вектора, х* - реальное состояние сис-
а
темы, х - получаемая оценка.
Результаты численных экспериментов приведены на рисунке для алгоритма оценки эмиссии пассивной примеси с применением различных увеличивающих множителей.
Я
Я
О «
и н о и о и
о <ц
8 13
а
и
и <ц
я
<ц
а о
□.□1 о.сгв 0.0118 0.0117 0.0)6 О.ООь 0Ш4 о.ооз
0.032 0.001 а
1 .2 3 4
номер шага по времени
и о ю 8
ш о
«
и &
и ар
в о
И
3 я и
!■
е «
о
а^аптийныл с ими и мультипликативный СМ зелоный <— мультипликативный а=1 1 - фиолетовый
ЁВ2 уННОТКЗКИЦВГП М тЧ:'НТЙ Л Я - ЬфЭСНЫЙ
1 :а- 3 4 5
номер шага по времени
а)
б)
Результаты численных экспериментов для алгоритма оценки эмиссии
пассивной примеси:
а) график среднеквадратического отклонения оценки эмиссии; б) график следа матрицы ковариаций ошибок оценок эмиссии
На рисунке различными цветами изображены результаты, полученные при использовании мультипликативного увеличивающего множителя (а = 1,1 -
фиолетовый, а = 1,04 - зеленый), результаты адаптивного подхода - синим, а результаты без использования увеличивающего множителя - красным. По горизонтальной оси отложен номер шага по времени. Все эксперименты проводились с ансамблями размером 150 элементов. Сетка по пространству была размером 100 узлов. Сетка по времени была размером 100 шагов, на графике представлены результаты первых шести шагов.
Самым действенным способом предотвратить расходимость алгоритма является адаптивный способ вычисления увеличивающего множителя. Среди графиков следа матрицы ковариаций ошибок самым высоким ростом следа отмечается график, полученный при реализации адаптивного подхода. На графике средне-квадратической ошибки оценки заметим, что данный вывод подтверждается самым быстрым падением значения отклонения при применении адаптивного увеличивающего множителя, что дает меньшую ошибку в получаемой оценке.
Заключение
В работе рассмотрены различные подходы к применению ансамблевого фильтра Калмана и ансамблевого сглаживания в задаче оценивания концентрации пассивной примеси. Проведен сравнительный анализ различных подходов к улучшению сходимости, таких как применение увеличивающего множителя и использование локализации. Выполнен сравнительный анализ разных подходов к вычислению коэффициента увеличивающего множителя. Сравнение результатов работы алгоритмов при различных значениях увеличивающего множителя показало, что наилучший способ - адаптивный. При использовании адаптивного способа вычисления увеличивающего множителя получается более точная оценка параметра системы.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Carrassi A., Bocquet M., Bertino L., Evensen G. Data assimilation in the geosciences: An overview of methods, issuers and perspectives // Wiley interdisciplinary reviews: Climate Change. 2018. V. 131, Issue5, e535, doi: 10.1002/wcc535.
2. Evensen, G., P.J. van Leeuwen An ensemble Kalman smoother for nonlinear dynamics // Monthly Weather Review. 2000. V. 128. P. 1852-1867.
3. Klimova, E. A suboptimal data assimilation algorithm based on the ensemble Kalman filter // Quarterly Journal of the Royal Meteorological Society. 2012. V. 138, P. 2079-2085.
4. Klimova E. G. Application of ensemble Kalman filter in environment data assimilation // IOP Conference Series: Earth and Environmental Science. - Vol. 211 (2018) 012049 .-doi:10.1088/1755-1315/211/1/012049.
5. Evensen, G. Data assimilation. The ensemble Kalman filter.Berlin Heideberg: Spriger-Verlag, 2009. 307 p.
6. Houtekamer, H.L. Zhang, F. Review of the ensemble Kalman filter for atmospheric data assimilation // Monthly Weather Review. 2016. V. 144. P. 4489-4532.
7. Jazwinski, A. H. Stochastic processes and filtering theory. New York: Academic Press, 1970. 376 p.
8. Feng L., P.I.Palmer, H.Bosch, and S.Dance Estimating surface CO2 fluxes from spaceborne CO2 dry air mole fraction observations using an ensemble Kalman filter // Atmospheric chemistry and physics. 2009. V. 9. P. 2619-2633.
9. Feng L. et al. Consistent regional fluxes of CH4 and CO2 inferred from GOSAT proxy XCH4:XCO2 retrievals, 2010-2014 // Atmospheric chemistry and physics. 2017. V. 17. P. 47814797.
10. Климова Е. Г. Стохастический ансамблевый фильтр Калмана с трансформацией ансамбля возмущений // Сибирский журнал вычислительной математики. - 2019. - Т. 22, № 1. С. 27-40.
© М. В. Платонова, Е. Г. Климова, 2019