Содержание к диссертации
Введение
Глава I. Рк-операторы и влияние процедуры регуляризации для интегральных уравнений I рода типа свертки 12
1. О некоторых обратных задачах математической физики 13
2. Регуляризирущие операторы для уравнений типа свертки 20
3. Асимптотические оценки уклонения регулярне ованного решения от точного 31
Глава II. Оптимальные Р-алгоритмы поиска приближенных решений 43
1. Оптимальный выбор параметра регуляризации 44
2. Оптимальные регуляризованные решения и фильтры, определяемые интегральными операторами свертки 53
3. Определение вероятностных характеристик сигнала и шума для оптимальных Рк-решений 63
4. Последовательное оценивание решений и статистическая регуляризация 68
Глава III. Р-операторы в обратной задаче гравиметрии для нескольких контактов 76
I. Оценка числа контактных поверхностей 77
2. Определение параметров геометрически подобных контактных поверхностей . 84
Заключение 89
Литература 91
Приложение 100
- Регуляризирущие операторы для уравнений типа свертки
- Оптимальные регуляризованные решения и фильтры, определяемые интегральными операторами свертки
- Определение вероятностных характеристик сигнала и шума для оптимальных Рк-решений
- Определение параметров геометрически подобных контактных поверхностей
Введение к работе
Диссертация посвящена разработке численных методов решения обратных задач математической физики,сводящихся к линейным интегральным уравнениям I рода, и созданию эффективных алгоритмов решения таких уравнений.
Актуальность исследования обусловлена тем,что рассматриваемые задачи возникают при моделировании важных для научно-технического прогресса проблем восстановления сигналов, искаженных приборами или окружающей средой, идентификации линейных систем,интерпретации геофизических наблюдений и зд.
Функциональная связь между характеристиками изучаемого объекта I и результатами наблюдений и* задается оператором А таким образом, что математическая модель об -ратной задачи описывается уравнением I рода
//j = и. , <"
где решение г и правая часть 1Л, принадлежат заданным метрическим пространствам:
Во многих случаях задача решения уравнения (I) является некорректно поставленной. Основы общей теории и методов решения некорректных задач были заложены в трудах А.Н.Тихонова [?4- 18 ], развиты в работах М.М.Лаврентьева, В.К.Иванова и многих других советских ученых. Подробная библиография приводится в монографиях \tt,i9,ьі, si, rs", Щ ЦtoyЩ
В частности, упомянутые задачи сводятся к линейным интегральным уравнениям Фредгольма I рода
Г*
А) - jW^)jfs)^S^^),bf^Jl; (2)
и их важному подклассу - линейным интегральным уравнениям I рода типа свертки
Аг н \М-$%Шь г и.#). (3)
w - О
Поиск решения \ уравнения (I) в случае, когда пра -вая часть и, известна приближенно, связан с рядом труд -ностей.Так, приближение к і , найденное по ІА клас -сическими численными методами (типа наименьших квадратов) , может как угодно сильно отличаться в заданной метрике про -странства У от точного решения, вследствие чего возникает проблема устойчивости приближенного решения относительно малых изменений исходных данных. Отсутствие в реальных задачах большой априорной информации для построения прибли -женных решений с оптимальными характеристиками (необходи -мой, например, для метода винеровской фильтрации) требует исследования возможностей уменьшения а также наиболее полного использования такой информации. Наконец, программная реализация автоматической обработки результатов наблюдений в условиях большого объема данных и многократного повторения процедуры вычислений при разных значениях параметров вызывает необходимость создания регуляризирующих алгоритмов для уравнений (1-3), эффективных по использованию вычислитель -ных ресурсов - памяти ЭВМ и времени счета.
Научная новизна и структура диссертационной работы
В работе предложены и исследованы новые регуляризирую-щие операторы к-го порядка для некорректно поставленных задач (1-3), изучены их свойства и построены эффективные устойчивые алгоритмы, реализованные в виде программ для ЭВМ .
Регуляризирущие операторы для уравнений типа свертки
В работе предложены и исследованы новые регуляризирую-щие операторы к-го порядка для некорректно поставленных задач (1-3), изучены их свойства и построены эффективные устойчивые алгоритмы, реализованные в виде программ для ЭВМ .
Глава I посвящена разработке специального класса регу-ляризирующих операторов , зависящих от натурального параметра "к" и называемых Рк-операторами, для интегральных уравнений типа свертки (3).
В первом параграфе вводятся необходимые для дальнейшего изложения определения и приводится ряд примеров, связанных с математической обработкой и интерпретацией данных физических измерений.
Во втором параграфе разработана новая конструтада Рк-операторов на основе метода интегральных преобразований В.Я.Арсенина [ ,??] , не требующая дополнительных ресурсов памяти ЭВМ. Ранее Рк-операторы исследовались В.А.Морозовым, П.Н.Заикиным, И.Ф.Дорофеевым Щ з.б/.НЧІ. В диссертации показана связь с этими подходами и вариационная интерпретация Рк-решений. Впервые изучено поведение невязки Рк-решения как функции параметра регуляризации оО в зависимости от порядка регуляризации "к" и порядка стабилизирующего функ -ционала. Выяснено, что разным "к" соответствуют непересекающиеся семейства асимптотически -близких стабилизирующих множителей.
В третьем параграфе получены новые асимптотические оценки уклонения Рк-решения от точного при об - О характеризующие влияние процедуры регуляризации. Рассмотрены два важных для приложений класса уравнений (типа I и II), определяемых асимптотическим поведением фурье-образов ядер [ . ?1 При к = 0 результат совпадает с оценками из [ ] , полученными другим способом, при о 4 с і» и дополнительных предположениях о гладкости точного решения увеличивается порядок по оС ( oL - о ) стремления к нулю уклонения приближенного решения от точного.
В главе II выполнены исследования оптимальных регуляри-зирующих операторов к-го порядка, связанные со статистическим подходом к описанию исходных данных для уравнений (1-3). В I анализируется случай, когда возмущение Ь правой части и, представляет реализацию стационарного случайного процесса со спектральной плотностью вида и имеет, как и всюду в дальнейшем, аддитивный характер Точное значение правой части Ц является детерминирован -ной функцией. Выведены новые асимптотические формулы (при сС +о ) для дисперсии влияния р на Рк-решение уравнений I и II типов. Для конечных порядков "к" определены, следуя [?, ?], значения параметра об , близкие к (С, f)-оптимальным, ко -торые приводят к наилучшим Рк-решениям с заданным стабили -зирующим множителем і . В 2 исследуются оптимальные Рк-операторы в винеров-ской постановке задачи о фильтрации стационарных процессов (5 /,13 Для стационарного процесса \Л с нормальным рас -пределением вероятностей построен Рк-оператор,задающий оптималыши линейный фильтр. При построении используются спектральные представления случайных процессов и метод регрессии [i3.so"i . Приведена оценка роста соответствующего опти -мального параметра регуляризации при возрастании "к". Получено р-оптимальное значение параметра регуляриза -ции для Рк-операторов при использовании простейших стабилизаторов р-го порядка для I типа ядер уравнений, имеющих точное решение со спектральной плотностью и исследованы свойства погрешности Рк-решения при заданном "к" как функции от параметра "р". 3 посвящен определению вероятностных характеристик Sh , CL шума и S? , і сигнала для априори заданных соотношений (4) и (6) по семейству" Рк-решений при разных значениях параметра регуляризации об согласно методу В.Я. Арсенина С5"" ]. Построен и исследован соответствующий функционал, зависящий от натурального параметра "к". Показано, что для эргодических процессов метод приводит к уменьшению априорной информации при построении оптимального ли -нейного фильтра, как и в [6 ] . Изменение порядка "к" регу -ляризации позволяет варьировать области определения неиз -вестных параметров и уточнять их оценки. В 4 рассмотрен вопрос использования Рк-алгоритмов для последовательного оценивания приближенных решений интегрального уравнения (2) путем перехода к дискретному аналогу задачи - системе линейных алгебраических уравнений. Правая часть задается в виде вектора наблюдений с нормально распределенной ошибкой, имеющей нулевое математическое ожидание и единичную матрицу ковариации. Исследование проводится с использованием байесовского подхода к идентификации вектора состояния линейной системы методом Е.Л.Жуковского и В.А.Морозова [4v4s\6.t. Показано, что каждый вектор Р: -решения алгебраической системы, j = 1,...,к , является нормальным решением отно -сительно предыдущего вектора Е 1 -решения. Этот факт служит основанием для гипотезы об априорной плотности вероятности вектора состояния при переходе к следующему порядку "к" ре-гуляризирующего оператора, а также для получения байесов -ских оценок вектора состояния. Установлено, что применение Рк-операторов для получе -ния оптимальной оценки ненаблюдаемой компоненты марковского процесса по наблюдаемой позволяет уменьшить число наблюде ний в (к+1) раз для совпадения в рамках критерия % с результатом фильтрации Калмана-Бьюси при к = 0.
Оптимальные регуляризованные решения и фильтры, определяемые интегральными операторами свертки
Учитывая монотонное стремление к нулю icW) при больших об , получаем существование точки локального максимума сіц, , расположенной правее d . Доказательство закончено.
Теорема I. Пусть случайные процессы сигнала J" и шума Ь удовлетворяют условиям, сформулированным в 2, и в частности выполняются соотношения (2.12) при некоторых значениях параметров. Тогда для каждого значения "к" по семейству Рк-операторов может быть найден оптимальный (или близкий к нему) линейный фильтр для Т из уравнения (2.1).
Доказательство. Составим таблицу Рк-реше-ний для уравнения (2.1) при некотором значении "к" с достаточно мелким шагом по параметру регуляризации оО . Поль -зуясь соотношением (I), можно найти значения функции %№ ) в заданном диапазоне изменения аргумента d, , определить точки локальных экстремумов d1 d . При этом на участке (С, Uл) можно определить (например, методом наименьших квадратов) параметры шума Ъ , %, , так как для функции определяющим здесь будет слагаемое ,$ .Поэтому интервал (о,1 ) можно назвать областью влияния характерис тик шума на Рк-решение. Аналогичным образом, интервал (U ol с определяющим слагаемым Ci Sr- можно назвать об ластью влияния характеристик сигнала. Рассмотрим изменение точки локального минимума в зависимости от порядка "к". Вследствие соотношения, полученного в утверждении I, для («&") можно написать Теперь правую часть равенства разложим в произведение двух множителей, один из которых - Й не зависит от "к" : Таким образом, в зависимости от соотношения между параметрами модели d и к+( точка локального минимума может сме -щаться влево или вправо и тем больше, чем больше "к".То есть, с помощью параметра "к" можно менять зоны влияния шума и сигнала, при этом на скорость стремления к нулю функции УїсМ при больших dj модификации существенного влияния не оказывают. Это обстоятельство может быть использовано для уточнения оценок неизвестных параметров S а, S- 4 и проверки результатов счета.
Итак, в рамках предполагаемой модели шума и сигнала мы можем определить значения соответствующих параметров и затем использовать их для оптимальной линейной фильтрации, напри -мер согласно формуле (2.7). Использование аппарата Рк-решений позволяет при этом существенно уменьшить априорную информацию об искомом решении и шуме. Теорема доказана.
Заметим, что в случае эргодических случайных процессов вычисление функции % (оЬ) может осуществляться более прос -тым способом (#1 ] по формуле
В предыдущих параграфах мы рассматривали оптимальные Рк-алгоритмы решения интегральных уравнений типа свертки в предположении, что нам доступна информация о спектральных плотностях шума и искомого решения. Сейчас же рассмотрим интегральное уравнение Фредгольма I рода[Дг 7 гЫг. л,??]
Рк-решение ищем в прежнем предположении аддитивности возмущения правой части уравнения (I) u(-i)-u,(i)+ fl) и некоррелированности с точным решением Т() , но с учетом другой априорной информации. Будем считать известными не спектральные плотности решения и шума, а их законы распределения вероятностей. Принимается, что J и ? представляют реализации гауссовских случайных процессов, причем g имеет нулевое математическое ожидание и единичную дисперсию , то есть у « ы(о,4).
При численном решении уравнения (I) необходимо провес -ти его разностную аппроксимацию, в результате которой мы получаем систему линейных алгебраических уравнений
Такой переход к дискретному аналогу задачи можно провести по-разному ц р /, однако важно отметить факт сходимости регуля-ризованных решений системы (2) к соответствующему решению (I).
Таким образом, вместо уравнения (I) мы будем рассматривать непосредственно систему (2). При этом статистический характер данных позволяет связать эту алгебраическую систему с такими физическими явлениями, когда соответствующие эксперименты или наблюдения могут быть повторены большое число раз при одинаковых условиях. Может быть использована следующая стохастическая модель С4ъ] . На вероятностном пространстве ( X, Л, Я) задается двумерная марковская цепь (jU, ) , і - О,... , которая управляется рекуррентными соотношениями так что IL). - наблюдаемая, а 3±, - ненаблюдаемая компо -ненты. Равенства понимаются в смысле стохастической эквивалентности случайных величин.
Обратимся к задаче нахождения оптимальной в среднеквад-ратическом смысле оценки fy компоненты 1 процесса (3). Известно Г 1 , что, если л/ - число проведенных изме -рений, то AV совпадает с регуляризов энным решением 3 №} системы линейных уравнений (2) со значением параметра регуляризации об равным У/-/ . При этом тихоновский слу -чайный функционал имеет вид
Определение вероятностных характеристик сигнала и шума для оптимальных Рк-решений
Так, при использовании метода невязки приходится решать уравнение вида (1.2.19) Если при этом используется какой-либо итерационный алгоритм, то на каждом шаге по необходимо вычислять значение функции невязки Т ЦЛ , "- )2V.- . С подобной проблемой мы встречаемся и при расчете высокочастотных характеристик сигнала и шума в главе II , когда ищется функция (2.3.1)
Наконец, и при решении линейных алгебраических систем уравнений (2.4.2) в самых разных ситуациях требуется находить Р-решения не для одного, а для целого ряда значений пара -метрам Ь.\\» Использование Рк-алгоритмов в таких однотипных задачах и позволяет получить одновременно несколько расчетных значений такой функции.
Пусть Ф(о(.) - заданная функция от параметра регуля -ризации, и ищется ряд ее значений. Вычисляем последовательность Рк-решений 3 /3/ 1 " ) fitfe) и по ней находим фр")( $ЛО, ,,. PK(cL) ДЛЯ ОДНОГО значения об . Согласно результатам для функций Р(и) , полученным в главах I, II , существует "к" значений параметра d \ о ,..., U таких , что РЫГ ) Pf-C ) і t -i,--, К- т При этом вычисление іф( г") ;ґ=І,...,к: \ требует К-У(#,//) операций (функция V (ЯУ) определена в главе I,стр.42), в то время как для вычисления Pf.(+), r l --у16 \ необходимо лишь V(#,/t/9+vcV операций. Ясно, что в некоторых задачах за счет этого может быть существенно уменьшен объем вычислений. Например, для поиска корня функции невязки f fa) можно использовать метод деления отрезка пополам. Тогда при -менение Рк-операторов позволяет на каждом шаге метода получать "к" дополнительных значений невязки, расположенных по убыванию. Таким образом, мы приближаемся к корню со стороны превышающих его значений используя Г ІЗ. 1 дополнительных подшагов (здесь Г 1 обозначает целую часть числа). В случае вычисления функции («О применение Рк-алгоритма позволяет в "к" раз увеличить получаемую информацию в виде найденных значений функции, что значительно об -легчает дальнейшее использование методов типа наименьших квадратов для определения параметров сигнала и шума. Итак, Рк-операторы являются эффективным средством при численном решении целого ряда задач. В практической работе за критерий выбора параметра "к" можно принять условие 1&С«0 - Фк-/Иіи " » где " -" устанавливается из дополнительных соображений о точности исходных данных. Геофизические методы исследования земных недр основаны на изучении физических полей, связанных со свойствами земных пород. Гравиразведка использует измерения гравитационного поля. Изменения (аномалии) силы тяжести вдоль земной поверхности на фоне среднего (нормального) значения, обусловленного " . полем Земли в целом, несут информацию о распределении соот -ветствующих физических параметров избыточных масс. Таким образом, выяснение свойств геологических объектов по косвенным наблюдениям является специфической особенностью геофизических исследований вообще и гравиметрии в частности. Задачи такого рода называются обратными. Если решения прямых задач построения локальных аномалий над залегающим телом получаются в виде интегралов, то реше -ния обратных задач сводятся к интегральным уравнениям I рода. Так как телам, относительно сильно отличающимся по плотности, могут отвечать близкие наблюдаемые аномалии, и поля измеряются всегда с некоторой погрешностью, то следует использовать устойчивые алгоритмы приближенного решения обратных задач. Итак, некорректность задач, стремление получить наиболее полную информацию из результатов эксперимента приводят к необходимости применения метода регуляризации с использова -нием средств современной вычислительной техники[І-Г. -ЇГ, . Существенную часть гравитационной разведки полезных ископаемых составляет построение геолого-геофизической модели исследуемой области и ее анализ на ЭВМ с помощью математи -ческих методов. Если учитывается лишь элементарная связь между локальными источниками и характеристиками наблюдаемого гравитационного поля, то модель может быть представлена уравнением I рода с необходимыми ограничениями на оператор Д , множество характеристик источников, и множество характеристик измеряемого поля V . Мы ограничимся двумерным случаем и будем пользоваться системой координат г. О/ с осью Oz. , направленной вниз. I. Начнем с простейшей модели для обратной задачи гра -виметрии, а затем усложним ее. Пусть функция З.Ы) X представляет границу раздела двух сред с различной плотностью масс f"= 6 w f и % - ши, соответственно. При этом контактная поверхность 1(}С) не отклоняется от некоторого горизонтального уровня А всюду кроме отрезка С 1 . Уклонение на отрезке создает на поверхности Земли аномалию напряжения силы тяжести & ? - -/ 2-12.=0 » где U- - по тенциал масс с плотностью р - - Требование постоянности плотностей является существенным, ибо известны при -меры структур с разной геометрией и различным распределением плотностей, создающих одинаковые аномалии.
Таким образом, зная потенциал на земной поверхности, мы приходим к задаче восстановления контактной границы раздела сред г(а} . Вопросам единственности решения обратной зада -чи в такой постановке посвящен ряд работ Остромогильского А.Х. (см., например,рчі ). Случай нескольких контактных поверх -ностей изучался в работах Гласко В.Б., Филатова Б.Г. [ 3- ,/3] и других.
Определение параметров геометрически подобных контактных поверхностей
Тогда обратная задача для нескольких контактов имеет единственное решение [ ъ] . Наблюдаемое поле ЯбО в этом случае можно представить как суперпозицию полей, вызванных V слоями В Региональный фон считается снятым из наблюденных данных и каждая составляющая &: ( ) в (5) возникает в результате уклонения контактной поверхности Ъ\ (уС) от горизонтального уровня It: . Для обеспечения устойчивости процесса интерпретации отрезок Tcv l , на котором осуществляются измерения поля И-60 » полагается, как обычно, в три раза превосходящим отрезок t 1 о Лї Ч 1- лі і 4- " " j и расположен -ным относительно него симметричным образом Г 6 ъ ] . 3. Пусть теперь решается задача нахождения контактных поверхностей в приведенной выше постановке с одним дополни -тельным требованием то есть формы поверхностей являются подобными некоторой од -ной функции fy} с неизвестными коэффициентами Qj . На практике в качестве такого эталона U/) могут использоваться дуги окружностей, парабол или других кривых в зависимости от конкретного района исследований,
В некоторых случаях Я (У) можно определить предвари -тельно разделяя наблюдаемое на поверхности z. - о поле W- ( ) на "полезную" составляющую - сигнал и шумовую составляющую Б Д« ,?о]. При этом, если представляет поле того контакта, который определяет главный вклад в и#) , а # ft) - влияние остальных контактов,включая и другие помехи. Методы выделения ft) основываются на каких-либо характерных признаках поля и чаще всего относятся к преобразованиям интегрального типа. В линейном случае это методы сглаживания, высших производных, или использующие решения граничных задач теории потенциала. Все они сводятся к линейным интегральным операторам типа свертки b ,?0 Wl. Если допустимо представление о контактной поверхности как о функции, являющейся реализацией стационарного эргодического процесса, то можно воспользоваться теорией оптимальной линейной фильтраііииІ .ЧЧї 57 ] Получив (V) , определяем эталонную форму ft) контактной поверхности при неизвестной глубине контакта согласно работе С? V4]. Если JZ ft) обусловлено первым слоем с известным (например, по данным бурения) значением глубины, то просто ищем Рк решение уравнения (2).
Таким образом, функцию з/W) можно считать известной и задача на данном этапе заключается в нахождении оценки числа л/ контактных поверхностей. Алгоритм определения глубин залегания слоев А; и коэффициентов $: будет рассмотрен в следующем параграфе. Сведем рассматриваемую проблему к алгебраической задаче. Для этого рассмотрим преобразование Фурье функции и.(уС)
Здесь фурье-образ и,(и ) с помощью численных методов (например,любым из алгоритмов быстрого преобразования Фурье) вычисляется по измеряемой функции (А (1С) . Фурье-образ Р (Ю) рассчитывается либо численно, когда з (\С) получается в результате решения интегрального уравнения, либо аналитически, когда Я 60 задается в виде простой параметрической кри -вой. Мы хотим разделить обе части равенства (7) на 1 ) , чтобы получить справа сумму экспоненциальных функций. Но , так как в измерениях присутствуют погрешности, рассматриваемые функции являются приближенно заданными, и мы находимся в рамках некорректных задач, то после деления можем получить результат, сильно отличающийся от истинного значения суммы в правой части. В такой ситуации можно использовать регуляризирующий алгоритм. Например, как это делается для уравнений типа свертки при переходе к фурье-образам [ . ] :
Параметр регуляризации ы, может быть найден любым из известных способов в соответствии с величиной погрешности для коэффициентов Q\ и показателей экспонент L: .Необходимо также учесть, что получаемая функция должна быть в рамках заданного уровня погрешности близка к действительной. Тогда, обозначая /« / -а г и проекцию l(u/) на действительную ось через Т?(Ґ) , получаем с учетом малости ар -гумента комплексной функции (о") соотношение