CC BY
39
ВЕСТНИК ПЕРМСКОГО УНИВЕРСИТЕТА
2011 Геология Вып. 4(13)
УДК 550.831.017
Использование функции локализации с целью определения параметров источника моногеничной аномалии силы тяжести
А.Ф. Шархимуллин, А.С. Долгаль
Пермский государственный национальный исследовательский университет, 614990, Пермь, ул. Букирева, 15 E-mail: dolgal@mi-perm.ru, art-perm@mail.ru
(Статья поступила в редакцию 11 сентября 2011 г.)
Описан эффективный алгоритм, предназначенный для выделения и локализации одиночного объекта по данным гравиметрических наблюдений. Построена репрезентативная выборка допустимых решений обратной задачи гравиметрии, позволяющая оценить вероятности обнаружения аномалиеобразующего объекта в тех или иных элементарных фрагментах изучаемого объема геологической среды. Результативность метода подтверждается модельными и практическими примерами. Намечены направления дальнейших исследований в данной области.
Ключевые слова: решение обратной задачи гравиметрии, «функция локализации», аномалиеобразующий объект, невязка полей, аномальная плотность.
кальной мощностей mz, на которые также заданы интервальные ограничения:
m\ < mx < mx m7 < m7 < m7 mm max, ™” соответст-
венно. Значения аномального гравитационного поля Ag заданы вдоль профиля (оси ОХ) на криволинейной поверхности Ag = Ag(x, z) и осложнены региональным фоном Ag<^ = Ax + B и помехами s. Задается также требуемая точность решения 5 (невязка наблюденного и модельного полей в евклидовой метрике).
Следует добавить, что с физикогеологических позиций количественная интерпретация моногеничной аномалии авторами рассматривается как задача локализации ее источника (т.е. как получение наиболее достоверных сведений о пространственном положении объекта, без выявления особенностей его границы). Именно поэтому в представленной реализации метода в качестве базовой конфигурации выбрано тело простой формы - че-
Введение
Задача определения параметров одиночных возмущающих объектов, создающих локальные гравитационные аномалии, до сих пор актуальна [1]. Ее решение широко используется на начальных этапах интерпретации материалов гравиметрической съемки и приобретает особую значимость при задании буровой скважины с целью подсечения конкретного геологического объекта.
Постановка задачи
Рассмотрим задачу определения параметров источника моногеничной аномалии силы тяжести - трехмерного возмущающего тела £Т, обладающего аномальной плотностью о: о™п < а < а^, которое может залегать в интервале глубин И: И^п < И < И^ и обладать некоторыми значениями горизонтальной т и верти-
© Шархимуллин А.Ф., Долгаль А.С., 2011
тырехугольная призма. Двухсторонние ограничения на величину а отвечают пе-троплотностным параметрам геологических тел в естественном залегании (точное значение а всегда остается неизвестным !).
Для решения ОЗГ используется следующий алгоритм:
1. Случайным выбором геометрических параметров генерируются «стартовая модель» - прямоугольная призма, параметры которой отвечают заданным априорным ограничениям глубины залегания, и размеры аномалиеобразующего объекта.
%
2. Определяются значения плотности о и коэффициенты фона А, В, минимизирующие невязку полей 5
3. Если 5 >35, то происходит переход к операции 1, т.е. генерируется новое начальное приближение; в противном случае происходит уточнение геометрических параметров призмы методом случайного поиска, а также вычисление ее плот-
%
ности о и коэффициентов фона А, В.
4. Если достигнуто значение невязки полей, которое не превышает заданной величины 5, то полученные результаты являются одним из допустимых решений обратной задачи О и записываются в файл.
В результате многократного выполнения данного алгоритма мы получаем набор независимых частных решений ОЗГ (т.е. векторы геометрических О], О2, ■■■, О и плотностных о], о2,..., о параметров модельных источников). Этот набор является репрезентативной выборкой, при большом числе N характеризующей структуру множества Q(5) допустимых решений ОЗГ.
Каждое частное решение представляется в виде совокупности элементарных объемов оа, образующих регулярное замощение некоторой части нижнего полупространства. Необходимо пояснить, что регулярным замощением плоскости называется представление этой плоскости некоторым числом правильных замкнутых многоугольников (элементов замощения)
соа, плотно прилегающих друг к другу по
целой стороне [3]. В данном случае элементами замощения являются прямоугольные призмы, а для конечноэлементного описания частных решений ОЗГ используется алгоритм В.В. Ломтадзе [4]. По аналогии с детерминистской оценкой достоверности решения ОЗГ [2] строится функция ц(оо), такая, что ц(аО) = 1, если оа е О; 1л{а0) = 0, если оа £ О; при этом получается N вариантов пространственных распределений ее значений (бинарных матриц), отвечающих всем имеющимся решениям ОЗГ. Затем проводится объединение всех этих решений, сводящееся к суммированию N значений Дюа), характеризующих одни и те же относительно малые области изучаемого объема геологической среды. Функция V отвечает пересечению допустимых решений ОЗГ, каждое из которых представлено объединением конечного числа элементов соа: О]П О2^Оз ... г^Оы.
Эту функцию г(оа) можно также рассматривать как частоту попадания элемента замощения соа во фрагмент моделируемой конфигурации О в серии решений ОЗГ, что позволяет естественным образом перейти к вероятностной оценке искомых геометрических параметров. Нормированная функция ь(оа) = v(0а)/N, , названная авторами «функцией локализации», характеризует вероятность Р присутствия аномалиеобразующего объекта в элементарном объеме среды оа, полученную в результате N равновозможных исходов эксперимента.
Примеры решения обратной задачи
Имеется объект, c формой правильного шестиугольника, расположенный в интервале глубин 6-9 км, его максимальная горизонтальная мощность составляет 4 км, избыточная плотность - 0.1 г/см3.
«Наблюденное» гравитационное поле, заданное в 51 точке профиля с шагом 500 м, имеет максимальную амплитуду 1.41 мГал и осложнено помехой, обладающей нормальным законом распределения с ну-
левым средним и среднеквадратическим отклонением (СКО) ±0.1 мГал. Для генерации помехи использовался алгоритм Бокса-Мюллера. Поверхность наблюдения горизонтальная, z = 0, региональный фон отсутствует.
Рис.1. Допустимые решения обратной задачи гравиметрии для правильного шестиугольника: 1 - аномалиеобразующий объект; 2 - график исходного поля, 3 - графики модельных полей, 4 - область Б2, 5 - область
Решение ОЗГ будем искать в классе четырехугольных призм, горизонтальная мощность которых может изменяться от 2 до 15 км, вертикальная мощность - от 0.5 до 5 км, предполагаемый интервал глубин залегания 2-12 км, избыточная плотность в интервале 0.05 - 0.3 г/см3. Выбрано пороговое значение невязки полей 8 = 0.15 мГал (т.к. помимо помех расхождение наблюденного и модельного полей будет обусловлено несоответствием модельных представлений о форме аномалиеобразующего объекта), в качестве элемента замощения используется вытянутый вдоль
оси OY параллелепипед, имеющий размеры 125x125 м в вертикальной плоскости.
Всего построено #=100 частных решений ОЗГ, графическое представление этих решений дано на рис. 1.
Синтез частных решений ОЗГ иллюстрирует функция локализации и контуры областей D1 и D2. Область D1 содержит всю совокупность полученных решений ОЗГ; область D2 объединяет элементарные объемы геологического пространства, внутри которых наиболее высока вероятность присутствия искомого объекта. Принцип оконтуривания области D2 близок к выбору значения доверительной вероятности при проверке статистических гипотез, т.е. несколько субъективен.
В данном случае максимальное значе-*
ние v = 0.8, что свидетельствует о том, что при имеющемся уровне помех и глубине залегания объекта мы не можем указать ни одного элемента та, который ему гарантированно принадлежит. Однако по пороговому значению v > 0.6 (что состав*
ляет 75% от v max) уверенно локализуется связная область D2, с вполне достаточной для практики точностью совпадающая с аномалиеобразующим объектом. Отме-^ ^ * - * тим, что данный критерий v >kv max, где
0.7<^<0.8, позволяет успешно оконтури-вать источники поля в целом ряде случаев, для которых также проводились подобные расчеты, но их рассмотрение выходит за рамки данной статьи.
Как очевидно, высокая степень совпадения полей и соответствие геометрических и физических параметров подобранной модели имеющимся априорным сведениям (заданным в виде интервальных оценок) не всегда влекут за собой геологическую содержательность результата, т.е. полученные в частном решении ОЗГ модели могут существенно отличаться от реального источника аномалии. О степени этих различий можно судить по размерам области D1. Синтез результатов серии допустимых решений задачи (целесообразно выбирать #=50-100), сводящийся к построению дискретной функции v*, характеризующей вероятность наличия ано-
мальных масс в элементарных объемах среды (оа, т.е. построению области D2, позволяет значительно повысить достоверность результата количественной интерпретации.
При использовании того же источника поля и той же априорной информации, что и в предыдущем примере, попытаемся оценить помехоустойчивость метода (рис. 2). Зададим три варианта нормальной помехи в «наблюденном» поле с разными амплитудами СКО - ±0.05, ±0.15, ±0.25 мГал; для каждого из этих вариантов выполним серию из #=50 решений ОЗГ. Полученные результаты приведены на рис. 3, который свидетельствует о том, что даже при высоком уровне помех, достигающем 18% от максимальной амплитуды поля, удается уверенно локализовать возмущающий объект (в данном случае при построении области D2 использовался кри-** терий v > 0.7 v max).
Следует обратить внимание на закономерное увеличение площади области D1 с ростом амплитуды помех и на эффекты
нарушения связности области Д2, также происходящие при этом. По аналогии с детерминистскими методами [5] представляется целесообразным меру те8{р2Д{) рассматривать как характеристику степени проявления 5-
эквивалентности для конкретной практической задачи.
Рассмотрим еще один модельный пример: гравитационная аномалия, зафиксированная на профиле длиной 10 км с постоянным шагом между измерениями 200 м, обусловлена объектом с плотностью 0.25 г/см3, с глубинами залегания кровли 250 м, подошвы - 1500 м.
Высотные отметки точек изменяются от 90 до 240 м, в поле присутствуют помеха (СКО = ±0.3 мГал) и региональный фон вида Лgф = 0.5х + 1.25. Решение задачи выполнялось в двух вариантах, при использовании различных априорных данных.
Рис.2. Результаты синтеза допустимых решений обратной задачи гравиметрии для правильного шестиугольника при различном уровне помех в «наблюденном» поле (А - ±0.05 мГал; Б -±0.15 мГал; В - ±0.25 мГал): 1 - аномалиеобразующий объект; 2 - график «наблюденного» поля; 3 - контур области D2, 4 — контур области D1
Как очевидно (рис. 3), результаты решения ОЗГ в варианте 2 обладают меньшей неоднозначность (е-
эквивалентностью). Следовательно, можно построить методику двухэтапного решения практических задач: на 1-м этапе выполняется уточнение априорных ограничений; на 2-м - проводится повторное решение ОЗГ с новыми, более жесткими, ограничениями на искомые параметры.
Необходимым элементом рациональной методики должно являться определение требуемой точности совпадения полей
5, т.к. уровень помех |£| в наблюденном
поле силы тяжести неизвестен. К тому же в реальных условиях определенный вклад в величину 5 будут вносить погрешности аппроксимации источника для используемого модельного класса тел. Целесообразно задаться несколькими допустимыми значениями невязки 5 : 51 < 52 < ••• <5к и начать вычисления при минимальном из них, т.е. при 5]. В случае, если после выполнения значительного числа (~10 ) циклов поиска квазирешение ОЗГ построить не удастся, следует переходить на следующее, большее пороговое значение невязки.
На рис. 4 приведены результаты решения ОЗГ по данным крупномасштабной гравиметрической съемки, выполненной над месторождением платино-медно-никелевых руд Норильск-1. Предполагается, что аномалия силы тяжести в основном обусловлена рудоносной интрузией базит-гипербазитового состава, ее избыточная плотность (по отношению к вмещающим породам трапповой формации) составляет от 0.2 до 0.3 г/см3, глубина залегания - 0.2 км < И < 1.5 км; размеры -2км < шх < 5км и 100м < тШ < 500м . Построение N = 50 решений ОЗГ при 5= 0.25
* *
мГал и их синтез при V = 0.7у тах позволяют уверенно локализовать в разрезе наиболее мощную часть интрузии.
Количественная интерпретация гравитационного поля по профилю скв. 2- скв. 4 проводилась, исходя из предположения о том, что отрицательная локальная ано-
малия обусловлена одним объектом, залегающим в интервале глубин от 5 до 300 м, обладающим аномальной плотностью от -0.02 до -0.1 г/см3. Размеры 2.5Б аномалиеобразующего объекта по горизонтали могут изменяться от 10 до 500 м, по вертикали - от 20 до 150 м, размер по простиранию - 300 м.
Рис. 3. Результаты решения обратной задачи гравиметрии при разной априорной информации: А - графики поля; Б - решение 1; В - решение 2. Графики полей: 1 - "наблюденного", 2 - "наблюденного" без помехи, 3 - фонового; 4 - локальной составляющей; 5 - аномалиеобразующий объект; 6 - контуры подобластей D2 (сплошная линия) и D1 (пунктирная линия)
На рис. 5 представлена интерпретация, которая проводилась путем построения функции локализации V, характеризующей вероятность Р присутствия аномалиеобразующего объекта в элементарном объеме геологической среды (оа, полученную в результате N равновозможных исходов эксперимента. Для этого выполнялось N = 30 независимых решений обратной задачи гравиразведки путем случайного поиска, а далее проводился синтез полученных результатов. Очевидно, что в пределах области с нулевыми значениями функции V =0 аномалиеобразующие массы будут отсутствовать, а область значений v=1 представляет собой фрагмент, гарантированно принадлежащий возмущающему объекту.
Дц. мГал
Рис. 4. Интерпретация гравитационного поля над месторождением платино-медно-никелевых руд Норильск-1: 1 - породы туфолавовой толщи; 2 - отложения тунгусской серии; 3 - силлы габбро-долеритов; 4 - рудоносная интрузия; 5 - график наблюденного поля; 6 - область 02 (сплошная линия), область В1(пунктирная линия)
Рис. 5. Результаты количественной интерпретации локальной гравитационной аномалии (А - съемка 2010 г., Б - съемка 2011 г.)
Опыт применения метода свидетельст-
**
вует, что критерий V >№ тах, где
0.7<^<0.8 позволяет успешно оконтури-вать источники поля в подавляющем большинстве случаев.
Результаты интерпретации, проведенной по результатам съемок 2010 и 2011 гг., представлены на рис. 5,А и 5,Б соответственно. Средняя невязка для наблюденного и модельного полей в серии из 30 решений составила ±0.019 мГал (2010 г.) и ±0.026 мГал (2011 г.). Сплошной линией оконтурена сравнительно широкая область, в пределах которой располагаются все допустимые решения обратной задачи. В пределах этой области изолинией 0.8 локализовано наиболее вероятное местоположение источника, пунктирная изолиния 0.99 отвечает фрагменту источника, фиксирующемуся при всех решениях задачи. Средняя величина аномальной плотности возмущающего объекта для 2010 г. - -0.087 г/см3, для 2011 г. - -0.099 г/см3. Наблюдается отчетливая тенденция к увеличению размеров зоны разуплотнения и приближению ее верхней границы к поверхности. Глубина этой границы в 2010 г. была порядка 50 м, в 2011 г. - уже порядка 15 м. При этом следует отметить, что местоположение нижней границы, находящейся на уровне около 100 м, практически не изменилось, т.е. отражающиеся в гравитационном поле объекты не связаны с глубинными процессами разуплотнения горных пород
Заключение
В настоящее время вычислительные возможности компьютеров возросли до такой степени, что появилась реальная возможность решения принципиально новых интерпретационных задач для многих геофизических методов, в т.ч. и для гравиразведки. Выявление закономерностей
Библиографический список
1. Бабаянц П.С., Блох Ю.И., Трусов А.А. Ло-
кальная интерпретация: от интерпретационных профилей к интерпретационным окнам // Матер. 32-й сессии Междунар.
при рассмотрении конечного множества допустимых решений ОЗГ является эффективным способом повышения качества получаемых геологических результатов. Техническая реализация построения такого множества для данного алгоритма (при числе вариантов решений N = 100) не представляет каких-либо затруднений и требует не более 5-7 мин машинного времени.
Практическое выполнение интерпретационных построений с применением функции локализации V может существенно повысить вероятность вскрытия искомых аномалиеобразующих объектов в заданных интервалах глубин поисковыми и разведочными скважинами, рекомендованными по гравиметрическим данным.
Данный подход к повышению достоверности количественной интерпретации можно попытаться развивать на более сложные распределения масс, т.е. на присутствие в моделируемом разрезе нескольких геоплотностных неоднородностей.
Функцию локализации V можно применить и для комплексной интерпретации данных гравиразведки и магниторазведки при интервальном представлении априорных данных о плотности а : ать < а <
атах и намагниченности J : /тт < J < /тах изучаемых объектов.
Работа выполнена при поддержке РФФИ (проект 10-05-96023_р_урал_а) и программы исследований ОНЗ РАН (проект 09-Т-5-1031).
Авторы выражают глубокую благодарность доктору физико-математических наук П.И. Балку за конструктивное обсуждение результатов исследований.
семинара им. Д. Г. Успенского. Пермь, 2007. С. 17-20.
2. Балк П.И., Долгаль А.С. Детерминированный подход к проблеме достоверности результатов интерпретации гравиметрических данных // Докл. Академии наук. 2010. Т. 431, № 1. С. 334-338.
3. Гольдшмидт В.И. Оптимизация процесса
количественной интерпретации данных гравиразведки. М.: Недра, 1984. 184 с.
4. Ломтадзе В.В. Программное обеспечение
обработки геофизических данных. Л.: Недра, 1982. 280 с.
5. Страхов В.Н. Общая схема и основные итоги развития теории и практики интерпретации потенциальных полей в СССР и России в ХХ веке // Развитие гравиметрии и магнитометрии в ХХ веке: тр. конф. (Москва, 23-26 сентября 1996 г.) М.: ОИФЗ РАН, 1997. С. 98-120.
The Localization Function for the Definition of Source Parameters of Gravity Monogenetic Anomaly
A.F. Sharhimullin, A.S. Dolgal
Perm State National Researching University, 614990, Perm, Bukirev st., 15 E-mail: dolgal@mi-perm.ru, art-perm@mail.ru
In the paper we describe an effective algorithm for dedication and localization of the single object on gravity data. The algorithm creates the representative sampling of the feasible solution of the gravity inverse problem. This allows estimating the location probability of the gravity source on different fragments of studied geological medium. The effectiveness of the method is confirmed on model and practical examples. We plan directions for subsequent researches in this field.
Key words: the solution of the gravity inverse problem, localization function, the object which create the anomaly, field residual, anomalous density.
Рецензент - доктор физико-математических наук П.И. Балк
