авторефераты диссертаций БЕСПЛАТНАЯ БИБЛИОТЕКА РОССИИ

КОНФЕРЕНЦИИ, КНИГИ, ПОСОБИЯ, НАУЧНЫЕ ИЗДАНИЯ

<< ГЛАВНАЯ
АГРОИНЖЕНЕРИЯ
АСТРОНОМИЯ
БЕЗОПАСНОСТЬ
БИОЛОГИЯ
ЗЕМЛЯ
ИНФОРМАТИКА
ИСКУССТВОВЕДЕНИЕ
ИСТОРИЯ
КУЛЬТУРОЛОГИЯ
МАШИНОСТРОЕНИЕ
МЕДИЦИНА
МЕТАЛЛУРГИЯ
МЕХАНИКА
ПЕДАГОГИКА
ПОЛИТИКА
ПРИБОРОСТРОЕНИЕ
ПРОДОВОЛЬСТВИЕ
ПСИХОЛОГИЯ
РАДИОТЕХНИКА
СЕЛЬСКОЕ ХОЗЯЙСТВО
СОЦИОЛОГИЯ
СТРОИТЕЛЬСТВО
ТЕХНИЧЕСКИЕ НАУКИ
ТРАНСПОРТ
ФАРМАЦЕВТИКА
ФИЗИКА
ФИЗИОЛОГИЯ
ФИЛОЛОГИЯ
ФИЛОСОФИЯ
ХИМИЯ
ЭКОНОМИКА
ЭЛЕКТРОТЕХНИКА
ЭНЕРГЕТИКА
ЮРИСПРУДЕНЦИЯ
ЯЗЫКОЗНАНИЕ
РАЗНОЕ
КОНТАКТЫ


Pages:   || 2 | 3 |
-- [ Страница 1 ] --

Государственный комитет Российской Федерации

По высшему образованию

Красноярский Государственный Университет

В.А. Кочнев

АДАПТИВНЫЕ МЕТОДЫ

РЕШЕНИЯ ОБРАТНЫХ ЗАДАЧ ГЕОФИЗИКИ

Учебное пособие

Красноярск

РОССИЙСКАЯ ФЕДЕРАЦИЯ

ГОСУДАРСТВЕННЫЙ КОМИТЕТ РОССИЙСКОЙ ФЕДЕРАЦИИ ПО

ВЫСШЕМУ ОБРАЗОВАНИЮ

КРАСНОЯРСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ

ВЫЧИСЛИТЕЛЬНЫЙ ЦЕНТР СО РАН (КРАСНОЯРСК) В.А. Кочнев АДАПТИВНЫЕ МЕТОДЫ РЕШЕНИЯ ОБРАТНЫХ ЗАДАЧ ГЕОФИЗИКИ Учебное пособие Красноярск 1993 2 УДК 550.834 ББК 26.2В К 758 Рецензенты: С.В. Гольдин- чл.-корр. РАН, профессор, зав.кафедрой НГУ А.М. Федотов- профессор, зав.кафедрой КПИ Кочнев В.А.

К 758 Адаптивные методы решения обратных задач геофизики учеб. пособие :

Краснояр. гос. ун-т;

Красноярск, 1993. 120 с.

Обосновывается и раскрывается суть адаптивного метода оценки одного параметра, решения линейных уравнений и обратных задач гравиметрии и магни тометрии. Проводится обзор литературы по адаптивному методу в теории управ ления, при прослеживании сейсмических волн;

при интерпретации сейсмических данных и при решении обратных задач геофизики. Имеются иллюстрации, при меры и упражнения.

Предназначена для студентов, спирантов и геофизиков.

Табл. 20 Ил. 11. Библиогр.: 57 назв.

180320000 Без объявл. ББК 26.2В К 6 Л 5(02) ©В.А. Кочнев, Красноярский Государственный Университет, ISBN 5-230-07988- Содержание Предисловие............................................................................................................................................ Глава 1. Введение в метод...................................................................................................................... Глава 2. Оценки одного параметра....................................................................................................... 2.1. Адаптивный вариант оценки среднего и дисперсии по неравноточным измерениям......... 2.2 Условия сходимости оценок...................................................................................................... 2.3 Упражнения............................................................................................................................... Глава 3. Адаптивный метод решения систем линейных Уравнений............................................... 3.1 Постановка задачи...................................................................................................................... 3. 2. Обоснование метода решения................................................................................................. 3.3. О сходимости метода................................................................................................................. 3.4 Примеры, иллюстрирующие свойства метода......................................................................... 3.5 Заключение.................................................................................................................................. 3.6 Упражнения................................................................................................................................. Глава 4. Адаптивный метод решения обратных задач гравиметрии............................................... 4. 1. Введение.................................................................................................................................... 4.2. Выбор модели и метода решения прямой задачи................................................................... 4.3 Особенности системы уравнений и обратной задачи............................................................. 4.4 Обоснование метода решения обратной задачи...................................................................... 4.6 Проверка точности решения прямой задачи............................................................................ 4.7. Проверка метода решения обратных задач на простых модельных примерах..

.................. 4.8. Решение обратной задачи с привлечением априорной информации................................... 4. 9. Погрешность определения параметров модели при различных горизонтальных размерах блоков................................................................................................................................................. 4.10. Проверка метода на сложных задачах................................................................................... 4.11. Результаты и свойства метода................................................................................................ 4.12. Некоторые рекомендации по методике решения прямых и обратных задач гравиметрии с применением пакета ADG-2.......................................................................................................... 4.12.1. Решение прямых задач..................................................................................................... 4.12.2. Решение обратных задач.................................................................................................. 4.12.3. Работа с избыточными и полными плотностями........................................................... 4.12.4. Редуцирование кривых dg................................................................................................ 4.13. Упражнения.............................................................................................................................. Глава 5. Решение прямых и обратных двумерных задач магнитометрии....................................... 5.1. Введение..................................................................................................................................... 5.2. Выбор модели среды................................................................................................................. 5.3. Выбор метода решения прямой задачи магнитометрии........................................................ 5.4. Разработка алгоритма решения обратной задачи магнитометрии........................................ 5.5. Пакет программ решения прямых и обратных задач магнитометрии (АDM-2 )................ 5.6. Исследование особенностей магнитных полей на моделях.................................................. 5.6.1. Необходимость модельных исследований....................................................................... 5.6.2. Проверка правильности работы программы решения прямой задачи......................... 5. 6. 3. Моделирование аномалий магнитных полей на простых объектах......................... 5.6.4. Исследование на более сложных моделях........................................................................ 5.7. Исследование устойчивости решения обратных задач по модельным данным.................. 5.7.1. Предисловие........................................................................................................................ 5.7.2. Исследование устойчивости обратных задач для однослойной модели среды............ 5.8. Решение обратных задач на реальных данных....................................................................... 5.9. Выводы........................................................................................................................................ 5.10. Рекомендации по моделированию и решению обратных задач с использованием пакета ADM-2................................................................................................................................................ 5.10.1. Рекомендации по моделированию магнитных полей.................................................... 5.10.2. Рекомендации по решению обратных задач.................................................................. 5.11. Упражнения.............................................................................................................................. 5.11.1. Моделирование аномалий от простых объектов........................................................... 5.11.2. Моделирование полей от сложных аномальных объектов........................................... Глава 6. Обзор литературы.................................................................................................................. 6.1. Адаптивные методы в теории управления.............................................................................. 6.2. Методы решения обратных задач в геофизике....................................................................... 6.3. Адаптивные методы в сейсморазведке.................................................................................... ЗАКЛЮЧЕНИЕ..................................................................................................................................... Предисловие Геофизические методы, основанные на измерении полей: магнитного- маг нитометрия, гравитационного- гравиметрия, электрического- электрометрия и уп ругих волн- сейсмометрия- получили широкое применения при исследовании строения планеты Земля как с целью познания, так и для поисков и разведки ме сторождений полезных ископаемых.

Для изучения тех или других объектов на поверхности Земли или внутри её проводятся специальные наблюдения естественных и искусственно создаваемых физических полей. Полученные данные предварительно обрабатываются и, затем, используются для суждения о свойствах изучаемого объекта.

Для того, чтобы получить представления о свойствах объекта, необходимо создать его модель, имеющую некоторые параметры. Так, можно предположить существование слоистой модели среды (объекта) с горизонтальными границами раздела слоев, обладающих различными физическими свойствами. Тогда модель будет включать в себя следующие параметры: глубины границ раздела слоев и свойства пород в каждом слое. Зная эти параметры модели и используя математи ческие зависимости между элементами модели и полем, мы можем вычислить теоретические значения поля для заданных условий его наблюдения. Процесс пе рехода от модели к полю называют моделированием поля ли решением прямой задачи. Переход от значения поля к параметрам модели среды - решением обрат ной задачи.

Одним из простейших вариантов решения обратной задачи является подбор такой модели, которая дала бы теоретическое поле, совпадающее или близкое к наблюдаемому.

Сопоставляя прямую и обратную задачи, необходимо отметить следующие их особенности.

Прямая задача, как правило, имеет единственное решение. Заданной модели при заданных условиях наблюдения соответствует единственное поле.

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

Прямые задачи являются, как правило, устойчивыми. Используя оговорку «как правило»- можно подразумевать, что в многообразии прямых задач может найтись и такая, которая окажется неустойчивой.

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

Для преодоления неединственности и неустойчивости используются методы регуляризации [41]. Основным способом регуляризации является поиск решения вблизи априорно заданного начального приближения параметров модели.

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

На основе адаптивного подхода созданы методы решения различных обрат ных задач сейсмометрии, гравиметрии и магнитометрии, которые положены в ос нову многих пакетов программ.

Учебное пособие предназначено как для студентов и аспирантов, изучаю щих спецкурс по адаптивным методам, так и для специалистов, использующих эти метод в практической интерпретации.

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

Для первого знакомства достаточно прочесть главу 1- «Введение в метод».

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

В главах 4-5 показаны возможности адаптивных методов при решении задач гравиметрии и магнитометрии. В процессе изучения этих глав предполагается ис пользование пакетов программ на персональных ЭВМ.

Обзор литературных источников и список литературы приведены в шестой главе.

Со стороны автора глубокая признательность всем, кто способствовал под готовке и улучшению пособия.

Глава 1. Введение в метод В основе метода должна быть, хотя бы одна, конструктивная идея.

Автор Вы приступили к изучению адаптивных методов и Вам хочется побыстрее узнать, что такое адаптивные методы и чем они отличаются от неадаптивных?

Прежде чем ответить на этот вопрос, начнем с простого примера.

Представим себе, что имеется некоторый постоянный, но неизвестный нам параметр x*, который можно определить или измерить с погрешностью. резуль тат i-го измерения обозначим через xi. модель процесса измерения примем xi = x * + i (1.1) Будем считать, что i – является случайно нормально распределенной поме хой [17, 22] с математическим ожиданием (М), равным 0 и с дисперсией (D=2).

При этих условиях оптимальной оценкой (имеющей наименьшую дисперсию) будет среднее, т.е.

1k 1k k k1 xi = x * + i (1.2) xk = k i =1k i = Эта хорошо известная оценка среднего найдена неадаптивным методом.

Представим себе, что мы провели еще одно k+1 измерение и решили узнать сред нее из k+1 измерений. Мы можем воспользоваться формулой (1.2), высчитав сно ва сумму и поделив ее на k+1, но можем уточнить среднее, построив следующую формулу ( ) (1.3) x k +1 = x k + x k +1 x k k + Убедитесь, что оценка среднего, полученная по формуле (1.3), точно совпа дает с оценкой 1 k + xi (1.4) x k +1 = k + 1 i =1k Формула (1.3) имеет следующие преимущества перед формулой (1.4):

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

2. не может возникнуть переполнения при больших k и x.

Оба эти преимущества имеют решающее значение при включении их в сис темы реального времени, в которых необходимо знание усредненных оценок ка кого-либо параметра в любой момент времени. Совместно с оценкой среднего можно оценивать и дисперсию, по следующей формуле (xk2+1 Dk ) (1.5) Dk +1 = Dk + k где x k +1 = x k +1 x k, т.е. разность, стоящая в скобках формулы (1.3).

Далее будем называть ее невязкой, а множитель, стоящий перед скобкой, корректирующим коэффициентом. В формуле (1.5) корректирующий коэффици 1 ент k +1 = отличается от такового в формуле (1.3), где k +1 = и это связано с k + k тем, что если для оценки параметров x* можно иметь одно измерение, приняв x1 = x1, то для оценки дисперсии необходимо, как минимум, иметь два значения.

Квадрат разности между ними и будет первой оценкой дисперсии. При больших k корректирующие коэффициенты в формулах (1.3) и (1.5) будут мало различимы, их можно брать равными.

Алгоритмы (1.3) и (1.5) являются адаптивными, поскольку в них оценка па раметров x и D ведется через уточнение предыдущих значений по невязке. Зна ние таких параметров, как мы увидим дальше, очень важно для построения опти мальных более сложных алгоритмов. Два простых примера оценки параметров показывают нам идею последовательного пополнения недостающих априорных сведений о модели в ходе поступления текущей информации. Накопление и не медленное использование текущей информации для устранения неопределенно сти из-за недостаточной априорной информации с целью оптимизации избранно го показателя качества и является наиболее характерной чертой адаптации [48].

Задачу оценки одного параметра более подробно обсудим во второй главе.

Там же покажем условия сходимости оценок.

Глава 2. Оценки одного параметра Все начинается с простого.

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

При однократном измерении погрешность получаемой оценки может опре деляться (оцениваться) точностью прибора или другими соображениями, связан ными с условиями измерения или эксперимента. Продолжая тему оценки одного параметра, будем считать, что измерения являются неравноточными и что для каждого из них нам известна оценка погрешности i. Тогда оптимальной оценкой (имеющей наименьшую дисперсию) будет k + x P i i x k +1 = i = k + P i i = где Pi = i Выбор веса обратно пропорционального дисперсии измерения дает оценку x с минимумом дисперсии [22]. Адаптивный вариант оценки, будет иметь вид (x ) Pk + x k = x k + k +1 x k +1 (2. 2) x k +1 = x k + k + Pk + Pk + k где Pk = Pi -вес среднего значения на k-ом шаге.

i = Перейдя от весов к 2 и используя выражения Pk +1 = k2+ Pk = k P k +1 = k + Получим k Pk + k +1 = = Pk + Pk +1 k2+1 + k Вес k+1 оценки x будет равен сумме весов P k +1 = P k + Pk +1. Записав последнее вы ражение через дисперсию, получим 1 1 = +, отсюда k2+ 2 (2.3а) k +1 k k2+1 = k2 (1 k +1 ) Значение k2+1 иногда может служить оценкой погрешности среднего, если k2 точны. Реально же такое бывает редко и поэтому возникает необходимость в по иске независимой оценки дисперсии, получаемой численно через невязки. Для оценки дисперсии соответственно получим Dk +1 = Dk + (2 x k +1 Dk ) (2.4) где - = если по аналогии с оценкой x ввести первый шаг, в котором задается априорное значение дисперсии D0 и погрешности дисперсии DD. В частном случае можем принять D0 = DD = 3 14 [17].

Нетрудно убедиться, что оценки x, получаемые по формулам 2.1 и 2.3, дают одинаковые значения при равноточных измерениях, т. е. при одинаковом i при дем к формулам 1..2 и 1.3.

Рассмотрим еще одно обобщение, предполагая, что ведется измерение зна чения функции U i, которое связано с x * следующим образом:

U k +1 = ax * + k +1 (2.5) Тогда для x k +1 будем иметь ( ) x k +1 = x k + k +1 U k +1 a x k (2.6) a k Где k +1 = k2+1 + a 2 k k - оценка дисперсии среднего на k-ом шаге.

Произведение в ( 2.6 ) a x k будем называть прогнозным значением U i, а раз ность - невязкой между фактическим и прогнозным значениями ( U ). Т.е. невязка на k+1 шаге: U = U k +1 a x k Для оценки дисперсии получим ( ) Dk +1 = Dk + k +1 2U a 2 Dk (2.7) что следует из свойства дисперсии DU=a2D. В частном случае при а=1 выражения 2.6 и 2.7 эквивалентны 2.3 и 2. 4.

На примере оценки одного параметра видны некоторые особенности адап тивных методов и главные из них:

1.Неизвестный параметр уточняется по невязке между очередным измерени ем и прогнозным значением измеряемого параметра 2.Переход от невязки к параметру происходит путем умножения невязки на кор ректирующий коэффициент k +1, зависящий в свою очередь от k2+1 и k и коэффи циента а (2.6).

3.Независимо от оценки k +1 осуществляется численная оценка дисперсии Dk+1.

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

2.2 Условия сходимости оценок А. Дворецкий предложил для исследования сходимости рассматривать любую про цедуру стохастической аппроксимации как обычный детерминированный метод, но с наложением на него случайной составляющей [42]. Следуя этому подходу, формулу 2. запишем в виде ( ) x k +1 = x * + k + k +1 x * + k +1 x * k (2.8) где k - ошибки среднего значения на k-ом шаге;

k +1 - ошибка к+1 измерения.

k Учитывая, что k = i i, получим i = k x k +1 = x * + (1 k +1 ) i i + k +1 k +1 (2.9) i = В этом выражении случайная составляющая сосредоточена в последних двух слагаемых. Обозначим каждое из слагаемых через ri. На этом примере проиллюстрируем условия сходимости, доказанные А. Дворецким. В [42] показано, что для сходимо сти процедур стохастической аппроксимации в среднеквадратическом и с веро ятностью единица, требуется выполнение двух условий;

- помеховая составляющая ( ri ) должна быть несмещенной:

M [ri ] = 0 (2. 10) - сумма дисперсий случайной составляющей должна быть конечной при любой беско нечной процедуре поиска:

M [r ] (2.11) i i = Для выполнения первого условия необходимо, чтобы второе и третье слагаемые выражения (2.9) стремились к 0, если k +1 - конечно и k +1 стремится к 0. Второе сла гаемое стремится к 0, если выполняется условие M [ i ] = 0. При его невыполнении оцен ка, среднего окажется смещенной.

Для анализа второго условия, следуя Дж. Уайлду [42] будем считать, что. суммар r. Если помехи образуют, последо ная ошибка, внесенная всеми помехами равна i i = вательность независимых случайных величин с дисперсией 2, то дисперсия суммар ной ошибки равна математическому ожиданию суммы квадратов отдельных ошибок D(ri ) = M ri 2 = ri 2 (2. 12) i =1 i = Учитывая сказанное, условие (2. 11) для (2. 9) будет D(rk ) = (1 k +1 ) 2 2 + k2+1 (2.13) Если при k и k 0, то второе слагаемое конечно и тоже стремится к0.

Для анализа первого слагаемого при 0 л +1 1 запишем.

k k (1 i )2 2 i2 2 i2 (2.14) i =1 i = Правая часть неравенства, будет конечной лишь в том случае, если i2.

Таким образом, для сходимости стохастической процедуры уточнения парамет ров необходимо выполнение трех условий:

() M ( i ) = 0, M 2 = 2 (2.15) k 0 при k (2.16) (2.17) i i = Из необходимости поиска правильного решения при сколь угодно удаленных начальных значениях возникает еще одно условие сходимости (2. 18) = i i = Условия (2.15) налагает ограничения на помеху, а условия (16) - (18) - на свойства алгоритма.

Как известно [42] условиям (16-18) в частности удовлетворяет гармониче ский ряд =.

k При k k 0, S1 = =, S 2 = k2 - конечна.

k =1 k k = При работе с компьютером физически невозможно требовать бесконечности, т. к. количество итераций уточнения всегда конечно. В этом отношении интерес но знать, чему же равны суммы S1 и S2 при конечных k для гармонического ряда.

Ниже приведены их значения k 10 100 1000 10000 20000 S1 2.93 5.10 7.48 9.79 10.48 10. S2 1.55 1.63 1.64 1.64 1.64 1. Как видно, S1 растет, хотя очень медленно, а сумма квадратов практически стабилизируется уже при первых 100 значениях.

Более подробно и основательно вопросы, связанные со сходимостью проце дур стохастической аппроксимации интересующийся найдет в 6 главе интерес ной и весьма доступной книги Уайлда “Методы поиска экстремума” [42].

2.3 Упражнения 1. Будет ли сходиться процедура стохастической аппроксимации (2.3).

При k =10;

1;

0.5;

0.1;

0. 01;

1/ к;

1/ к2;

1: k ;

к: (к+1 ). Почему?

Каковы свойства таких алгоритмов?

2. Будет ли сходиться алгоритм (2.6)? Почему? При каких сходимость будет больше: при =10;

1 ;

0. 1?

3. Из 10 случайных чисел i постройте для некоторого x* выборку по мо дели xi = x * + i. Рассчитайте х и D, используя k k +1 и k +1 = = k +1 k +1 + k Постройте графики xk, xk и Dk. Сопоставьте со стандартными.

* 4. дополнительно к заданию 3 на 7-ой точке сымитируйте ураганную по меху, внеся число x7 = x * + 10 D. Включите в алгоритм перед вычислени ем очередного значения проверку: если 2 x 9 Dk, то уточнение на k-том шаге не производится. Получите результат при ураганной помехе и при её исключении.

5. Напишите программу и исследуйте процедуру для оценки амплитуды и дисперсии реальной части спектра для разных w, для выборки, форми рующейся по следующей модели xi = a cos w0 t i, 2 k где t i = 0;

;

;

....;

ww w Используя формулы 2. 1 и 2. 7, рассчитайте значения x k и Dk для w = 0.1w0 ;

0.2w0 ;

....;

0.9w0 ;

w0 ;

1.1w0 ;

....;

2w0. Постройте и проинтерпретируйте гра фики x k = f ( w) и Dk = f ( w) 6. То же, что и в 5-ом пункте, сделайте для функции xi = a cos w0 t i + a cos1.2w0 t i 7. То же, что и в 5-ом сделайте, для выделения реальной и мнимой части 2 k для функции xi = a cos w0 t i + a sin w0 t i при t = 0;

, постройте и ;

;

...;

2w 2w 2w проинтерпретируйте графики x( w) и D x ( w) - для реальной и y ( w) и D y ( w) - для мнимой компонент.

2 k 8. То же, что и в 7-ом пункте, посчитайте при t = 0;

. Построй ;

;

...;

4w 4w 4w те и объясните отличие графиков x( w) и D( w), полученных в 7 и 8 пунк тах;

Выберите правильный алгоритм.

9. То же, что и в 7-ом пункте, посчитайте для функции xi = a cos w + t i.

Постройте графики x( w) и D x ( w), y ( w) и D y ( w) и сделайте анализ.

Глава 3. Адаптивный метод решения систем линейных Урав нений "Дорогу осилит идущий".

Народная мудрость.

Для изучения адаптивного метода данная глава является ключевой. Это объ ясняется тем, что, освоив суть метода на хорошо известной математической за даче, нетрудно разобраться в методах решения конкретных геофизических задач.

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

3.1 Постановка задачи Предположим, имеется система линейных уравнений [a ] [x ] = [U ] ij j i (3.1) i = 1, m, j = 1, n где. [aij ]- матрица коэффициентов системы.

[x j ] - вектор неизвестных.

[U i ] - вектор- столбец, являющийся результатом физических изме рений или некоторых оценок.

Допустим, что известны:

1. Средняя квадратическая погрешность каждого i-того измерения, дающая [ ] вектор Ui [] 2. Вектор начальных приближений неизвестных - x 0 j.

[] 3, Вектор погрешностей неизвестных x j По поводу этих допущений необходимо заметить, что при решении конкрет ных задач, связанных с оценкой параметров модели, указанные сведения обычно имеются. Мерой погрешности правых частей в первом приближении может быть точность прибора. Например, при определении времени прихода отраженной волны погрешность сильно варьирует и определяется как функция от соотношения сигнал-помеха на конкретном участке трассы [13]. При таком подходе каждое уравнение системы (1) становится неравноточным [14], что бо лее адекватно физической сути задачи, чем подход, который вытекает из реше ния равноточных алгебраических уравнений. Следовательно, на поиск путей оценки значений вектору U при предлагаемой постановке задачи должно обра щаться серьезное внимание. Приступая к оценке параметров физической мо дели, исследователь обычно знает о приближенных значениях параметров и мо жет задавать начальное приближение вектора неизвестных. На минимизации не которого функционала с учетом этого вектора и построен принцип регуляриза ции по Н. А. Тихонову [41]. Из физической сути задачи следует, что сте пень достоверности начальных значений параметров известна не одинаково точ но. И это очень важно учесть при решении, для характеристики достоверности и [] введен вектор погрешностей неизвестных x j. В качестве первого приближе ния значений x j могут быть взяты возможные пределы изменения неизвестно го.

В самой неопределенной ситуации значения вектора неизвестных могут быть приняты равными 0, а значения x j одинаковые, равные какому-либо большему числу, которое будет характеризовать возможную наибольшую раз ницу между истинными x *j, и начальными значениями.

x = max x0 j x *j (3.1) Определив условия задачи, перейдем к конструированию алгоритмической части процедуры решения систем линейных уравнений.

3. 2. Обоснование метода решения В соответствии с адаптивным подходом (14) уточнение неизвестных ведется по каждому алгебраическому уравнению раздельно. Очередность поступления урав нений может быть нормальной или определяться величиной невязки, как в методе релак саций [24]. Запишем i--тое уравнение системы (1) (3.2) ai1 x 1 +... + aij x j +... + ain x n = U i Возьмем в качестве прогнозных значений неизвестных x j начальные значения x0 j, т.е. x j = x0 j. Подставим прогнозные значения в левую часть уравнения, получим некото рое значение U i, которое назовем прогнозным.

(3.3) U i = aU x 1 +... + aij x j +... + a n x n Найдем разность между фактическим значением (правая часть) и прогнозным и будем называть её невязкой (3.4) U i U i U i Величина невязки в общем случае обусловлена следующими причинами:

1. Отклонениями прогнозных значений x j от истинных x *j x j = x * x j j 2. Ошибками измерения U i = U i* + i 3. Неточностями коэффициентов [aij ] или, другими словами, неадекватностью модели физической, из которой получены измерения, модели математиче ской, из которой получены прогнозные значения. Неадекватность обычно выявляется и устраняется на этапе отладки алгоритма. Дальнейший разго вор будем вести предполагая, что система уравнений адекватна и что доми нирующими в создании невязок являются первые два фактора.

Следуя [27], представим невязку в виде dU i dU i dU i x n + i (3.5) U i = x1 +... + x j +... + dx1 dx j dx n что в линейной системе дает U i = ai1 x1 +... + aij x j +... + ain x n + i (3.6) Обозначив слагаемые через z j, получим n U i = z j, где z 0 = i j = Задача уточнения n неизвестных по одной невязке, в общем случае, имеет мно жество решений. Все слагаемые являются случайными величинами, распределен ными, по нормальному закону, что следует из введенных ранее предположений и свойств случайных величин [17]. Найдем совместную плотность вероятности рас пределения параметров (от которых зависит U i.) в n+1 - мерном пространстве.

z n (z0.....zn ) = exp j 2 z j =0 2 z j j (3.7) n где z j = x j aij ;

z0 = Ui ;

z0 = U z j 2 22 2 j = (3.7) Значения z j выберем таким образом, чтобы плотность вероятности была мак симальной (метод максимума апостериорной плотности вероятности) [17].

z2 max ln P( z 0,...., z n ) = ln + j 2 z2 2 z j j откуда, исключив первую сумму констант и поменяв знаки, получим новый вид целевой функции n U i z j n z + min O(z 0,..., z n ) = j (3.8) U z j = i j Если принять U = const = ;

= const = ;

= 2 2 2 2 i, Ui U zj z z то локальная целевая функция будет иметь вид min O(z1,..., z j,..., z n ) = U i z j + z 2, n n j j =1 j = что по форме соответствует регуляризующему функционалу А.Н. Тихонова [41 с. 117], т.е. в анализируемом методе имеет место локальная регуляризация, и параметру регуля ризации для каждого неизвестного и в каждом уравнении определяются величинами z2 и U. Продифференцировав по каждой из неизвестных, получим систему из n i уравнений с n неизвестными. Она имеет следующий вид:

1 1 1.. 1 z1 U i.

...

1...... =.

. (3.9).....

. 1. j...

..

n z n U i 1. 1... z2j Где j = 2 + U i Решая систему, получим z (3.10) z j = U i j n zj j = Таким образом, невязка разбрасывается пропорционально дисперсии неизвестного.

В знаменатель входит супца дисперсий неизвестных и дисперсии правой части. При n U z2 неизвестные практически уточняться не будут. При обратном неравен i j j = стве невязка будет полностью распределена на уточнение параметров.

Получив оценки zj, перейдем непосредственно к оценкам неизвестных х, полу чим:

a ij x j U i = x (jk ) + ijk +1) U i (3.11) ( k +1) =x + (k ) ( x j j n + a 2 2 2( k ) Ui ij xj j = Если принять = 0, а все xi = 1, то 2 Ui aij (3.12) x (jk +1) = x k + U i j n a ij j = А это известный математикам метод Качмажа [2 1]. Таким образом создаваемый адаптивный метод является обобщением известного итерационного метода.

Нам удалось построить процедуру уточнения неизвестных по одному уравнению и она схожа с процедурой оценки одного параметра ( формулы 2.2 и 2.6 ). Следовательно, мы можем использовать условия сходимости, полученные ранее. Главное из ник k должно стремиться к 0 при k, стремящемся к бесконечности. Дисперсии уточняемых параметров будем находить по формуле D (j k +1) = D (j k ) + aij ijk +1) (x 2 D (j k ) ) (3.13) ( j Анализируя видим, что дисперсия может увеличиваться (при положительной раз ности т. к. произведение а на всегда положительно). Следовательно, она не может использоваться как вычисления корректируемых коэффициентов. Для этих целей будем вычислять оценку 2, через сумму весов, т.е.

j ( k +1) (k ) (3.14) + p ( k +1) =p p ( k +1) (k ) Где p - веса оценок параметра на k+1 и k-шаге, иp - вес уточняющей добавки ( x ).

( k +1) p Обобщим приведенный в главе 2 подход для изменения параметров, измерен ных не непосредственно, а оцениваемых через функциональные связи. В этом случай оценка x ( k +1) зависит от точности измерения, от слагаемых и коэффициентов, входящих в уравнение. Допустим, измеряется U :

U = ax + by + U (3.15) Тогда (U by U ) x= a Полагая, что a и b известны точно, для оценки ошибки x запишем ( b y ) x = U (3.16) a Дисперсия, оценки х будет равна математическому ожиданию квадрата выраже ния правой части уравнения (3.12) ( ) + b 2 y 2 U b y 2 U (3.16а) = x a Выражение окажется близким к 0, если U b 2 y и b 0, и будет большим, если 2 U b 2 y и b 0. Во избежание этих крайних ситуаций третье слагаемое в скобках 2 отбросим, т.е. будем постулировать независимость ошибки измерения и ошибки случайной составляющей второго слагаемого. С учетом сказанного получим ( ) + b 2 y 2 U (3.17) = x a Учитывая, что P = и подставив полученную оценку в (3. 14), найдем a 1 (3.18) = + x ( k +1) k2( k ) U + by y (k ) 2 2 Преобразовав формулу, получим 2( k ) a2 x = (3.19) 2 ( k +1) 2( k ) U + a 2 x ( k ) + b 2 y ( k ) x x 2 2 Аля многомерного случая формула для j-го неизвестного в i-том уравнении имеет вид aij 2 ( k ) 2( k ) 2( k +1) = j j (3.20) j n U i + aij j 2( k ) 2 j = Или, учитывая (2.11), найдем 2( k +1) = 2 ( k ) (1 aij (j k +1) ) (3.20а) j j Здесь вывод формулы (3.17) сделан методом индукции из общих подходов метода наименьших квадратов.

В модифицированном варианте формула имеет вид aij 2 ( k ) ( ) 2( k ) 2( k +1) = j 1 = j 1 j j (3.21) ( k +1) 2( k ) j n U i2 + U i + aij 2 ( k ) 2 j j = при о 1.

Для нелинейного случая следует иметь в виду, что U i (3.22) aij = x j при x *j (k ) для всех j от 1 до n.

Таким образом, в соответствии с поставленной задачей мы обосновали алго ритм итерационного уточнения параметров в одном уравнении. Он включает формулы: 3, 4, 11 и 21.

При переходе к следующему уравнению ( от k+1 к k+2 ) в качестве априор ных будут использоваться уточненные значения параметров и их дисперсии. Так будет продолжаться до тех пор, пока не будут исчерпаны все уравнения системы.

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

D ( L ) D ( L 1) (3.23) где - малое наперед заданное число, a D ( L ) -дисперсия невязок, вычисляе мая, по формуле n U i2 Pi = n, i = (3.24) D ( L) Pi i = при P = или по более удобным рекуррентным формулам (2.2 или 2.4 ). В предель ном случае, когда X X *, из (3. 24) следует n i Pi = i =1m (3.25) D (L) Pi i = т. е, D(l ) будет характеризовать средневесовую дисперсию помехи, которая от номера итерации не зависит. Следовательно, будет выполнено условие (3.23).

U = const) за критерий оста При равноточных значениях вектора U ( при i нова может быть принято неравенство:

D ( L) U (3.26) Аля наглядности представим алгоритм адаптивного метода решения сис тем линейных алгебраических уравнений в виде блок-схемы (рис. 1).

В центре схемы видим блок памяти (БП) для хранения матрицы коэффициен тов, значений правых частей, неизвестных и их дисперсий и других вспомога тельных переменных. Он соединен с вычислительными блоками двойными ли ниями, указывающими передачу информации. В блохах А0, А1, A2 устанавливаются начальные значения констант и организуются циклы по итера циям и уравнениям. Блоки A3 – А6 реализуют вычислительные функции, указан ные в них, а блоки А7-А8 переключают алгоритм на очередные уравнение или итерацию и на окончание.

3.3. О сходимости метода Созданный итерационный метод по классификации Ф.А. Фаддеева и В. Н.

Фаддеевой (Вычислительные методу линейной алгебры-М. : Физматгиз 1963) является не стационарным (циклическим), т. к. матрица Н(L) в уравнении x ( L ) = x ( L 1) + H ( L )UAx ( L 1) (3.27) существенно зависит от номера итерации (L). Поэтому доказательство общей сходимости построить не удается, да это и не имеет особого смысла, т. к. в реаль ных задачах при огромных матрицам исследовать сходимость метода для кон кретной системы - задача более трудоемкая, чем ее решение и анализ на сходи мость по уменьшению невязки. Поэтому в дальнейшем дли анализе сходимости будем использовать подходы, разработанные в методах стохастической аппрок симации [42]. В методах стохастической аппроксимации число независимых на блюдений предполагается неограниченным и за счет этого, при определенных ог раничениях на свойства помех, возможно получение оценок со сколь угодно вы сокой точностью. В наших условиях мы вынуждены компенсировать недостаток независимых наблюдений многократным использованием одной и той же выбор ки с тем, чтобы максимально ослабить влияние априорной информации.

Условия сходимости с двумя неизвестными исследовались в работе [14]. Приведем заключительную формулу критерия сходимости на очеред ном k-том шаге 1 + k k2 ok 1 1 + k2 k W= 0, 1 + k ok 2 1 + k2 k (3.28) y y y* a где л = k ;

k = 2 ;

ok = k x xk x * bk На k-том шаге алгоритм приблизится к решению, если будет выполнено приве денное неравенство.

В частном случае при = 1;

W =, т, е. имеет место сходимость независимо от соотношения других параметров. Этот частный случай соответствует методу Качмажа [21].

3.4 Примеры, иллюстрирующие свойства метода Приведем простые задачи, допускающие проверку и наглядное представле ние результата.

Пример 1. От начального приближения x0=0,5;

y0=0.3, используя различ ные итерационные методы, найдем решение следующей системы x y = 0.5 x + 2 y = Точное решение системы x * = 0;

y * = 1.

Отклонение x = x0 x * = 0.5;

y = y 0 y * = 2, т.е. из (7.1) o=4. Рассмотрим сна чала наиболее благоприятные условия, положив априорные средние квад ратические ошибки (АСКО) равными x = 0. y = Используя адаптивный метод за 4 итерации приближаемся к точному решению (табл. 3.1 и рис. 3.2) Таблица 3. x y xk U k Номер ите- Номер Номер ша yk рации уравнения га k (k ) (k ) 0 0 - 0,5 3 0,5 1 1 1 1,5 0,59 1,59 0,48 0, 2 2 -0,88 0,69 1,17 0,47 0, 2 1 3 -0,52 0,2 1,2 0,114 0, 2 4 -0,3 0,24 1,06 0,11 0, 3 1 5 -0,18 0,07 1,07 0,03 0, 2 6 -0,1 0,08 1,02 0,03 0, 4 1 7 -0,06 0,02 1,02 0,006 0, 2 8 -0,03 0,02 1,01 0,006 0, o Этот же пример решим положив x = y = 1 = 1 и = 4. Как видно из рис.

0 3.2 такое же приближение, удается получить за пять итераций.

И, наконец, предположим, что априорная информация достоверности o неверна, приняв x = 2;

y = 0.5 = 0.25;

= 16. На рис. 3.2 видим, что для по 0 лучения того же результата необходимо на две итерации больше, чем в первом варианте и на одну по сравнению со вторым.

Для сопоставления на рис. 3.2 приведены траектории результатов, получен ных с помощью методов покоординатного спуска и Качмажа. В данном при мере большая скорость сходимости оказалась у адаптивного метода (при пер вом и втором варианте АСКО), а самая медленная - у метода Качмажа.

Пример 2. Начиная от точки x0=y0=2, найдем решение численными ите рационными методами системы x 5y = x + 5y = Следует заметить, что решить эту систему методами покоординатного спуска без поворота системы координат нельзя. Они Не дают решения, т.к. не сходятся. Метод Качмажа (рис. 3) сходится на каждом шаге, но по мере при ближения к решению скорость сходимости постепенно уменьшается. Что каса ется адаптивного метода (рис. 3.3 и табл. 3.2), то его траектория становится бо лее «целеустремленной» по мере приближения к решению. Это объясняется тем, что из-за больших коэффициентов перед y, 2y сильно уменьшается и не вязка в большей степени идет на уточнение x.

Таблица 3. x (k ) y (k ) xk U k номер ите- номер 2 yk рации уравнения 2 0 0 2,3 0, 1 0, 1 0, 2,06 -0, -4, 2 0, 0, 1,7 0, -4, 1 0, 2 0, 1,2 -0, -3, 2 0, 0, 0,83 0, -2, 1 0, 3 0, 0,44 -0, -1, 2 0, 0, 0,16 0, -0, 1 0, 4 0, 0,03 -0, -0, 2 0, 0, 0,006 0, -0, 1 0, 5 0, 6 2 10 3 0, 2 0, 0, Пример 3. Начиная от точки x0 = 0.5;

y 0 = 3 с использованием адаптивного метода найдем решение несовместной системы трех уравнений x y = 0.5 x + 2 y = 0.333 x + y = 2. Имеющих различные погрешности правых частей.

U = U = 0. 1 U = 0. На рис.3.4 и табл.3.3 видно, что траектория решения приближается к точке пересечения 1 и2 линий, а уравнения 3 оказывает влияние только при первом (с его участием) уточнении (отрезок от 2 к 3 шагу), пока x и y являются больши ми.

Таблица 3. x (k ) y (k ) xk U k Номер ите- Номер урав- 2 yk рации нения 0, 0 0, 0, 2, 1, 0 1 -1, 0, 0, 1, 1, 1 2 -1, 0, 0, 1, 1, 3 0, 0, 0, 1, 0, 1 -1, 0, 0, 1, 0, 2 -0, 0, 0, 1, 0, 3 0, 0, 0, 0, 0, 1 -0, 0, 0, 1, 0, 2 -0, 0, 0, 1, 0, 3 1, 0, 0, 1, 0, -0, 0, 0, 1, 0, -0, 4 0, 0, 1, 0, 1, Пример 4. Используя адаптивный метод, найдем решение системы из примера 3, приняв U 1 = U 2 = U 3 = 0.5.

Как видно из рис.3.5 и табл. 3.4 решение постепенно приближается к центру треугольника.

Пример 5. Используя адаптивный метод, решим нелинейную систему 0. t2 + = 1. v t2 + = 1. v при следующих начальных данных:

v0 = 2.2;

t 0 = 0.05;

v 0 = 0.2;

U = 0. t 0 = 0.95;

Точное решение t = 1, v = Результаты уточнения для = 0 приведены в таблице 3. Таблица 3. x (k ) y (k ) xk U i Номер ите- Номер урав- 2 yk рации нения 0, 0, 0, 2, 1, -1, 0 0, 0, 1, 1, - 0, 0, 1, 1, 0, 0, 0, 1, 0, - 0, 0, 1, 1, -0, 2 0, 0, 1, 1, 0, 0, 0, 1, 0, -0, 0, 0, 1, 0, -0, 0, 0, 1, 1, 0, 0, 0, 1, 0, -0, 0, 0, 1, 0, -0, 0, 0, 1, 0, 0, Таблица 3. tk vk tk U k Номер ите- Номер урав vk рации нения 0 0,95 2,200 0,050 0, 1 0,051 1,001 2,196 0,001 0, 2 0,016 1,008 1,983 0,001 0, 1 -0,001 1,000 1,983 0,0008 0, 2 -0,002 1,000 1,997 0,0007 0, Пример 6. (Взят из [42, c.263]). От начального приближения x=10 найти корень для функции y = x 3 + x + k, где k - случайное число, которое предлагается выби рать ± c, бросая монету. Не меняя сути задания, в качестве помехи используем гармоническую функцию k = c sin(0.5k ).

На этом примере сопоставим три метода: Ньютона, стохастической аппроксима ций и адаптивной. При отсутствии помех, т.е. при с=0, метод Ньютона и адаптивный при решении этой задачи практически совпадают, отличие результатов имеет место при наличии помех, что видим из таблицы 3. 6, в которой сведены результаты поиска корня при с=2. Адаптивный метод испытывался при двух вариантах задания АСКО начально го приближения:

x20 = 10 и x 0 = Таблица 3. Адаптивный метод к Метод Ньютона Вариант 1 Вариант x y k y k y k xk xk xk x 9, 6, 6, 6, 9, 4, 4, 4, 9, 2, 2, 2, 6, 1, 27, 1, 1, 27, 1, 27, 4, 0, 10, 6, 0, 10, 0, 10, 0, 0, 1, 0, 0, 1, 0, 1, 0, 0, -1, 0, 1, -1, 1, -1, 0, 0, 1, 0, 0, 2, 0, 3, 0, 0, 2, 0, 0, 2, -0, 3, 0, 0, 0, 0, 0, 1, -0, -0, 0, 0, -1, 0, 0, -1, 1, -2, 0, 0, 0, 0, 0, 0, 1, 9, 0, -0, 2, 0, -0, 2, 0, 4, 14 0,26 0,025 -0,03 -0,02 0,16 -0,02 -0,017 0, 15 -0,19 2,0 -2,0 0,25 0,15 -2,0 0,25 0, Из таблицы видно, что методы Ньютона и адаптивный (при двух вариантах дисперсии x20 = 10 и x20 = 100 ) первоначально одинаково быстро приближаются к решению, но, начиная с 7-го шага, видны отличия, которые отчетливо проявляются на 11 и 15 шагах. В результатах адаптивного метода флюктуации оценок корня, вы званные помехой, постепенно уменьшаются из-за постепенного уменьшения x2, а в результатах метода Ньютона они остаются на одинаковом уровне ( ± 2 ).

Сходимость метода стохастической аппроксимации [42] сильно зависит от выбора шага. Если в этом примере взять k =, то процесс уточнения k будет быстро расходиться. Уже на первом шаге будет получено x1 500, а на 0. втором x 2 40 10 6. Поэтому, по рекомендации [42], возьмем k = и дополни k 0. тельно k =.

k Результаты счета сведем в таблицу 3. 0.01 0. МСА k = МСА k = k k k yk + k yk + k xk xk 1 1012 -0,12 1012 - 2 -0,12 -0,12 -1084 0, 3 -2,12 -0,11 -1,18 0, 4 -0,11 -0,11 0,83 0, 5 1,89 -0,11 2,8 0, 6 -0,11 -0,11 0,8 0, 7 -2,1 -0,11 -1,19 0, Из таблицы видно, что результаты метода, стохастической аппроксимации практически сильно зависят от первого или второго шага, а дальше они уточняются очень мало, Приведенные выше примеры показывают нам, что адаптивный метод по зволяет:

а) ускорить сходимость решения задачи не только за счет учета достовер ности начальных приближений, но и за счет собственных адаптивных свойств алгоритма (примеры l, 2) б) решать несовместимые неравноточные системы (примеры 3 и 4), т. е, системы, не имеющие решения;

в) решить нелинейные системы (примеры 5 и 6).

3.5 Заключение В начале главы, мы предположили, что кроме системы уравнений нам из вестны:

1.Средняя квадратическая погрешность каждого i-того измерения, дающая вектор [ Ui ] 2.Вектор начальных приближений неизвестных - [x0 j ] 3.Вектор погрешностей неизвестных [ xj ] В такой постановке мы обосновали новый метод итерационного уточнения параметров в одном уравнении. Он включает формулы 3, 4, 11, 21.

При переходе к следующему уравнению ( от k+1 к k+2 :) В качестве априор ных будут использоваться уточненные значение параметров и их дисперсии. Так будет продолжаться до тех пор, пока не будут исчерпаны все уравнения системы.

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

Анализируя, мы увидели следующие свойства метода:

1. Метод относится к классу итерационных.

2. Он обобщает метод Качмажа.

3. Метод имеет гибкую локальную регуляризацию: параметры регуляризации для каждого неизвестного и в каждом уравнении определяются величи нами z2 и U. Приведенные выше примеры показали, что адаптивный ме тод позволяет: ускорить сходимость решения задач, решать как совместные, так и несовместные, как равноточные, так и неравноточные, как линейные, так и нелинейные системы.

Для дальнейшего изучения возможностей метода полезны приведенные ниже упражнения.

3.6 Упражнения 1. Выберите систему двух уравнений с двумя неизвестными. Задайте начальные приближения неизвестных и x и y. Решите систему адаптивным методом, по строив таблицу решения (см. пример 1).

На координатную плоскость XY нанесите уравнения и траекторию решения ме тодов: адаптивного, Качмажа и покоординатного спуска (см. рис 3.2). Сопоставь те скорость сходимости методов. Смените значения x и y и повторите решение, дайте интерпретацию полученных результатов.

2.Выберите систему нелинейных уравнений с двумя неизвестными. Решите адаптивным методом и методом Ньютона. Постройте таблицы и графики, дай те интерпретацию.

3.То же что и в третьем примере, но в одно из уравнений введите помеху cos k, где k = 0,1,2,...,20. Соответственно примите U = 1.

4.Используя программу решения систем линейных уравнений (ADASLU) на персональных ЭВМ, решите выбранную вами систему линейных уравнений при различных параметрах x и y, распечатайте невязки и результаты. Проин терпретируйте результаты.

5.То же, что и 4, но введите в уравнение помехи: anois = 0.1;

1 при U = 0,1, U = 1.

Выпишите СКН и дайте объяснение уровню невязок (см. формулу 3. 25) 6. Выберите и решите недоопределенную систему с различными начальными данными, опишите результат.

7. Выберите и решите переопределенную систему при U = 0.01. Сделайте опи сание результата.

8. То же что и в 7, но с U = 1.

Глава 4. Адаптивный метод решения обратных задач грави метрии С гравиметрического метода берет свое начало геофизика Автор 4. 1. Введение Обратные задачи гравиметрии (особенно трехмерные) относятся к классу наиболее трудных [9] в силу нескольких причин.

1. Отсутствие единственности решений. Одному, и тому же гравитационно му полю, в общем случае, может быть подобран не один вариант модели среды.

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

3. Трудности решения алгебраических систем с большим числом уравнений и неизвестных (103 и более).

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

Эта задача решалась многими исследователями [2,4,19,34,36-40], но наиболее близкими к нашему направлению являются методы сеток, предложенные Г. Г.


Булахом и А. А.Юньклвым [4]. В дальнейшем для решения систем и учета огра ничений на плотность С. В. Шалаевым было предложено использовать методы линейного программирования [49]. Позднее этот подход был продолжен в рабо тах В. В. Ломтадзе [19]. Использование методов линейного программирования не всегда оказывается конструктивным. Его неконструктивность выявляется при не обходимости ввода жестких ограничений. Например, при отрицательных анома лиях необходимо решать "зеркальную" задачу, т. е. переводить и аномалии, и из быточные плотности из отрицательных в положительные. А как быть, когда есть и отрицательные, и положительные? Кроме того, трудности в методе возникают и тогда, когда число неизвестных велико.

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

Используем для решения обратных задач гравиметрии изучаемый нами адап тивный метод. Впервые он был применен для решения кинематических задач сейсморазведки [11-16]. Накопленный опыт позволил создать методику разра ботки и исследования алгоритма конкретной обратной задачи [13]. В соответст вии с методикой необходимо пройти последовательно следующие этапы:

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

Намеченный порядок разработки и используем для последующего описания.

4.2. Выбор модели и метода решения прямой задачи За основу примем слоистую модель среды, которая является достаточно общей при изучении осадочных толщ. Она же может использоваться и при моде лировании большого числа интрузивных тел. Каждый слой по оси X разбива ется на прямоугольные параллелепипеды, вертикальные размеры которых оп ределяются мощностью слоя, а горизонтальные - заданными величинами dx и dy, соответствующими шагу прямоугольной сетки наблюдения или скани рования гравитационного поля. Точки наблюдения будем располагать на по верхности над центрами параллелепипедов.

Проследим логику вывода формулы расчета гравитационного поля для пря моугольного параллелепипеда. В соответствии с законом Ньютона имеем Gm1 m (4.1) F= r где m1 и m2 -точечные массы тел, r 2 - расстояние между ними, см G - гравитационная постоянная, равная 6.67 10 г / с Если F поделить на m1, то получим вместо силы напряженность гравитацион ного поля, которую обозначим через g. Введя вместо m2 – переменную M по лучим GM (4.2) g= r Напряженность - векторная, величина, имеющая то же направление что и сила F. Следуя [23], введем радиус вектор R, соединяющий точку наблюдения и ис точник. Направление по радиусу к источнику считаем положительным, тогда вектор g GM R (4.3) g= R Для расчета напряженности поля объемного тела, с плотностью вещества, рассмотрим элементарный объем dv, а затей интегрируя по всему телу, получим Rdv g = G (4.4) vR Рассмотрим вертикальную составляющую гравитационного поля в декарто вой системе. Формулу для нее получим из (4. 4) (z z 0 )dv g z = G (4.5) vR где z и z0- координаты точки источника и наблюдения.

Преобразуем формулу (4.5), учитывая, что в гравиметрии интерпретируется не все поле, а его аномальная составляющая g, вызванная аномальной или из быточной плотностью.

( z z0 ) g = G (4.6) dv vR Предположим, что точка наблюдения расположена в начале координат, ис точник напряженности представлен в виде параллелепипеда, ограниченного сверху горизонтальной плоскостью z1, а снизу z 2, по оси X плоскостями x1 и x 2, а по оси Y соответственно y1 и y 2. Переходя от интегрирования по объему к ин тегрированию по осям декартовой системы координат, получим вклад в поле g от элементарного параллелепипеда x2 y 2 z z dx dy dz g = G (4.7) (x ) + y2 + z / x1 y1 z Проигнорировав, получим [2, 9, 36, 39] g = G yR1 + xR2 + zR3 | x12 | y12 | z12 (4.8) xyz где R1 = ln ( x + R ) R2 = ln ( z + R ) xy R3 = arctg zR ( ) 1/ R = x2 + y2 + z Для двухмерного случая при y будем иметь x g = G + x ln( R z ) + 2 z arctg | x12 | z12 (4.9) xz z Вводя множество параллелепипедов (которые будем называть, для кратко сти, блоками или призмами) и, обозначая избыточную плотность в каждом блоке через ijk, а гравиметрический эффект от блока с = 1 через cijk, полу чим следующее уравнение для расчета g от блока в точке наблюдения с координатами х0,у0,z I J K g ( x0, y 0, z 0 ) = ijk cijk (x0, y 0, z 0 ) (4.10) i =1 j =1 k = Для каждой точки наблюдения получим своё уравнение.

4.3 Особенности системы уравнений и обратной задачи Проиллюстрируем особенности систему уравнений и обратной задачи на простом двухмерном варианте с одним слоем. В этом случае, g в i--той точке наблюдения от J элементарных блоков будет J g = cij j (4. 11) j = Взяв совокупности из I наблюдений, получим систему из I уравнений с J не известными. c1J 1 g c.....

. ciJ j = g i ci1 (4.12).....

c I 1. c IJ J g I При I= J система будет иметь единственное, но не всегда устойчивое реше ние, что будет проиллюстрировано далее. В частном случае, когда поверхность наблюдения и границы слоя будут горизонтальны, то матрица коэффициентов обладает следующими свойствами:

1. диагональные элементу матрицы равны и симметричны относительно глав ной диагонали.

2. Элементы главной диагонали будут наибольшими.

3. Значения остальных элементов будут убывать по мере удаления от главной диагонали.

Все эти свойства четко объясняются из физических предпосылок. На главной диагонали находятся коэффициенты cii равные g ii при избыточной плотности, равной единице. Они отражают гравиметрический эффект от блоков, располо женных под точкой наблюдения. Вспомним, что точки наблюдения находятся на поверхности над центрами блоков (4.2). По мере удаления блоков от точки на блюдения их влияние на гравиметрическое поле будет уменьшаться.

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

Приведем для иллюстрации коэффициенты для трех слоев, кровли которых располагаются на глубинах 0, 500 и 1000 м, при горизонтальной величине бло ков 500 м.

Таблица 4. c13 c15 c16 c17 c c11 c12 c N слоя интервал 1 0-500 11,6 2,6 0,8 0,4 0,2 0,1 0,1 0, 2 500-1000 4,43 3,08 1,6 0,89 0,55 0,37 0,26 0, 3 1000-2000 4,57 4,06 3,06 2,19 1,58 1,16 0,88 0, В случае, когда поверхность наблюдения или границы слоев не горизон тальны, упомянутые свойства коэффициентов не сохраняются. Поэтому необхо димо вычислять их значения для всех точек наблюдения, что приводит к увеличе нию числа коэффициентов в I раз (где I число точек). Для двухмерной задачи, включающей сто точек по профилю и пять слоев, потребуется хранить 100 100 5 = 50000 коэффициентов. Для трехмерной задачи их число увеличивает ся как минимум на порядок.

4.4 Обоснование метода решения обратной задачи Предположим, что априорно известны начальные приближения избыточных плотностей j в каждом блоке. В частном случае они могут, быть заданы константой, равной нулю.

Будем считать, что (j0 ) = * + (j 0) (4.13) j где *j - точное, но неизвестное нам значение избыточной плотности, а (j0 ) - погрешность начального приближения.

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

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

Предположим также, что правые части, т. е. значения g i, осложнены поме хой g i = g i* + i (4.14) и что нам известен стандарт ошибок правой части i, в частном случае он может быть задан константой 0, равной погрешности съемки.

В предлагаемом методе [11-16] система уравнений не решается в целом, а осуществляется многократное, по всей системе, последовательное уточнение не известных по каждому уравнению отдельно, как это осуществляется в некото рых пошаговых итерационных алгоритмах [24, 42].

Предположим, что неизвестные уточняются на k+l-ом шаге l-ой итерации, по i-ому уравнению, которые связаны следующим соотношением: k = i + I (l 1).

Подставим в него сформировавшиеся на k-ом шаге значения неизвестных (k ) :

j J = cij (j k ) (k ) (4.15) g i j = (k ) Значение g i будем называть прогнозным значением гравитационного поля в точке i.

Найдем разность между фактическим и прогнозным значениями гравитационного поля (4. 16) g i( k ) = g i g i( k ) Наличие ненулевой невязки в общем случае обусловлено следующими причи нами:

• отклонением прогнозных значений избыточных плотностей (k ) от j истинных, т. е. j являются не нулевыми;

• ошибками измерения i ;

• неадекватностью математической и физической модели.

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


Итак, предполагая, что основная причина невязки обусловлена первой и второй причинами, представим невязку в виде суммы J g i = cij j + i (4.17) j = Детально повторив рассуждения, приведенные в главе 3, придем к следующим оценкам значений j cij ( 2 ) (k ) j (4.18) ( k +1) = g (k ) j i + c ( ) 2 (k ) 2 0 ij j j где k- номер шага, а i- номер уравнения.

Из приведенной формулы видим, что в большей степени уточняются избыточ ные плотности тех блоков, у которых большие коэффициенты и большие диспер сии погрешностей, отсюда уточненные избыточные плотности в i-том уравнении на k-том шаге будут равны:

(j k +1) = (j k ) + (j k +1) (4.19) В адаптивном методе на каждом шаге изменяется (уменьшается) значение параметра 2 в соответствии с обоснованием, приведенным в главе j примем следующую формулу, пересчета дисперсий ( ) = ( ) (1 ), ( k +1) 2 (k ) 2 (k ) j j j (c ) ( k 1) (k ) = ij j (4.20) (k ) где (g ) + + (c ) j (k ) 2 ( k 1) i 0 ij j при 0 Из выражений (17-20) видим, что значению параметра, уточненному в большей степени, будет соответствовать и большее уменьшение 2.

j После уточнения оценок неизвестных и изменения дисперсий переходим на очередное уравнение и процесс уточнения неизвестных повторяем, используя формулы (15, 16, 18-20). Продолжая до последнего уравнения, пройдем первую итерацию уточнения. Из-за того, что процесс уточнения идет приближенно, необ ходимо повторное уточнение на нескольких итерациях до тем пор, пока не будет выполнено одно из двух условий: k = k max или ( ) 1I g i k 0 (4.21) I i = Где 0 - средняя квадратическая погрешность g.

Обоснованные выше алгоритмы были положены в основу пакетов программ решения прямых и обратных двух- и трехмерных задач гравиметрии. Дадим им краткую характеристику.

4.1 Пакеты программ решения прямых и обратных задач гравиметрии Пакеты созданы для персональных ЭВМ и запрограммированы на языке Пас каль. Авторы разработки В. А. Кочнев, В. И. Хвостенко. Двухмерный пакет (ADG-2) по зволяет решать прямые задачи с любым числом слоев. В обратных – число слоев огра ничено пятью, а число блоков в каждом слое - сто. Время счета зависит от числа слоев и блоков и колеблется от минут до десятков минут.

В пакета предусмотрено несколько режимов работы:

• INV=1 - решается прямая задача;

• INV=2 - решается прямая, а затем обратная задачи (неизвестными являются плотности);

• INV=3 - решается только обратная задача (неизвестными являются плотности);

• INV=4 - решаются прямая и обратные задачи (неизвестными являются гра ницы слоев);

• INV=5 - решается обратная задача (неизвестными являются границы разде ла слоев).

Для связи с сейсморазведкой рассчитываются скорости, нулевые времена и статические поправки. Для этого предусмотрена возможность задания коэффици ентов линейной связи между скоростью и плотностью.

Трехмерный пакет ADG-3 может решать задачу с пятью слоями. Размерность матрицы блоков в каждом 'слое ограничена 100x100. Режимы работы те же, что и в двухмерном пакете.

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

4.6 Проверка точности решения прямой задачи Проверим первоначально правильность решения прямой задачи. В качестве теста используем расчет аномального эффекта над центром горизонтального слоя.

Используя формулу (4.9), возьмем предельный случай, приняв x1 =, x 2 =, z1 = 0, z 2 = h, получим:

g = 2Gh Проверяя при = 1, построим таблицу значений g в зависимости от h и метода.

Таблица 4. h(м) Аналитическ. ADG-2 ADG-3 Палеточный g в мГал 20, 20, 20, 41, 41, 41, 83, 82, 83, 167, 163, 167, Из таблицы видно, что за счет ограничения горизонтальных размеров слоя (по км от точки наблюдения) с увеличением толщины сдоя погрешность увеличивает ся при расчете g пакетом ADG-2 и, особенно, при расчете g с помощью палет ки. Погрешность "палеточных" значений возрастает с увеличением мощности слоя, что объясняется ограничением горизонтального размера палетки. В пакете ADG-3 влияние краевых эффектов учитывается путем продолжения крайних бло ков в бесконечность по горизонтали. Поэтому при любых горизонтально слои стых моделях прямая задача этим пакетом решается точно.

Используя пакет, численно оценим погрешность расчета g, из-за неточно сти аппроксимации криволинейной границы прямоугольными призмами. Для этой цели в качестве криволинейной поверхности возьмем синусоиду. Будем по ложительный полупериод синусоиды аппроксимировать 65, 33, 17, 9 и 3 блоками таким образом, чтобы центральная точка находилась на максимуме синусоиды (рис. 4.1). Отношения погрешностей, вычисленных в этой точке, при различной детальности аппроксимации и при разных уровнях границы приведены в таблице 4.3.

min расстоя- Число блоков ния от по 65 33 47 9 5 верхности Относительная погрешность g 0 0 0,01 0,03 0, 1000 0 0,01 0,04 0, Как и следовало ожидать, наибольшая относительная погрешность возникает при меньшем числе точек. При трех точках погрешность достигает 10 - 12%, Она.

значительно уменьшается при уменьшении шага в 2 раза, т. е. зависимость по грешности примерно пропорциональна квадрату размера блока или шага между точками. Здесь, следует иметь в виду следующие обстоятельства:

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

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

2. Погрешность получена для слоя, меняющего мощность от 0 до 1000 м. Ес ли же слой может быть разложен по мощности на два подслоя: один из них с по стоянной (hconst), а второй переменной толщиной (hvar ), то погрешность вычисле ния g от постоянного подслоя будет равна 0. Относительная погрешность от переменного подслоя вычислена нами выше (Еу). Тогда суммарная относитель ная погрешность ЕS будет меньше. Она может быть оценена по формуле EV hvar ES = hvar + hconst Таким образом, реальная относительная погрешность при большой мощно сти слоя и малой амплитуде границ будет значительно меньшей, чем та, которая приведена в расчетной таблице.

4.7. Проверка метода решения обратных задач на простых модель ных примерах.

Проведем серию экспериментов по решению обратных задач на двух простых наглядных моделях. Первую модель примем однослойной с горизонтальными границами. Кровля, она же поверхность наблюдения, располагается на высоте – 500 м, а подошва на 0. Расчленим слой на 21 блок с шагом 500 м. Примем 5 бло ков находящихся в центре, аномальными с избыточной плотностью 0,5 г/см3.

Вторую модель произведем от первой, поместив этот же аномальный слой на бо лее глубокий уровень от О до 500 м, а слой от -500 до 0 сделаем с нулевой избы точной плотностью. Решив прямую задачу для этих двух моделей получим две кривые g (рис. 4. 2). Как видно, они значительно отличаются величиной мак симумов и крутизной.

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

Эксперимент 1. Решая обратную задачу при этих условиях, получено прак тически точное решение. Средняя квадратическая ошибка (СКО) полученных оценок плотности составляет 0,01. Ниже приведены средние квадратические не вязки (СКН) между фактическими и прогнозными g на каждой итерации:

номер итерации 1 2 3 4 5 6 СКН 1.66 0.60 0.24 0.10 0.04 0.02 0. Эксперименты 2-5 посвятим проверке устойчивости метода к помехам. Для этого будем вводить после решения прямой задачи случайные помехи с нулевым мате матическим ожиданием и заданным уровнем помехи -. Результаты счета после 5 итераций увидим ниже:

Уровень шума 0 0,1 0,2 0,4 0, СКН 0,01 0,04 0,14 0,34 0, СКО 0,00 0,01 0,02 0,04 0, К эксперименту 6-10 проведем для проверки метода при решении обратной зада чи в нелинейном варианте, т.е. известной является плотность, а неизвестной явля ется граница раздела двух слоев. Начальное приближение границы примем рав ным 500 м. Погрешность начального приближения для всех блоков равна 200 м. В результате решения при разном уровне шумов на 5-ой итерации получим сле дующие оценки (время счета 18 с).

уровень шума 0 0.1 0.2 0.4 0. СКН 0. 01 0. 05 0. 12 0. 30 0. СКО 0.11 0.38 5. 20 18.6 27. Анализируя результат, видим, что СКН практически точно отображает уро вень помех. СКО увеличивается, достигая 27 м при 0 = 0. Проведем подобные эксперименты для того же аномального слоя, располо женного на глубине 500 м от поверхности наблюдения. В экспериментах 11- получим следующие результаты после 7-ой итерации:

уровень шума 0 0.1 0. 1 0. 4 0. СКН 0. 02 0.10 0. 16 0. 30 0. СКО 0.07 0. 08 0.08 0. 08 9. Анализируя таблицу, видим, что средние квадратические отклонения от мо дели во всех экспериментах, независимо от уровня помех, близко к 0. 1 г/см3.

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

4.8. Решение обратной задачи с привлечением априорной информа ции Для ответа на поставленный вопрос проведем эксперименты на этой же модели с незашумленной g. При решении контактной задачи, вводя дополни тельную априорную информацию.

В эксперименте 16 априорные глубины границ зададим равные 750 м, по грешность наблюдения - 0. 01 мГал, погрешность положения границы для всех точек примем равной 1000 м. На 11 итерации получим СКН 0.05, а СКО резуль тата - 47 м.

Предположим, что мы уточнили модель, а затем в центре аномалии пробу рили скважину и уточнили глубину верхней границы. Приняв ее равной 500м (абсолютная отметка 250 м), а погрешность в этой точке зададим равной нулю. В результате на 5 итерации получим СКН равной 0. 05 мГал, а СКО результата - м.

В эксперименте 17 уточним положение границы еще в двух точках - 8 и 13.

В результате на 3 итерации получим СКН равной 0.00, а СКО – 2.3 м.

Таким образом, исследованиями показана возможность уточнения ре зультатов путем вовлечения априорной информации в некоторых точках изучае мой модели среды.

Дополнительные эксперименты показали, что не все точки оказываются в равной степени априорно информативными. Для проверки этого утверждения проведем серию экспериментов с тремя априорными точками. Причем, средняя будет на ходиться в центре на 11 блоке, а две другие располагаться симметрично относи тельно центральной. Ниже приведены результаты экспериментов.

номер эксперимента 18 19 20 положение априорных точек 10, 11, 12 9,11 и 13 8,11 и 14 7, 11 и номер итерации 20 4 20 СКН 0,2 0 0, 1 0. СКО 10, 9 2, 3 12, 9 25, Результат наиболее информативного 19 эксперимента приведен на рис. 4.3.

Исследуя причину разного влияния априорных данных на результат, прихо дим к выводу, что наиболее информативными оказываются те точки, в которых первоначальный результат наиболее сильно отличается от модельного. Для под тверждения этого приведем таблицу оценок, полученных в эксперименте 17.

номер блока 7 8 9 10 11 12 13 14 результат 750 707 323 217 250 222 312 713 модель 750 750 250 250 250 250 250 750 разность 0 -43 73 -33 0 -28 62 -37 - Из приведенных данных видно, что именно точки 9 и 13 имеют наиболь шие, а точки 7 и 15 наименьшие отклонения от модели. Это подтверждает извест ное в теории планирования эксперимента правило - новый эксперимент по изуче нию надо ставить там, где дисперсия оцениваемого параметра наибольшая. И ес ли в процессе его постановки удается получить большое отличие - между резуль татом предыдущей интерпретации и новой, то это говорит об удаче эксперимента.

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

Далее рассмотрим погрешность оценки параметров модели при различных го ризонтальных размерах блоков.

4. 9. Погрешность определения параметров модели при различных горизонтальных размерах блоков Введем параметр, характеризующий отношение горизонтальных размеров блока dx к средней глубине hср т. е.:

dx ср = hср Проведем эксперименты на той же модели при разных dx, a, следовательно, и.

Номер экспе 16 22 23 24 римента 250 333 500 666 dx ср 0.5 0.6 1 1.33 СКН 0.04 0.04 0.03 0.02 0. СКО 57.8 51 24.5 13 2. % 11.5 10 5 2.6 0. Как видно, увеличение (за счет увеличения dx) приводит к уменьшению среднего квадратичного отклонения результативных значений глубин от мо дельных. Относительная погрешность составляет около 5% при =1. Эту оценку можно принять за граничную, характеризующую разрешающую способность гравиметрического метода.

Проведем подобные эксперименты, изменяя параметр за счет изменения положения аномального блока относительно поверхности наблюдения (при dx=50).

N эксперимента 24 25 26 27 0,5 0,66 1 1,33 СКН 0,05 0,03 0,03 0,02 0, СКО 68 59 27 8,5 1, % 5,4 2,8 0, Результаты оказались близкими к тем, что получены в предыдущих экспе риментах. Таким образом, эксперименты, проведенные на простых моделях, позволяют сделать выводы:

1. о помехоустойчивости методов решения линейной и нелинейной (контактной) обратных гравиметрических задач;

2. о зависимости погрешности результата от соотношения горизонтальных раз меров блоков и средней глубины границы;

3. о возможности дополнительного учета априорных сведений, полученных из данных бурения или других источников.

4.10. Проверка метода на сложных задачах.

В заключение приведем примеры решения обратных задач, иллюстри рующие некоторые возможности адаптивного метода и пакета ADG- Взяв из [9, с.400] сложный пример обратной контактной задачи, построим близкую к нему модель и решим ее при тех же условиях, т.е. приняв начальное положение уточняемой границы на глубине 11000 м. В результате решения об ратной задачи по данным без помех получим СКН 0. 01 мГал, а СКО глубины 62. м. При гауссовских помехах с нулевым математическим ожиданием и с 0 =0, мГал получим СКН – 0,77 мГал, СКО - 285 м. Относительная погрешность опре деления глубин ср составила 1% в варианте без помех и 5%-с помехами. Резуль тат с помехой видим на рис. 4.4. На верхнем графике - модельная и подобранная g. На втором графике –модель, а на третьем - средние скорости пробега волн до нижней границы, т.е. от поверхности до 11000 м.

В пакете предусмотрена возможность перехода от плотностной подели к скоростной, по формулам вида V = V0 + b( 0 ). В данном случае принято. Пере ход от плотностной модели к скоростной необходим при совместной интерпре тации гравиметрических и сейсмических данных. Более полно эти возможности пакета представим, вернувшись к рис. 4.3, где показаны средние скорости (Vср ), статические поправки () и нулевые времена (t0).

В заключение приведем результаты решения реальной обратной задачи в сложной модели среды, полученные сотрудниками гравиметрической экспедиции #3 с применением пакета ADG-2 с учетом данных интерпретации первых вступ лений преломленных сейсмических волн (рис. 4.5). Главная цель совместной об работки всех геофизических данных - построить модель разреза и, главным обра зом, модель верхней части разреза (ВЧР). На верхнем графике приведены исход ная и подобранная кривые g. Они практически совпали. На среднем графике модель среды, а на нижнем - график нулевых времен, т.е. времен пробега сейсми ческих волн по вертикали от поверхности до нижней границы и обратно.

Метод адаптивного решения обратной задачи гравиметрии нашел примене ние как при оценке параметров ВЧР, так и при оценке параметров рудных тел, оценке положения фундамента. Опыт, накопленный при решении модельных и реальных задач, позволил сформулировать некоторые рекомендации, приводи мые в разделе 4.10.

4.11. Результаты и свойства метода В главе прослежен путь разработки и исследовании адаптивного метода ре шения обратных задач гравиметрии.

Для персональных ЭВМ IBM AT разработаны пакеты решения прямых и обратных задач.

За основу принята слоистая модель среды. Каждый слой по осям X и Y раз бивается на прямоугольные параллелепипеды, вертикальные размеры которых определяются мощностью слоя, а горизонтальные - заданными величинами dx и dy.

Модельные исследования погрешности расчета g показали, что при трех блоках, перекрывающих аномалию (положительный полупериод синусоиды) она может достигать 10 – 12% и значительно уменьшается (до 3-4) при уменьшении шага в 2 раза, т.е. зависимость погрешности примерно пропорциональна квадра ту блока или шага между точками. Реальная относительная погрешность слоя и малой амплитуде границ будет значительно меньшей.

Эксперименты, проведенные на моделях, позволяют сделать выводы:

1. о помехоустойчивости методов решения линейной и нелинейной (кон тактной) обратных гравиметрических задач;

2. о зависимости погрешности результата от соотношения горизонтальных размеров блоков и средней глубины границы;

3. о возможности дополнительного учета априорных сведений, полученных из данных бурения или других источников;

4. относительная погрешность определения глубин ср на сложной модели со ставил 1 % - в варианте без помех и 5% - с помехами.

4.12. Некоторые рекомендации по методике решения прямых и об ратных задач гравиметрии с применением пакета ADG- 4.12.1. Решение прямых задач Главные трудности здесь связаны с ограниченностью задаваемой модели, в то время как реальная среда, влияющая на гравитационное поле, является сущест венно большей и более детальной, чем задаваемая модель.

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

Но, даже увеличив этот параметр до 10, мы можем столкнуться с конкретной за дачей, где понадобится слоев больше чем 10. Преодолеть эти затруднения можно с помощью процедуры укрупнения слоев, суть которой в том, что несколько сло ев объединяются в один с усредненной плотностью. Описанный подход позволяет переходить от моделей детальных к моделям обобщенным, которые только и име ет смысл рассматривать при решении обратных задач, связанных с большой мощностью осадочной толщи.

В составе пакета имеется процедура «Интегрировать dg», которая дает воз можность подсчитывать суммарное поле от большого числа слоев, не создавая укрупненных слоев и не подбирая усредненных плотностей, если они нас не ин тересуют.

Более детально об этом сказано в документации к пакетам.

4.12.2. Решение обратных задач При решении обратных задач главное внимание необходимо обратить на вы бор обрабатываемого профиля. Он должен проходить вкрест простирания ано мального поля, Если это условие не выполнено, то решение следует считать при ближенным. Опыт показывает, что погрешности могут доходить в этом случае до 60% и более. Они особенно велики, если обрабатываемый профиль идет вдоль длинной оси вытянутой аномалии.



Pages:   || 2 | 3 |
 





 
© 2013 www.libed.ru - «Бесплатная библиотека научно-практических конференций»

Материалы этого сайта размещены для ознакомления, все права принадлежат их авторам.
Если Вы не согласны с тем, что Ваш материал размещён на этом сайте, пожалуйста, напишите нам, мы в течении 1-2 рабочих дней удалим его.