Построение разрезов методом аналитического продолжения геофизических полей

Размещено на http://www.allbest.ru/

Построение разрезов методом аналитического продолжения геофизических полей

Ермохин К.М.

Методам аналитического продолжения геофизических полей с профиля измерений в нижнее полупространство посвящено много работ ведущих геофизиков: Б.А. Андреева, И.Г. Клушина, В.Н. Страхова. Явление распадения поля в окрестности особых точек констатированное в [4], оказалось, скорее, не информативным признаком, а препятствием для практического применения метода (см. рис. 2).

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

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

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

Единственной относительно удачной попыткой обойти сложности является метод полного градиента В.М. Березкина. Однако его трудно признать удовлетворительным ввиду невозможности определения необходимого числа членов рядов Фурье при суммировании [2]. Используемые им ряды в нижнем полупространстве, как нетрудно заметить, расходятся. Это приводит к прямой зависимости результата от количества использованных членов ряда, т.е. к фактической неопределенности.

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

Частным случаем задачи аналитического продолжения является метод особых точек [6], в котором объектом поиска являются нули знаменателя, а о числителе речь не идет вообще, что, очевидно, не достаточно для представления об общей структуре поля в нижнем полупространстве.

Предлагаемый нами метод основан на аппроксимации исследуемых функций цепными дробями. Мы назвали его - CFCM (Continued Fraction Continuation Method).

Основания для применения аналитического продолжения к полям, применяемым в геофизике

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

Аналитическим продолжением (АП) функции, заданной на отрезке, называется:

1. аналитическая функция, совпадающая с заданной на этом отрезке.

2. аналитическая функция, вещественная часть которой совпадает с заданной на отрезке.

Мы в этой работе будем пользоваться вторым определением.

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

Аналитически продолжить с отрезка задания (профиля измерения) можно любую непрерывно дифференцируемую функцию. Это означает, что решение задачи АП всегда существует [4].

АП заданной функции единственно [4].

В классе мероморфных функций АП асимптотически устойчиво [1]. (мероморфная функция - функция, у которой все полюсы имеют конечную кратность, упрощенно можно сказать, что это рациональная дробь с бесконечным числителем и знаменателем).

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

Решение задачи

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

I. Рассмотрим измеренную по профилю функцию , и представим ее рядом Фурье по многочленам Чебышёва первого рода (предварительно приведя аргумент к отрезку [-1,1]:

,

где - многочлены Чебышёва первого рода, - коэффициенты.

Поскольку функция реально может быть измерена только на конечном отрезке, представление ее рядом, а не интегралом, является вынужденным, но вполне естественным. Аппроксимация синус - рядом Фурье, использованная в методе Березкина, с общей точки зрения является неудачной. Такая аппроксимация подразумевает периодическое продолжение измеренной функции на всю вещественную ось, что, очевидно, не соответствует действительности и приводит к непредсказуемым искажениям (особенно на краях профиля). Кроме того, выравнивание значений на концах путем вычитания линейного тренда без последующего восстановления, вносит значительные искажения в спектр функции.

Коэффициенты ряда Фурье - Чебышёва можно вычислить по формулам:

, (1)

Для упрощения вычислений и дальнейших преобразований, используем связь многочленов Чебышёва с тригонометрическими функциями. Сделаем замену переменной в интеграле: . Тогда [5]:

,

Коэффициенты Чебышёва есть коэффициенты Фурье функции F()=f(cos()) при четном продолжении f() на отрезок [-,0]. Эти коэффициенты вычисляются прямым интегрированием без применения специальных формул посредством интерполяции функции F() периодическим кубическим сплайном. Такая аппроксимация позволяет получить достаточно большое количество коэффициентов (точнее, 2N-1 коэффициент, где N- число точек измерения), которые убывают по величине быстрее, чем . Аппроксимация сплайном является более естественной, чем применение дискретного преобразования Фурье (ДПФ), т.к. аналитичность исследуемой функции подразумевает ее бесконечную гладкость, а при ДПФ происходит ее замена на кусочно постоянную, что изначально создает логическое противоречие.

Здесь интересно отметить, что сравнение рядов Фурье и Чебышёва позволяет понять причину эффекта Гиббса (сильные колебания частичных сумм ряда Фурье на краях отрезка). Как видно из (1), концевые точки функции входят в интеграл с бесконечным весом, что при численном интегрировании и порождает неустойчивость.

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

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

Представление функции в виде ряда Чебышёва имеет целью восстановить мешающий вычислениям линейный тренд без ущерба для скорости убывания коэффициентов и без формального периодического продолжения на всю вещественную ось, а сводится, фактически, к нелинейному перемасштабированию комплексной плоскости (чего нет в методе Березкина).

Прямая подстановка в ряд Чебышёва комплексной переменной вместо вещественной с целью аналитического продолжения с отрезка невозможна в силу того, что ряд этот расходится в нижней полуплоскости (y<0) комплексной плоскости, за исключением отрезка , со скоростью геометрической прогрессии со знаменателем большим единицы. Однако если этот ряд представить в виде степенного, к нему можно будет применить представление цепной дробью и в таком виде просуммировать.

Формально можно привести ряд Чебышёва к виду степенного собрав одинаковые степени в выражении , но эта процедура некорректна ввиду неустойчивости, порожденной быстрым ростом коэффициентов в выражениях для .

II. С целью приведения ряда Чебышёва к виду степенного ряда, перепишем ряд функции в комплексной форме на основе определения полиномов Чебышёва [5]:

Обозначим (знак минус перед корнем дает продолжение вниз (y<0), плюс - вверх (y>0)) и заменим вещественную переменную комплексной

,

где y означает глубину/высоту. Соответствие взаимно-однозначное.

Отрезок [-1,1] перейдет при этом в нижнюю единичную полуокружность на комплексной плоскости Y т.к. x изменяется от -1 до 1, при изменяющемся от 0 до - и . При имеем .

Ряд Чебышёва для тогда перепишется в виде:

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

(2),

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

III. Для продолжения вниз заменим ряд соответствующей цепной дробью, а точнее С-дробью общего вида:

где - целые положительные числа, - вещественные коэффициенты.

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

Существует много способов представления функций, заданных формальным степенным рядом, цепными дробями. Наиболее эффективные алгоритмы разработаны для правильных С-дробей, J-дробей, g-дробей, S-дробей Стилтьеса, Р-дробей Магнуса и Т-дробей Трона [1,2]. Области сходимости ряда и соответствующей ему цепной дроби могут принципиально различаться. В частности, таким образом можно вычислять суммы расходящихся формальных степенных рядов, порожденных изучаемой функцией. Как показали наши исследования, наиболее подходящей математической конструкцией для продолжения геофизических полей оказалась С-дробь общего вида.

Если последовательно брать в качестве приближения подходящие дроби (т.е. учитывать все большее количество членов), то формула (3) эквивалентна представлению F(Y) в виде дробно-рациональной функции посредством аппроксимаций Падэ:

где P_L(Y) и Q_M(Y) - полиномы степени L и M от Y.

Близость аппроксимаций Падэ, каковыми и являются цепные дроби, к представляемой ими функции при L+M описывается теоремой Бейкера [1]:

“Пусть последовательность аппроксимаций Падэ функции , аналитической в начале координат, и последовательность равномерно ограничена в односвязной области D, содержащей начало координат. Если , то равномерно в каждой области R, компактно принадлежащей D, и функция аналитически продолжается в R ”.

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

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

Другое, кардинально более эффективное, решение было предложено в 1960г. Рутисхаузером [1]. Коэффициенты цепной дроби можно последовательно вычислить на основе построения QD-таблицы по правилам ромба. Первый столбец этой таблицы представляет собой последовательность отношений Это может воспрепятствовать её построению при нулевых или попарно несоизмеримых по величине коэффициентах ряда Чебышёва, что при практических вычислениях не редкость, а правило.

Для построения С-дроби, соответствующей ряду (2) мы применили предельно эффективный, надежный и устойчивый метод В.И. Висковатова [1,7], предложенный им в 1805г. Он позволяет по заданным коэффициентам формального степенного ряда построить последовательность коэффициентов соответствующей непрерывной дроби . Приведем этот алгоритм в авторском варианте с небольшими упрощающими изменениями и заменой архаичных обозначений.

Пусть исходная функция задана степенным рядом по переменной , тогда:

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

Таким образом мы получаем С - дробь общего вида:

где все , - целые.

IV. Получаемая посредством этого алгоритма дробь (4) сходится, но не монотонно и достаточно медленно. Для ускорения сходимости применим следующий прием: для заданного Y вычисляются значения последовательных подходящих дробей C-дроби _, по рекуррентным формулам прямого вычисления числителей (A_i) и знаменателей (B_i) [3]:

Вычисленные таким образом комплексные числапроецируются на сферу Римана с радиусом 1 путем соединения соответствующих им точек на комплексной плоскости с ее северным полюсом. Точки пересечения соответствующих отрезков со сферой отождествляются с точками комплексной плоскости. При этом северный полюс соответствует бесконечности, а южный - нулю. Сходимость понимается в смысле метрики на сфере Римана. Это означает, что в качестве предельной точки берется средняя . По завершении процесса точка Р обратно проецируется на комплексную плоскость и ее проекция F принимается в качестве предельной для . Описанный процесс поясняет рис. 1.

Размещено на http://www.allbest.ru/

Рис. 1. Схема соответствия точек комплексной плоскости и сферы

Примечания к решению

Информация об аналитичности продолжаемой функции является не дополнительной, а главной информацией о множестве N ее измеренных значений. Аппроксимация Падэ это, по существу, экстраполяция степенного ряда за интервал задания коэффициентов. Практика расчетов (не только наших, но и авторов монографии [1]) показала, что информация о глубинных неоднородностях содержится в дальних знаках коэффициентов Фурье-Чебышёва и недоступна без определения высокочастотных коэффициентов, поэтому первостепенной задачей является предельно точное их вычисление.

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

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

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

Примеры применения метода CFC

Перейдем теперь к примерам программной реализации описанного выше метода.

Рис. 2. Аналитическое продолжение гравитационного поля точечной массы в нижнее полупространство (здесь и далее в модельных примерах единицы измерения поля - условные)

На рис. 2 слева представлен результат аналитического продолжения поля “аппроксимационным методом” В.Н. Страхова [4], из которого видно, что при приближении к источнику поля возникает интенсивный колебательный процесс, приводящий к “эффекту распадения поля”. Это означает, что “аппроксимационный метод” позволяет определить (и то - весьма приближенно) лишь верхнюю кромку аномалиеобразующего объекта. Заметим также что поле , в отличие от , практически не измеримо, поскольку реально ось гравиметра может быть установлена только по полю силы тяжести, а направление оси Z - лишь математическая фикция. Выбор объясняется тем, что удовлетворяет уравнению Лапласа, а - нет и, соответственно, к нему “аппроксимационный метод” вообще не применим.

Справа показан результат аналитического продолжения гравитационного поля методом CFC. Индекс “ext” подчеркивает, что аналитически продолжается внешнее поле, реально измеренное на поверхности. Это означает, что внутренняя структура поля (внутри объекта) остается недоступной и может интерпретироваться только на основании внешних проявлений.

Графики поля на глубинах 400, 800, 1200 и 1500м (источник на глубине 1000м) свидетельствуют об отсутствии “эффекта распадения поля”. В самой правой части рисунка показан график поля по вертикальной оси. Выше источника функция положительная, ниже - отрицательная, что абсолютно соответствует поведению поля при измерении вдоль вертикальной оси. Характер изменения продолженного поля с глубиной отражает процесс измерения на различных глубинах.

Внизу - разрезы в изолиниях собственно продолженного поля и, ниже его, модуля . - это вещественная часть комплексной аналитической функции, которая совпадает на поверхности с измеренной (в данном случае - рассчитанной). Модуль функции можно интерпретировать как плотность энергии гравитационного поля в разрезе.

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

В примерах расчетов, приводимых далее, показан модуль функции |F(Y)|, достаточно хорошо отражающий распределение неоднородностей в разрезе в большинстве случаев. Физический смысл этой функции - плотность энергии продолженного поля. В сложных для интерпретации случаях привлекается информация о вещественной и мнимой частях комплексной функции F(Y).

Рис. 3. Аналитическое продолжение гравитационного поля от двух шаров, находящихся один под другим

На рисунке 3 представлен график гравитационного поля g по центральному профилю от двух шаров (3D модель) с радиусами 50 и 100м на глубине 200 и 400м с избыточной плотностью 1 (белые контуры), находящихся один под другим при шаге измерения 10м. Очевидно, что визуально предположить по графику наличие двух объектов невозможно. Внизу показано распределение поля в разрезе.

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

Рис. 4. Аналитическое продолжение гравитационного поля от трех шаров на глубине 200, 300 и 400м радиусом 100, 150 и 200м одинаковой избыточной плотности при отсутствии помехи (вверху) и наличии помехи (средний уровень 15%) в виде множества мелких шаров в верней части разреза на глубине 50м (внизу)

Рис. 4 не только подтверждает теоретический вывод [1] об асимптотической устойчивости метода к помехам, но и дает основания к более сильному утверждению о глобальной устойчивости к помехам геологического типа.

Рис. 5. Аналитическое продолжение магнитного поля изогнутого пласта при вертикальном (вверху) и горизонтальном (внизу) намагничивании

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

Рис. 6. Аналитическое продолжение поля вызванной поляризации (кажущейся поляризуемости в %) при реальных измерениях по профилю. Красный пунктир - разведочная скважина вскрывшая рудное тело в интервале 176-246м). Синий пунктир - пустые скважины, пробуренные до глубины 200м под локальные максимумы

Рис. 6 демонстрирует возможности метода CFC в условиях реальных измерений при низком уровне аномалий и высокой геологической зашумленности поля. (Технические помехи, как показали контрольные полевые измерения, практически отсутствовали).

Рис. 7. Комплексный анализ результатов аналитического продолжения поля концентрации пропана, магнитного и гравитационного поля над вскрытым нефтегазоконденсатным месторождением (контур месторождения - черный пунктир)

Рис. 7 дает возможность сопоставить результаты аналитического продолжения полей различной природы с целью выявления перспективных областей для разведочного бурения. В правой части разрезов (пикеты 40 000-50 000) отмечаются ярко выраженные особенности трех полей, что свидетельствует о наличии аномалиеобразующего объекта. На пикетах 15 000 и 25 000 объекты выделяются только по полю концентрации пропана, в магнитном и гравитационном полях особенностей не наблюдается, что говорит о низкой вероятности нахождения там залежей.

Заключение

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

Область применения метода не ограничивается потенциальными полями. Метод опробован на задачах электроразведки (поле сопротивлений-_к и вызванной поляризации-_к), магниторазведки (поле Т), гравиразведки (поле g), полей концентрации газов и может быть распространен на другие геофизические методы. Развитие метода может состоять в его модификации для сейсмических полей, МПП, МТЗ, ЧЗ и т.п., разработке 3D - версии алгоритма, учете рельефа местности.

Преимущества метода:

Не требует “априорной” информации о направлении и величине вектора намагниченности, строении окружающего пространства, положении источников поля (естественного или искусственного) и др.

Устойчив к геологическим помехам.

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

Объективно не зависит от выбора модели среды.

Универсален.

Недостатки:

Разработанная на сей день 2D - реализация алгоритма требует “разумного” выбора интерпретационных профилей во избежание геологического истолкования “теней” от близлежащих объектов по сторонам от изучаемого профиля, что не позволяет полностью автоматизировать процесс. Простое “склеивание” 2D разрезов по параллельным профилям измерений не дает адекватной трехмерной картины поля в нижнем полупространстве.

Большой пробел в области инженерного математического образования в ХХ веке - забвение теории цепных дробей. Этот древний раздел математики ведет свое начало со времен Евклида, достиг расцвета в трудах Эйлера, Стилтьеса, Чебышёва, Адамара, Маркова, но ныне практически забыт. Однако, успехи геофизической теории, достигнутые на основе линейных математических методов, почти исчерпаны, а практика требует двигаться дальше. Предложенный подход представляется шагом в новом направлении.

Литература

геофизический поле степенной дробь

1. Дж. Бейкер, мл, П. Грейвс-Моррис. Аппроксимации Падэ. М. Мир, 1986.

2. В.М. Березкин, М.А. Киричек, А.А. Кунарев. Применение геофизических методов разведки для прямых поисков месторождений нефти и газа. М. Недра, 1978.

3. У. Джоунс, В. Трон. Непрерывные дроби. М. Мир, 1985.

4. В.Н. Страхов. Аналитическое продолжение потенциальных полей. Гравиразведка: Справочник геофизика / Под ред. Е.А. Мудрецовой, К.Е. Веселова. М. Недра, 1990.

5. П.К. Суетин. Классические ортогональные многочлены. М. Наука, 1979.

6. Г.А. Трошков, А.А. Грознова. Математические методы интерпретации магнитных аномалий. М. Недра, 1985.

7. Viskovatov B. De la methode generale pour reduire toutes sortes des quantitees en fraction continues. Memoires de l' Academie Imperiale des Sciences de St. Petersburg, 1, 1805.

8. K.M. Ermokhine. Analytical continuation of geophysical fields into the area of anomaly sources by the Continued fraction method (CFCM). Vienne, EAGE-2006, abstr. P324.

9. Ермохин К.М. Аналитическое продолжение геофизических полей в область источников аномалий с помощью цепных дробей. Материалы 34 семинара им. Д.Г. Успенского “Вопросы теории и практики интерпретации гравитационных, магнитных и электрических полей”. М. ИФЗ РАН, 2007.

10. Жданова Л.А., Ермохин К.М. Результаты применения метода продолжения полей непрерывной дробью (CFCM) при геофизических исследованиях в Карело-Кольском регионе. Материалы 34 семинара им. Д.Г. Успенского “Вопросы теории и практики интерпретации гравитационных, магнитных и электрических полей”. М. ИФЗ РАН, 2007.

11. Ермохин К.М. Аналитическое продолжение геофизических полей в область источников аномалий методом аппроксимации цепными дробями. Геофизика (ЕАГО), 1.2007.

12. Ермохин К.М., Жданова Л.А. Аналитическое продолжение геофизических полей, как отражение глубинного геологического строения. РАН, Петрозаводск, материалы XIV международной конференции “Связь поверхностных структур земной коры с глубинными”, ч. 1, 2008.

13. Ермохин К.М. Аналитическое продолжение геофизических полей методом цепных дробей. Записки Горного института, т. 183, СПб, 2009.

14. Ермохин К.М., Жданова Л.А. Эффективный метод аналитического продолжения модельных и практических геофизических полей в область источников. Материалы 37 семинара им. Д.Г. Успенского “Вопросы теории и практики интерпретации гравитационных, магнитных и электрических полей”. М. ИФЗ РАН, 2010.

Размещено на Allbest.ru