УДК 004.9:504:519.6 И.В. КОВАЛЕЦ
АССИМИЛЯЦИЯ МЕТЕОРОЛОГИЧЕСКИХ ДАННЫХ В СИСТЕМАХ ПОДДЕРЖКИ ПРИНЯТИЯ РЕШЕНИЙ________________________________________________________________________________
Abstract: The procedure for the assimilation of the meteorological measurements in the making decision support systems (DSS), dealing with the atmospheric dispersion of the pollutions is described. The results of calculations are mpared with the measurement data in the field experiments. The effectiveness of the developed algorithms and their applicability for the use in the real-time DSSs are shown. Further developments are discussed.
Key words: data assimilation, atmospheric dispersion models, meteorological pre-processor.
Анотація: Розроблений та представлений алгоритм асиміляції даних для систем підтримки прийняття рішень (СПР), які пов'язані з розповсюдженням атмосферних забруднень в оточуючому середовищі. Результати обчислень порівнюються з вимірюваннями у натурних експериментах. Підтверджена ефективність розроблених алгоритмів та їх придатність для використання у СПР, що працюють у реальному часі.
Ключові слова: асиміляція даних, моделі атмосферного перенесення, метеорологічний препроцесор.
Аннотация: Разработан и представлен алгоритм ассимиляции метеорологических данных для систем поддержки принятия решений (СПР), имеющих дело с распространением атмосферных загрязнений в окружающей среде. Результаты вычислений сравниваются с измерительными данными натурных экспериментов. Подтверждена эффективность разработанных алгоритмов и их пригодность для использования в СПР, работающих в режиме реального времени.
Ключевые слова: ассимиляция данных, модели атмосферного переноса, метеорологический
препроцессор.
1. Введение
Модели атмосферного переноса загрязнений (МАП) играют ключевую роль в последовательности моделей систем поддержки принятия решений (СПР), используемых в случае выбросов загрязняющих веществ в атмосферу. Одна из таких систем - система поддержки принятия решений в реальном времени (оп-Ппэ) при авариях на атомных станциях РОДОС [1]. СПР обычно подключены к данным метеорологических измерений и/или к данным численного прогноза погоды (ЧПП). Измерения несут информацию о действительном локальном состоянии атмосферы, в то время как данные ЧПП охватывают значительно большую (по сравнению с данными измерений) пространственно-временную область. Интерфейсом между входными метеорологическими данными, поступающими в СПР и МАП, обычно служат метеорологические препроцессоры (МПП), которые специальным образом интерполируют данные ЧПП или данные измерений на сетку МАП. При этом возникает задача согласования данных ЧПП и данных измерений или задача ассимиляции данных. В моделях ЧПП используется для этой цели четырехмерная (пространственно-временная) ассимиляция (4ДА), например, Калмановская фильтрация. Оперативные СПР не имеют непосредственного доступа к моделям ЧПП, а получают данные ЧПП извне несколько (от 1 до 4) раз в сутки. Таким образом, время анализа ЧПП может отстоять на несколько часов от текущего момента времени. Как следствие, могут возникнуть сравнительно большие отличия между результатами МАП, полученными с использованием только данных метеорологических измерений, и результатами, полученными с использованием данных ЧПП.
Такие отличия наблюдались, в частности, в практике оперативного использования системы РОДОС [1], ранее внедренной в Украине [3]. На рис. 1 представлены результаты использования системы РОДОС во время учений по внеобъектному ядерному аварийному реагированию,
проведенными Национальной атомной энергогенерирующей компанией (НАЭК Энергоатом) в августе 2002 г. на Запорожской АЭС (ЗАЭС). Во время этих учений с помощью системы РОДОС рассчитывалось распространение радионуклидов вследствие гипотетической аварийной ситуации на ЗАЭС. Было проведено два варианта расчетов: в первом варианте входными
метеорологическими данными для системы РОДОС служили данные измерительной станции,
расположенной вблизи источника выброса. Во втором варианте были использованы данные ЧПП.
Как видно из рис. 1, между результатами
обоих вариантов расчетов имеется
существенное отличие. Очевидно, что
метеорологической модели в этом
случае не удалось правильно
рассчитать локальную скорость ветра в
точке выброса. Как следует из рис. 1в,
вектор ветра вблизи источника выброса
резко меняет свое направление
(вследствие влияния на ветровой поток
водоема), так что его направление
чувствительно к расстоянию от
источника. Поэтому использование информации об измерениях ветра вблизи источника для
расчета поля ветра является чрезвычайно важным для точности расчетов.
Одновременное использование данных ЧПП и измерений в диагностической
метеорологической модели, какой является МПП, - задача трехмерной ассимиляции данных (3ДА). Ранее методы 3ДА были разработаны и использовались для вычисления начальных условий
моделей ЧПП, но не использовались в МПП МАП. Поскольку МПП используются для вычисления
на более мелкой пространственно-временной сетке, чем модели ЧПП, с целью более точного описания атмосферного пограничного слоя, основной целью настоящей работы является разработка подходящих методов 3ДА для использования в оперативной версии СПР РОДОС. Будет изучено влияние разработанных методик 3ДА на качество результатов путем сравнения результатов, полученных МПП, с использованием разработанных методик и данных натурных экспериментов по атмосферному переносу загрязняющих веществ.
В данной статье разработаны и представлены методы 3ДА. Валидация этих методов осуществлена путем сравнения с данными метеорологических наблюдений, проведенных во время натурного эксперимента по атмосферному переносу ЕТЕХ [4].
2. Метод решения
Методы 3ДА, представленные в данной работе, были разработаны для метеорологического препроцессора системы РОДОС [2]. МПП системы РОДОС рассчитывает двумерные и трехмерные поля всех метеорологических переменных, необходимых для МАП системы РОДОС. Это
Рис. 1. Пример вычислений
системы РОДОС во время учений по радиационной безопасности на ЗАЭС 22.08.2002: а)
интегральные концентрации 1311 в воздухе, рассчитанные с
использованием только
измерений ветра в точке выброса; б) рассчитанные на основании данных ЧПП; в) линии тока ветра, рассчитанные моделью ЧПП для момента, соответствующего началу выброса
трехмерные поля модуля и направления скорости, температура, влажность, давление и коэффициенты диффузии, двумерные поля осадков, облачности, излучения, турбулентного потока тепла, категории устойчивости, динамическая скорость, масштаба Монина-Обухова, высота слоя перемешивания (пограничного слоя). С использованием предложенных методов 3ДА расчет упомянутых метеорологических переменных осуществляется следующим образом:
1. Расчет первого приближения полей ветра и температуры на расчетной сетке МПП путем
1/ г2 интерполяции соответствующих значений с расчетной сетки модели ЧПП.
2. Ассимиляция измерений скалярных переменных (таких, как температура у поверхности, облачность, интенсивность излучения, осадки).
3. Ассимиляция измерений скорости и направления ветра.
4. Вычисление оставшихся метеорологических переменных, таких, как категория устойчивости, масштаб Монина-Обухова, динамическая скорость, на основе существующих методик в текущей версии МПП системы РОДОС.
Ассимиляция измерений скалярных переменных проводится с использованием процедуры Крессмановского типа объективного анализа, предложенной в [5] и называемой "Итерациями к оптимальному решению" (ИОР). Согласно математической статистике, если есть два вектора наблюдений Г0 и Гв величины Г , измеренных в одном и том же наборе точек 0 < к < К +1,
с К X К ковариационными матрицами О = ({/0 -/' 11)), в = (^/'в -)•(/£-)),
0 < г, 1 < К (здесь квадратные скобки означают осреднение), оптимальной оценкой ГА вектора Г в точках к , 0 < к < К +1 будет [5]:
Если определить поле первого приближения данной метеорологической переменной как Гв с постоянной среднеквадратической ошибкой <гв и вектор измерений той же величины Г0 с постоянной среднеквадратической ошибкой <г0, то, как показано в [5], следующая итерационная процедура сходится в точках измерений к соотношениям, определенным (1):
(1)
к
если элементы ковариационных матриц
и
введенных выше, могут быть
представлены как Ьк = '№{ггк )&В, °к =&01, где I - единичная матрица. В (3) индекс "г" означает номер узла расчетной сетки МПП; индекс "к" - точку наблюдения; м{гш)- весовую функцию, определенную, например, в [5] как
где Я0 - радиус влияния.
Метод ИОР (2) - (4) используется для ассимиляции наблюдений температуры, облачности, интенсивности излучения и осадков с параметрами (У20 и Я0 , определенными в [6, 7].
Ассимиляция наблюдений скорости производится в два этапа: вертикальная экстраполяция наблюдения на высоте 10 м значений скорости на вертикальные уровни МПП и ассимиляция экстраполированных значений скорости (которые затем интерпретируются как “измерения”). Экстраполяция производится в нижнем 200 м слое атмосферы, где обычно применима теория Монина-Обухова, с использованием экстраполяционных формул для модуля и направления горизонтальной составляющей скорости, приведенных в [6]. Эти формулы учитывают отклонение направления горизонтальной составляющей средней скорости за счет влияния силы Кориолиса.
Ассимиляция измеренных значений скорости может производиться как с использованием описанного ИОР метода, так и с использованием метода оптимальной интерполяции. В первом случае относительная ошибка определяется с помощью соотношений, предложенных в [6, 7], которые описывают рост ошибки экстраполированных значений скорости (используемых в качестве псевдоизмерений в процедуре ДА) с ростом высоты экстраполяции. При этом радиус влияния и его зависимость от типа рельефа (равнина, холмистая местность, горы) определены в [7]. Для случая равнины радиус влияния Я = 150 км.
В методе оптимальной интерполяции требуется задать формулы для ковариационных функций и и V компонент скорости. Известно, что для двумерного изотропного бездивергентного поля скорости соответствующие корреляционные функции могут быть приближенно описаны с помощью соотношений [5]:
тии (г-, г1) = ехр(-г / я0)
1 — яп
т (гг, г1 ) = ехр(-г / я0)
т (гг , г1 ) = ехр(-г / я0)
1 - (Уг - У1 )
(5.1)
У
(5.2)
_Т_ (Хг - Х1 )(Уг - У1 )
я г2
(5.3)
Тогда оценки иА {гг) и vA {гг), и и V компонент скорости в данной точке гг находятся
из
системы уравнений:
и.
V
1 (гг ) - и В (гг ) = Кии {И0 - И в )+ ^ {Х 0 - V в ) ; !(г- ) - V в (гг ) = ^и {И0 - И в )+ &0 - Хв ),
(6)
где ив {г(.),vb {г(.) - поля первого приближения х и у компонент горизонтальной составляющей скорости, определенные на вычислительной сетке МПП на 1-м шаге, и0 ,у0 - векторы длины К результатов наблюдений (где К - суммарное число точек измерений скорости), ие,уе - векторы результатов первого приближения, интерполированные с вычислительной сетки МПП в точки
г
измерений, векторы весов Щ_ии,Щ_т длины К определяются из следующей системы уравнений [5]:
Ь:+со к+[£+со, к=с(г); кВ,+со к+\ы+с. =с., (г,);
к:,+со К+[£+со.к, = св„(г);
Ісі,+со к,+Ь+со к=с" (г). (7)
В (7) сі,,сі,сі, - это ковариационные матрицы размера ошибок поля первого приближения в точках измерений размера К XК , например, сі, =е е Г) . Здесь еив - вектор длины К
ошибок поля первого приближения в точках измерений. Матрицы с,.,С?с,, определяются
аналогично. Векторы Свии,Свп,Свт длины К состоят из ковариаций ошибок поля первого приближения в точках измерений и ошибок поля первого приближения в точках анализа /:
сВии (г) = ЄЄ1(Г-)) . Здесь еив (г) - ошибка поля первого приближения в точке анализа /.
Векторы Свии (/",),С„ (г,) определяются аналогично. Система уравнений (7) решается численно для каждого узла / МПП.
Ковариационная функция СВи (/,/) (индексы " и", " В" использованы для примера)
определяется с помощью корреляционной функции тии (//) соотношением
св: (г,, г,) = т (г,, т, )аВ (г №В (г), (8)
где аВ(г,),&В(/) - среднеквадратические отклонения (в данном случае "и" - компоненты скорости поля первого приближения) от истинного значения в точках /, г. соответственно. Поэтому для
определения ковариационных матриц в системе уравнений (7) необходимо определение
соответствующих значений среднеквадратических отклонений, а также корреляционных функций ошибок измерений, тогда как корреляционные функции ошибок полей первого приближения уже определены соотношениями (5.1, 5.2, 5.3).
Общепринятой в метеорологической практике [8] является гипотеза об отсутствии
корреляций ошибок измерений, поэтому т°ии (Г,/) = №°п (/, /і ) = ^ї . Дополнительно можно
предположить, что Ци° (/,/) = 0. Предположим также, что выполнено равенство среднеквадратических ошибок наблюдений обеих компонент скорости на данном вертикальном уровне МПП а(/) =а°(/) = а°(/). Напомним, что “наблюдениями” считаются и значения, интерполированные или экстраполированные с более низких измерительных уровней.
Среднеквадратическая ошибка поля первого приближения предполагается постоянной для данного вертикального уровня, на котором производится ассимиляция: аВ(/) = аВ(/) = ав. После
этих предположений система уравнений (7) может быть представлена в виде
B I в B
Л,т + 1кии + Лш,к«у = Лии(г,);
^В_ +1 }1иу = Л*у (г,);
^уу
и I и и
Л + Ікуи + Л куу = лиу (г,);
^ши -*— =иу— ±иу^1!
В
Л кш +
=иу —
В і в
Луу + 1куу =Луу(г,), (9)
^=УУ
где матрицы и векторы, входящие в (10), определяются следующим образом:
(л* )и = тВи (г ? гі) (аВ / (аО \),1 £ к ?1 £ К; (л*. (г) )к = тій (г ? г) (а / (аО )к),1 £ к £ К?
(10)
Ч) У(У0 ^
В В в в
где ^- количество наблюдений. Определения матриц Л ,Л и векторов Л ,Л аналогичны.
= ЫУ =УУ —иу—уу
Таким образом, система уравнений (9) содержит, кроме корреляционных функций, только относительные ошибки полей первого приближения и наблюдений. Относительная ошибка прогнозов и наблюдений компонент приземной скорости может быть взята из литературы [8] и
равна приблизительно У2/у0 »10. Однако в процедуре ассимиляции используются не только
“истинные” измерения скорости, но и значения, экстраполированные с нижележащих уровней измерений на уровни МПП. Ясно, что для таких значений их ошибка будет определяться не только неточностью наблюдений, но и ошибкой, вносимой экстраполяцией. Для параметризации этой ошибки был использован подход, предложенный авторами МПП ОДЬМБТ [9], в котором поле ветра Га, определенное в узлах МПП, рассчитывалось как взвешенная сумма полей первого приближения Гв и поля, полученного путем интерполяции/экстраполяции измерений во все узлы МПП: Га = Ж0 • {0 +(1 - Ж0 )• {в .
Отметим, что использование взвешенной суммы существенно отличается от предложенного в данной работе метода ассимиляции данных. Для определения Г0 в этом методе
необходимо экстраполировать измерения не только по вертикали, но и по горизонтали, что может приводить к большим неточностям. Однако между весовым коэффициентом, используемым в [9], и
относительной ошибкой Уд/У20 может быть установлено соответствие. Отметим, что в случае одноточечных измерений случайных величин /0 и /в одной и той же величины / и двух различных среднеквадратических ошибок измерений у0 и ув , "оптимальной" оценкой /А величины / будет [5]: fA =((у0/ув)2 0 + fB)/((у0/ув) 2 +1). Из этого следует оценка:
аВ
ао=w0/(1-). (11)
^Cherbourg {VaIognes|
pBA0^LLl . AMCAL^;Vt! '
LOU ARGA]j_ « KERPERT;
Г Mill*
- 'R<pSTRENEN
• • У
PLOURAY • -
НоУмаЫ-HffiMcFp —
^MONTERFIL
;РГОЕКМЕМ
Domain of calculations
© Center point, sodar location a. Observation points (used in D.A. procedure)^,
• Observation points (for comparison)
■i- NWP nodes
Grid 5x5 km of met. preprocessor •
ARGENT AN . . . Alencore-
|N|
A
Рис. 2. Вычислительная область МПП с точками
метеорологических измерений в БД эксперимента ЕТЕХ: треугольники - точки, измерения в которых использовались в ДА; кружки - точки, измерения в которых использовались для сравнения; крестики - узлы модели ЧПП; сплошные линии - сетка МПП
С учетом (11) соотношения, предложенного в [9] для вычисления
весового коэффициента
Wn
учитывающего возрастания
погрешности за счет вертикальной экстраполяции измерений, могут быть использованы для
определения относительной ошибки у2в/у2 . Некоторые модификации
этих соотношений в деталях описаны в [7].
Следует отметить, что в разработанном методе измерения следующих величин могут быть ассимилированы значения модуля и направления скорости, температуры, облачности, суммарной
интенсивности излучения и осадков. Все остальные метеорологические величины, такие, как высота
пограничного слоя, категория устойчивости, масштаб Монина-Обухова, для которых обычно прямые измерения не проводятся, подвергаются косвенному влиянию ассимиляции измерений первой группы переменных благодаря их взаимозависимости через параметризации. Используемые параметризации детально описаны в [6].
3. Сравнения с результатами измерений в эксперименте ЕТЕХ
Рассчитанные метеорологические переменные сравнивались с измерениями, проведенными в натурном эксперименте ЕТЕХ [4]. База данных измерений ЕТЕХ детально описана в [10] и открыта к доступу в Интернете. Она содержит синоптические измерения наземных станций (средние значения температуры, ветра, облачности и осадков за каждые три часа), собранные во время экспериментов с распространением трасера (2 периода в течение трех дней каждый, начиная с 23 октября 1995 г. и 14 ноября 1995 г. соответственно). Измерения проводились в квадрате 400x400 км вокруг точки выброса (Монтерфиль, Франция, рис. 2). Область, изображенная на рис. 2, совпадала с вычислительной областью МПП. Сетка МПП имела горизонтальный пространственный шаг 5x5 км (сплошные линии на рис. 2). Вертикальное разрешение сетки было переменным. В нижнем 1 км слое атмосферы было расположено 4 вычислительных уровня МПП, а выше в 10 км слое - остальные шесть уровней. Данные восьми наземных метеорологических станций, показанных треугольниками на рис. 2, использовались в процедуре ассимиляции данных. Данные 71 наземных станций, отмеченные на рис. 2 кружочками, использовались для сравнения расчетов с наблюдениями. Также в базе данных эксперимента ЕТЕХ хранятся данные анализа
метеорологической модели ЧПП ECMWF [10], то есть метеорологических полей, рассчитанных метеорологической моделью с учетом уже имеющихся на момент произведения анализа измерений. Эти метеорологические поля определены на сетке, обозначенной крестиками на рис. 2. От ECMWF были получены также данные ЧПП, произведенного для указанных выше периодов, которые рассчитывались на той же сетке, что и анализы. Были проведены четыре варианта расчетов: первый - используя данные только ЧПП; второй - с использованием ЧПП и с ассимиляцией данных измерений; третий - с использованием данных только анализа ECMWF; четвертый - с использованием данных анализа и с ассимиляцией данных измерений. Повторная ассимиляция данных измерений в четвертом случае необходима в связи с интерполяцией результатов анализа с грубой сетки метеорологической модели на мелкую сетку МПП [9].
В базе данных ЕТЕХ также находятся данные измерений анемометра и содара, расположенных в центре области (Монтерфиль), которыми измерялись турбулентные потоки тепла и вертикального потока импульса (динамической скорости), а также вертикальные профили скорости и направления ветра, которые использовались для сравнений.
Четыре основных параметра характеризуют близость рассчитанных и наблюдаемых полей ветра. Это среднеквадратическое и систематическое отклонения модуля скорости ветра (гтэи, Ывви), а также среднеквадратическое и систематическое отклонения направления ветра. Отметим, что ветер в области направлен приблизительно в одну сторону в один момент времени, но в различные моменты времени средние по области направления ветра могут существенно отличаться. Поэтому определение систематического отклонения для направления ветра было несколько модифицировано:
rmsf ■
Ык N
к і
'(Ык • Мі), f = и, ё ;
(12)
biasu
Ык Мі
II (< - ик)
к і
'(Ык • Мі), Ыа8ё :
Ык
і
Мі
I (и - и)
'(Ык • N),
где "Ык" означает количество временных горизонтов, для которых доступны данные ЧПП и измерений, а "№" - количество измерительных станций, данные которых использовались для сравнения. Для первого эксперимента ЕТЕХ Ык = 9 , для второго - Ык = 13 , № = 71 для обоих случаев.
Таблица 1. Статистические характеристики ошибок рассчитанных метеорологических полей в экспериментах ЕТЕХ1 и ЕТЕХ2 с использованием и без ассимиляции данных
к
Переменная Первое приближение -ЧПП Ассимиляция Первое приближение -анализ Ассимиляция
гтэи, м/с, ЕТЕХ1 2,57 2,44 3,19 2,52
Ывэи, м/с, ЕТЕХ1 0,67 -0,08 1,06 0,25
гтэб, град, ЕТЕХ1 33 33 31 31
Ывэб, град, ЕТЕХ1 10,3 5,5 9,55 4,26
гтэи, м/с, ЕТЕХ2 2,21 1,94 — —
Ывэи, м/с, ЕТЕХ2 0,5 -0,18 — —
гтэб, град, ЕТЕХ2 56 48 - -
Ывэб, град, ЕТЕХ2 18,7 5,9 - -
Таблица 2. Результаты 4ДА, полученные различными метеорологическими моделями с использованием четырехмерной ассимиляции данных [13]
Модель, шаг сетки Без ассимиляции гтви (м/с), biasd (град) С ассимиляцией rmsu (м/с), biasd (град)
RAMS, dx = 1,32 км 3,63, - 0,93, -
RAMS, dx = 0,33 км 2,58, - 2,23, -
MM5, dx = 4 км 2,63, -8,4 1,72, -9,2
Величины гтзы,гтзй,Ыазы,Ыазй, рассчитанные по полю первого
приближения и с использованием ассимиляции данных, представлены для каждого из экспериментов ЕТЕХ в табл. 1. В первом варианте поле первого
приближения рассчитывалось на основе данных ЧПП. Во втором - на основе данных анализов метеорологической модели. В обоих случаях сравнение показывает определенное улучшение полей ветра,
рассчитанных с использованием ассимиляции данных. Для случая первого эксперимента ЕТЕХ значения гтзы и гтзй почти не изменяются под воздействием ассимиляции данных, тогда как для второго эксперимента ЕТЕХ улучшение параметров гтзы, гтзй под воздействием ассимиляции данных очевидно. Это может быть объяснено тем фактом, что, согласно [4], первый эксперимент ЕТЕХ характеризуется погодными условиями, более точно описываемыми параметризациями ЧПП модели, чем погодные условия второго эксперимента. Систематические отклонения Ыазы и Ыазй в большей степени подвергаются воздействию ассимиляции данных и улучшаются во всех представленных в табл. 1 случаях.
_1-----------,----------------------------------- 0 1 і —
О 100000 О 100000
х, т х, т
Рис. 3. Поля ветра на высоте 10 м, рассчитанные МПП для 25.10.95 г., 00:00иТС: а) поле первого
приближения, рассчитанного на основе ЧПП; б) поле ветра, рассчитанное с ассимиляцией и радиусом влияния Я0 = 50 км; в) наблюдаемые векторы скорости; г) поле ветра, рассчитанное с ассимиляцией и радиусом влияния Я0 = 150 кт
Особо следует отметить, что когда поле первого приближения в эксперименте ЕТЕХ1
рассчитывалось на основании данных анализа, среднеквадратические
отклонения абсолютной величины скорости гтэи существенно улучшались под влиянием
ассимиляции данных. Это говорит о том, что даже если ассимиляция данных производилась моделью ЧПП, при интерполяции ее результатов на более мелкую сетку необходимо
23.5 24 24.5 25 25.5 26 26.5
повторно производить ассимиляцию,
время, дни
поскольку мелкомасштабные
компоненты движений фильтруются
моделью с грубым пространственным
разрешением.
Сравним результаты 3ДА,
представленные в табл. 1, с
результатами 4ДА из [11],
представленными в табл. 2. Видно, что
качество улучшения статистических
характеристик полей ветра в случае
использования 3ДА согласуется с
улучшением тех же характеристик при время, дни ’ ’ г г г
использовании 4ДА, которые
Рис. 4. Зависимость от времени: а) скорости используются в прогностических м°делях.
трения; б) потока тепла. Точки - измерения; Поскольку 3ДА методы гораздо более сплошные линии - расчет с ассимиляцией;
пунктирные - без ассимиляции эффективны в смысле вычислительного
времени, чем 4ДА методы, эти результаты являются хорошим обоснованием развития 3ДА методов для их использования в СПР.
Примеры метеорологических полей, рассчитанных МПП для 25 октября 1995 г., 00:00 ч, показаны на рис. 3. Как видно из рисунка, поле первого приближения, рассчитанное с использованием только данных ЧПП, гораздо более гладко, чем измеренное. Это естественно, поскольку метеорологическая модель ЧПП фильтрует коротковолновые возмущения за счет грубого шага вычислительной сетки (рис. 2). Однако использование ассимиляции данных, как видно из рис. 3, в некоторой степени исправляет распределение ветра, хотя для того, чтобы добиться более полной согласованности рассчитанного и измеренного полей ветра, необходимо объединение представленных в настоящей работе статистических методов ассимиляции данных с вариационными методами, описанными в [12].
На рис. 4 представлены сравнения рассчитанных МПП и измеренных анемометром, расположенным в центре области (Монтерфиль, рис. 2), значений
вертикального турбулентного потока импульса
(скорости трения и* = (и'м>)) и потока тепла
Рис. 5. Вертикальные профили скорости (а1, б1) и направления ветра (а2, б2) для различных моментов времени: рассчитанные МПП с ассимиляцией (■, сплошная линия); на основе только ЧПП (▲, пунктирная линия); измеренные содаром (♦)
д = {и'Т') во время первого эксперимента
ЕТЕХ. Как было отмечено в предыдущем разделе, ассимиляция данных измерений скорости, температуры и облачности может лишь косвенно влиять на значения величин скорости трения и турбулентного потока тепла.
Как видно из рис. 4, использование ассимиляции данных существенно улучшило значения скорости трения. В то же время значения турбулентного потока тепла
остаются почти неизменными в результате ассимиляции. Это связано с тем, что уже в первом приближении значения потока тепла достаточно близки к измеренным.
На рис. 5 показано сравнение вертикальных профилей скорости и направления ветра, рассчитанных МПП с измерениями содаром. Как видно из представленных на рисунке результатов, вертикальные профили скорости, рассчитанные с использованием ассимиляции, лучше согласуются с измеренными, чем вертикальные профили первого приближения. Это оправдывает вертикальную экстраполяцию измеренных приземных значений скорости и направления ветра на уровни МПП и использование экстраполированных значений в качестве наблюдений. Поскольку экстраполяция приземных значений скорости производится в нижних 200 м атмосферы, улучшения можно добиться только в слое соответствующей толщины, что и видно на рис. 5. Имеющиеся расхождения на рис. 5 в направлении ветра объясняются расхождением экспериментальных данных между собой: в некоторых случаях направление ветра, измеренное содаром на 30°, отличается от направления ветра близлежащей метеорологической станции за счет пульсации метеорологических полей и различного времени осреднения измерений содара (30 мин.) и синоптических станций (3 часа).
Особо следует отметить, что в настоящей работе никакие измерения, используемые для сравнений с расчетами, не использовались в процедуре ассимиляции. В случае использования измерений содаром вертикальных профилей скорости в процедуре ассимиляции можно было бы добиться еще более лучшего согласования рассчитанных вертикальных профилей и измеренных в более глубоком (чем 200м) слое атмосферы.
4. Выводы
В данной работе представлена процедура ассимиляции данных, разработанная для метеорологического препроцессора системы поддержки принятия решений при авариях, связанных с выбросами загрязняющих веществ в атмосферу. Решена задача одновременного использования данных ЧПП и данных метеорологических измерений на местности для расчета метеорологических полей на расчетной сетке МПП. Для этой цели были адаптированы методы трехмерной ассимиляции данных (оптимальная интерполяция, [5, 8]), ранее разработанные для инициализации метеорологических моделей ЧПП. Новым в данном подходе является сочетание методов 3ДА с метеорологическими параметризациями пограничного слоя атмосферы, уменьшение числа независимых параметров задачи (требуются только относительные, но не абсолютные ошибки метеорологических величин) за счет некоторых предположений о структуре ошибок, сочетание представленных в статье методик 3ДА с ранее разработанными методиками для определения “весового коэффициента” в МПП CALMET [9]. Однако наиболее важным результатом работы является убедительное подтверждение эффективности разработанных методик 3ДА за счет сравнения результатов расчетов, проведенных с использованием и без использования ассимиляции данных с метеорологическими измерениями в экспериментах ETEX [4]. Показано, что уровень улучшения статистических характеристик ошибок полей ветра за счет применения 3ДА в МПП приблизительно совпадает с уровнем улучшения за счет применения 4ДА в метеорологических моделях ЧПП [11]. Как подтверждено результатами расчетов, использование 3ДА позволяет в некоторой степени “реставрировать” мелкомасштабные компоненты движений, которые фильтруются метеорологической моделью ЧПП за счет ассимиляции данных измерений. Особенно интересным является тот факт, что даже в случае, когда первое приближение определялось из данных анализа метеорологической модели (т.е. ассимиляция данных измерений уже имела место), повторная ассимиляция данных все равно приводит к улучшению результатов расчетов, что связано с различным горизонтальным разрешением метеорологической модели ЧПП и МПП. Это говорит о том, что возможности метода не исчерпаны. В дальнейшем следует сочетать разработанные в этой работе статистические методы 3ДА с вариационными методами для адекватного описания влияния сложной подстилающей поверхности на метеорологические поля. В настоящее время представленные в статье методики внедряются в СППР РОДОС в рамках проекта РП6 Европейской Комиссии EURANOS.
Благодарности
Представленная работа полностью финансировалась за счет фанта для молодых ученых из третьих стран программы EURATOM Европейской Комиссии. Автор благодарит всех сотрудников лаборатории моделирования окружающей среды Национального научного центра “DEMOKRITOS” (Греция) и особенно докторов Спироса Андронопулоса, Александроса Венецаноса, профессора Джона Бартиса за обсуждение представленных в работе результатов и ценные рекомендации.
СПИСОК ЛИТЕРАТУРЫ
1. Raskob W., Ehrhardt J. The RODOS System: Decision Support For Nuclear Off-Site Emergency Management In Europe, RODOS report GEN_TN99_02. - 1999. - 112 p. http://www.rodos.fzk.de.
2. Andronopoulos S. et al. FILMAKER: a computer model fro pre-processing of the weather data. RODOS report WG2_TN98_08. - 1998. - 82 p. http://www.rodos.fzk.de.
3. Zheleznyak M., Kovalets I., Sorokin Y., Dvorzhak A., Kushchan A., Bogorad S., Shlyahtun N. Implemetation in Ukraine of the Rodos system // Proc. of International Symposium on “Off-site Nuclear Emergency Management -Capabilities and Challenges”. - Salzburg, Austria. - 2003.
4. Dop van H., Addis R. et al. ETEX: a European tracer experiment; observations, dispersion modelling and emergency response // Atmospheric Environment. - 1998. - Vol. 32. - P. 4089 - 4094.
5. Daley R. Atmospheric Data Analysis. - Cambridge University Press. - 1991. - 621 p.
6. Kovalets I. Introduction of data assimilation techniques in the meteorological pre-processor of the RODOS system, Final Report, EUROATOM grant № FIKR-CT-2002-50024. - 2003.
7. Kovalets I., Andronopoulos S., Gounaris N., Kushchan A. Introduction of data assimilation procedures in the meteorological pre-processor of atmospheric dispersion models used in emergency response systems // Atmospheric Environment. - 2004. - Vol. 38/3. - P. 457 - 467.
8. Гандин Л. Статистические методы интерпретации метеорологических данных. - Ленинград: Гидрометеоиздат, 1976. - 221 с.
9. Scire J. A users guide for the CALMET meteorological model // Earth Tech, Inc. - 2000. http://www.src.com/calpuff/calpuff1.htm.
10. Straume A.G. et al. Meteorological observations collected during the European Tracer Experiment (ETEX). -Report. - Joint Research Center, Italy, 1997. http://rem.jrc.cec.eu.int/etex/.
11. Seaman N. Meteorological modelling for air-quality assessments // Atmospheric Environment. - 2000. - Vol. 34. P. 2231 - 2259.
12. Fisher B.E.A. et al. COST Action 710, Final report, Harmonisation of the pre-processing of meteorological data for atmospheric dispersion models. - L-2985 European Commission, Luxembourg, 1998. - EUR 18195 EN (ISBN 92-828-3302-X).