Диссертация (Электромагнитные геометрические зондирования с донными косами при поисках углеводородов на мелководье), страница 11

Вместе с тем числопакетов, которые доведены до уровня, позволяющего использовать их для обработкипроизводственных данных, относительно невелико.Большая часть таких программных пакетов основаны на конечно-разностном методе, вкотором система уравнений Максвелла аппроксимируется системой конечно- разностныхуравнений на пространственной сетке (Спичак, 1983; Druskin and Knizhnerman, 1994; Controlledsource electromagnetic inversion for resource exploration. Oldenburg [et al.], 2005; Commer andNewman, 2004; 3D inversion of marine CSEM data using a fast finite-difference time-domain forwardcode and approximate Hessian-based optimization. Zach [et al.], 2008; 2.5D forward and inversemodeling for interpreting low-frequency electromagnetic measurement. Abubakar [et al.], 2008;Mackie, Watts and Rodi, 2007; и др.).К этой же группе относиться пакет Otze, разработанный C.Scholl и N.Yavich (Yavich andScholl, 2012), который был использован в настоящей работе для выполнения 2.5D и 3Dинверсии (см.

ниже).Достаточно активно в практике морской электроразведки применяются программы,использующие метод конечных элементов. Список программ, которые применялись длярешения практических задач морской геоэлектрики с контролируемым источником, включает,но не ограничивается: (Li and Key, 2007; Решение трёхмерных нестационарных задачимпульсной электроразведки. Иванов [и др.],2007; Nechaev, Shurina and Botchev, 2008;Сравнение методов решения трёхмерных задач становления поля с использованиемаппроксимаций в частотной и временной областях. Персова [и др.], 2013; и др).52Программы, основанные на методе интегральных уравнений, применяются относительнонечасто, однако известны впечатляющие примеры инверсии реальных морских измерений.Неполный список включает: (Dmitriev and Barashkov, 2012; 3D inversion of towed streamer EMdata, 2011; и др.).

К этой же группе относиться программа А.А. Петрова (Петров, 1992), котораябыла использована в настоящей работе для 3D моделирования.Для 2.5D и 3D моделирования и инверсии использовался пакет Otze, реализованныйC. Scholl и Н.Б. Явичем (Yavich and Scholl, 2012). Прямая 2D задача, которая ввидутрёхмерности источника часто обозначаемая как 2.5D, решается методом конечных разностей вчастотной области. Алгоритм минимизации основан на Тихоновской регуляризации (Тихонов иАрсенин, 1979) с выбором параметра регуляризации на каждой итерации.

Кластерное времябыло предоставлено Fugro Electromagnetic Italy Srl.Ниже приведено конспективное описание 2D модулей пакета, которые активноиспользовались автором.Решениепрямойзадачизаключаетсявпоискекомпонентнапряжённостейэлектрического и магнитного поля:E  ( E x ( x, y, z), E y ( x, y, z), E z ( x, y , z )),(2.39)H  ( H x ( x, y, z ), H y ( x, y, z), H z ( x, y, z ))для заданной геологической модели УЭС:00  xy ( x, y , z )00  xy ( x, y, z ).00 z ( x, y, z ) Предполагается(2.40)трансверсальная анизотропия УЭС. Компоненты электромагнитногополя удовлетворяют уравнениям Максвелла в частотной области, которые могут записаны ввиде:rotE  i 0 H  0,1rotH   E  J ,где  – круговая частота, J – объёмная плотность сторонних токов,(2.41)i – мнимая единица, 0  4 10 7 Гн/м – магнитная проницаемость вакуума.

Единственность, а также физическийсмысл решений системы (2.41) обеспечивается краевыми условиями на бесконечности:E( x, y, z)  0 и H ( x, y, z)  0 при x  y  z  .(2.42)Для 2D моделей трудоёмкость вычислений сокращается выполнением преобразованиеФурье в направлении однородности модели (ось y ):53Eˆ ( x, k y , z )  E ( x, y, z )e ik y ydy ,Hˆ ( x, k y , z )  H ( x, y , z ) e ik y ydy ,(2.43)где k y – пространственная частота. Исключение зависимых неизвестных из Фурье-образовуравнений (2.41) приводит к следующей системе для Eˆ y ( x, k y , z), Hˆ y ( x, k y , z) (Electromagneticmethods in applied geophysics, 2006):   1 Hˆ y    1 Hˆ y   1 Eˆ y    1 Eˆ y 1 ˆiky2 K 2 x    xy E y x  K xy2  xy x  z  K z2  z z xzzK z  xy   Jˆ    Jˆ   Jˆ y  ik y   x2    z2  , x  K xy  z  K z   1 Hˆ y    1 Hˆ y  ik y    1 Eˆ y    1 Eˆ y  ˆ  2Hy x  K z x  z  K xy2 z   0  x  K z2 z  z  K xy2 x   Jˆ xz  K xy2(2.44)   Jˆ z  , x  K 2  z здесь Jˆ  ( Jˆ x ( x, k y , z), Jˆ y ( x, k y , z), Jˆ z ( x, k y , z)) – преобразование Фурье плотности стороннихтоков,K xy2  k y2 i 0, xyK z2  k y2 i 0.z(2.45)Найдя решение (2.44) для каждой из интересующих пространственных частот (возможноприближённое), и выполняя обратное преобразование Фурье, получают решения исходнойпрямой задачи,1E ( x, y, z ) 2ik y y Eˆ ( x, k y , z )e dy,1H ( x, y , z ) 2 Hˆ ( x, ky, z )eik y ydy.

(2.46)Для системы уравнений (2.44) краевые условия примут вид,Eˆ y ( x, k y , z)  0 и Hˆ y ( x, k y , z)  0 при x  z  ,(2.47)Однако такие условия затруднительно использовать в вычислительном методе. Поэтомув рассмотрение вводят достаточно большую прямоугольную расчётную область V , и краевыеусловия (2.47) заменяют аналогичными на её границе:Eˆ y ( x, k y , z)  0 и Hˆ y ( x, k y , z)  0 на границеV .(2.48)Замена условий (2.47) на (2.48) вносит некоторую ошибку, так как решение задачи(2.44),(2.47) заведомо отличается от задачи (2.44),(2.48).

Однако эту ошибку можно сделать54сколь угодно малой в силу экспоненциального затухания ЭМ полей в геологических породах исолёной воде (Ваньян, 1997).Вычислительная часть решения прямой задачи (2.47),(2.48) основывается на конечноразностном (КР) подходе (Самарский, 1977). Расчётную область V покрывают линиями сетки,параллельными координатным осям. Характеристики используемойКР сетки определяютточность и эффективность решения прямой задачи.

Шаг КР сетки диктуется длинной волны,удалённостью от источника, размерами источника. Это приводит к КР сеткам, неоднородным вгоризонтальном и вертикальном направлениях. Рисунок 2.13 иллюстрирует типичнуюиспользуемую сетку.Рис. 2.13. Конечно-разностная сетка в части расчётной области (источник находиться в началекоординат)Зависимость шага КР сетки от длины волны влечёт за собой зависимость от частотызондирования.

Поэтому для зондирований, использующих несколько частот, необходимо либогенерировать отдельные КР сетки, либо генерировать одну, удовлетворяющую требованиямнаибольшей и наименьшей частоты (шаг сетки пропорционален длине волны от наибольшейчастоты, а размер расчётной области пропорционален длине волны от наименьшей частоты).В конечно-разностном подходе поиск полей Ê y и Ĥ y , определённых в области V ,заменяется поиском полей, определённых лишь конечном числе точек, называемыми степенямисвободы:i  N 1, k  K 1Eˆ y   E 1 1  i  2 k  2  i 1,k 1,i N ,k KHˆ y  H ik i 1, k 1(2.49)Де факто, стандартом является использование степеней свободы, изображённых нарисунке 2.14, предложенных в работе (Yee, 1996).55Рис.

2.14 Степени свободы для дискретизации уравнений (5), электрическое поле •, магнитное поле ×Такое расположение позволяет избегать интерполяции на этапе дискретизацииуравнений.Затем дифференциальные операторы системы (2.44) заменяют конечно-разностнымианалогами, что в результате приводит к системе линейных алгебраических уравнений. Вопрососреднения коэффициентов УЭС, необходимого при построении этой системы, рассмотренниже.Её размер составляет от тысяч до сотен тысяч неизвестных, a матрица имеетразряженную структуру.

Эллиптичность задачи (2.44) и краевые условия (2.48) влекутневырожденность матрицы. Перечисленные свойствапозволяют находить решение системыкак точными, так и итерационными методами. Однако, с практической точки зрения,итерационные методы применяют лишь в случаях, когда уже найдено приемлемое начальноеприближение (например, с близкой частоты или геологической модели), так как в противномслучае потребуется существенное количество итераций метода. Это связано с отсутствиемалгоритмов эффективного предобусловливания системы (2.44).Конечно-разностное моделирование требует определения УЭС в ячейках КР сетки.Поэтому для сложных геологических моделей (включающих, например, рельеф дна,выклинивания и пр.) зачастую рекомендуют прибегать к конечно-элементному моделированию,имеющему существенное преимущество в геометрической гибкости расчётных сеток.

Но это необязательно. Например, работа (A finite difference scheme for elliptic equations with roughcoefficients using a Cartesian grid nonconforming to interfaces. Moskow [et al.], 1999), предлагаетнесколько алгоритмов осреднения коэффициентов дифференциального уравнения по КР сетке.Методика осреднения, использованная в рамках данного проекта, основана на принципесовпадения сопротивлений. Пусть для некоторой прямоугольной ячейки, заданной накоординатной плоскости [ x1 , x2 ] [ z1 , z2 ] , требуется определить осреднённое горизонтальное56hудельное сопротивление  xy на основе удельного сопротивления геологической модели xy ( x, z) .