Математическое моделирование. Оптимальное управление Вестник Нижегородского университета им. Н.И. Лобачевского, 2011, № 3 (2), с. 108-114
УДК 517.956.4:517.983.54
ДВОЙСТВЕННАЯ РЕГУЛЯРИЗАЦИЯ В ОБРАТНОЙ ЗАДАЧЕ ФИНАЛЬНОГО НАБЛЮДЕНИЯ ДЛЯ ПАРАБОЛИЧЕСКОГО УРАВНЕНИЯ
© 2011 г. Ф.А. Кутерин
Институт прикладной физики РАН, Н. Новгород
Поступила в редакцию 09.07.2010
Рассматривается обратная задача нахождения начального условия и распределенной правой части операторного параболического уравнения по приближенно известному в финальный момент времени решению. Обсуждается применение алгоритма двойственной регуляризации и итеративного двойственного метода для решения этой некорректной задачи. Приводятся результаты численных расчетов, демонстрирующих поведение алгоритма.
Ключевые слова: некорректные задачи, двойственная регуляризация, параболические уравнения, обратные задачи, дифференциальные уравнения с частными производными.
В работах [1-4] был предложен и развит метод (алгоритм) двойственной регуляризации для исследования и решения задач математического программирования [4], оптимального управления [1, 3], а также обратных задач [2]. Этот метод заключается в непосредственном решении на основе метода регуляризации Тихонова задачи, являющейся двойственной к исходной оптимизационной задаче, и, как следствие, одновременном решении исходной (прямой) задачи.
Хорошо известно, что двойственные численные алгоритмы являются одними из самых популярных и эффективных при решении оптимизационных задач с ограничениями [5, 6]. По этой причине представляется естественным, что развитие методов регуляризации на основе теории двойственности позволило дополнить существующие методы регуляризации оптимизационных задач (см., например, [5]) эффективными методами, в которых самым существенным и непосредственным образом используется классическая идея «снятия» ограничений, заложенная в принципе Лагранжа.
Целью настоящей работы является численная иллюстрация особенностей метода двойственной регуляризации применительно к обратной задаче финального наблюдения для операторного параболического уравнения по определению правой части уравнения и начального состояния. Рассматриваемая формально в работе достаточно абстрактная обратная задача содержит в себе в качестве частных случаев различные конкретные обратные задачи финального наблюдения для классического уравнения
теплопроводности. В частности, в заключительной части настоящей работы приводятся результаты численного решения обратной задачи финального наблюдения по определению начального состояния в третьей краевой задаче для уравнения теплопроводности. Посредством численного решения и анализа этой обратной задачи иллюстрируется тот факт, что метод двойственной регуляризации сохраняет работоспособность как при наличии ошибки задания исходных данных, так и в отсутствие седловой точки функционала Лагранжа оптимизационной задачи, эквивалентной исходной обратной задаче. Последними обстоятельствами обсуждаемый двойственный алгоритм принципиально отличается от классического двойственного алгоритма Удзавы [6, 7], заключающегося, как известно, в непосредственном решении задачи, двойственной к исходной оптимизационной, на основе градиентного метода. Этот весьма популярный алгоритм широко используется при решении различных прикладных задач (см., в частности, [8, 9]), вопросы его сходимости рассматривались, например, в работах [8, 9]. В то же время, к его существенным недостаткам относятся, как отмечено в [4], неустойчивость по возмущению исходных данных и необходимость существования седловой точки функции Лагранжа при доказательстве сходимости [8, 9].
1. Постановка задачи
Пусть Н - гильбертово пространство, и1, и2 с Н - выпуклые замкнутые множества, =
= (ие L2(0,T;H): u(t) е U при почти всех (п.в.) te [0,T]}, D2 = U2 с H. Множество D = Dj x D2 назовем множеством допустимых управлений.
Рассмотрим задачу Коши для операторного параболического уравнения [10]
y'(t') + Ay(t) = u(t), 0 < t < T, j(0) = v. (1)
Будем предполагать, что оператор A является энергетическим расширением линейного неограниченного симметричного положительно определенного оператора с областью определения, плотной в вещественном сепарабельном гильбертовом пространстве H. Это значит, что A е L(V ^ V*), где V - энергетическое гильбертово пространство, V* - пространство, сопряженное к V, D(A) = V - область определения и R( A) = V * - область значений оператора A . Будем предполагать также, что имеют место * *
вложения V с H; H с V , причем эти вложения плотные и непрерывные. Решением задачи Коши (1), соответствующим паре (u, v) е D, называется регулярная обобщенная функция у(-) = у[%\(:) е L2 (0, T; V), имеющая регулярную производную y'(t) e L2(0,T;V*), удовлетворяющая уравнению (1) при п.в. t е (0,T) и начальному условию Необходимые подробности, связанные с понятием обобщенного решения задачи (1) и определения пространств L2(0,T;V) и L2(0,T;V*) см. в [10].
Обсуждаемая обратная задача финального наблюдения заключается в нахождении пары п = (u, v) е D с минимальной нормой |У22 +
II IIL (0J;#)
+ IV H по известному в момент времени T наблюдению q е H. Предполагаем при этом, что решение обратной задачи существует, обозначим его через п0. Таким образом, имеем равенство y[n° ] = q.
Поставленная обратная задача эквивалентна задаче оптимального управления с ограничением типа равенства
/0 (п) ^ inf, I (п) = q, (2)
hi <§, II 11я
§ > 0.
(5)
где п е D , q е H и
Г / Ч II II2 II II2
0 (п) = ||М||l2(0,T;H) +1 HIH ’
Il(n) = y[n\(T).
С учетом погрешности финального наблюдения имеем возмущенную задачу оптимального управления
I0(п) ^ inf, /j(n) = q8 = q + h8, пе D, q8 e H.
(6)
2. Необходимые теоретические сведения
Для задачи (1) при сделанных выше предположениях об операторе А справедлива теорема существования и единственности решения [10].
Теорема 1. Пусть А е Ь(¥ ^ V*) - энергетическое расширение линейного неограниченного симметричного положительно определенного оператора с областью определения, плотной в гильбертовом пространстве Н, и вложение пространств V с Н компактно. Тогда для любых и е £2(0,Т;У*) и V е Н существует единственное решение задачи (1).
Определим функционал Лагранжа задачи (6)
(n, X) = /0(п) + (Х, I1(n) - q8),
пе D, X е H
и двойственную задачу Хе H.
К (Х) = min L° (п,X) ^ sup,
(7)
(8)
(3)
(4)
Будем считать, что оператор А исходной задачи (1) задан точно, а финальное наблюдение -с ошибкой Нъ, такой, что
Заметим, что операция шт в определении целевой функции двойственной задачи законна, так как функционал является сильно выпуклым на выпуклом замкнутом множестве Б с с I? (0,Т; Н) х Н. При этом минимум достигается для любого q в единственной точке л8 [X] = ащшт 1?д (п, X), Хе Н.
Лемма 1. Производная Фреше функционала в точке Хе Н равна
д¥д5 (Х) = у[п5[Х]](Г) - ^8.
Эта производная удовлетворяет условию Липшица.
Доказательство этой леммы полностью аналогично доказательству леммы 2.6 в [2].
3. Метод двойственной регуляризации
Метод двойственной регуляризации [1-4] применительно к решению задачи (2) заключается в максимизации регуляризованного по Тихонову функционала
(9)
Этот функционал является сильно вогнутым и обладает в силу леммы 1 равномерно непрерывным (липшицевым) градиентом. Последнее обстоятельство позволяет применить для максимизации функционала (9) тот или иной градиентный метод, например, хорошо известный метод наискорейшего подъема.
Обозначим через А8’а единственную в Н точку, дающую безусловный максимум функционалу Яъ,а.
Теорема 2. При условии согласования ошибки задания исходных данных 5 с параметром регуляризации а
8/а(8) ^ 0, 8^ 0 регуляризованные элементы л8[А,8,а(8) ] для любого ц сильно сходятся при 8 ^ 0 к решению
невозмущенной задачи п°.
Доказательство проводится аналогично доказательству в [2].
Следует отметить, что, как и в [2], утверждение теоремы справедливо как в случае, когда функционал V0 достигает максимума на Н, так и в случае, когда этого нет.
4. Итеративная двойственная регуляризация
Для каждого фиксированного а может быть найдено решение Х5,а задачи
Я8’“ (А) ^ тах, Хе Н. (10)
Будем применять для решения двойственной задачи (10) процедуру итеративной регуляризации [5]. Поясним кратко смысл этой процедуры.
Предположим, мы задали последовательность параметров ак, согласованных с последовательностью ошибок наблюдения 8к. В этом
случае с возрастанием номера к задачу (10) при а = ак придется решать со все более высокой
точностью, поэтому число итераций пк метода
градиентного подъема может оказаться весьма большим при больших к. Это обстоятельство может создать трудности при практической реализации метода двойственной регуляризации. Поэтому для практического решения обратной задачи (2) предлагается использовать процедуру итеративной двойственной регуляризации [2-4].
Принцип итеративной регуляризации градиентного двойственного метода решения задачи (10) заключается в том, что при каждом фиксированном к делается один шаг метода градиентного подъема с шаговым множителем в = вк
при а = ак и затем делается переход при а=ак+1-
Начальную точку Х0 можно выбрать произвольным образом из гильбертова пространства
Н. Чтобы полученная таким образом процедура была устойчива по отношению к ошибкам исходных данных, необходимо, чтобы последовательность параметров регуляризации ак и
последовательность шаговых множителей рк были определенным образом согласованы с последовательностью ошибок в исходных данных 8 к.
Итак, пусть последовательность Лс,
к = 1,2, к, в соответствии с методом итеративной двойственной регуляризации конструируется по правилу
X*+1 = X* + рк(дУ*к (X*)- 2акX*), то
(11)
к = 0,1,к; X0 єН,
где числовые последовательности 8к,ак,Рк , к = 0,1, к, удовлетворяют условиям
8* ^ 0, ак > 0, рк > 0,
Нш(8к +а* + р* )-0,
а
-< С
Кк+1 а*+і- а* I _
Ііт . з
к^“ (а*) Рк 8
ак
к1 -0, Ііш-%-0,
(12)
к (ак)
Іі^-^-3 -0, ХакРк - +да-к (ак) к=і
Тк
Теорема 3. При условии, что элементы X при к = 1,2, ... находятся из рекуррентной формулы (11) и выполняются условия согласования (12), имеет место предельное соотношение
Ііт
к
^ к Г'і к п 0
п к [X ] -п
= 0,
где п0 - решение исходной обратной задачи (2).
Замечание. В качестве одного из возможных примеров указанных выше числовых последовательностей можно взять, в частности,
ак = к-1/6, Рк = к-3/5.
При практической реализации схемы двойственной регуляризации в случае, когда исходные данные задаются с определенной фиксированной погрешностью 8 > 0, может быть использовано следующее условие останова итеративного процесса (11) (см., например, [4]). При фиксированном уровне погрешности 5 > 0 итерации продолжаются до такого наибольшего номера к = к (8), при котором выполняются неравенства
5к >5, к = 1,2,к,к(8). (13)
Теорема 4. Пусть X0 е Н, а элементы А при к = 1,2, к находятся из рекуррентной формулы (11) и выполняются условия согласования (12). Тогда верно предельное соотношение
Иш
5^-ГО
[Хк(5)] -:
где П 5к [Хк (5) ] - результат к (8) итераций
8 г
цесса (11) при дУя к (X) = дУ8 (X).
про-
5. Численный эксперимент
Рассмотрим пример, демонстрирующий работу и особенности метода итеративной двойственной регуляризации [2-4].
В качестве модельной рассмотрим обратную задачу нахождения начального условия у(х) в третьей краевой задаче для уравнения теплопроводности
У, - Ухх = 0,
_у(х,0) = х),
Ух(0, * ) = Ж 0> (14)
.Ух (ь> * ) = -У(ь * )>
0 < * < Т,
Ь = 1, Т = 0.001
из множества допустимых управлений В = = Д = (у(х)еЬ2(0,Ь): у(х) е [-2,2] при хе[0,Г]}.
Чтобы воспользоваться описанной выше процедурой, поясним, как задача (14) записывается в виде абстрактной задачи Коши (1). Следуя схеме, приведенной в [10, гл. 2, §2, с. 66], введем оператор
Л/ = — f''(х), 0 < х < Ь, с областью определения
ДА) = {/ = /(х) е С2[0,Ь]:
-/'(0) + /(0) = 0, /'(Ь) + /(Ь) = 0},
плотной в Н = Ь (0, Ь).
Оператор А - линейный, симметричный, неограниченный и положительно определенный. Энергетическое пространство V, получаемое пополнением 0(А) в метрике (А/, /)1/2,
здесь совпадает с пространством Н 1(0, Ь). Энергетическое расширение
А -.V ^ V * = (Н 1(0, Ь))*
исходного оператора А определяется равенством
(Л/, ^ = / (0) я (0) + / (Ь) я (Ь) +
I
+ { /'(х) я'(х)Лх, У/, я е Н1 (0, Ь).
о
Это позволяет рассматривать задачу (14) как частный случай задачи (1) и применять для решения обратной задачи финального наблюдения описанные выше метод двойственной регуляризации и итеративный двойственный метод.
Численные эксперименты проводились по замкнутой схеме. Это означает, что для некоторой функции у0, «точного решения» обратной
задачи, решалась задача (14) при V = у0 , и функция ут = у[у0 ](•, Т) рассматривалась как финальное наблюдение. Далее функция ут возмущалась с целью моделирования «погрешности измерения». «Приближенно известное» финальное наблюдение ц5 вычислялось по формуле
8 . >8 Я = ут + п .
Далее решалась обратная задача нахождения V при финальном наблюдении ц. При этом использовался метод итеративной двойственной регуляризации (11) с описанным выше правилом останова (13). В качестве последовательностей, удовлетворяющих условиям согласования (12), использовались последовательности ак = а0к-1/6, вк = Ро*-3/5 при а0 = 0.0000001, р0 = 2000. Значения констант а0 и р0 были
подобраны опытным путем таким образом, чтобы метод сходился «достаточно быстро». Полученное решение обратной задачи V сравнивалось с у0, кроме того, оценивалось значение
нормы соответствующей двойственной переменной X.
Задача (14) и сопряженная задача решались с помощью разностной схемы Кранка-Никол-сона, которая, как хорошо известно, устойчива при любых шагах по времени и имеет второй порядок точности. Соответственно, вместо функций из гильбертовых пространств рассматривались сеточные представления этих функций. Для обеспечения достаточной точности решения использовалось равномерное разбиение по пространственной переменной с помощью 256 точек.
А К
Ошибка к в пространстве сеточных приближений задавалась следующим образом. Бра-
Л1
лась сеточная функция й со случайными значениями, распределенными равномерно на отрезке [-1,1], вычислялась ее норма й1 , а также
0
п
«собственно погрешность» по формуле
/8 = б(йу| |/гЦ).
В качестве нормы в пространстве сеточных приближений рассматривался сеточный аналог нормы в Ь2
начальной функции взята гладкая функция V
(1)
„(1) _ X2 (X - У)2
224 У4
(15)
«Погрешность наблюдения» равнялась 1%,
то есть 5 = 0.01 останова (13).
Использовалось условие
Таблица 1
Уменьшение невязки с возрастанием точности решения задачи
Число итераций, N 11^ N1 |К - ^]
10 0.555782 0.00554065
100 1.00724 0.000299552
1000 1.79222 0.000202258
10000 3.18836 0.000137671
100000 5.67136 9.37373е-05
IIй!)2 ^1 /к -1),
где п - число точек разбиения, а и1 - значения
функции й в соответствующих точках.
При проведении численных экспериментов ставились три цели:
1. Продемонстрировать сходимость алгоритма итеративной двойственной регуляризации в случае, когда у функции Лагранжа (7) задачи (2)-(4) есть седловая точка.
2. Показать поведение алгоритма в ситуации, когда нет седловой точки у функции Лагранжа (7).
3. Продемонстрировать неустойчивость алгоритма (11), (12) в отсутствие регуляризирую-щего слагаемого при наличии ошибки в исходных данных.
Для достижения перечисленных целей проводились три серии численных экспериментов. Ниже при описании их результатов в каждом из трех случаев приводятся графики начальной функции у0, результата восстановления методом итеративной двойственной регуляризации у[Ам ] и график функции Xм (двойственной
переменной) на последнем Ж-м шаге итеративного процесса. Кроме того, подсчитывается норма двойственной переменной Xм и невязка
1К - ^ N1.
В первой серии экспериментов в качестве
Рис. 1. Двойственная переменная на 10-м шаге
Рис. 2. Начальная функция и восстановленная начальная функция, 10 итераций
Уменьшение нормы разности
к} - VI
[*■ N ]||
говорит о сходимости метода итеративной двойственной регуляризации. Следует отметить низкую скорость сходимости итерационного двойственного метода, обусловленную априорным выбором последовательностей шаговых и регуляризующих множителей а к и рк.
Во второй серии экспериментов использовалась разрывная начальная функция
\ -1, х < а, х > Ь,
1, а < х < Ь,
у(2) = У0 -
(16)
принимающая значения внутри допустимого множества [-2,2]. В этой ситуации, как показано в [3], функционал Лагранжа (7) задачи (2)-(4) не имеет седловой точки. Это значит, что норма двойственной переменной X по мере работы алгоритма должна неограниченно увеличиваться. Чтобы наблюдать этот эффект «в чистом виде» (независимо от вклада погрешности наблюдения) «погрешность наблюдения» 5 была положена равной нулю.
Анализируя приведенные ниже результаты, можно наблюдать, что с ростом числа итераций происходит приближение решения к начальной функции и одновременное увеличение нормы X к.
0
Таблица 2
Рост нормы двойственной переменной
Число итераций, N 1* Л ||У0 - N1
100000 647.621 0.0299701
1000000 1764.2 0.00961037
10000000 3601.9 0.00121545
Рис. 3. Двойственная переменная на 10000-м шаге
Рис. 4. Начальная и восстановленная функции на 10000-м шаге
Рис. 5. Двойственная переменная на 100000-м шаге
1.5 1
0.5 0
-0.5 -1 -1.5
0 0.2 0.4 0.6 0.8 1
Рис. 6. Начальная и восстановленная функции на 100000-м шаге
В третьей серии экспериментов, как и в первой, в качестве начальной функции была взята гладкая функция V® (15), но все элементы а к,
отвечающие за регуляризацию в формуле (11), были положены равными нулю. В то же время, «финальное наблюдение» возмущалось равномерно распределенной случайной ошибкой к8, как было описано выше. Уровень ошибки был равен 3% (§ = 0.0з|| д 0||). Результаты численного
эксперимента представлены на следующих двух графиках и в таблице.
Таблица 3
Неустойчивость итерационного двойственного метода без регуляризации
Число итераций, N ^ N Уо - ЧЛ N ]
10000 299.144 0.0923157
100000 1054.59 0.127038
1000000 4742.26 0.222071
4000 3000 2000 1000 0
-1000 -2000 -3000
0 0.2 0.4 0.6 0.8 1
Рис. 7. Двойственная переменная на 100000-м шаге
--V — V-
Рис. 8. Начальная и восстановленная функции на 100000-м шаге
По мере работы алгоритма как норма двойственной переменной, так и невязка увеличивались. И, как можно заметить из графиков начальной и восстановленной функций, приближенное решение обратной задачи «не тяготеет» к ее точному решению, что говорит о неустойчивой работе итеративного двойственного метода без регуляризации и о ее существенности в задаче с не точно известными начальными данными.
Результаты проведенных экспериментов показали эффективность итеративного двойственного метода [1-4] в задаче финального наблюдения для одномерного уравнения теплопроводности и необходимость регуляризации при решении подобного рода задач.
В заключение настоящей работы выражаю благодарность своему научному руководителю профессору М.И. Сумину за постановку задачи и внимание к работе.
Работа подготовлена при финансовой поддержке Российского фонда фундаментальных исследований (проект 09-01-97019-р_поволжье), а также аналити-
ческой ведомственной целевой программы «Развитие научного потенциала высшей школы (2009-2010 годы)» Минобрнауки РФ (проект 2.1.1/3927).
Список литературы
1. Сумин М.И. Оптимальное управление параболическими уравнениями: двойственные численные методы, регуляризация // Распределенные системы: оптимизация и приложения в экономике и науках об окружающей среде: Сб. докладов к Международной конференции (Екатеринбург, 30 мая-2 июня 2000 г.). Екатеринбург: Изд-во Ин-та математики и механики УрО РАН, 2000. C. 66-69.
2. Сумин М.И. Регуляризованный градиентный двойственный метод решения обратной задачи финального наблюдения для параболического уравнения // Ж. вычисл. матем. и матем. физ. 2004. Т. 44. № 11. C. 2001-2019.
3. Сумин М.И. Регуляризованный двойственный алгоритм в задачах оптимального управления для распределенных систем // Вестник Нижегородского университета. Серия «Математическое моделирование и оптимальное управление». 2006. Вып. 2(31). С. 82-102.
4. Сумин М.И. Регуляризация в линейновыпуклой задаче математического программирования на основе теории двойственности // Ж. вычисл. матем. и матем. физ. 2007. Т. 47. № 4. С. 602-625.
5. Васильев Ф.П. Методы оптимизации. М.: Факториал-Пресс, 2002. 824 с.
6. Мину М. Математическое программирование. Теория и алгоритмы. М.: Наука, 1990. 488 с.
7. Эрроу К.Дж., Гурвиц Л., Удзава Х. Исследования по линейному и нелинейному программированию. М.: ИЛ, 1962. 336 с. [Англ. оригинал: Arrow K.J., Hurwicz L., Uzawa H. Studies in linear and nonlinear programming. Stanford University Press, 1958.]
8. Экланд И., Темам Р. Выпуклый анализ и вариационные проблемы. М.: Мир, 1979. 400 с.
9. Гловински Р., Лионс Ж.-Л., Тремольер Р. Численное исследование вариационных неравенств. М.: Мир, 1979. 576 с.
10. Осипов Ю.С., Васильев Ф.П., Потапов М.М. Основы метода динамической регуляризации. М.: Изд-во Московского университета, 1990.
DUAL REGULARIZATION IN THE INVERSE PROBLEM WITH FINAL OBSERVATION FOR
PARABOLIC EQUATION
F.A. Kuterin
An inverse problem is considered of finding the initial conditions and the distributed right-hand side of a parabolic equation by the approximation known at the final time. The application of the dual regularization algorithm and of the iterative dual method for the solution of this ill-posed problem is discussed. Numerical calculation results demonstrating the algorithm behavior are presented.
Keywords: ill-posed problems, dual regularization, parabolic equations, inverse problems, partial differential equations.