Оценка погрешности численного метода решения одной обратной задачи Текст научной статьи по специальности «Математика»

УДК 517.927.2

ОЦЕНКА ПОГРЕШНОСТИ ЧИСЛЕННОГО МЕТОДА РЕШЕНИЯ ОДНОЙ ОБРАТНОЙ ЗАДАЧИ

Рассмотрен линейный дифференциальный оператор и система краевых условий, задаваемая линенйыми в пространстве п раз непрерывно дифференцируемых функций линейно-независимыми функционалами. Функция Грина для краевой задачи, определенной этим оператором и упомянутыми функционалами, строится как решение интегрального уравнения Фредгольма II рода, параметры которого определяются функцией Грина вспомогательной задачи. Предложенный метод обращения дает возможность эффективно решить как прямую (т.е. задачу нахождения решения), так и обратную (т.е. задачу нахождения правой части уравнения по экспериментально полученному решению) задачи. Обсуждены особенности численной реализации метода и возможности оценки точности полученных решений.

Ключевые слова: кравевая задача, интегральные уравнения, функция Грина.

Введение

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

Многие подобные задачи могут быть сформулированы как линейные краевые задачи:

( L[x] = x(n) + Pn-ix(n l) +-+ Pix' + pox = f(t),

\ Uj(x) = a,j = 1,2,...,n,

где Pi(t),f(t) - непрерывные на [a,b] функции, aj - числа, Uj(x) - линейные, линейно-

Краевая задача (1) — это прямая задача нахождения х(Ь) при заданной f (Ь). Задачу нахождения правой части f (Ь) по экспериментально измеренной функции Х(Ь) будем называть обратной задачей.

1. Метод обращения

Одним из наиболее эффективных методов решения обратной задачи является метод обращения дифференциального оператора (1) с помощью функции Грина, использующий хорошо известное соотношение для решения полуоднородной краевой задачи:

ь

В.И. Заляпин, Ю.С. Попенко, Е.В. Харитонова

независимые функционалы.1

(2)

a

где С(Ь,г) - функция Грина этой задачи.

хБез ограничения общности можно считать, что а^ = 0.

Содержательной частью описываемого метода является построение функции Грина G(t,T), осуществляемое в два этапа - сперва строится функция Грина G(t,T) вспомогательной задачи

х(П = f (t), (3)

Uj (x) = 0, j = 1, 2,... ,n,

которая может быть найдена непосредственно по определению с использованием фундаментальной системы решений соответствующего однородного уравнения.. На втором этапе численно отыскивается функция Грина основной задачи G(t,T), оказывающаяся решением интегрального уравнения Фредгольма II-го рода:

b

G(t,s) — G(t,s) = j G(t,T)V(T,s)dT. (4)

a

Детали и необходимые подробности можно найти в [2, 3].

На основании этих соображений был построен алгоритм решения обратной задачи теории, и написана компьютерная программа с использованием пакета Mathematica 8.0. 2

Предложенный подход поддается распараллеливанию, которое было реализовано [4] с использованием языка Си++ и стандарта MPI-2 (Message Passing Interface). Программа для высокопроизводительного вычислительного кластера «СКИФ Аврора> суперкомпью-терного центра ЮУрГУ (8832 вычислительных ядра, 117,64 триллиона операций в секунду) находится в стадии разработки. Подробности в [4].

2. Регуляризация

Важным вопросом, возникающим в контекстсте описанной выше процедуры, является вопрос о точности получаемого решения. Численная реализация метода предполагает последовательное решение двух интегральных уравнений - уравнения (4) для нахождения функции Грина С(Ь,г) и уравнения (2) для нахождения функции f (Ь). И если первая из задач особых проблем в плане точности получаемого численного решения перед исследователем не ставит - уравнение (4) является уравнением Фредгольма II рода, то вторая относится к классу задач, решения которых неустойчивы к малым изменениям исходных данных. Поскольку правая часть рассматриваемого уравнения - экспериментально определяемая функция х(Ь), а ядро - результат применения численных процедур решения уравнения (4), то для нахождения решений уравнения (2) необходимо так построить процедуру решения, чтобы погрешности в определении ядра и правой части не сильно портили ее решение -рассматриваемую задачу необходимо регуляризовать.

Мы использовали метод регуляризации А.Н. Тихонова [5], сводящий задачу решения уравнения (2) к задаче минимизации функционала

Ма = || Аf - Х||2 + а ■ О2, (5)

где - О2 стабилизатор, взятый в виде

ь

О2 = ^ р(Ь)/2(Ь) йЬ, р(Ь) > 0,

а

а - параметр регуляризации, который определяется из условия минимума невязки и зависит от точности задания правой части уравнения (2).

2 Свидетельство о государственной регистрации № 2012619481.

Напомним, что псевдорешением операторного уравнения А/ = х называется решение /, минимизирующее невязку

х = а^т£ ЦА/ — х||.

Уравнение может иметь несколько псевдорешений. Обозначим через Q совокупность всех его псевдорешений, и пусть ф некоторый элемент из Q. Нормальным относительно ф решением или ф - нормальным решением операторного уравнения называется псевдорешение /о с минимальной нормой ||/ — ф||, то есть такое, что

Ш — ф1 = ^ V — ф1.

1

Если ф = 0, то /о называется просто нормальным решением.

Рассмотрим операторное уравнение

А/ = х, / е Г, х е X,

X и Г - гильбертовы пространства, А - линейный вполне непрерывный оператор из Г в X. Если вместо точных х и А известны их приближения х и А, такие, что

1х — хЦх < 5,

А — А

< е,

то фактически решается уравнение

А-/ = х, / е Г, х е X.

Известно [5], что задача нахождения элемента /а, на котором функционал (5) достигает минимального значения, имеет единственное решение. Это решение является нормальным решением, то есть при точных А и х из всех точных решений операторного уравнения в методе Тихонова автоматически выбирается нормальное решение. Если же 5 = 0 и/или е = 0, то метод дает решение /а, близкое к нормальному. Известные на сегодня оценки точности получены сравнением численного решения /а и нормального решения /о

Аа/ =

3. Оценивание точности по Винокурову В.А.

Наиболее удачными (с теоретической точки зрения) на этом пути нам представляются оценки, полученные В.А. Винокуровым [6].

Подробнее: пусть (С = А* А, а (С + - псевдообратный оператор. Тогда

А а

II А/а | < — + 7----------- 11/0 I , (6)

21/5 (а + |#|

где А = 5 + е I/1|^ .

В некоторых ситуациях оказывается полезной также оценка относительной ошибки решения. Она имеет вид:

| А/а

где

Е(а) =

А

< Е(а), (7)

(5оЬп + еоЬп) + Ра

л/а ра + 1

2

причем

р

С *

Л*

2 л — А

, лоЬп — и и , єоіи

Ж

При А ^ 0 ||А/а|| ^ 0, если только выполнено а(А) = 0(А2).

Одновременно, полученные В.А. Винокуровым оценки дают возможность получить оптимальное значение параметра регуляризации ар из соотношения

а

||Л||(Лоіп + Єоіи)

4р

(1 + ра)3.

(8)

Для нахождения решения уравнения (8) можно воспользоваться итерационными процедурами либо приближенным соотношением

аорі

||Л||(Лоіп + ^оьП

4р

1 І 4 ( ІІЛІІ(Лоіп + £оы)

+ з р\ 4р

(9)

которое справедливо, если только

р

||Л||(Лоіп + єоЬп) 4р

<< 1.

Эффективные с теоретической точки зрения соотношения (6) - (9) мало пригодны для практической оценки точности получаемых решений.

4. Эмпирическая оценка погрешности 4.1. Дискретизация

Зададим на промежутке [а; Ь] сетку а < Ьо < Ь\ < ... < Ьп = I. Положим /(Ь) Un(t),

где

п.

.и — £ пит

і=0

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

0 , Ь < Ьг-1 и Ь > ^г+1 ,

1 1%_1 и_1 < ь < и,

\і(і) — <

Иі Иі+1

І — Іі+1

Іі < І ^ Іі+і.

- Ьг+1

При этом функционал Ма(/) заменится функцией Ма(ип(Ь)) = Ма(ио,Щ, ■■■ип), где и'П значения функции ип(Ь) в узлах:

2 / п \ 2

Ма —

С ь ( г ь и \ 2 г ь ( и \

/ I / С(Ь,г)^п'ПХі(г)йт — Ж(і) I (М + а / р ^^^п'ПЛі(т)1 йт.

■’а Vа і=і У -’а \і=і У

Обозначив через 3%(£) интеграл ^ С(і,т)Лі(т)йт, получим:

Ма —

/ ... ЗоЗп ёо п п ^ ^ ^2 — 2 у пи^і(і) ж + ж

(пп...пп) йі+

V ЗиЗо ... Зи е е п = у

є

ь

її (РЛ2) (т -. її (РЛ0Ли) йт ' п

— ^ Міі пипи

- її (РЛпЛ0) йт ... її (рЛП) йт _ V пП ) і,і=0

+а (пП ... пП)

Необходимое условие минимума функционала Ма запишется в виде

-2^2 ЖпП+Цх"2

і=0

— > МізпП — 2Ш3 — 0, 8 — 0,1, 2, ...п.

диП ^

5 г=0

В итоге получаем систему п + 1 линейных уравнений

п

Мг5иП = 2Шз, в = 0, п,

г=о

которая и подлежит решению.

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

Ма(ип(Ь)) = Ма(и0,и1, ■■■ип) ^ ш1п

прямым методом Хука - Дживса.

4.2. Оценивание точности решения уравнения (4)

Ядро уравнения (2) определялось в результате решения интегрального уравнения (4), являющегося при каждом фиксированном значении переменной а < Ь < Ь уравнением Фред-гольма второго рода относительно функции С(Ь,т) = Ф*(т):

ь

Ф*(в) - С(Ь,в) = ( Фг(т)У(т,в)1т,

ядро которого У(в,т) определяется функцией Грина С(Ь,т) вспомогательной задачи. Положив при каждом фиксированном значении а < Ь < Ь

а = то < Т1 <...<Тп = Ь, Ф* = Ф(Ь,тг), С\ = &(Ь,Тг),

пп

Ф*(т) = 5>&(т)М,Т) = £ С*Хг(т)■

получим

і=0

і=0

Обозначая

V- — Лі(,в)У(в, ті)(в, і,і — 1,2, ...,п,

а

и подставляя в (4), приходим к дискретному аналогу рассматриваемого уравнения

и

(Лі (ті) — VI) — С, і — 0,1,...,п,

і=0

ь

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

(I - V)Ф* = Сг,

где

/ Ф? N / С? \ / 1 1 °1 V?1 . . V?” \

Фг = Ф1 , С г = С} I 1 = V0 1 -VI1 . . V”

V ФП V СП V V? VI . . 1 - V” )

Поскольку задача решения уравнения Фредгольма второго рода корректна, и исходное интегральное уравнение (4) однозначно разрешимо, то для каждого значения N соответствующая дискретная задача однозначно разрешима. При этом конечномерные решения сходятся к точному решению уравнения.

Пусть СN(Ь,т) - решение уравнения (4), полученное при п = N. Из вышеизложенного следует, что точность нахождения ядра уравнения (2) из уравнения (4) может быть оценена разностью двух соседних приближений для достаточно больших значений N:

АС = \\См+1 — См У,

что дает возможность оценить величину е - точность задания оператора А, поскольку:

АА = || А - А\\ =АС < е.

4.3. Оценивание точности решения уравнения (2)

Основные соображения, положенные в основу метода оценки точности решения уравнения (2), состоят в следующем: параллельно с задачей решения уравнения (2) относительно неизвестной функции / решается по тому же алгоритму модельное уравнение

ь

! С(г,т)1'(т)йт = ХтосС^)- (10)

а

Пусть - /тос - известное точное решение уравнения (10) для точно известной правой части хтос, /тос - решение, полученное регуляризацией. Точность метода оценивается величиной

АМ/ = ||/тоС /тос||.

В случае неточного задания ядра и/или правой части уравнения (2) поступаем следующим образом: возмутим точное решение уравнения (10), положив /регг(^ = /тос + £(£)> и построим возмущенную правую часть для уравнения (10), положив Хтос = хтос + ц(Ь), где П(Ь) дается соотношением

ь

п(г) = I С(г,т)£(т)(1т.

а

Регуляризованное решение уравнения (10) с правой частью Хтос обозначим через /регг и сравним его с /тоС, положив Аз/ = \\/тоС - 1Рен\\.

Модельная задача может быть построена многими способами. Например, если есть априорная информация о структуре и поведении отыскиваемого решения / (£), то в качестве хтос(^) следует взять функцию

ь

Хтос(^) = J С(ь,т )д(т )йт,

а

где g(t) - аналог (в каком-то смысле) функции f (t).

Если априорная информация отсутствует, то на грубой сетке (п не очень велико) отыскивается регуляризованное решение уравнения (2) и в качестве аналога g(t) функции f берется сглаженное (например - сплайнами невысокого порядка) регуляризованное решение.

Описанная выше эмпирическая методика оценивания точности показала свою эффективность при решении различных прикладных задач.

Литература

1. Грановский, В.А. Динамические измерения: Основы метрологического обеспечения / В.А. Грановский. - Л.: Энергоатомиздат. Ленингр. отд-ние, 1984.

2. Zalyapin, V.I. Inverse Problems of the Measurements Theory / V.I. Zalyapin,

H.V. Kharitonova, S.V. Ermakov // Inverse problems, Design and Optimization Symposium.

- Miami, Florida, U.S.A., 2007. - P. 91-96.

3. Асфандиярова, Ю.С. Метод интегральных уравнений построения функции Грина / Ю.С. Асфандиярова, В.И. Заляпин, Е.В. Харитонова // Вестник ЮУрГУ. Серия: Математическое моделирование и программирование. - 2012. - № 27 (286), вып. 13. - С. 16-23.

4. Асфандиярова, Ю.С. Численный анализ обратной задачи теории измерений / Ю.С. Асфандиярова // Тр. 53-й науч. конф. МФТИ «Современные проблемы фундаментальных и прикладных наук>. - М., 2010. - Ч. VII, т. 3. - С. 6-7.

5. Тихонов, А.Н. Методы решения некорректных задач / А.Н. Тихонов, В.Я. Арсенин. -М.: Наука, 1986.

6. Винокуров, В.А. О погрешности решения линейных операторных уравнений / В.А. Винокуров // Журнал вычислительной математики и математической физики. - 1970. -Т. 10, № 4. - С. 830-839.

7. Menikhes, L.D. On the Convergence of Approximations of the Regularization Method and the Tikhonov Regularization Method of n-th Order / Menikhes L.D., Tanana V.P. // J. Inv. and Ill-Posed Problems. - 1998. - V. 6, №3. - P. 241-262.

Владимир Ильич Заляпин, кандидат физико-математических наук, профессор, кафедра «Математический анализ>, Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), [email protected].

Юлия Сагитовна Попенко, ассистент, кафедра «Математический анализ>, ЮжноУральский государственный университет (г. Челябинск, Российская Федерация), [email protected].

Елена Владимировна Харитонова, кандидат физико-математических наук, доцент, кафедра «Математический анализ>, Южно-Уральский государственный университет (г. Челябинск, Российская Федерация), [email protected].

Bulletin of the South Ural State University. Series «Mathematical Modelling, Programming & Computer Software>,

2013, vol. 6, no. 3, pp. 51-58.

MSC 34B27

Error Estimate of Numerical Method for Solving an Inverse Problem

V.I. Zalyapin, South Ural State University, Chelyabinsk, Russian Federation,

[email protected],

Yu.S. Popenko, South Ural State University, Chelyabinsk, Russian Federation,

[email protected],

Ye.V. Kharitonova, South Ural State University, Chelyabinsk, Russian Federation,

[email protected]

Linear differential operator and the system of boundary conditions were considered. The boundary conditions are linear and linear independent functionals. The Green functions for the boundary problem defined by this operator and the functionals was build as a solution of the Fredholm integral equation of the second kind. Characteristics of the Fredholm equation was defined by the Green function of the auxiliary problem. The suggested method enables to solve both direct (the problem of finding solutions) and inverse (the problem of finding the right part of the equation from the experimentally obtained solution) problems. The characteristics of the numerical implementation of the method and the possibility of assessing the accuracy of the solutions were discussed.

Keywords: boundary problem, integral equations, Green’s function.

References

1. Granovsky V.A. Dynamic Measurements: The Basics Metrological Support. Leningrad, 1984. (in Russian)

2. Zalyapin V.I., Kharitonova H.V., Ermakov S.V. Inverse Problems of the Measurements Theory. Inverse Problems, Design and Optimization Symposium,. Miami, Florida, U.S.A., 2007, pp. 91-96.

3. Asfandiyarova Yu.S., Zalyapin V.I., Kharitonova Ye.V. The Method of the Integral Equations to Construct the Green’s Function. Bulletin of the South Ural State University. Series «Mathematical Modelling, Programming & Computer Software», 2012, no. 27 (286), issue 13, pp. 16-23. (in Russian)

4. Asfandiyarova Yu.S. Numerical Analysis of the Inverse Problem of Measurement. Proceedings of the 53-rd Conference MIPT «Modern Problems of Fundamental and Applied Science». Part VII, 2010, vol. 3, pp. 6-7.

5. Tikhonov A.N., Arsenin V.Ya. Metody resheniya nekkorektnykh zadach [Methods for Solving Ill-Posed Problems]. Мoscow, Nauka, 1986.

6. Vinokurov V.A. On the Error of Solutions of Linear Operator Equations. Zhurnal vychislitel’noy matematiki i matematicheskoy fiziki [Computational Mathematics and Mathematical Physics], 1970, vol. 10, no. 4, pp. 830-839. (in Russian)

7. Menikhes L.D., Tanana V.P. On the Convergence of Approximations of the Regularization Method and the Tikhonov Regularization Method of n-th Order. J. Inv. and Ill-Posed Problems. 1998, vol. 6, no. 3, pp. 241-262.

Поступила в редакцию 15 мая 2013 г.