Практикум по электроразведке (5 курс) Моделирование мт-поля в двумерной среде методом конечных разностей Практическое задание

МГУ им. М.В. Ломоносова, Геологический Факультет

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

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

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

1. О методах решения прямых задач электроразведки

Рис. 1.1. Классификация методов

Для простых моделей геоэлектрического разреза прямые задачи различных методов электроразведки имеют аналитическое решение. Электромагнитное поле в сложных геоэлектрических разрезах изучается с помощью методов физического и математического моделирования.

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

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

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

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

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

,

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

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

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

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

Рис. 1.2. Классификация методов математического моделирования.
Методы интегральных уравнений (МИУ) основаны на применении аппарата интегральных уравнений Фредгольма. Их характерной особенностью является то, что поле рассматривается как сумма нормального и аномального полей. Физически идея МИУ заключается в замене неоднородностей некоторыми вторичными источниками, расположенными либо в объеме (метод объемных интегральных уравнений - МОИУ), либо на поверхности границ раздела (метод поверхностных интегральных уравнений - МПИУ). Поэтому сеткой покрывается лишь область неоднородности – в МОИУ, или вовсе лишь ее границы – в МПИУ. Вследствие этого достоинством программ, основанных на МИУ, является относительно низкий объем требуемой памяти и небольшое время счета. С точки зрения разработчика программы важным достоинством метода является также возможность контролировать работу программы на промежуточных этапах. Недостатком МИУ является необходимость разделения неоднородностей и нормального разреза, вследствие чего трудно работать с моделями, выходящими за рамки горизонтально-слоистой структуры нормального разреза. Поэтому, в частности, при моделировании переменных полей МИУ в настоящее время применяется редко.

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

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

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

2. Магнитотеллурическое поле в двумерно-неоднородной среде
Основной особенностью магнитотеллурического поля в двумерно-неоднородной (2D) среде является его разделение на две независимые части, называемые поляризациями или модами. Рассмотрим эту особенность подробно.

П
Рис. 2.1. 2D среда.

редположим, что свойства среды меняются только вдоль осей x и z (рис. 2.1). Ось y в таком случае направлена перпендикулярно плоскости x0z и вдоль нее свойства среды не меняются. Поэтому ось y называется осью однородности среды.

Напомним, что два первых уравнения Максвелла для комплексных амплитуд напряженностей электрического и магнитного полей записываются в виде:

Здесь - круговая частота колебаний поля, - магнитная проницаемость вакуума. При этом среда в окрестности точки полагается однородной с электропроводностью . Кроме того, полагается, что отсутствуют токи смещения и сторонние токи.

В декартовой системе координат эти уравнения принимают вид:

Поскольку среда вдоль оси y не меняется, все производные по y равны нулю:

Переходя к скалярной записи первого уравнения, получим:

(1.1)

Аналогично, второе уравнение в скалярной форме примет вид системы:

(1.2)

Теперь выпишем в одну систему второе уравнение системы (1.1), а также первое и третье уравнения системы (1.2):

(1.3)

Составим также систему, включающую первое и третье уравнения системы (1.1) и второе уравнение системы (1.2):

(1.4)

Легко заметить, что система (1.3) содержит только компоненты , и , а система (1.4) - только компоненты , и .

Таким образом, электромагнитное поле распадается на две независимые части. Первая из них, описываемая системой (1.3), называется E-поляризацией, поскольку в этом случае электрическое поле поляризовано вдоль оси y (т.е. имеет только одну компоненту, которая параллельна оси y). Аналогично, вторая часть поля описывается системой (1.4), которая определяет H-поляризацию, так как содержит поляризованное вдоль оси y магнитное поле. В английском языке E-поляризации отвечает термин TE-mode, а H-поляризации – термин TM-mode. Поэтому в литературе E-поляризацию часто называют TE-модой, а H-поляризацию – TM-модой.

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

Рассмотрим сначала случай E-поляризации. Для этого выразим компоненты магнитного поля из второго и третьего уравнений системы (1.3):

(1.5)

Теперь подставим эти соотношения в первое уравнение рассматриваемой системы:

Это выражение нетрудно привести к виду:

Где квадрат волнового числа . Полученное уравнение называется двумерным уравнением Гельмгольца. Если найти из этого уравнения и далее по формулам (1.5) рассчитать и , то задача определения E-поляризованного поля будет решена полностью. Далее можно рассчитать продольный импеданс на земной поверхности и построить амплитудные и фазовые продольные кривые МТЗ.

Чтобы получить дифференциальное уравнение для H-поляризации, выразим из первого и второго уравнений системы (1.4) компоненты электрического поля:

(1.6)

Подставив эти соотношения в третье уравнение системы (1.4), получим:

В данном случае выносить электропроводность из-под знака дифференцирования нельзя, поскольку в общем случае она меняется вдоль осей x и z. Деля левую и правую части последнего выражения на , получим уравнение:

Определив из этого уравнения и рассчитав и по формулам (1.6), мы полностью решим задачу для H-поляризации. Остается лишь рассчитать поперечный импеданс на земной поверхности и построить амплитудные и фазовые поперечные кривые МТЗ.

Основное различие E-поляризованного и H-поляризованного полей заключается в том, что в первом случае аномалии имеют индукционную, а во втором – гальваническую природу. В E-поляризации электрическое поле поляризовано вдоль структур, поэтому токи не пересекают границ областей различного сопротивления (рис. 2.2 а). Аномалии проявляются за счет избыточных токов, индуцируемых в проводящих неоднородностях или низкого числа токов, индуцируемых в непроводящих неоднородностях. В то же время в H-поляризации токи текут вкрест структур (рис 2.2 б). При этом они пересекают границы неоднородностей, и аномалии проявляются за счет гальванического затекания в них токов (для низкоомных неоднородностей) или обтекания их токами (для высокоомных неоднородностей).

Рис. 2.2. Распределение токов в 2D среде: а) при Е-поляризации, б) при Н-поляризации.

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

3. Двумерное моделирование МТ-поля методом конечных разностей
Метод конечных разностей наиболее широко используется при моделировании магнитотеллурического поля в двумерной среде. В частности, на нем основана и программа, с которой требуется работать при выполнении данной задачи. Поэтому рассмотрим этот метод применительно к решению прямой 2D задачи МТЗ подробнее.

Рис. 3.1. 2D модель, применяемая в методе конечных разностей.

Пусть изучаемая среда (рис. 3.1) состоит из аномальной прямоугольной области (, ) с произвольным распределением электропроводности и из трех нормальных областей, в которых электропроводность меняется только по : области (, ) с электропроводностью (левый нормальный разрез), области (, ) с электропроводностью (правый нормальный разрез) и области (, ) с электропроводностью (нижний нормальный разрез). Среда контактирует с непроводящей атмосферой (). Модель возбуждается плоской E-поляризованной или H-поляризованной электромагнитной волной, падающей на поверхность среды. Требуется рассчитать значения продольного и поперечного импедансов на поверхности области (при ).

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

Постановка краевой задачи. Как было показано в предыдущем параграфе, компонента Е-поляризованного поля может быть получена в результате решения уравнения:

(3.1)

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

Классическая теория описывает три способа постановки краевых задач:

Краевая задача Дирихле. На границе задаются значения искомой функции.
Краевая задача Неймана. На границе задаются значения производной искомой функции по нормали к границе.
Третья краевая задача. На границе задаются значения определенной линейной комбинации искомой функции и ее производной по нормали к границе.

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

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

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

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

Конечно-разностная аппроксимация краевой задачи. Для конечно-разностной аппроксимации уравнения (3.1) в области моделирования покроем эту область прямоугольной неравномерной сеткой (рис. 3.2):

Введем вспомогательную сетку , составив ее из узлов, помещенных в центры ячеек сетки :

Представим y-компоненту электрического поля как дискретную функцию на сетке , а волновое число – как дискретную функцию на сетке .

Для того, чтобы понизить порядок производных в уравнении (3.1), обычно конечно-разностно аппроксимируется не само уравнение (3.1), а интегральное тождество, получающееся в результате интегрирования этого уравнения по элементарной ячейке :

и последующего применения к левой части полученного выражения формулы Грина, связывающей интеграл по поверхности с криволинейным интегралом по контуру :

, (3.2)

где - единичный вектор нормали (рис. 3.2).

Итак, выразим интегралы в (3.2) приближенно в терминах дискретных функций и :

;

.

Подставляя эти формулы в (3.2) и вводя для краткости обозначения, получим уравнение:

, (3.3)

причем может принимать значения от до , а - от до . В этом уравнении - комплексные, а , , и - действительные константы:

; ; ; ;

; ; ;

; ;

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

Решение системы линейных уравнений. Записывая уравнения вида (3.3) для каждого внутреннего узла сетки , получим систему из разностных уравнений. Каждое уравнение этой системы связывает значения электрического поля в пяти соседних узлах.

Система уравнений вида (3.3) и соответствующих краевых условий может быть записана в матричном виде:

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

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

Для решения системы уравнений применяют как прямые, так и итерационные методы. Напомним, что в теории численных методов прямыми называют такие методы, которые приводят к точному решению системы за определенное конечное число действий, а итерационными – такие методы, в основе которых лежит последовательное уточнение приближенного решения. Соответственно, первые более надежны, а вторые - менее требовательны к объему оперативной памяти и, в большинстве случаев, позволяют получать результат за меньшее время. В настоящее время количество программ, использующих прямые, и количество программ, использующих итерационные методы, примерно одинаково.

В результате решения системы линейных уравнений получаются значения компоненты в узлах сетки .

Расчет компонент поля, не участвующих в разностной схеме. На последнем этапе рассчитываются компоненты магнитного поля и . Для этого используются соотношения (1.5), причем производные горизонтального электрического поля обычно рассчитываются разностным способом. Далее может быть рассчитан продольный импеданс и продольные кривые кажущегося сопротивления.

Решение рассматриваемой задачи для случая Н-поляризации может быть получено аналогичным образом.

Контроль точности моделирования. На всех рассмотренных этапах конечно-разностного решения 2D задачи МТЗ возникают определенные погрешности.

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

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

4. Особенности работы с программой К. Новожинского

Для моделирования МТ-поля в 2D среде в рамках настоящей работы необходимо использовать программу К. Новожинского (Институт Геофизики Польской Академии Наук).

Основные особенности алгоритма данной программы таковы:

При постановке краевой задачи используются граничные условия Дирихле. На левой и правой границах значения компонент поля получаются в результате расчета одномерных задач для левого и правого разрезов, а на верхней и нижней границах – в результате линейной интерполяции.
При конечно-разностной аппроксимации краевой задачи для понижения порядка производных, входящих в уравнения поля, применяется формула Грина.
При решении системы линейных алгебраических уравнений применяется быстрый прямой метод - метод Банашевича-Холецкого.
При расчете компонент поля, не участвующих в разностной схеме, в случае Е-поляризации применяется разностная трехточечная, а в случае Н-поляризации –двухточечная схема.

Программа позволяет решать не только прямую, но и обратную задачу МТЗ. При выполнении данной работы мы будем решать только прямую задачу. Поэтому вопросы, связанные с использованием программы К. Новожинского для решения обратной задачи, в рамках данной разработки рассматриваться не будут.

Рис 4.1 А. Меню File. Рис 4.1 Б. Меню Data.

Рис 4.1 В. Меню Calculate. Рис 4.1 Г. Меню Results.

Для удобства ввода данных и анализа результатов моделирования расчетные программы снабжают графическими пользовательскими интерфейсами, называемыми также графическими оболочками. Наличие графического интерфейса позволяет уменьшить число ошибок ввода-вывода данных, повысить производительность работы с программой и превратить эту работу в приятное и увлекательное занятие. У программы К. Новожинского также имеется графическая оболочка – программа IGF_MT2D. Она рассчитана на 32-битные операционные системы семейства Microsoft Windows. Рассмотрим вкратце особенности работы с программой К. Новожинского и графической оболочкой к ней.

В начале работы необходимо запустить программу IGF_MT2D. Меню этой программы представлено на рис. 4.1. Вызов команд управления файлами данных ("Создать", "Открыть", "Сохранить") осуществляется из меню "File". Ввод данных производится с помощью четырех диалоговых окон, которые могут быть вызваны при выборе соответствующего пункта меню "Data". Для запуска расчетной программы необходимо затем выбрать один из пунктов меню "Calculate". Полученные на этом этапе результаты расчета можно анализировать в виде кривых зондирования или графиков профилирования и сохранять в удобном формате (меню "Results"). Для доступа к инструкции (на английском языке) и информации о программе имеется также меню "Help".

Ввод данных. Для создания нового файла данных следует выбрать пункт "New" меню "File". После этого одно за другим появятся четыре диалоговых окна ввода данных, которые могут также быть вызваны из меню "Data".

Рис 4.2. Диалоговое окно Grid.

Диалоговое окно Grid (рис. 4.2) предназначено для ввода размеров сетки. В нем задаются число ячеек вдоль оси z в атмосфере (Number of meshes along Z above the surface), вдоль оси z Земле (Number of meshes along Z under the surface) и вдоль оси x (Number of meshes along X). Для расчетов в данной задаче размеры сетки могут быть порядка 60x30 ячеек.

Для каждого из этих чисел предусмотрена возможность задания числа частей, на которые будут разбиваться ячейки в соответствующем направлении (Number of parts to divide mesh). Технически это сводится к увеличению дробности сетки (например, если задать значения, равные двум, сетка просто станет в два раза более частой). Напомним, что этот прием позволяет изучать влияние погрешностей конечно-разностной аппроксимации краевой задачи и численного расчета компонент поля, не участвующих в разностной схеме. Кроме того, он позволяет ослабить эффект превращения наклонных границ в лестнице-подобные ломаные, возникающие за счет того, что среда разбивается на прямоугольные ячейки. Конечно, при этом сильно возрастает время счета. Отметим также, что при интерпретации данных, сводящейся обычно к многократному решению прямой и обратной задачи с некоторыми промежуточными исправлениями параметров, делить ячейки на части для более-менее грамотно составленной модели не рекомендуется. Однако есть смысл сделать это для итоговой модели, предназначенной, например, для публикации. В рамках данной задачи необходимо выполнить расчеты для чисел частей, на которые нужно делить ячейки, равных 1 (ячейки не делятся на части) и 2 (ячейки делятся на две части). Совпадение результатов расчетов для этих двух случаев будет свидетельствовать о малости величины влияния указанных видов погрешностей.

Для выполнения граничных условий необходимо, чтобы аномальное поле вблизи границ было пренебрежимо малым. Для этого границы области моделирования нужно отнести далеко от неоднородностей. С этой целью вблизи границ приходится задавать ячейки большого размера с фоновыми сопротивлениями, увеличивая размеры этих ячеек с геометрическим шагом (обычно равным 2) по мере приближения к границам. Программа предусматривает возможность автоматического добавления слева и справа таких дополнительных ячеек, обеспечивающих удаленность боковых границ области моделирования, примерно равную некоторому подобранному опытным путем значению. Для этого взводится флажок "Additional large meshes at left and at right". Эта возможность призвана помочь неопытным пользователям снизить вероятность ошибок при создании модели и упростить этот процесс. Однако в большинстве случаев квалифицированный пользователь лучше справится с решением этой задачи вручную. В рамках настоящей работы требуется отключить процедуру автоматического создания таких ячеек. Кроме того, следует выполнить расчеты для двух разных расстояний до боковых границ области моделирования, различающихся в 3-5 раз. Совпадение результатов этих расчетов покажет, что погрешность постановки краевой задачи невелика.

Рис 4.3. Диалоговое окно Meshes.
В диалоговом окне Meshes (рис. 4.3) задаются размеры ячеек: вдоль оси z в атмосфере (Mesh sizes along Z above the surface), вдоль оси z в Земле (Mesh sizes along Z under the surface) и вдоль оси x (Mesh sizes along X). Размеры всех ячеек в соответствующем направлении могут быть сделаны равными при взводе флажка “All meshes of one size”, однако в данной работе делать этого нельзя.

Вертикальный размер верхнего ряда ячеек в Земле должен быть в несколько раз меньше глубины проникновения поля на самом маленьком периоде. По мере роста глубины размеры ячеек обычно увеличивают. Глубина до нижней границы области моделирования должна превышать глубину исследования на самом большом периоде (в противном случае могут не выполниться соответствующие граничные условия). Размеры ячеек в атмосфере лучше задавать таким образом, чтобы они увеличивались снизу вверх с геометрическим шагом, равным 2, от значения, равного высоте верхней ячейки в Земле, до достижения высоты, примерно равной глубине до нижней границы области моделирования.

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

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

Рис 4.4. Диалоговое окно Model.
В диалоговом окне Model (рис. 4.4) на заданной ранее сетке рисуется модель среды. Сетка подписана номерами узлов и (в скобках) километрами. Для лучшего восприятия лишь некоторые границы между ячейками показаны черными линиями (например, каждая пятая). Текущая позиция отмечена, жирной точкой, а ее координаты приведены в верхней части окна. Изменение текущей позиции осуществляется с помощью мыши или клавиш перемещения курсора. После нажатия клавиши Insert в текущей позиции создается вершина тела. Изменяя текущую позицию и нажимая Insert пользователь создает контур тела (нажатие в процессе этого Delete уничтожает контур). После того, как контур оказывается замкнутым, ограниченному им телу предлагается присвоить сопротивление. Новые тела могут рисоваться поверх существующих. Если какое-либо тело оказывается полностью перекрыто более новыми телами, оно автоматически уничтожается. Изменение положений вершин тел осуществляется с помощью клавиши перемещения курсора при нажатой клавише Shift. Изменение сопротивлений тел осуществляется после щелчка над соответствующим телом правой кнопкой мыши. Сопротивления тел показаны цветом, причем в правой части окна приведена поясняющая колонка. Для изменения сопротивления тела нужно щелкнуть в его пределах правой кнопкой мыши и в появившемся диалоговом окне ввести новое значение.

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

Рис 4.5. Диалоговое окно Curves.
В диалоговом окне Curves (рис. 4.5) водятся полевые данные, необходимые для инверсии. В рамках данной работы нужно заполнить лишь первый столбец таблицы, т.е. ввести периоды, для которых будет производиться расчет. Периоды должны изменяться с геометрическим шагом, покрывая весь необходимый диапазон.

Выполнение расчетов. После ввода данных необходимо сохранить их в файле с расширением “dat”. После этого станут активными пункты меню “Calculate” и можно будет вызвать процедуру расчета прямой задачи (пункт “Forward”). При этом запустится программа К. Новожинского, которая читает файл данных, производит вычисления и записывает результаты в файл с расширением “res”. По окончании расчета становятся доступными пункты меню “Results” и можно проводить анализ результатов.

Рис 4.6. Главное окно программы IGF_MT2D.
Просмотр результатов. Полевые данные и результаты расчета могут быть представлены либо в виде кривых зондирования (после выбора пункта меню Curves) - как на рис. 4.6, либо в виде графиков профилирования (после выбора пункта меню Profiles) - как на рис. 4.7. После выбора пункта меню Save результаты расчета, представленные в данный момент на экране, будут сохранены в файле формата программы Grapher фирмы Golden Software.

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

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

Более подробная информация о работе с программой содержится в ее инструкции.

Рис 4.7. Главное окно программы IGF_MT2D.

5. Контрольные вопросы

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

6. Литература

М.С. Жданов. «Электроразведка». Москва, Недра, 1986.
М.С. Жданов, И.М. Варенцов, Н.Г. Голубев, В.А. Крылов. «Методы моделирования электромагнитных полей (Материалы международного проекта COMEMI)». Москва, Наука, 1990.
К. Новожинский, П.Ю. Пушкарёв. «Анализ эффективности программ для двумерной инверсии магнитотеллурических данных». Физика Земли, 2001, № 6, с. 72 - 85.