CC BY
0
2023 Геология Том 23, № 1
УДК 550.312
Оценка параметров источника гравитационной аномалии методом роя частиц
А.С. Долгальа' ь, Р.Н. Петросянь с, Н.В. Рыжовь, с
а Горный институт УрО РАН
614007, Пермь, ул. Сибирская, 78а. Е-шаИ: dolgal@mi-perm.ru
ь Пермский государственный национальный исследовательский университет
614068, Пермь, ул. Букирева, 15. Е-шаП: petrosyan@antikarst.ru
с ООО «Противокарстовая и береговая защита»
614068, Пермь, ул. Данщина, 4а. Е-mail: ryzhov@antikarst.ru
(статья поступила в редакцию 29 декабря 2023 г.)
Рассматривается алгоритм решения нелинейной обратной задачи гравиразведки методом роя частиц для двухмерной горизонтальной призмы. После сравнительно небольшого числа итераций проведено определение размеров и глубины залегания возмущающего объекта. Использованы три варианта корректировки скоростей перемещения роя в процессе оптимизации. Результаты интерпретации представлены в виде ге-оплотностной модели, отвечающей минимуму невязки наблюденного и теоретического полей, а также вертикальной карты изолиний функции локализации. Выполнена количественная интерпретация силы тяжести Березовской аномальной зоны (Пермский край). Сделан вывод о высокой перспективности применения роевого интеллекта для оценки параметров источников гравитационных аномалий.
Ключевые слова: гравиразведка, аномалия, эволюционная оптимизация, метод роя частиц, источник поля, функция локализации.
БОТ: 10.17072^и^ео1.23.1.85
Введение
В настоящее время методы искусственного интеллекта (ИИ) являются одним из приоритетных направлений развития науки в России и в мире. В английском языке термин artificial intelligence, в русском переводе -ИИ, означает лишь искусственное умение рассуждать разумно, а не интеллект (intellect) (Гаврилова, Хорошевский, 2000). Единого определения для ИИ не существует. Согласно А.В. Смолину (2007), «интеллектуальной называется система, способная целеустремленно, в зависимости от состояния информационных входов, изменять не только параметры функционирования, но и сам способ своего поведения». Академик В.Н. Страхов (1988) еще в конце ХХ в. писал, что в будущем необходимо «создание основ специализированного искусственного интеллекта как базовой позиции в становлении автоматической компьютерной интерпретации». Следует отметить, что геофизика является хорошо структурированной предметной областью, в которой уже существуют апробированные надежные методы, и некоторые
задачи, тесно связанные с ИИ, в ней уже успешно решаются. В качестве примера можно привести алгоритмы классификации многомерных данных (Никитин, Петров, 2008) и распознавания образов (Калинин, 2011).
Метод роя частиц (МРЧ, англ. particle swarm optimization, PSO) рассматривается в теории ИИ как метод эволюционной оптимизации (Саймон, 2020) и позволяет моделировать поведение децентрализованной самоорганизующейся системы. Первоначально этот метод был предложен для описания движения птичьих стай с целью выяснения законов, управляющих синхронным движением и быстрой сменой направления полета. С этих позиций стая птиц рассматривается как система роевого интеллекта (англ. swarm intelligence), состоящая из множества частиц, локально взаимодействующих между собой и окружающей средой (Eberhart, Kennedy, 1995).
Текущее состояние частицы характеризуется координатами в пространстве решений, а также вектором скорости перемещения. Оба этих параметра выбираются случайным образом на этапе инициализации. Кроме того,
© Долгаль А.С., Петросян Р.Н., Рыжов Н.В., 2024
85
каждая частица хранит координаты лучшего из найденных ей решений, а также лучшее из пройденных всеми частицами решений -этим имитируется мгновенный обмен информацией между птицами. В итоге локальные и в некоторой степени случайные взаимодействия приводят к возникновению интеллектуального группового поведения, неконтролируемого отдельными особями (Казакова, 2022).
В настоящее время МРЧ широко применяется в задачах машинного обучения, параметрической и структурной оптимизации в области проектирования, в областях робототехники, биохимии, биомеханики и др. Иранские геофизики успешно использовали МРЧ для моделирования поверхности фундамента нефтегазоносного осадочного бассейна по аномальному гравитационному полю в провинции Голестан (Ьош, МеЬгат^, 2022). МРЧ также достаточно эффективно использовался для моделирования аномалий высших производных гравитационного потенциала, связанных с дизъюнктивными нарушениями (ЕШшБет е! а1., 2023). В отечественной литературе пока имеется единственный пример применения МРЧ при инверсии данных гравиразведки (Долгаль, 2023).
Основной теоретической предпосылкой для решения обратных задач гравиразведки (ОЗГ) с использованием МРЧ является существование множества равновероятных (в общем случае) допустимых решений, отвечающих априорным представлениям о геометрических параметрах и эффективной плотности источников поля, а также обеспечивающих требуемое значение невязки наблюденного и модельного полей е0 (Долгаль, 2018). Рассмотрим возможность решения нелинейной обратной 2Б-задачи гравиразведки рудного типа для простейшей ситуации (одиночный изолированный объект с известной эффективной плотностью) с использованием МРЧ. В частности, таким объектом может быть интрузивное тело или рудная залежь, которые требуется локализовать в разрезе.
Алгоритм метода роя частиц и результаты его применения
Будем использовать аппроксимацию ано-малиеобразующего геологического объекта с
известной эффективной плотностью а = const прямоугольной призмой с бесконечными размерами по простиранию, находящейся внутри области поиска решения обратной задачи D. Введем обозначения для координат {x0, zc} центра тяжести C призмы и размеров d, h ее горизонтальных и вертикальных граней соответственно. Таким образом, призма полностью характеризуется четырехмерным вектором параметров р = [х0, z0, d, h} и в дальнейшем рассматривается в качестве отдельной частицы роя (Долгаль, 2023). Рой состоит из ^-частиц со случайными параметрами pj Е R4, i = 1, JV, равномерно распределенными в подобласти D0 с D. Каждой частице роя соответствуют модельное гравитационное поле ДдМ°д и некоторое значение функционала F 2 = ||Д,днабл —
Д^Мод lli2.
Состояние всего роя характеризуется наилучшим состоянием его частиц, таким образом ОЗГ сводится к минимизации функционала min F2(p), которое осуществляется в
pED
течение некоторого времени, разбитого на к = 1,М интервалов. Скорость частицы v Е R10 в текущий момент времени к+1 определяется суммой трех слагаемых:
v, (к + 1) = avi (к) + £Ui®[L (к) — р (к)] + YV2®[G(k) — р (к)], (1)
т.е. инерционной, когнитивной и социальной компонент. Инерционная компонента зависит от скорости частицы в предыдущий момент времени Vj(к), когнитивная отражает стремление частицы к собственному наилучшему положению (локальному рекорду) L(к), социальная показывает ее движение к лучшему положению роя (глобальному рекорду) G(к). Символ (В> означает прямое (покомпонентное) произведение векторов. Длины случайных векторов U равномерно распределены в интервале [0,1], a, fi, у - коэффициенты, которые являются постоянными величинами или зависят от времени к. Схема алгоритма МРЧ показана на рис. 1. Этот алгоритм был реализован авторами в виде программы PDBR_ROI, написанной на языке Delphi 7.0 (рис. 2).
Рис. 1. Блок-схема алгоритма МРЧ
Рис. 2. Главное окно программы РБЬЯКОТ
Рассмотрим синтетический пример: на профиле длиной I = 50 км зафиксирована аномалия гравитационного поля Дднабл с амплитудой до 25 мГал, обусловленная двухмерной горизонтальной пятиугольной призмой, обладающей эффективной плотностью 0,25 г/см3. Шаг Дх между точками составляет 1 км. В области поиска решения обратной задачи Б размером 50x25 км (разрезе) разместим N = 100 призм с плотностью 0,25 г/см3 (частиц).
На рис. 3 представлены исходные данные и начальное положение используемого роя частиц. Динамику поиска наилучшего положения частицами роя при числе итераций М = 40 отражает рис. 4. После выполнения М = 40 итераций была достигнута средняя величина невязки наблюденного и модельного полей ^2 = 0,79 мГал для всего роя частиц. Для наилучшего решения 0(40) величина невязки составила 0,34 мГал, т.е. около 1,5 % от максимальной амплитуды аномалии.
-10-
-15-
-20-
Z, км
Рис. 3. График аномалии силы тяжести поля и источник поля (черный контур с заливкой). Синие точки - центры призм С для начального положения роя
Рис. 4. Пространственное расположение центров призм С после выполнения 5 (а), 15 (б), 30 (в) и 40 (г) итераций
Устойчивость процесса оптимизации в МРЧ обеспечивает ввод ограничений на максимальные значения скорости частиц V,<Утах. Выбор величины \тах столь же важен, как выбор параметра регуляризации в классических методах решения обратной задачи гравиразведки, где решаются плохо обусловленные системы алгебраических уравнений. В общем случае для рассматриваемой задачи целесообразно опираться на соотношение Утах < 0,5Дх для параметров (х0, г0] и условие 0,9ру < ру+1 < 1,1ру, ] = 3,4 для граней призмы.
Используются 3 различных варианта изменения скоростей во времени, отвечающие выражению (1). В первом варианте используются постоянные коэффициенты <х= 0,7298, р = у = 1,4962. В втором варианте а = ^(0,4 - 0,9)/М + 0,9, £ = ^(0,4945 - 1,4945)/М + 1,4945, у = ^1,4945 - 0,4945)/М + 0,4945, т.е. величина Р линейно уменьшается с ростом итерации, а величина у линейно возрастает. Это позволяет сначала искать лучшее положение каждой частицы, со временем смещая акцент на лучшее положение роя. В третьем случае используется коэффициент сужения д = 0,5714, представляющий собой множитель в
правой части (1), при постоянных значения коэффициентов к= 1, р = у = 2,05:
Vj (к + 1) = /¿{avi (к) + £U!®[L(к) -р(Ю] + 7U2®[G(fc) -p(fc)]}. (2)
Каждый из вариантов определяет собственный ход процесса оптимизации роя частиц и может иметь преимущества при решении конкретной задачи. Однако, согласно теореме NFL (no free lunth theorem), «все алгоритмы оптимизации работают одинаково хорошо при усреднении по всем возможным задачам» (Саймон, 2020), поэтому выбор варианта может быть только экспериментальным. Применительно к рассматриваемому синтетическому примеру вариант 3 заметно проигрывает по отношениям к вариантам 1 и 2 (рис. 5).
МРЧ может оказаться весьма эффективным инструментом для построения репрезентативного подмножества допустимых решений ОЗГ. Можно задаться некоторым малым пороговым значением невязки F2n0p и выбрать п частиц роя, для которых величина F2 < F2n0p. Каждую из этих частиц можно описать связной совокупностью конечных элементов ш в форме квадратов, образующих замощение области D.
Рис. 5. Характеристики итерационного процесса различных вариантов решения ОЗГ: значения Е2 отвечают глобальному рекорду О(р). Номер варианта указан на графиках
Обозначим через п* число всех призм, в каждой из которых содержится элемент ш. Тогда р = п*/п - оценка вероятности ситуации, когда элемент ш является фрагментом неизвестного носителя масс (Балк, Долгаль, 2020).
Функцию пространственных координат = р(х,г) с областью определения [0,1] назовем функцией локализации, характеризующей структуру полученного множества допустимых решений ОЗГ. В качестве новой математической формы представления
результатов интерпретации будем использовать одну или несколько «наилучших» результативных моделей и вертикальную карту функции локализации 'ф(й)). На рис. 6 показаны результаты решения ОЗГ в варианте 1, отвечающие значению ^2пор = 0,4 мГал, числу частиц п = 28 и размеру элементов замощения 100x100 м. Как очевидно, качество выполненных интерпретационных построений высокое и является вполне достаточным для заложения заверочной скважины.
Рис. 6. Результаты решения обратной задачи гравиразведки: а — призма, отвечающая минимальному значению невязки полей Р2 = 0,34 мГал (синий пунктир); б — карта изолиний функции локализации Черный контур — аномалиеобразующий объект
Приведем результаты применения МРЧ при моделировании источника Березовской аномальной зоны в пределах Юрюзано-Сыл-венской депрессии на юго-востоке Пермского края. Выполнена интерпретация данных сред-немасштабной гравиметрической съемки по профилю АБ в предположении, что избыточная плотность аномального объекта равна 0.45 г/см3. Соотношение между протяженностью (свыше 100 км) и шириной (около 20 км) аномальной зоны допускает использование 2D-модели источников (рис. 7 а).
а
Размерность роя составила N = 100 частиц, их перемещение происходило в течение М = 50 интервалов времени. Величина средней невязки для роя после завершения процесса эволюционной оптимизации составила 0,70 мГал. Можно предположить, что аномалия силы тяжести обусловлена интрузивным телом ультраосновного состава мощностью ~1 км, обладающим высокой избыточной плотностью и низкой намагниченностью, кровля которого располагается на глубине около 2 км от дневной поверхности (рис. 7 б).
б
0 км 10 км 20 км 30 км
А
: 1 2 I
Рис. 7. Карта изоаномал гравитационного поля (а) и результаты интерпретации Березовской аномальной зоны (б): 1 - интерпретационный профиль; 2 - график наблюденного поля; 3 - модель источника гравитационной аномалии, обеспечивающая минимум невязки полей Е2 = 0,69 мГал
Заключение
Полученные результаты подтверждают высокую перспективность применения МРЧ в области прикладной геофизики для моделирования источников аномалий геопотенциальных полей. Разработанный алгоритм обладает высокой скоростью сходимости, что делает возможным получение приемлемого по качеству решения ОЗГ даже при небольшом числе итераций. Нужно отметить, что для получения близких результатов с помощью метода Монте-Карло (Долгаль, Петросян, 2021) требуется примерно в 100 раз больше пробных решений прямой задачи гравиметрии.
Особый интерес представляет синтез МРЧ и аддитивных технологий количественной
интерпретации геопотенциальных полей, один из примеров которого представлен выше. Дальнейшим направлением исследований в этой области является разработка модификации МРЧ для многоугольных призм, нескольких аномалиеобразующих 2Б-тел и для решения 3Б ОЗГ. Также может оказаться эффективным применение других методов роевого интеллекта - муравьиного и пчелиного алгоритмов, метода серых волков, гравитационного поиска и др.
Разработанная программа РББК_К01 может использоваться студентами при выполнении лабораторных работ, что позволит им познакомиться с методикой имитационного моделирования и закрепить основные понятия из области теории решения некорректных
задач математической физики. Важной особенностью программы как средства обучения является уникальный характер каждого вычислительного цикла, позволяющий судить о самостоятельном получении имеющихся результатов расчета. Программа PDBR_ROI также позволяет осуществлять решение ряда практических задач и может эксплуатироваться в производственных организациях.
Исследование выполнено при финансовой поддержке Минобрнауки России в рамках государственного задания (регистрационный номер НИОКТР 1022040500598-4-1.5.6).
Библиографический список
Eberhart R., Kennedy J. Particle swarm optimization. Proceedings of the IEEE International Conference on Neural Networks, 4, IEEE, 1995. pp. 19421948. DOI: 10.1109/ICNN.1995.488968.
ElhusseinM., Abo-EzzE.R., Gomaa O.A. et al. Two-sided fault parameter assessment through the utilization of a particle swarm optimizer to different gravity horizontal gradients-orders with application to various cases studies. Environ Earth Sci 82, 502 (2023). DOI: 10.1007/s12665-023-11081-0
Loni S., Mehramuz M. Comparison of improved particle swarm optimization with Marquardt Algorithm for simulation of sedimentary basin with parabolic density contrast using gravity data // J. Ind. Ge-ophys. Union. 2022. 26 (1). pp. 1-15. DOI: 10.1186/2193-1801-2-462
Балк П.И., Долгаль А.С. Аддитивные технологии решения обратных задач гравиразведки и магниторазведки. М., Научный мир. 2020. 455 с.
Гаврилова Т.А., Хорошевский В.Ф. Базы знаний интеллектуальных систем. СПб: Питер, 2000. 384 с.
Долгаль А.С., Петросян Р.Н. Решение обратной задачи гравиразведки для 2D призматических тел методом статистических испытаний // Вестник Пермского университета. Геология. 2021. Т. 20. № 4. С. 334-343. БО1: 10.17072Zpsu.geoL 20.4.334
Долгаль А.С. Элементы искусственного интеллекта в монтажном методе интерпретации гравитационных аномалий // Стратегия и процессы освоения георесурсов: сборник научных трудов. Вып. 16. Пермь, 2018. С. 146-149. БО1: 10.7242Zgdsp.2018.16.38
Долгаль А.С. Локализация источника аномалии силы тяжести методом роя частиц // Теория и практика разведочной и промысловой геофизики: сборник научных трудов, Пермь, 23-24 ноября 2023 года. Пермь: Пермский государственный национальный исследовательский университет, 2023. С. 98-104.
Казакова Е.М. Краткий обзор методов оптимизации на основе роя частиц // Вестник КРАУНЦ. Физико-математические науки. 2022. Т. 39. № 2. С. 150-174. DOI: 10.26117/2079-6641-2022-39-2150-174
Калинин Д.Ф. Информационно-статистический прогноз полезных ископаемых. СПб: ФГУНПП «Геологоразведка», 2011. 163 с.
Никитин А.А., Петров А.В. Теоретические основы обработки геофизической информации. М: ГНЦ РФ ВНИИгеосистем, 2008. 116 с.
Саймон Дэн. Алгоритмы эволюционной оптимизации / пер. с англ. А.В. Логунова. М.: ДМК Пресс, 2020. 1002 с.
Смолин А.В. Введение в искусственный интеллект: конспект лекций. 2-е изд., перераб. М.: ФИЗ-МАТЛИТ, 2007. 264 с.
Страхов В.Н. Что делать? (о развитии гравиметрии и магнитометрии в России в начале 21 века). М.: ОИФЗ РАН. 1988. 24 с.
Estimation of Source Parameters of Gravitational Anomaly by Particle Swarm Optimization
A.S. Dolgala b, R.N. Petrosyanb c, N.V. Ryzhovb c
a Mining Institute UB RAS
78a Siberian Str., Perm 614007, Russia. E-mail: dolgal@mi-perm.ru b Perm State University
15 Bukirev Str., Perm 614068, Russia. E-mail: petrosyan@antikarst.ru c LLC Antikarst and coastal protection
4a Danshchina Str., Perm 614068. Russia. E-mail: ryzhov@antikarst.ru
An algorithm for solving a nonlinear inverse problem of gravity exploration using the particle swarm optimization for a two-dimensional horizontal prism is considered. After a relatively small number of iterations, the size and depth of the disturbing object were determined. Three options were used to adjust the speed of movement of the
swarm during the optimization process. The interpretation results are presented in the form of a geodensity model that corresponds to the minimum discrepancy between the observed and theoretical fields, as well as a vertical isoline map of the localization function. A quantitative interpretation of the gravity force of the Berezovskiy anomalous zone (Perm region) was carried out. It is concluded that the use of swarm intelligence is highly promising for assessing the parameters of gravitational anomalies sources.
Key words: gravity prospecting; anomaly; evolutionary optimization; particle swarm optimization; field source; localization function.
References
Eberhart R., Kennedy J. 1995. Particle swarm optimization. In: Proceedings of the IEEE International Conference on Neural Networks, 4, IEEE, pp. 19421948. doi: 10.1109/ICNN.1995.488968.
Elhussein M., Abo-Ezz E.R., Gomaa O.A. et al. 2023. Two-sided fault parameter assessment through the utilization of a particle swarm optimizer to different gravity horizontal gradients-orders with application to various cases studies. Environ Earth Sci 82:502. doi: 10.1007/s12665-023-11081-0
Loni S., Mehramuz M. 2022. Comparison of improved particle swarm optimization with Marquardt Algorithm for simulation of sedimentary basin with parabolic density contrast using gravity data. J. Ind. Geophys. Union. 26(1): 1-15. doi: 10.1186/21931801-2-462
Balk P.I, Dolgal A.S. 2020. Additivnye tekhnologii resheniya obratnykh zadach gravi-razvedki i magnitorazvedki. [Additive methods for solving inverse problems of gravity and magnetic prospecting]. Moskva, Nauchnyy mir. p. 455 (in Russian).
Gavrilova T.A., Khoroshevskiy V.F. 2000. Bazy znaniy intellektualnykh sistem [Intelligent systems knowledge bases]. St. Petersburg, Piter, p. 384 (in Russian).
Dolgal A.S., Petrosyan R.N. 2021. Resheniye obratnoy zadachi gravirazvedki dlya 2D prizmatich-eskikh tel metodom statisticheskikh ispytaniy [Solution of the inverse problem of gravity exploration for 2D prismatic bodies by the statistical tests method]. Vestnik Permskogo universiteta. Geologiya. 20(4): 334-343 (in Russian). doi: 10.17072/psu.geol. 20.4.334
Dolgal A.S. 2018. Elementy iskusstvennogo intel-lekta v montazhnom metode interpretatsii gravita-
tsionnykh anomaliy [Elements of artificial intelligence in the installation method for interpreting gravitational anomalies]. In: Strategiya i protsessy osvoy-eniya georesursov. Perm, p. 146-149 (in Russian). doi: 10.7242/gdsp.2018.16.38
Dolgal A.S. 2023. Lokalizatsiya istochnika anom-alii sily tyazhesti metodom roya chastits [Localization of the source of gravity anomalies by the particle swarm method]. In: Teoriya i praktika razvedochnoy i promyslovoy geofiziki. Perm, Perm State University, p. 98-104 (in Russian).
Kazakova Ye.M. 2022. Kratkiy obzor metodov op-timizatsii na osnove roya chastits [Short overview of particle swarm optimization methods]. Vestnik KRAUNTS. Fiziko-matematicheskiye nauki. 39(2): 150-174 (in Russian). doi: 10.26117/20796641-2022-39-2-150-174
Kalinin D.F. 2011. Informatsionno-statisticheskiy prognoz poleznykh iskopayemykh [Information and statistical forecast of mineral resources]. St. Petersburg, Geologorazvedka, p. 163 (in Russian).
Nikitin A.A., Petrov A.V. 2008. Teoreticheskiye osnovy obrabotki geofizicheskoy informatsii [Theoretical foundations of geophysical information processing] Moskva, VNIIgeosistem, p. 116 (in Russian).
Simon Dan. 2013. Evolutionary Optimization Algorithms. John Wiley & Sons. p. 784.
SmolinA.V. 2007. Vvedeniye v iskusstvennyy intellect [Introduction to artificial intelligence]. Moskva, FIZMATLIT, p. 264 (in Russian).
Strakhov V.N. 1988. Chto delat? (o razvitii gravi-metrii i magnitometrii v Rossii v nachale 21 veka) [What to do? (on the development of gravimetry and magnetometry in Russia at the beginning of the 21st century)]. Moskva, IEP RAS. p. 24. (in Russian)
