CC BY
12
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ МОРСКИХ СИСТЕМ
УДК 551.46(262.5)
Идентификация мощности источника загрязнения в Казантипском заливе на основе применения вариационного алгоритма
© 2015 В.С. Кочергин, С.В. Кочергин
Морской гидрофизический институт РАН, Севастополь, Россия E-mail: vskocher@gmail.com
Поступила в редакцию 09.06.2014 г. После доработки 20.06.2014 г.
Рассматривается модель переноса пассивной примеси в Азовском море. На ее основе реализован вариационный алгоритм идентификации мощности источника загрязнения, в том числе переменной по пространству. На тестовом примере показана работоспособность алгоритма поиска оптимального распределения по пространству мощности источника, согласованного с данными измерений. Тестовые расчеты проведены для акватории Казантипского залива при восточном ветровом воздействии.
Ключевые слова: идентификация мощности источника, вариационный алгоритм, функционал невязок, поле концентрации, усвоение данных измерений, Азовское море.
Введение. При изучении динамики распространения примесей необходимо использование как современных математических моделей [1], так и методов усвоения данных измерений [2, 3], которые позволяют идентифицировать входные параметры модели. Алгоритмы усвоения данных измерений основаны, как привило, на минимизации квадратичного функционала качества прогноза, характеризующего отклонения модельного решения от данных измерений. При этом модель переноса пассивной примеси выступает в качестве ограничений на вариации входных параметров. В работе [4] всесторонне рассмотрен вариационный алгоритм идентификации мощности источника для двумерной модели, показана его работоспособность при наличии данных измерений на периферии пятна загрязнения в случае действия точечного источника постоянной и переменной по времени мощности. В настоящей работе такой подход применен для трехмерной модели переноса пассивной примеси в Азовском море. Рассматривается задача идентификации переменной по пространству мощности постоянно действующего источника.
Вариационный алгоритм усвоения данных измерений. Рассмотрим модель переноса пассивной примеси в а -координатах
dDC dDUC dDVC dWC д . dDC д . dDC д K DC
-+-+-+-= — AH-+ — AH-+---(1)
dt dx dy да dx dx dy dy да D да
с условиями на боковых границах
Г :DC = 0, (2)
дn
краевыми условиями на поверхности и на дне
а = 0: ^ = Qss(x - xc У - Ус) ' а = -1: ^ = QBS(x - x0,y - y0) , (3) да да
или
а = 0: да = Qs(x,y) , а = -1: ^ = QB{x,y) (4)
да да
и начальными данными
C (x, у,а,0) = Co (x, y, а), (5)
где t - время; x0, y0- координаты точечного источника; D - динамическая глубина; C - концентрация примеси; C0 - начальная концентрация примеси; U,V,W - компоненты поля скорости; AH и K - коэффициенты горизонтальной и вертикальной турбулентной диффузии соответственно; n - нормаль к боковой границе.
В условиях (3) QS, QB = const, а в (4) QS (x, y), QB (x, y) - переменные мощности источника на поверхности и на дне соответственно.
Задача усвоения данных измерений Сизм состоит в минимизации квадратичного функционала
I0 = 2 (P(C - СИзм ), P(C - СИзм DMt , (6)
где Mt = M х[0,Г ]; M - область интегрирования; P - оператор расширения нулями функций невязок, заданных на множестве точек измерений; скалярное произведение определяется стандартным способом. Минимизация (6) с краевыми условиями (4) эквивалентна поиску экстремума следующего функционала:
I = 10
(dDC dDUC dDVC dWC д . dDC
--Ah
v
дt дx дy да дx дx
д . дDC д К дС ^ (дС Л (, „ \ --Ан-----,С I +1—,С I + С -С0,С + (7)
ду н ду да D да )М { дп )Г V ]'=0 1'=0'м
+[£ - е. Xу) с*)о+{£ -еX У) С'
Записывая вариацию функционала (7) и интегрируя по частям с учетом краевых условий и аналога уравнения неразрывности в а -координатах
дD дБи дDV дW Л
-+-+-+-= 0, (8)
д1 дх ду да
получим
SI = (х, у), C % + (х, у), C, (9)
где C: - множители Лагранжа, которые выбираются из решения следующей сопряженной задачи:
дDC: дDUC дDVC* дWC: ^ д л дC
--------^ан~
дt дх ду да дх дх
_ D—а. да-± £ дс! = - р(с - с,„,),
н ^ ^ ^ V изм / '
ду ду да D да
(10)
Г : — = 0, а = 0: — = 0, а = -1: — = 0, (11)
дп да да
t = Т: е = 0. (12)
В случае, когда данные измерений имеются на конечный момент времени Т , в (10) задаем правую часть, равной нулю, а при t = Т в (12) используется условие
t = Т: С: = Р(Сизм - С). (13)
Из условия стационарности функционала (7) 81 = 0 и определения градиента функционала имеем
т
vqSI = |С*(х,у,0,^ , (14)
0 т
= |С*(х, у,-1, ^ . (15)
0
Аналогично для краевых условий (3) можно получить
т
Уе/ = 1С * (х0, У0,0, М, (16)
т
= |С*(х0,У0,-1,t^ . (17)
0
Далее осуществляется итерационный спуск в направлении соответствующего градиента функционала.
Результаты численных экспериментов. Численные эксперименты проводились с использованием модели из работы [1] для акватории Азовского моря. Для тестирования вариационного алгоритма идентификации мощности
МОРСКОЙ ГИДРОФИЗИЧЕСКИЙ ЖУРНАЛ № 2 2015 81
источника в Казантипском заливе был проведен расчет на установление модельного поля течений под воздействием постоянного ветра восточного направления со скоростью 10 м/с. В результате моделирования также было получено пространственное распределение коэффициентов Ан и К . Поля скоростей и коэффициентов турбулентной диффузии использовались в качестве входной информации при интегрировании модели переноса пассивной примеси на срок 10 сут.
Рис. 1. Область загрязнения О и абсолютные значения мощности источника QB
На рис. 1 показана область загрязнения О в Казантипском заливе, в которой задается переменное по пространству значение мощности источника на
дне QB . При этом максимальное значение QBlax = _1 задается в центральной части области. На рис. 2 представлен результат моделирования распространения пассивной примеси при восточном ветровом воздействии и указанном переменном по пространству потоке вещества на дне Казантипского залива. Рис. 3 характеризует падение нормированного значения функционала качества прогноза в случае усвоения полной информации о поле концентрации, при этом в процессе итерации происходит идентификация пространственной структуры QB .
Рис. 2. Поле концентрации загрязнения при восточном ветровом воздействии и переменной по пространству мощности источника
1П о
Рис. 3. Падение нормированного значения функционала качества прогноза в зависимости от номера итерации п при переменной по пространству мощности источника Qв
На рис. 4, 5 представлены абсолютные значения мощности источника QB на первой и двадцатой итерациях работы алгоритма идентификации. Видно, что в процессе итераций существенно улучшается соответствие найденной мощности QB (рис. 5) первоначально заданной (рис. 1). При уменьшении общего объема информации о пространственной структуре поля концентрации сходимость итерационного процесса замедляется. Существенно ускоряется процесс идентификации в случае, если QB в указанной области является постоянной величиной, тогда для достижения минимума функционала при работе данной процедуры достаточно одной итерации. Отметим, что QB = const в области Q можно оценить и другим способом, например методом линеаризации или на основе метода сопряженных уравнений.
45.5
45
Азовское море
10.25 0.2 0.15 0 1 0.05 0
Черное море
36 36.5 в. д.
Рис. 4. Абсолютное значение мощности источника QB на первой итерации вариационного алгоритма идентификации
Используемый алгоритм позволяет восстанавливать пространственную структуру потока вещества. В данной работе без ограничения общности рассматривается величина QB . Аналогично решается задача для QS . В статье [4] показана работоспособность алгоритма при идентификации переменной по времени мощности Q{t). Главное требование при этом - достаточное ко-
личество информации для сходимости итерационного процесса. Понятно, что Сизм может не принадлежать одному моменту времени Т , а относиться к различным моментам времени t е[0, Т]. В соответствии с работой [5], модель в данном алгоритме выступает в роли пространственно-временного интерпо-лянта, поэтому главное также, чтобы общего количества информации было достаточно для сходимости итерационного процесса. Наиболее информативными будут точки измерений в области максимальных значений поля концентрации и в районах, обладающих существенной информационной связью с областью О. Эти точки можно определить по функциям влияния на основе решения соответствующих сопряженных задач.
45.5
45
Азовское море
Черное море
10.8 I Об
10.4
36
36.5
в .д.
Рис. 5. Абсолютное значение мощности источника QB на двадцатой итерации вариационного алгоритма идентификации
Моделирование распространения загрязнения от постоянно действующего источника единичной мощности проводилось при различном ветровом воздействии. Программный код предусматривает задание источника как на поверхности моря QS = 1, так и на дне QB = -1. Рассмотрим случай QB = -1 в Казантипском заливе при восточном ветре. При таком ветровом воздействии происходит достаточно интенсивное «проветривание» залива, и примесь распространяется в северо-западном направлении (рис. 6). При ассимиляции в качестве данных измерений используется информация с периферии области загрязнения (левее вертикальной линии). На рис. 6, 7 приведены шка-
лы значений полей концентрации, нормированных на соответствующие максимальные величины.
с.ш.
45.5
45-
Черное море
Азовское море
1
0.8 0 6 | 0.4 0.2 0
36
36.5
в.д.
Рис. 6. Поле концентрации загрязнения при восточном ветровом воздействии и постоянно действующем точечном источнике
Рис. 7. Функция влияния (решение сопряженной задачи)
86 МОРСКОЙ ГИДРОФИЗИЧЕСКИЙ ЖУРНАЛ № 2 2015
Рис. 8. Падение нормированного значения функционала качества прогноза в зависимости от номера итерации п при постоянной мощности источника QB
В результате решения сопряженной задачи (10) - (12) при первой итерации имеем распределение С * на момент времени ^ = 0 (рис. 7), которое характеризует влияние невязок прогноза на мощность QB в точке (х0, у0). В процессе итераций (рис. 8) происходит падение нормированного значения функционала качества прогноза и восстанавливается известное значение QB . Результаты численных экспериментов показали, что сходимость итерационного процесса зависит от количества усваиваемой информации. В случае ассимиляции всей информации о моделируемом поле на конечный момент времени для достижения минимума функционала требуется одна итерация. Наибольшей информативностью обладают точки, расположенные ближе к источнику загрязнения.
В статье [4] показано, что для переменной по времени мощности источника Q(t) необходима информация обо всем поле концентраций на конечный момент времени. Учитывая возможность распределения данных измерений по пространству и по времени, можно также утверждать, что для более точ-
ной идентификации Q(t ) необходимо располагать точки измерений в области, прилегающей к источнику загрязнения.
Заключение. В целом проведенные численные эксперименты показали надежную работу вариационного алгоритма идентификации мощности источника загрязнения применительно к модели переноса пассивной примеси в Азовском море. Результаты могут быть использованы для решения различных задач экологической направленности при изучении воздействия источников загрязнения антропогенного характера в акваториях Азовского и Черного морей.
СПИСОК ЛИТЕРАТУРЫ
1. Иванов В.А., Фомин В.В. Математическое моделирование динамических процессов в зоне море - суша. - Севастополь: ЭКОСИ-Гидрофизика, 2008. - 363 с.
2. Кочергин С.В., Кочергин В.С., Фомин В.В. Определение концентрации пассивной примеси в Азовском море на основе решения серии сопряженных задач // Экологическая безопасность прибрежной и шельфовой зон и комплексное использование ресурсов шельфа. - Севастополь: ЭКОСИ-Гидрофизика, 2012. - Вып. 26. - Т. 2. - С. 112 - 118.
3. Marchuk G.I., Penenko V.V. Application of optimization methods to the problem of mathematical simulation of atmospheric processes and environment // Modelling and Optimization of Complex Systems / Ed. G.I. Marchuk. - Proc. оf the IFIP-TC7 Working conf. - NewYork: Springer, 1978. - P. 240 - 252.
4. Кочергин В.С., Кочергин С.В. Использование вариационных принципов и решения сопряженной задачи при идентификации входных параметров модели переноса пассивной примеси // Экологическая безопасность прибрежной и шельфовой зон и комплексное использование ресурсов шельфа. - Севастополь: ЭКОСИ-Гидрофизика, 2010. -Вып. 22. - С. 240 - 244.
5. Пененко В.В. Методы численного моделирования атмосферных процессов. - Л.: Гидро-метеоиздат, 1981. - 350 с.
Identification of a pollution source power in the Kazantip Bay using the variation algorithm
V.S. Kochergin, S.V. Kochergin
Marine Hydrophysical Institute, Russian Academy of Sciences, Sevastopol, Russia e-mail: vskocher@gmail.com
Model of passive impurity transport in the Azov Sea is considered. It serves a basis for implementing the variation algorithm aimed at identifying the pollution source power including a variable over space. The test example shows efficiency of the algorithm for searching optimal distribution of the power source in space agreed with measurement data. Test simulations are performed for the Kazantip Bay area at the eastern wind impact.
Keywords: identification of power source, variation algorithm, discrepancy functional, concentration field, assimilation of measurement data, Sea of Azov.
