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

Вычислительные технологии

Том 13, Специальный выпуск 3, 2008

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

А. В. Пененко Институт вычислительной математики и математической геофизики СО РАН, Новосибирск, Россия е-таП: aleks@ommgp.sscc.ru

This work addresses a variational approach applied to the localization of the point-wise sources of atmospheric substances. Our modification is based on an incorporation of both data of concentration measurements of the pointwise substances and a certain a priori data on their capacities. Using a priori data and the weight functions, the domains which have to contain sources are selected. The mentioned weight functions are constructed with the help of the solutions of adjoint problems formulated for a model describing advection-diffusion of atmospheric substances. Localization is sought as an intersection of the selected domains.

Введение

Задача об обнаружении источников примесей в атмосфере рассматривалась, например, в работах [1-6] Наиболее теоретически изученным в данной теме является вопрос построения связи между данными наблюдений и неизвестными источниками с помощью сопряженных задач и весовых функций (см., например, [1-3]). Эта связь позволяет свести задачу об обнаружении источников к задаче о восстановлении функции по ее проекциям. Сложность полученной задачи состоит в том, что множество функций, на которое производится проектирование функции-источника, не является полным, а следовательно, задача, вообще говоря, имеет неединственное решение. Одним из подходов, предлагаемых в [3, 4], является подход, при котором решение задачи ищется как линейная комбинация весовых функций. Однако в случае точечных источников вид весовых функций и вид источников различаются. Следовательно, для обнаружения таких источников требуется предложить способ интерпретации получаемых распределенных решений, позволяющий на их основе получить решения в виде систем точечных источников. Подход, развиваемый в данной работе, состоит в том, чтобы, используя априорную информацию о мощности источника, выделить на носителях весовых функций области, в которых могут располагаться источники. Тогда локализация источника ищется как пересечение таких областей.

* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (грант № 07-05-00673) и контракта Европейской комиссии (грант № 013427).

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

1. Задача локализации источника пассивной примеси в атмосфере

1.1. Постановка задачи

Предположим, что имеется набор из £ вещественных значений {1%|£ € 5} — результатов измерений в разных точках одного и того же пространственно-временного поля концентрации некоторой пассивной примеси, представленного вещественной функцией р, определенной на некоторой сеточной области:

Шт = {(Кв],Пк^и)1г = 1..Мх,] = 1..Му,к = 1..Мг,и = 1...М} .

Сеточная область шт получается в результате введения на шаровом слое

От = о х [о,т] = {(\,е,п,Ш е [Лш,\е] ,в е [в3,ви] е [1, о] ,г е [о,т]},

равномерной сетки с шагами АЛ, АО, Ац, т. Здесь Л,О,п — долгота, широта и высота соответственно, причем п = 0 отвечает верхней границе атмосферы, а п = 1 — нижней. Кроме того, для каждого из рассмотренных шагов введем величины бИ,1 = бНМК = АН/2 и бН2 = ... = = АН. Определим скалярное произведение сеточных функций:

(и, V) = ^ и%к]бъ,

г],к,п

где Б]к — элемент метрики на шаровом слое. Пространство сеточных функций на шт, снабженное скалярным произведением, становится пространством Ь2(шт). Используя это скалярное произведение и вводя дискретные дельта-функции б(х — х%) и б (г — ):

1/в% ,х = х%, — ^_) = (1/бг,г = , о, х = х%, ( %) \о, г = ,

запишем выражение

1% =

где х% и — это соответственно время и место проведения измерения

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

Ьр = Q.

Далее, предполагая, что источники были точечными и располагались на поверхности земли, мы можем записать предположение, что поле наземных источников Q имеет структуру

Мд

Q = ^2 дь б(х — хь )б(г — и),

1=1

где Мд — число точечных источников; дь — мощность источника с номером /; (хь— его координаты.

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

б(х — х%) =

1.2. Метод локализации

По определению сопряженного оператора L* к данному линейному оператору [1], имеем

(L<p, Ф€> = (<р^*Ф<>,

откуда на основании равенств (1),(2) получим набор соотношений на функцию источника:

(q, ф€ > = е s.

Система функций |£ е S} в общем случае не является полной, поэтому задача о нахождении Q по и может не иметь единственного решения, а следовательно, при интерпретации данных будем руководствоваться критерием экономичности описания системы источников, достаточной для объяснения данных измерений (потому что если, например, все измерения производятся в один момент времени, то самым простым решением будет предположить, что все источники находились непосредственно в точках измерений, однако такая интерпретация данных может быть неэкономичной). Заметим, что обычно интерпретацию данных измерения производят в терминах весовой функции [1, 2, 4], т.е. пытаются представить неизвестный источник в виде линейной комбинации весовых функций. Мы же в данной работе попытаемся интерпретировать данные измерений в терминах набора точечных источников.

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

= ((x, t) е шт|Ф(Х, t) > а} .

На основании знания о структуре поля источников можно сделать следующие утверждения.

Лемма 1. Пусть I — некоторый результат измерения и Ф — соответствующая весовая функция. Если известно, что при

(Q, i>z < q

внутри области Z нет точечных источников, то тогда, внутри области S^, где ai = I¡q, нет источников. Здесь (., .>z обозначает взятие скалярного произведения, по Z 1 1

Условие утверждения выполняется, если, например,

q < min [qi (6(. - xi)6(. - U), 1>],

leNq

q

Лемма 2. Пусть ск2 > . . , I > 0, тогда в области гарантированно имеется

(Q, 1>

источник.

Следовательно, множество ОФ = Бф2 \Бхф1 гарантированно содержит источник. Просуммировав индикаторы этих множеств, получим функцию локализации:

Ф = ^ 1п(ОФ

ч-

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

2. Построение системы весовых функций

2.1. Дифференциальная модель

Уточним оператор модели. В качестве оператора Ь : Ь2(шт) ^ Ь2(шт) рассмотрим дискретизацию модели адвекции-диффузии пассивной примеси в атмосфере, определенной От

дф

— + щгд.(\р> — сИу(^^гас1</?) = 0, (х, € От,

д^ дф

К---й„р = 0, х е дП/дПсимв, К---й„р = <^,х е дПСиМв, = 0.

дп дп

Здесь дО — граница области О; дОскмо — приземная граница области О; и — скорость ветра; К — коэффициент диффузии. Коэффициенты модели вычисляются как функции от заданных скоростей ветра, температуры и приземного давления. По вертикали введена гибридная система координат.

2.2. Дискретизация модели

Дискретизация модели производится методом конечных объемов с использованием теоремы Остроградского—Гаусса. При этом основные агрегаты модели дискретизуются соотношениями типа

дф д (тф) дт ф — фп-1

др\П (\тиМц - ^г+г) + (фг - ¥>г-1 У

ти

а соэ 9 д\) I acosвj5\i

тд дф

п

^ дФ \

цк

та еов 9д\ х а еов 9д\

[Кхт\{ ,1-—— - [Кхт^_ 1-——-

г+2 асоявуААг г 2 асоэ^-ааг-!

а еов вj8Хг

\ )

jk

т

разностной схемы уравнение неразрывности; а — радиус Земли; [.]+, [.]- — операции взятия положительной и отрицательной частей соответственно. В результате дискретизации получается дискретная модель вида

п

т" — т"-1 Оп

тп~1 ^-^-+ Г>га + Сп<рп = тпН- 4

т —8ц 1

Vп = V" + V" + V", Сп = Сп + С" + С",

т0 = о,

где Н — оператор продолжения поверхностных источников на всю область; Vп —

Сп

разностная схема расщепляется по пространственным координатам:

,ПП— д _ -П—1 глп

тп~^-£-+ + = тпН—~-—;

т -8ц1

ш™"1^--+ + = 0;

шп-1^га-^га 3 + уп^ + ¿^п = 0 ^

т

2.3. Весовые функции

Соответствующая (4) — (6) сопряженная модель выглядит следующим образом:

т™"1^- + [С + С':г Г = тпУ--

т т 8£„

тт

Здесь [.]* — операция взятия сопряженного оператора относительно заданного скалярного произведения. Заметим, что в (3)-(6) п = 1, Л^, а в (7)-(9) п = Л^, 1.

Лемма 3. Если функция ф удовлетворяет уравнениям (7)-(9) и условиям

кт Ф{т+1)'1 = ^

а функция т ~ уравнениям, (3) — (6) то они связаны, соотношением,

тй = (Я, ф> ,

где

т 2

Фга = --[шгаЯ] >га"з,

г, к — координаты, точек наблюдений.

49

А. В. Паиаико

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

В приводимом ниже численном эксперименте задаются источники примеси и с номощыо решения прямой задачи вычисляются предполагаемые данные измерения. Метеорологические данные, необходимые для вычисления коэффициентов модели, предоставлены Датским королевским метеорологическим институтом, они соответствуют периоду 23 октября 1995 г. — 27 октября 1995 г. (86 ч). Шаг но времени 1 ч. Рассматриваемая пространственная область — Европа. Число шагов сетки но долготе и широте 100. Число шагов но высоте 40. В качестве системы измерений взята система, применявшаяся в эксперименте etexl_vl. 1.960505, насчитывающая 168 наблюдательных комплексов.

а б

Поло приземной концентрации примеси в зависимости от времени (а);

функция локализации

На протяжении 10 ч действует источник постоянной мощности, через 10 ч после окончания его работы всей измерительной системой производятся измерения концентрации примеси, и па основе этих данных восстанавливается расположение источника. На рисунке а приведена концентрация примеси па поверхности земли в зависимости от времени, па рисунке б — получающаяся функция локализации. Из ее графика видно, что кроме точек расположения действительного источника функция локализации достигает больших значений еще в ряде точек. Для удаления подобных "фантомных" источников необходимо привлекать дополнительные данные.

Заключение

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

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

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

[11 марчук Г.и. Сопряженные уравнения и анализ сложных систем. м, Наука, 1992.

[2] Pudykiewicz J.A. Application of adjoint tracer transport equation for evaluating source parameters // Atmos. Environ. 1998. Vol. 32, N 11. P. 3039-3050.

[3] Penenko V., Baklanov A., Tsvetova E. Methods of sensitivity theory and inverse modelling for estimation of source parameters // FGCS. 2002. Vol. 18. P. 661-671.

[4] Issartel J.P. Rebuilding source of linear tracers after atmospheric concentration measurements // Atmos. Chem. Phys. Discuss. 2003. N 3. P. 3173-3203.

[5] Пененко В.В. Методы обратного моделирования и оценки экологических рисков от антропогенных воздействий // Обозрение прикладной и промышленной математики. 2003. Т. 10, № 1. С. 26-38.

[6] Пененко В.В. Методы численного моделирования атмосферных процессов. Л.: Гидроме-теоиздат, 1981.

Поступила в редакцию 21 февраля 2008 г.