CC BY
21Научный вестник НГТУ. - 2012. - № 3(48)
УДК 519.63
SD-инверсия данных вызванной поляризации на основе конечно-элементного моделирования*
М.Г. ПЕРСОВА, Ю.Г. СОЛОВЕЙЧИК, М.Г. ТОКАРЕВА, Е.И. СИМОН, М.В. АБРАМОВ, У.А. САФРОНОВА
Предлагается подход к выполнению BD-инверсий данных электроразведки методом вызванной поляризации с использованием для решения прямой задачи метода конечных элементов. Приводится описание математического аппарата решения обратной задачи и примеры численных экспериментов.
Ключевые слова: BD-инверсии данных вызванной поляризации, BD-моделирование, метод конечных элементов. ВВЕДЕНИЕ
Основными подходами к решению обратной задачи электроразведки при выполнении геофизических исследований по методу вызванной поляризации (ВП) являются подходы, основанные на моделировании электрического поля в одномерной среде. При этом в их основе лежит решение прямой задачи, в которой эффекты вызванной поляризации учитываются с использованием феноменологической модели Cole-Cole, основанной на кажущейся зависимости удельной проводимости среды от частоты. Кроме того, использование одномерных моделей среды при обработке данных, полученных при изучении существенно неоднородных трехмерных сред, может приводить к довольно серьезным ошибкам в определении изучаемых параметров среды и, как следствие, к ошибкам в геологическом прогнозе. В этой связи очевидна необходимость создания программно-математических средств для обработки данных по методу ВП, основанных на решении трехмерной обратной задачи.
В настоящее время в зарубежных программных разработках используются подходы, в которых поле ВП нелинейно зависит от параметров ВП (с использованием модели Cole-Cole или зависимости проводимости среды от поляризуемости на постоянном токе), что приводит к необходимости многократного решения (и при этом на достаточно большом количестве частот) набора векторных гармонических или, в лучшем случае, скалярных эллиптических задач [1-3]. Это влечет за собой очень высокую вычислительную сложность соответствующих алгоритмов инверсии. В отличие от них предлагаемые в данной работе алгоритмы 3D-инверсии будут основаны на однократном решении набора вычислительно не слишком трудоемких трехмерных эллиптических задач, причем без какого-либо снижения адекватности моделируемых процессов реальным физическим процессам, что позволит существенно повысить вычислительную эффективность и может способствовать скорейшему внедрению таких подходов к выполнению 3D-инверсий в практику геологоразведочных работ, проводимых по методу ВП.
* Получена 30 марта 2012 г.
Работа выполнена при поддержке ФЦП «Исследования и разработки по приоритетным направлениям развития научно-технологического комплекса России» на 2007-2013 гг., ФЦП
«Научные и научно-педагогические кадры инновационной России» на 2009-2013 гг. и гранта Президента Российской Федерации для государственной поддержки молодых российских ученых - докторов наук (№ грантаМД-1925.2011.5).
1. МАТЕМАТИЧЕСКИЙ АППАРАТ РЕШЕНИЯ ТРЕХМЕРНОЙ ОБРАТНОЙ ЗАДАЧИ
Для описания процессов ВП, возбуждаемых гальваническими источниками (горизонтальной электрической линией, либо вертикальной электрической линией, либо круговым электрическим диполем) [4], будем использовать математическую модель вида [5-6]
^у(с(х, у, z)grad V/P ) = div(а(х, у, z)р(т(х, у, z), t)сЁрс ),
1/1Р ТЭТ-Г Г РС
где V - потенциал поля ВП; Ё - напряженность поля постоянного тока, протекающего
через проводящую среду с удельной проводимостью с; а - коэффициент поляризуемости
(начальная поляризация) среды; Р - функция спада ВП (виды которой приведены, например, в
работе [7]); т - параметр функции спада.
Приведем кратко математический аппарат, на котором основана предлагаемая процедура 3D-инверсии для восстановления параметров поляризуемости. Параметры ВП будут определяться на основании минимизации функционала
п р / \2 Ф(а,т) = 1(Ь )(и (ь) - и, (Ь)) , (1)
,=1/=1
где и, (t) - сигналы ВП, полученные в , -м приемнике (, = 1...п , п - количество приемников для всех положений приемно-генераторной установки) в виде разности потенциалов на концах приемной линии М^ Ш, (t) - теоретические значения сигнала ВП на тех же приемных линиях; ю, (t) - некоторые веса, t/ - моменты времени, в которые выполнялась регистрация поля
ВП (/ = 1...р, р - количество времен), а а = (а1,...,ак)т и т = (т1,...,тк)т- искомые векторы параметров начальной поляризации и спада ВП, компонентами которых являются значения соответствующих параметров (а и т ) в ячейках 0.к исследуемой области.
Как было показано в работах [5-6], значения Ш, (t) могут быть вычислены в виде
и, (t)=1 а^, (2)
к=1
где Vik - значения сигнала ВП в , -м приемнике при единичной поляризации элементарной подобласти 0.к , на которые разбит изучаемый объем и в каждой из которых ищутся свои значения начальной поляризации ак и параметры спада ВП в виде постоянной времени Тк в функции спада Рк (t) = Р(тк,t) .
Таким образом, функционал (1) нужно минимизировать по параметрам ак и Тк. При этом минимизацию можно выполнять последовательно - сначала по ак, затем по Тк.
Входными данными для инверсии является набор сигналов ВП в определенном наборе приемников для нескольких положений генераторной линии. В качестве исходных данных задаются также параметры вмещающей горизонтально-слоистой среды с заданными значениями удельной проводимостью и параметров вызванной поляризации. В каждой ячейке 0.к задается начальное значение искомого параметра поляризации, совпадающее с его значением для вмещающей среды.
На первом этапе работы алгоритма выполнения инверсии с использованием специальной конечноэлементной процедуры, предложенной в работах [5-6] и реализующей расчет трехмерного поля ВП на основе технологии выделения поля, для всех ячеек и положений приемно-
генераторной установки вычисляются значения . Эта процедура может быть очень эффективно распараллелена (вычисление поля для каждой ячейки и отдельного положения генераторной линии может выполняться независимо), и поэтому даже для очень большого количества ячеек и достаточно большого количества расстановок генераторной линии процедура расчета трехмерного поля, порождаемого одной поляризующейся ячейкой с единичной поляризацией, при наличии соответствующего количества ядер в вычислительной системе может быть выполнена достаточно быстро. При этом процедура вычисления выполняется только один раз, и полученные значения могут быть использованы как для нахождения стартового распределения параметров ак, так и в итерационном процессе поочередного уточнения пара-метр°в {ак } и {хк }.
Стартовое распределение параметра поляризуемости в виде вектора а (с компонентами ак, являющимися значениями параметра поляризуемости в ячейках 0.к) ищется из решения системы линейных алгебраических уравнений (СЛАУ)
Ва = Ь, (3)
где компоненты матрицы В и вектора правой части Ь определяются соотношениями
В^ = 11Ю (Ч ) • {Ц V • Р* ('/ К, (4)
г=1/=1
Ь =1 I ю (ч )Л (ч V • й, (г1 ). (5)
г=1/=1
В качестве весов ю, (t/) будем брать квадраты величин, обратные значениям сигнала ВП в г -м приемнике в момент времени t/ для геоэлектрической модели вмещающей среды. Такой выбор весов позволяет правильно учесть масштаб изменения принимаемого сигнала по времени.
Кроме того, довольно часто (особенно, когда используется ячеистая структура с мелкими ячейками) значения ак могут находиться отрицательными, что не соответствует реальному процессу. Устранить такие значения можно введением в функционал (1) соответствующих параметров регуляризации у к :
" Р \2 К ! \2
Ф(а, т) = Цюг (г1 )(и, (г1)-и, (Ь )) + £у к («к-а0) , (6)
г=1/=1 к=1
где а0 - начальное распределение параметров поляризуемости (в виде значений в ячейках О к).
Получаемая в результате минимизации функционала (6) СЛАУ будет иметь вид
(В + у1 )а = Ь + у/а0 , (7)
где I - единичная матрица, у - вектор параметров регуляризации, а значения компонент матрицы В и вектора правой части Ь вычисляются по формулам (4) и (5).
После нахождения стартового распределения параметров а к запускается итерационный процесс, в котором поочередно ищутся значения %к и новые приближения параметров ак.
Значения хк ищутся с использованием для минимизации функционала (1) линеаризации теоретических значений сигнала ВП по параметрам хк, а новые приближения параметров ак определяются затем по описанному выше алгоритму с измененными в ячейках параметрами хк, полученными на предыдущей итерации. При этом минимизация функционала невязки по значениям хк выполняется по схеме, изложенной в [8], с заменой значений удельной проводимости в ячейках ок на значения в них параметра функции спада хк.
2. РЕЗУЛЬТАТЫ РАБОТЫ ПРОЦЕДУРЫ 3D-ИНВЕРСИИ НА СИНТЕТИЧЕСКИХ ДАННЫХ
Для исследования описанной процедуры 3D-инверсии рассмотрим геоэлектрическую модель для случая, когда поисковые объекты отличны от вмещающей среды только по параметру начальной поляризации а .
Пусть проводящее полупространство с удельным сопротивлением 10 Ом-м и параметрами поляризации = 0.01 (1 %) и х^ = 0.03 с содержит два объекта повышенной поляризации на разных глубинах так, как это показано на рис. 1. При этом «объект 1» с параметром поляризации а1 = 0.2 (20 %) имеет размеры 1 х 2 х 0.1 км3 и расположен на глубине от 400 до
500 м, а «объект 2» с параметром поляризации а2 = 0.1(10 %) имеет размеры 2 х 1 х 0.1 км3 и расположен на глубине от 100 до 200 м. Приемно-генераторная установка представляет собой горизонтальную электрическую линию длиной 500 м и приемные линии длиной 300 м, расположенные по разные стороны от генераторной линии на ее оси и на двух параллельных профилях. Приемно-генераторная установка перемещается по четырем профилям, на каждом из которых выполняется по семь ее расстановок. Эта система наблюдений и первое положение приемно-генераторной установки показано на рис. 1, а.
В
-2000 -1000 0
3000 4000 5000
2000 3000
а б
Рис. 1. План (а) и разрез (б) геоэлектрической модели, содержащей на разной глубине два объекта повышенной поляризации
При выполнении 3D-инверсии сигналов ВП было исследовано три типа разбиения исследуемой области на ячейки. Эти разбиения представлены на рис. 2. Первое разбиение содержало ячейки с размерами 1 х 1 х 0.1 км3, второе - с размерами 0.5 х 0.5 х 0.1 км3, а третье - с размерами 0.25 х 0.25 х 0.1 км3. При этом все разбиения охватывали область [-2; 5] х [-0.5; 3.5] х
х [-0.5; -0.1] км3, т. е. были построены таким образом, что границы обоих объектов геоэлек-
/ / / Г Г 7-1-J-1-1— 1 тщтпш
///JMJ-UU-
т -44-
в
Рис. 2. Разбиение исследуемой области на ячейки с размерами: а - 1 х 1 х 0.1 км; б - 0.5 х 0.5 х 0.1 км; в - 0.25 х 0.25 х 0.1 км
Рис. 3-4 иллюстрируют результаты работы алгоритма поиска распределения параметров поляризуемости (ак) для каждого из разбиений. Числами представлены значения ак • 100 % в соответствующей ячейке, пунктирными линиями показаны реальные границы объектов геоэлектрической модели.
Из этих рисунков видно, что при первом, самом грубом разбиении исследуемой области на ячейки параметры объектов определяются практически точно, а параметры вмещающей среды (т. е. поляризуемость в ячейках, не совпадающих с местоположением объектов) искажаются довольно незначительно. С уменьшением же размера ячеек параметры объектов подбираются все более неоднородными, и появляется все больше ячеек с ложным повышением параметра поляризуемости.
i 0.94 0.92 1.1 0.95 0.0 1 a 0.00 .4 1.4 0.51 1.2 1.7 0.27
1.1 0.99 1.2 0.88 JO 10 1 § 0.68 .2 0.046 1.6~" 0.27 0.086
1
1.1 0.98 0.95 1.1"^"- -053-____ .i...... 0.91 a 0.55 .2 1 0.57" 0.0 1.3
_ J_ __ _ 1.1 0.94" -859-., 1.1 ■ V ■ 0.41 1.4" 1--- 0.75 0.67
,_, ,........ ._, .......
-idoo " Ó юЬо 20b0 зоЬо 40 Ьо 50bt -idoo Ó юЬо 20 bo зоЬо 40 bo sobo
AIPha0.001 0.5 i 2 4 g*** а 8 10 15 20 Ib JJ Alpha0 DQ1 0 5 2 4 8 10 15 6 JJ -33- 30 40 E
1.3 0.15 0.69 2.1 0.59 0.44 2.1 зоро 0.79 .6 0.74 0.23 1.2 0.79 1.1
1.6 0.18 0.46 0.56 0.59 0.58 a 0.6 .0 1 -l... r 23 1 0.4"" 1:5-- 2.2 2.1
1.9 0.15 0.84 1.6 " 0.73 0.47 _ 0.42 .9 1 i ____L.. 0.48" 1.6 1.4
ют
£L84.. i 2.3 0.057"'"" -05Э-., 1.7 1.7 ■ V ■ ,2 0.059 1.7"" 856- 0.42 0.41
-idoo ■ Ó юЬо 20b0 зоЬо 40 Ьо 50b[ -idoo Ó юЬо 20 bo зоЬо 40 bo sobo
AlphaQ Q01 Q 5 1 2 4 8 10 15 20 25 30 Jfl A|pha0.001 0.5 2 4 t В 10 15 U 25 эй 4Й e
в
г
Рис. 3. Распределение параметра начальной поляризации по ячейкам в четырех слоях для первого разбиения: а - слой 100-200 м; б - 200-300 м; в - 300-400 м; г - 400-500 м
□ 8 0.5 1.4 □ 9 0.8 1 □ 1 1 2 1 16 1 9 2.7 09 1.3 1 8 0.2 1.2 23 04 12 02 0.2 0.9 2 0.8 1 1 0.81
□ 9 1.2 1.1 1 11 оя , 07 1 1 0.9 06 07 0й 0 25 1.2 09 0 3 1 2 04 16 08 2 0 8 1 7 06 08 04 3.2
1.1 0.7 0.0 11 4.2" "02- _1_._1 1.2 1 1 12 12 12 ц 1 0 9 2.5 0.7 0"б" 11Э--2Л. - 0.1 0.5 0.0 0.7 0.5 0.23
1 1.3 1.1 0.9 1.3 0 9 0.7 1.3 1.2 9.4 "6.2- 101 N 0.0 0.2 0.3 0.9 0.3 0.2 2 0.6 0.2 0.3 -2,1. 0.95
1 1
1 3 0.8 0.8 09 1 (№ и.9 _1_ 3 0 9 0.7 03 03 0.4 0 52 0.6 04 2.5 07 04 П-&--1-_ -Л4 0.8 2.1 о / 0.5 ОЬ 2.4
1 1 1.3 1.1 1 1 0.9 09 08 06 11 0.9 0" 9 "0 7- - 0-9 0 86 0.9 00 04 05 04 03 18 22 0 7 0 9 19 -1.5 1.5
1 1 0.9 0.9 12 1""" "1-1- 1.1 .1 4 1 1.1 07 07 0.8 1 3 1.2 06 1.7 0.3 1""" ' 4 - 11' -.0 0 0.8 0.8 07 0.9 07 0.18
П " Г71П 09 09 09 09 1 0" 9 "Пй- - 0-9 1 . ^ 1*7 -0.7__ ал. 12 1 1 1 1 1 ? 0 9 Тз "15- -1-4 095
-1(1)0 ь 1 иЬ 2О0О зоОо 4000 50Ь0 А|р|1а -11)0 0 й юЬо 2000 3000 400 5060
^Р^ОООЮй 1 2 10 15 И тг~? 001 0.5 1 2 4 6 8 10 \ь 20 2^ 5Г~ ЖП
а б
0.6 1.2 1.5 0.9 0.6 0.6 1.3 0.4 3.3 0.8 0.6 0.4 0.9 0.95 0.7 1.6 1.1 0.2 0.4 0.8 1 0.6 1.4 0.7 0.6 0.4 0.8 0.87
0.9 0.6 1.4 0.7 0.В 0.6 0.7 1.2 0.0 0.0 0.9 0.7 0.5 0.54 У 0.7 0.3 1.2 0.6 0.7 3.4 0.8 1.4 1.4 0.2 0.3 0.7 0.0 0.35
0 0.7 0.9 0.0 0.7 -0,-7, ■0.0, .0.4 0.4 0.5 0.2 0.6 0.2 0.39 0 1 3.2 0.6 О.Ё! 24 " 7:41 - -0Л, ^0.6 4.8 0.6 0.1 0.7 0.2 0.21
к 1.4 0.6 1.7 0.5 0.5 0.7 0.5 0.5 0.4 0.3 0.7 "0.5- • 0,7. 0.29 н 1.2 0.4 0.7 2 71 20 27 0.2 0.' 4.6 0.5 0.7 '05- 0.75
4.
О 0.8 0.7 0.5 0.6 5.1 "0:7-- ■0.7-, .0.1 0.7 1.2 0.2 0.6 0.7 0.023 о 1.8 2.5 0.4 2.а 1 14 22 - йб, ^0.2 1.7 0.5 0.3 0.7 0.1 0.73
0.5 0.5 1.7 0.6 4.8 0.7 0.2 0.7 0.4 1.2 0.7 "Ч.З-ОЛ. 2.8 \2__ 2.1 .16.. . .0^8. 0.7 ..17. J(.. Щ ..1.3
0.5 0.1 0.8 0.3 {.Г "от- €Л_ .0 3 0.6 0.6 0.8 1.1 0.9 0.46 1.4 3 0.9 2.7 0.9" - -2,5. .0.8 2.8 0.8 0.9 1.1 1.1 1.8
. 1-4. 0-7 1.4 1 ? 1.1 0.4 09 п"п "0.6- - од 1.2 П.Соя. -53. 07 П8 0.4 1.3 1 ? 1.4 "о". 9 "1--4- 1.2
-1000 А1рИа0 001 0 5 6 того 2000 300 0 4000 5000 -1с1оо Й 1000 2000 зоЬо 401 0 50Ь0
2 У 10 1Ь Л) 25 30 4У к А1рЬааш] 0 5 2 4 Б 8 10 15 20 Ж 411 к
в г
Рис. 4. Распределение параметра начальной поляризации по ячейкам в четырех слоях для второго разбиения: а - слой 100-200 м; б - 200-300 м; в - 300-400 м; г - 400-500 м
Так, при втором разбиении (рис. 4) в ячейках, попадающих в первый объект, значения поляризуемости были получены в интервале от 7.4 % до 27 %, а в ячейках, соответствующих положению второго объекта, значения поляризуемости были получены в интервале от 6.2 % до 12 %. Кроме того, во втором и третьем слоях (рис. 4, б, в), где в исходной модели вообще отсутствуют аномально поляризующиеся объекты, были получены ложные повышенные значения поляризуемости, которые на практике могут быть проинтерпретированы как поисковые объекты (хотя на самом деле их нет). При третьем разбиении картина получается еще более «пестрой», поисковые объекты - гораздо более «размытыми».
Таким образом, использование более подробных разбиений в данном случае не только не способствовало улучшению подбора геометрии поисковых объектов, но и даже наоборот, несколько ухудшило результат инверсии.
Немного улучшить качество инверсии с использованием сеток с мелкими ячейками можно, существенно увеличив плотность съемки. Результаты подбора для увеличенной вдвое плотности измерений для второго разбиения приведены на рис. 5. При этом, если для второго разбиения улучшение подбора параметра начальной поляризуемости за счет увеличения плотности съемки было довольно значимым, то для третьего (самого мелкого) разбиения это улучшение уже стало гораздо менее заметным.
Кроме того, при проведении практических работ увеличение плотности съемки, во-первых, не всегда возможно, а во-вторых, приводит к существенному увеличению затрат на их проведение. При этом уже достаточно хороший результат можно получить и при невысокой плотности съемки, если использовать крупные ячейки, - во всяком случае можно с достаточно высокой точностью определить местоположение основных аномалий и не получить ложных.
Вычислительные затраты при использовании первого разбиения (на крупные ячейки) составили 6 мин. (112 ячеек), при использовании второго разбиения - 32 мин (448 ячеек), при использовании третьего разбиения - порядка 4.5 часов (1792 ячейки). Расчеты выполнялись на 28 ядрах вычислительного кластера с процессорами 3 ГГц Intel Xeon.
а б
Рис. 5. Распределение параметра начальной поляризации по ячейкам в слоях с объектами для второго разбиения при более плотной сети наблюдений: а - слой 100-200 м; б - 400-500 м
Таким образом, использование разбиений с относительно крупными ячейками позволило получить более точные результаты при существенно меньших вычислительных затратах.
Для более точного подбора границ объектов с аномальной поляризуемостью можно использовать подход, предложенный в работе [8] для решения задачи восстановления проводимости среды. В этом подходе после первого этапа, где с помощью ячеистой структуры ищется стартовое распределение параметров, выполняется второй этап, на котором уточняются параметры локальных неоднородностей, сформированных из пробных объектов с близкими значениями искомых параметров.
Рассмотрим теперь ситуацию, когда у объектов исходной модели аномальными являются не только значения параметра начальной поляризации, но и значения параметра функции спада. У «объекта 1» значение параметра функции спада т зададим равным 0.5 с, а у «объекта 2» - равным 0.15 с. При этом фоновое значение (т. е. значение во вмещающей среде) параметра функции спада оставим, как и прежде, равным 0.03 с.
На рис. 6-7 приведены результаты подбора параметров начальной поляризации а и функции спада Р с помощью предложенного в данной работе алгоритма инверсии при использовании первого разбиения (с самыми крупными ячейками).
Для ускорения этой процедуры можно при подборе распределения параметра а использовать только часть измеренной в каждом приемнике кривой, которая определяется, например, процентом изменения принимаемого сигнала по времени. При подборе же распределения т (при фиксированном распределении а ) можно использовать всю регистрируемую в приемнике кривую. Так, в рассмотренном нами случае при подборе распределения параметра а использовался временной интервал до 10 мс, который определялся 15-20 % уровнем изменения кривой во времени. В этом случае, как правило, достаточно одной-двух итераций при последовательном подборе распределений а и т .
Отметим, что если при подборе обоих параметров использовать измеренные значения во всем временном интервале, то это не только приводит к многократным итерациям последовательного подбора по а и по т , но и может заметно ухудшить качество конечного результата по сравнению с ситуацией, когда для подбора распределения параметра а используются данные лишь на относительно ранних временах, до начала существенного спада ВП.
Вычислительные затраты на подбор распределений обоих параметров (с ограничением времен при подборе а ), включая решение прямой задачи на 28 ядрах вычислительного кластера, составили порядка 7 мин. При использовании же полной кривой при подборе а потребовалось порядка 19 итераций по а и т , что привело к увеличению вычислительных затрат до 13 мин., и результат при этом получился несколько хуже.
0.99 0.92 0.87 1 0.9!
0.96 1.1 0.9
-1ЙЮ " б юЬо 2оЬо " зоЬо ' доЬо ' 5оЬо 5 А|РИа1е-00Я.5 1 2 -5-ТП-15 2Й 25 30 4А М 1е-Ю02
1.9 0.74
0.095 1.9 Э.-67---- 0.024 0.072
0.28 1.4
0Л2 ._. 1.4
0.15 1.3
1.-2--... 0.57 0.49
ЛШо ' Ь юЬо 2оЬо ' зоЬо ' 4оЬо ' 5оЬо > 1е-005).5 Г 2 4 6 3 10 15 20 2Ь Зи е-Ю02
0.83 0.097 0.19 1.4 0.77 0.51 0.7
2 0.33 0.48 0.2 " " ~ -0.-93- . 0.16 0.066
2.5 0.0006 0.24 0.36 "" " ~ ~ ""6.73- - 23 0.53
■ °,78 ■ 0.32 2 2.2
-1000 0 1000 2000 3000 4000 5000 X
Л|Рпа1е-00Ж1.5 1 2 4 40 60 1 е+002
"1
0.37 2.3
----- 0.39 1.4
.9--.. 0.23 0.016
Ш ' 3 ' юЬо ' 2000 " зо'оо " 40в0 ' 5оЬо X А^^в-оовл 1 2 —5—та—ЗА 25 ЗА 55 1е-юо2
б
а
в
г
Рис. 6. Распределение параметра начальной поляризации по ячейкам в четырех слоях для первого
разбиения:
а - слой 100-200 м; б - 200-300 м; в - 300-400 м; г - 400-500 м
0.031 0 03 0 03 0.020 0.029 0 021 0 020
0.032 0.031 0.017 0.043
0.032 0 032 О 029 0.041
0.14 0.14
0 026 0 031
0 031 0 029
О 0 01 0 03 0 05
Ъ ' Ш ' 2Ш0 ' зШо ' До ' 5оЬо 01502 шетттиг!
0.024 0.03 0.033 0.057 0.031 0.039 0.009
0.025 0.0065 в_В 0^3, 0.013 0.003 0.03 " ""-"Э.-026,____0.001 0.0014
0.014 0.0044 0.02" ___0.013 0.02
■ ■ 0 038 0.023 0.-О30_____ 0.028 0.042
■ ■ ■ ■ ■ ■■ ■
О 0 01 0 03 0.05
007 0.1 0.15 02 Я5П
оЛ, о!,
0 043 0.011 0 001 0 012 0.027 0 011 0 0059
0.0012 0 0011 0 017 ""--9.047-.__ 0 0011 0 001
0 037 0.0024 0 0015 0 036 ---0.O26-.__ 0 023 0 021
0 02 0 01 -"--0.043.__ 0 032 0 025
-|Ло ' б ' ш5о ' 2оЬо ' зоЬо ' доЬо ' 5060 X
Та" 0 0.01 0 03 0 05 0 07 0.1 0 15 0.2 О 'А 6.3 0 3^ 0.4 0 4 6.^ 1
0 625 0 03 0 622 0.025 0 637 0.0017
0030 0003, 041 | II' - 1" 1_. _ 0.03 00
О 022 0 043
0 44 1 0.0052 -" -6044- _ _ _ 0.09 0 067
_ Д.017.___ 0 021 0.079
0.018 0 0078
-1Л0 ' 3 ' ШБо ' 2оЬо ' зоЬо ' доЬо ' 5оЬо х
Тгш 0 0 01 0 63 0 05 0 67 0 1 0.15 0 2 0.25 0 3 0 35 6 4 0 45 6 5 1
в г
Рис. 7. Распределение параметра функции спада по ячейкам в четырех слоях для первого разбиения: а - слой 100-200 м; б - 20-300 м; в - 300-400 м; г - 400-500 м
Заметим также, что дробление ячеек при подборе т , как и при подборе а , не улучшает, а даже несколько ухудшает ситуацию: картина становится более «пестрой», параметры объектов становятся более неоднородными и отклоняющимися от истинных, а во вмещающей среде появляются существенные ложные аномалии.
ЗАКЛЮЧЕНИЕ
Разработаны и опробованы алгоритмы автоматической 3D-инверсии данных вызванной поляризации, основанные на конечноэлементном решении прямой трехмерной задачи. В отличие от существующих аналогов предлагаемые алгоритмы 3D-инверсии основаны на однократном решении трехмерной эллиптической задачи, что позволяет существенно повысить их вычислительную эффективность.
В результате проведенных исследований установлено, что, несмотря на использование различных регуляризаций, дробление ячеистой структуры может не только не улучшать конечный результат инверсии, но и вносить значимые искажения в распределение искомых параметров ВП по сравнению с использованием структур с более крупными ячейками. Поэтому повышения точности подбора параметров поляризуемости имеет смысл добиваться не дроблением ячеистой структуры, а выполнением второго этапа инверсии, на котором будут уточняться параметры поляризуемости и искаться размеры крупных объектов, стартовое положение и параметры которых формируются из пробных объектов (ячеек) с близкими аномальными значениями параметров поляризуемости.
Использование же для проведения инверсий ячеистых структур с мелкими ячейками требует, как минимум, более плотных сетей наблюдений, что, во-первых, приводит к удорожанию полевых работ (а также к резкому увеличению вычислительных затрат при выполнении самих инверсий), а во-вторых, может дать улучшение при подборе геометрии поисковых объектов только до определенных, относительно небольших глубин. При этом разработанный про-
б
а
граммно-математический аппарат в зависимости от требований к глубинности предстоящих геофизических исследований методом ВП позволяет связать необходимую плотность системы наблюдений с уровнем подробности ячеистой структуры, которая будет использована при проведении инверсий практических данных.
СПИСОК ЛИТЕРАТУРЫ
[1] M. Karaoulis, A. Revil, D. D. Werkema, B. J. Minsley, W. F. Woodruff and A. Kemna. Time-lapse three-dimensional inversion of complex conductivity data using an active time constrained (ATC) approach. Geophys. J. Int. (2011) 187, 237-251.
[2] Michael Zhdanov. Generalized effective-medium theory of induced polarization. GEOPHYSICS,VOL. 73, NO. 5 (SEPTEMBER-OCTOBER 2008). - P. F197-F211.
[3] Michael Commer, Gregory A. Newman, Kenneth H. Williams, and Susan S. Hubbard. Three-dimensional inversion of EM coupling contaminated spectral induced polarization data. SEG Denver 2010 Annual Meeting. - P. 624-628.
[4] Могилатов В.С. Индуктивный, смешанный и гальванический источники в электроразведке становлением поля // Изв. РАН. Сер. Физика Земли. - 1997.- № 12. - С. 42-51.
[5] Моисеев В.С., Рояк М.Э., Соловейчик Ю.Г. Математическое моделирование процессов вызванной поляризации в сложно построенных средах для токовой линии с заземленными электродами // Сибирский журнал индустриальной математики. - 1999. - Т. 2. - № 1. - С. 79-94.
[6] Персова М.Г., Соловейчик Ю.Г., Тригубович Г.М.Компьютерное моделирование геоэлектромагнитных полей в трехмерных средах методом конечных элементов // Физика Земли - 2011. - Т. 47. - № 2. - С. 3-14.
[7] Моисеев В.С. Метод вызванной поляризации при поисках нефтеперспективных площадей. - Новосибирск: Наука, 2002. - 136 с.
[8] Персова М.Г., Соловейчик Ю.Г., Вагин Д.В., Токарева М.Г., Чернышев А.В. О подходе к выполнению 3D-инверсий данных зондирований нестационарным электромагнитным полем на основе конечноэлементного моделирования // Научный вестник НГТУ. - 2011. - № 2(43). - С. 97-106.
Персова Марина Геннадьевна, доктор технических наук, доцент, профессор кафедры прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование электромагнитных полей в задачах геоэлектрики и электромеханики. Имеет более 90 публикаций, в том числе 2 монографии. E-mail: persova@fpm.ami.nstu.ru.
Соловейчик Юрий Григорьевич, доктор технических наук, профессор, заведующий кафедрой прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование электромагнитных и тепловых полей. Имеет более 100 публикаций, в том числе 2 монографии. E-mail: kpmt@fpm.ami.nstu.ru.
Токарева Марина Георгиевна, кандидат технических наук, доцент, доцент кафедры прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование электромагнитных полей в задачах геоэлектрики. Имеет 23 публикации. E-mail: tokarevamg@mail.ru.
Симон Евгения Игоревна, аспирант кафедры прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - методы инверсии в задачах геоэлектрики. Имеет 8 публикаций. E-mail: simona_@mail.ru.
Абрамов Михаил Владимирович, кандидат технических наук, доцент кафедры прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование электромагнитных полей в задачах геоэлектрики. Имеет 15 публикаций. Email: mike_abramov@inbox.ru.
Сафронова Ульяна Анатольевна, магистрант второго года факультета прикладной математики Новосибирского государственного технического университета. Основное направление научных исследований - конечно-элементное моделирование электромагнитных полей в задачах геоэлектрики. E-mail: safronova.u.a@mail.ru.
M.G. Persova, Yu.G. Soloveichik, M.G. Tokareva, Eu.I. Simon, M.V. Abramov, U.A. Safronova
3D-Inversion of Induced Polarization Data by means of Finite Element Modeling
In this paper we suggest an approach to the execution of 3D-inversions of data obtained by induced polarization survey with the help of the finite element method for solving a direct problem. The mathematical apparatus for solving a direct problem and examples of computational experiments are described.
Key words: 3D-inversions of induced polarization data, 3D-modeling, finite element method.
