ФИЗИКА ЗЕМЛИ, 2004, № 7, с. 71-77
УДК 550.8
ОДИН НЕЛОКАЛЬНЫЙ МЕТОД МИНИМИЗАЦИИ ФУНКЦИОНАЛОВ, ВОЗНИКАЮЩИХ В ЗАДАЧАХ ГРАВИМЕТРИИ И МАГНИТОМЕТРИИ
© 2004 г. Е. Г. Булах, Т. Н. Кишман-Лаванова
Институт геофизики им. С И. Субботина НАН Украины, г. Киев Поступила в редакцию 01.10.2003 г.
Метод подбора в обратных задачах гравиметрии и магнитометрии сводится к минимизации многопараметрических функционалов. Часто эти функционалы могут быть не одноэкстремальными, тогда решение задачи существенным образом зависит от начально построенной геологической модели. Описан метод, который позволяет в многопараметрическом пространстве вдоль выбранных направлений проследить изменение целевой функции и установить начальные значения параметров той модели, которая обусловила аномальное поле. Далее может успешно работать какой-либо локальный метод минимизации.
ВВЕДЕНИЕ
Геологическая интерпретация полевого материала представляет собой один из важнейших этапов любого геофизического метода разведки. В конечном итоге она определяет успех геологоразведочных работ. В практике интерпретации гравиметрических и магнитометрических задач используют метод подбора. Можно разделить взгляд В.Н. Страхова о том, что подбор является не единственным методом решения интерпретационных задач. Возможно, алгебраический подход станет ведущим при решении обратных задач гравиметрии и магнитометрии [Страхов, 2002]. Для этого потребуется еще некоторое время. На данном этапе дальнейшее совершенствование метода подбора остается важной вехой как в теории, так и в практике интерпретации сложных гравиметрических и магнитометрических полей. Для использования метода подбора необходимо составить начальную геологическую модель ано-малиеобразующих объектов. Здесь сформулировано понятие - класс геологических моделей.
В методе подбора выделяется два подхода. Первый - реализует метод простого моделирования. Интерпретатор работает только с программными средствами, которые позволяют в выбранном модельном классе решать прямые задачи. Сопоставляя исходное поле и теоретическое, интерпретатор сам вносит изменения в численные значения параметров модели, и вычислительный цикл начинается заново. Этот подход хорошо зарекомендовал себя. Сошлемся только на две работы - С.С. Красовского [1981] и В.И. Андреева [1992].
Второй подход - метод математического моделирования. Для сопоставления исходного и теоретического полей составляют специально построенные функционалы. Они зависят от численных
значений параметров модели. Функционалы минимизируются методами математического анализа. Получают новую геологическую модель, которая порождает аномальное поле, близкое к исходному [Страхов, 1964; Булах, 1973; Старостенко, 1978; Кобрунов, 1995; Блох, 1998 и др.].
Сосредоточим внимание на втором подходе -методе математического моделирования. Пусть в точках координатной оси абсцисс задано исходное поле аномалии силы тяжести. Требуется решить обратную задачу гравиметрии - найти те источники, которые обусловили исходное поле. В этом поле зафиксируем п точек и выделим первый массив данных
X) = (I), I = 1, 2,..., п . (1)
Если для решения интерпретационной задачи выбран метод подбора, то должна быть построена геологическая модель начального приближения, так фиксируется определенный модельный класс. Геологическая среда в этом классе характеризуется совокупностью параметров
Р = Р(р 1, Р2, ..., Рт) = р(Р]), ] = 1, 2,., т. (2)
В соотношении (2) записана в установленном порядке последовательность величин. Обычно такую последовательность называют кортежем или вектором. Если выбрана модель аномалиеоб-разующих источников и зафиксированы численные значения составляющих вектора (2), то может быть вычислено теоретическое поле.
А^,(х) = А^,(I, Р), I = 1, 2,., п. (3)
Здесь подчеркнуто, что это поле зависит от численных значений параметров модели аномалие-образующих источников. Сопоставив исходное и теоретическое поля, получим массив невязок
е(I) = е(I, Р) = Ая„(I) - Ая,(I, Р), I = 1, 2, ..., п.
(4)
Требуется найти такой вектор Р*, который минимизирует массив (4). Найденный вектор будет принят за результат решения обратной задачи.
Для минимизации невязок (4) составляют серию функционалов. Например, поля (1) и (3) можно сопоставлять в метрике Ь2, тогда записывают
^ = ^(Р) = (I) - Ая,(I, Р)]2;
1 = 1
(5)
Р = { Рь } = 1, 2, ..., т }.
Поиск минимума можно вести различными способами. Обычно их разбивают на три группы [Гельфанд, 1966]. К первой группе относят методы слепого поиска. Делается просмотр целевой функции (5) в случайно выбранных точках. Иногда последовательность таких точек выбирается по установленному правилу. Часто результаты предшествующих расчетов не используются в фиксированном вычислительном цикле.
Во вторую группу входят многочисленные методы локального поиска решения задачи. Обращаются к модели начального приближения. Параметры вектора (2) получают численные значения. Будем говорить, что для минимизации функционала (5) выбрана начальная точка. Изучаются свойства функционала в окрестности этой точки. Алгоритм минимизации должен использовать эти свойства функционала и осуществить переход к другой точке, где функция цели (5) имеет меньшее значение.
К локальным методам минимизации функционалов относятся различные модификации градиентного спуска к решению. Функционал может быть не одноэкстремальным, тогда решение задачи существенным образом зависит от выбора начального приближения. Метод градиентного спуска является устойчивым. Сошлемся только на работу В.И. Старостенко и С.М. Оганесян [2001]. Локальные методы не исключают зацикливания поиска на второстепенном мелком минимуме. Если целевая функция имеет овражный характер, то требуется перейти к специальным алгоритмическим решениям.
Третья группа методов поиска минимума многопараметрического функционала называется методами нелокального поиска. Рабочая точка в пространстве параметров перемещается не обязательно по непрерывной траектории. Здесь, как правило, не используют свойств функционала в фиксированной точке. Часто нелокальный поиск осуществляется при сочетании случайного метода с локальным. Один такой метод был изложен ранее [Булах, 1983]. Часто минимизируемый
функционал имеет овражное строение, и сходимость решения сильно замедляется. Метод позволяет алгоритмически найти такое направление, которое согласуется с осью овражной структуры. Сходимость решения значительно улучшается.
О СТРУКТУРЕ МИНИМИЗИРУЕМЫХ ФУНКЦИОНАЛОВ
Этот вопрос рассматривался ранее, но опустить его в настоящей работе или просто сослаться на статью [Булах, 1999] - это значит нарушить последовательность, полноту и целостность изложения.
В общем случае функционал (5), который требуется минимизировать, определяется т параметрами. Число т может быть большим, тогда сделать суждения о структуре функционала и точках минимума затруднительно. Если т = 2, то двухмерное поле может быть легко изображено графически. Пусть геологическая модель состоит из двух прямых уступов (рис. 1) и описывается такими параметрами:
Р = {(а, Н1, Н2, й), ] =1; 2}. (6)
В каждом уступе массы имеют постоянную избыточную плотность а. Верхняя и нижняя грани -это горизонтальные плоскости г = Н1 и г = Н2. Собственно уступ формируется вертикальной гранью х = й.
Первый пример. Пусть модель имеет такие параметры (рис. 1а):
Р1 = {(0.2; 0.2; 1; 3); (0.2; 1; 2; 1)}.
Поскольку элементы вектора Р - запись (6) - получили численные значения, можно решить прямую задачу и записать массив значений
Аяп(Х1) = Аяп(I), I = 1, 2, ..., п . (7)
В нашем примере п = 21. Поле вычислено в точках оси абсцисс на отрезке [-9; 11].
Далее будем считать, что это исходное поле аномалии силы тяжести. В отношении аномалие-образующих масс положим, что они описываются вектором (6). Пусть в этой записи переменными будут только два параметра - й1 и й2. Численные значения всех других параметров известны.
Будем задавать численные значения двух параметров й1 и й2 и вычислять теоретическое поле
АХ1) = Ая,(I, йь й2), I = 1, 2,..., п. (8)
Здесь подчеркнуто, что теоретическое поле определено параметрами й1 и й2.
Исходное поле, записанное ранее формулой (7), и теоретическое (формула (8)) сопоставляются между собой. Пусть такое сопоставление сделано в метрике Ь2, тогда получим функционал -формула (5): ^ = Ё(й1, й2); п = 21; т = 2. Это функ-
п
¿2, км (а)
модель 1
2
г, км
¿2, км (б)
-6 -4 -2 0
-2 -1 0 -1
х, км
-2 -1 0 1 2 3 4 5 6
0 2 4
-2 7 /Я\
\ "V
-4 А ч \
-6
¿1, км
-6 -4 -2 0 2 4 6 х, км
I_I_I_I_I_I_1_
модель 2
г, км
J_I_I_I_I_1_
J_I_I_I_I_I_1_
эквивалентная 2 модель
г, км
Рис. 1. (а) - Первая геологическая модель и изолинии функции Г = Г^; ¿2); (б) - вторая геологическая модель и изолинии функции Г = 1; ¿2).
ция двух переменных. Она может быть изображена в изолиниях своих значений - рис. 1а. Функция с одним экстремумом Г(3; 1) = 0. Она имеет тенденцию к овражному строению. В качестве начального приближения можно выбрать любую точку (¿1, ¿2), и из нее итерационный спуск в направлении антиградиента позволит получить достаточно точное решение задачи.
Второй пример. Пусть теперь модель имеет такие параметры:
Р2 = {(0.2; 0.2; 1; 3); (0.2; 1; 2; -3)}.
Она изображена на рис. 16. Как и в первом примере, модель позволяет решить прямую задачу и получить поле аномалии силы тяжести. Это поле примем в качестве исходного - формула (7). Оно определено в точках оси абсцисс на отрезке [-9; 11]и получено 21 значение.
Аномалиеобразующее тело описывается вектором (6). Здесь только два параметра й1 и й2 могут изменять свои численные значения. Будем фиксировать эти значения и вычислять теоретическое поле - запись (8). Для сопоставления исходного поля - формула (7) и теоретического -формула (8) построим функционал (5). Он зависит от параметров ¿1 и ¿2. Функция двух переменных может быть изображена в изолиниях своих значений - рис. 16. Это функция с двумя экстремумами. Существует помимо основной или исходной модели эквивалентная (или квазиэквивалент-
ная) модель. Она порождает аномальное поле, близкое к исходному. В рассматри
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.