Регулярный алгоритм аппроксимации негладких решений для интегральных уравнений Фредгольма первого рода Текст научной статьи по специальности «Математика»

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

Том 15, № 2, 2010

Регулярный алгоритм аппроксимации негладких решений для интегральных уравнений Фредгольма

первого рода*

В. В. Васин, Т. И. Серёжникова Институт математики и механики УрО РАН, Екатеринбург, Россия e-mail: vasin@imm.uran.ru, sti@imm.uran.ru

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

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

Введение

Вычислительная практика решения интегральных уравнений первого рода убедительно показала, что применение в тихоновской регуляризации стабилизаторов в форме вариации различных типов позволяет значительно улучшить качество аппроксимации разрывных решений по сравнению с классической квадратичной регуляризацией, когда используется гильбертова норма L^h (см. [1-4]). При применении вариации функции в качестве стабилизатора удается обосновать кусочно-равномерную сходимость регуляризованных решений, что, по-видимому, и объясняет преимущество такого рода регуляризующих алгоритмов при восстановлении негладких решений.

Однако при численной реализации этих методов, т. е. при построении регуляризо-ванного семейства приближенных решений, приходится иметь дело с задачей негладкой минимизации. Здесь наметились два подхода: а) негладкий функционал предварительно аппроксимируется семейством дифференцируемых функционалов, а затем применяются традиционные методы гладкого анализа, например, методы градиентного или ньютоновского типа [1, 2]; б) для регуляризованной негладкой выпуклой задачи минимизации тихоновского функционала применяются субградиентные процессы [3, 4].

В работе [5] был предложен другой параметрический класс недифференцируемых стабилизирующих функционалов на основе нормы пространства Липшица и обоснована равномерная сходимость тихоновских аппроксимаций к непрерывному необязательно дифференцируемому решению исходного уравнения. Оказалось, что при подходящем выборе управляющих параметров регулярные алгоритмы, построенные в рамках упомянутой регуляризации с привлечением ргох-метода, вполне пригодны для г,осп а ног,-

* Работа выполнена при финансовой поддержке РФФИ (грант № 09-01-00053) и Междисциплинарного проекта УрО РАН. © ИВТ СО РАН, 2010.

ления негладких (разрывных) решений интегральных уравнений Фредгольма первого рода.

Заметим, что идея о возможности использования нормы пространства Липшица в качестве стабилизатора была высказана в обзорной статье [6].

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

1. Теоремы сходимости приближенных решений

Пусть П — компакт в Яп. Определим множество функций, удовлетворяющих на П условию Липшица

3 d> 0: |п(х^ — п(х2)| < d |я1 — х2|м V х\,х2 € П, (1)

/ п \ 1/2

где 0 < л < 1, |х1 — х2| = ^ = |х1г — , Обозначим через Нм = НМ[П] множе-

ство функций, удовлетворяющих на П условию (1) с фиксированным параметром По аналогии со случаем функции одной переменной определим норму

...... |п(х1) — п(х2) |

и(х) № = тах ф) + вир —:-:-,

Х6П Х1,Х2&П |х1 — х2|^

относительно которой Нм является банаховым пространством (см, [7]), Рассмотрим линейное уравнение

Ап = / (2)

с оператором А, действующим из пространства Нм в пространство непрерывных функций С(П), и единственным решением П € ННепрерывность обратного оператора А-1 не предполагается, поэтому (2) относится к классу существенно некорректных задач. При приближенно заданной правой части /, || / — /1| < 8, для построения регуляризо-ванного семейства приближенных решений используем метод Тихонова в виде [5, 7]

ш1п(||Ап — /||с[п] + а||п||ям : п € Н^[П]}. (3)

Теорема 1. Пусть А — линейный непрерывный оператор в С(П), Кег(А) = {0},

|| / — /б || < 8. Тогда задача имеет решение па, возможно неединственное, и при связи,

8

параметров а(8) —> 0, 0, 8 —> 0 имеет место равномерная сходим,ость

а(8)

б^о

кроме того,

б^о

Иш ||па(б) — П||с(п) = 0, Иш ||па(б)||Ям = ||П||Ям.

Замечание. Если вместо липшицевой нормы в (3) использовать стабилизатор

= |Н|Як + ^1Ы|£р, 1<р<оо, (4)

то решение па задачи (3) единственно.

Используя метод Рптца, можно построить аппроксимации регуляризованного приближенного решения ua элементами конечномерного пространства кусочно-линейных функций. Предположим, что компакт П допускает разбиение симплекеами n-ro порядка с максимальным диаметром h = maxAm. Обозначим через Uh конечномерное подпро-

m

странство кусочно-линейных (линейных на каждом симплексе) функций.

Теорема 2. Пусть выполнены предположения теоремы, 1. Тогда конечномерная задача

min{||Au - fs||с[п] + «ОЫЦрП + ||и^||ям] : uh G Uh} (5)

имеет единственное решение и^, для, которого выполнены соотношения,

lim ||< - иа||с[П] = 0, lim ||<||ям[П] = ||иа||ям[П],

h ^U h ^U

где иа — решение, полученное методом Тихонова со стабилизатором, (4)-

Как показали численные эксперименты, выполненные для интегральных уравнений Фредгольма первого рода с негладкими (разрывными) решениями, при малых значениях параметра а субградиентные методы работают неудовлетворительно, в результате чего не удается восстановить, например, хорошо выраженный разрыв или излом, С дру-

а

решения при данном методе регуляризации,

В значительной степени эти трудности удается преодолеть, если для устойчивой аппроксимации экстремального элемента иа использовать ргох-метод (итерированный вариант метода Тихонова), а именно: если обозначить целевой функционал в (3) через Фа(и), то итерационный процесс для аппроксимации иа строится следующим образом:

uk = argmin {Фа(и) + в ||и - uk-1||H : и G U}, (6)

где в > 0, || ■ ||h — некоторая гильбертова норма,

В отличие от выпуклого функционала Фа(и) функционал Фа'в = Фа(и) + в 11и — v||H является сильно выпуклым, что гарантирует устойчивый счет при нахождении uk субградиентным методом. По-видимому, это обстоятельство объясняет факт повышения качества аппроксимации негладкого решения, однако за счет увеличения трудоемкости алгоритма,

2. Алгоритм решения интегрального уравнения Фредгольма первого рода и численные эксперименты

Рассмотрим одномерное интегральное уравнение Фредгольма первого рода

1

Au = У K(x,y)u(x)dx = f (y), 0 < y < 1 (7)

U

K(x, y), f(y)

тально алгоритм, предложенный в предыдущем разделе. Прежде всего заметим, что конечномерная аппроксимация (см (5)), описанная в теореме 2, не ведет к полной дискретизации регуляризованной задачи, поэтому используем иную, более удобную для численной реализации схему [7].

Зададим равномерную сетку по х и по у с шагом к = 1 /М и аппроксимируем интеграл формулой правых прямоугольников.

Введем обозначения: п(х^) = п, /б(уг) = /г. Тогда дискретный аналог регуляризо-ванной задачи со стабилизатором (4) принимает вид

N

min < max ^ hK(yh Xj)uj — fi| +

+а

i<j<N j=1

1 N |u - u-I

X—^ 7 I I 2 I I i j

- > fiiUj + max Uj +max--—

2 j= jl i<j<N' jl i=j |x, - xjI1

{Uj }N e RN} . (8)

Теорема 3. Экстремальная задача, (8) имеет единственное решение Пм = (п1, п2,... и^), и последовательность кусочно-линейных восполнений пм (х), построенных по вектору пм, равномерно сходится к решению задачи,

-1,

min j || Au — f ||C + а

+ INI*"

: u e H1

Теперь дискретный аналог ргох-метода (6) принимает форму

12

uk ^^min {ФN(u) + ß ||u — uk-1||R : u e RN} = ФN k, (9)

где ФN — целевая функция в задаче (8),

Так как целевая функция в (9) является сильно выпуклой, то ргох-отображение

(см, [8]), определяющее итерационный процесс (9), принадлежит классу псевдосжима-

ющих отображений, и, следовательно, ргох-метод (9) сходится к решению задачи (8)

(см, [9, 10]), Поскольку выпуклая функция субднфференцнруема, то при фиксировап-uk-1 uk

Наиболее простой вариант субградиентного метода имеет вид [11]

Uk,u+1 = uk,u_ и = 0,1,2, ...,nk, (10)

где vk,V e ЗФ^к (uk,V), _ целевая функция в задаче (9), 5Ф означает субдифференциал функционала Ф. Как известно, в общем случае для выпуклой функции при

те те

Ak > 0, £Ak = оо, £ Ak < оо установлена сходимость итераций по функционалу

i=0 i=0

[11], но поскольку Ф^^ ~ сильно выпуклая функция, то имеет место сходимость по аргументу, т.е. lim ||uk,V — uk|| = 0.

Если известна некоторая оценка ФN, k для оптимального значения ФN k в задаче (9),

т, е,

ФN,k < ФN)k + £k, £k > 0,

то целесообразно перейти от задачи (9) к решению выпуклого неравенства

ФЙ (u) — ФN, k < 0. (11)

Теперь для приближенного решения неравенства (11) применим более эффективный

[12]

+ vk, v

uk>v+1 = - А (Ф%%(ик>») - Ф„,к) (12)

где vk'v £ (uk'v), 0 < Л < 2, v = 0,1, 2,..., nk, а величина Фможет уточняться

в процессе счета.

Важно отметить, что итерации uk ргох-метода (9) могут вычисляться субградиентными методами (10), (12) с некоторой погрешностью с сохранением сходимости, а именно, справедлива следующая теорема, обоснование которой следует из результатов [13, 14].

Теорема 4. Пусть Т — ргох-отображение, т. е.

T : v ^ arg min (Ф^(u) + в ||u - v||2 : u £ }. Пусть ргох-метод реализуется с погрешностью

ж

II/+1 - T(zk)||< 7fc, z0 £ , ]T7fc <

i=0

Тогда, lim ||zk — UN|| = 0, где UN — решение задачи, (8).

Замечание. Если известна априорная информация о решении уравнения (2) вида u £ Q, где Q — выпуклое замкнутое множество, заданное, например, системой линейных или выпуклых неравенств, то вместо ргох-метода можно использовать итерационный процесс [14]

uk+1 = Pq(T (uk)),

где T — ргох-отображение, Pq — фейеровекое отображение, отвечающее за ограничение u £ Q.

Таким образом, все описанные этапы алгоритма, построенного на основе тихоновской регуляризации со стабилизатором в форме нормы Липшица, ргох-метода, дискретной аппроксимации квадратурным методом и субградиентного процесса негладкой оптимизации, являются обоснованными, следовательно, мы имеем дело с регуляризую-щим алгоритмом,

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

ф(и) = тах т—--=7"- (13)

0<ij<N |x — x

Обозначим через J(u) множество пар индексов i0,j0, для которых в (13) достигается максимум. Определим вектор

11

-—, 0,..., 0, —--

/у . _ /у ^ /у . _ /у . ^

Кад ^Jo I КЗД Jo I

Vioh = 0,...,о,--—, о,...,о,---—, о,...,о

где ненулевые элементы занимают позиции с номерами го и тогда

М = с^„иеЛ«), если - > 0, = ^ со{Ми(—М)}, если — и70 = 0,

-М = Со{-г^}г0,^0€^(и), еСЛИ «¿0 — «70 < 0,

где co(M} — выпуклая оболочка множества M, u = (u0, u2,..., uN),

Для проведения численных экспериментов рассматривается решение интегрального уравнения Фредгольма первого рода [15]

которое моделирует ситуацию, когда необходимо найти продолженное на глубину х = Н гравитационное поле и(х) по измеренному на земной поверхности полю / (у) (плоский случай). Были проведены расчеты по восстановлению трех модельных решений. Использовалось нулевое начальное приближение. В формуле (8)для нормы Липшица полагалось ^ = 0.001, в (9) задавалось в = 300. Количество итераций по формулам (12) принималось равным 60. Для расчета вариантов с точной правой частью для задания параметра регуляризации использовались а = 10-4, 10-5, а для расчетов с возмущенной правой частью — а = 10-3, 10-4. Понятно, что с ростом величины параметра Н, 0 < Н < 2, процесс восстановления приближенного решения задачи (14) усложняется, поскольку свойства матрицы ухудшаются. Для моделей 1 и 2 полагалось Н = 2, для модели 3 — Н =1.

Для каждой модели рассчитывался вариант с точной правой частью и вариант с возмущенной правой частью, величина возмущения полагалась ||/ — /||с = 5 < 0.0108.

Точная правая часть для каждой модели рассчитывалась подстановкой точного решения и в сеточный вариант операторного уравнения Аи = /.

(14)

1 1 о о

0 20 40 60 80 100 120 х

а

1 1 0 о

0 20 40 60 80 100 120 х

б

0 20 40 60 80 100 120 *

в

Рис. 1. Модель 1. а — точное решение; б — восстановленное решение при точной правой части (Д1 = 1.7 х 10-2, Д2 = 1.0 х 10-6); в — восстановленное решение при возмущенной правой части (Д1 = 2.1 х 10-2, Д2 = 8.4 х 10-3)

Точное решение и(ж) (0 < ж, и < 2) для каждой модели рассчитывалось на равномерной сетке по ж, состоящей из 129 узлов с шагом Н = 2/128 ~ 0.016, по следующим формулам:

Модель 1 « =

Модель 2 « = <

Модель 3 « = <

1/32 х г,

г = 0,..., 64,

1/32 х (128 - г), г = 65,..., 128.

' 0, г = 0,..., 42, 2, г = 43,..., 85, 0, г = 86,..., 128.

0, г = 0, ... , 15,

1/8 х (г - 16), г = 16,..., 32,

2 - 1/32(г - 32), г = 33,..., 64,

, «128-г, г = 65, ..., 128.

В каждом варианте расчетов для полученного приближения вычислялись относительная погрешность по решению Д1 и относительная погрешность по невязке Д2:

Д1 =

- х\\ь2 \\х\\ь2

Д2 =

1|Ах - /||ь2

II/II

¿2

И -, 1.6" 1.20.80.40

20 40 60 80 100 120 х

1.61.20.80.4-

0 20 40 60 80 100 120 х

б

0 20 40 60 80 100 120 х

в

Рис. 2. Модель 2. а — точное решение; б — восстановленное решение при точной правой части (Д1 = 7.6 х 10-4, Д2 = 6.9 х 10-4); в — восстановленное решение при возмущенной правой части (Д1 = 3.9 х 10-2, Д2 = 2.7 х 10-2)

Рис. .3. Модель .3. а — точное решение; б — восстановленное решение при точной правой части (Д1 = 1.6 х 10-4, Д2 = 9.5 х 10-5); в — восстановленное решение при возмущенной правой части (Д1 = 7.1 х 10-3, Д2 = 5.1 х 10-3)

На рис. 1-3 приведены графики восстановленных решений как для точных, так и для возмущенных данных для трех модельных решений. На всех рисунках горизонтальная ось — номера узлов равномерной сетки по аргументу х.

Анализ численного моделирования показывает, что предложенный алгоритм достаточно хорошо восстанавливает решение с изломом (разрыв в производной) и разрывом первого рода самой функции, сохраняя структуру исходной модели. Несколько хуже алгоритм работает для более сложной модели 3. По-видимому, в данном случае для повышения качества решения необходимо использовать вместо и0 = 0 начальное приближение, учитывающее некоторые особенности решения, а также более тщательно выбирать управляющие параметры а,

Заключение

В работе предложен регулярный алгоритм, который позволяет восстанавливать негладкие решения интегральных уравнений первого рода. Алгоритм основан на тихоновской регуляризации с недифференцируемым стабилизатором (в форме нормы пространства Липшица) с привлечением ргох-метода и субградиентных процессов для решения задачи негладкой минимизации. Численные эксперименты подтверждают, что использование ргох-метода существенно повышает качество приближенного решения, поскольку дает возможность применять субградиентные процессы для сильно выпуклой функции, что обеспечивает устойчивость счета.

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

[1] Леонов A.C. Кусочно-равномерная регуляризация двумерных некорректных задач с разрывными решениями. Численный анализ // Журн. вычисл. математики и мат. физики. 1999. Т. 39, № 12. С. 1939-1944.

[2] Vogel C.R. Computational Methods for Inverse Problems. Philadelphia: SIAM, 2002. 183 p.

[3] Васин B.B., Серёжникова Т.И. Двухэтапный метод аппроксимации негладких решений и восстановление зашумленного изображения // Автоматика и телемеханика. 2004. № 2. С. 126-135.

[4] Vasin V.V., Korotkii M.A. Tikhonov regularization with non-differentiable stabilizing functionals // J. Inv. Ill-Posed Problems. 2007. Vol. 15, No. 8. P. 853-865.

[5] Васин B.B. Устойчивая аппроксимация негладких решений некорректно поставленных задач // Докл. РАН. 2005. Т. 402, № 5. С. 586-589.

[6] Тихонов А.Н., Васильев Ф.П. Методы решения некорректных экстремальных задач // Banach Center Publ. 1976. Vol. 3. P. 297-342.

[7] Васин B.B. Аппроксимация негладких решений линейных некорректных задач // Тр. Ин-та математики и механики УрО РАН. 2006. Т. 12, № 1. С. 64-77.

[8] Moreau J.J. Proximité et dualité dans un espace Hilbertien // Bull. Soc. Math. France. 1965. Vol. 93, No. 2. P. 273-299.

[9] Martinet B. Determination approachee d'un point fixe d'ne applications pseudo-contractante // C. R. Acad. Sei. Paris. 1972. Vol. 274. P. 163-165.

[10] Васин В.В., Агеев А.Л. Некорректные задачи с априорной информацией. Екатеринбург: УНФ "Наука", 1993.

[11] Поляк Б.Т. Введение в оптимизацию. М.: Наука, 1983.

[12] Васин В.В., Еремин И.И. Операторы и итерационные процессы фейеровского типа. Теория и приложения. Москва; Ижевск: НИЦ РХД, 2005.

[13] Rockape llar R.T. Monotone operators and the proximal point algorithm // SIAM J. Control and Optimizat. 1976. Vol. 14, No. 5. P. 871-898.

[14] Васин B.B. Проксимальный алгоритм с проектированием в задачах выпуклого программирования. Свердловск, 1982 (Препр. Ин-та математики и механики УНЦ АН СССР).

[15] Лаврентьев М.М., Романов В.Г., Шишатский С.П. Некорректные задачи математической физики и анализа. М.: Наука, 1980.

Поступила в редакцию 30 апреля 2009 г., с доработки — 20 ноября 2009 г.