автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.16, диссертация на тему:Математическое моделирование трехмерных электромагнитных полей в двумернонеоднородных средах

ДАЛЬНЕВОСТОЧНОЕ ОТДЕЛЕНИЕ РОССИЙСКОЙ АКАДЕМИИ НАУК ВЫЧИСЛИТЕЛЬНЫЙ ЦЕНТР

На правах рукописи

УДК 519.63

Пересветов Владимир Викторович

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ТРЕХМЕРНЫХ ЭЛЕКТРОМАГНИТНЫХ ПОЛЕЙ В ДВУМЕРНО- НЕОДНОРОДНЫХ СРЕДАХ

05.13.16 - применение вычислительной техники, математического моделирования и математических методов в научных исследованиях

Диссертация на соискание ученой степени кандидата физико- математических наук

Научный руководитель: профессор, доктор физико-математических наук С. И. Смагин

Хабаровск 1998

Оглавление

Введение 5

1 2.5-Б ЗАДАЧА МОДЕЛИРОВАНИЯ ЭЛЕКТРОМАГНИТНЫХ ПОЛЕЙ 24

1.1 Описание параметров модели....................................24

1.2 Уравнения электромагнитного поля............................25

1.2.1 Уравнения Максвелла для первичного и вторичного полей......................................................26

1.2.2 Уравнения электромагнитного поля в пространствах (кх,у,г) и (к х) я)..................29

1.3 Электромагнитное поле в частных случаях источника и структуры среды..................................................34

1.3.1 Уравнения поля для частных значений кх............34

1.3.2 2-Б задача моделирования с плоским источником . . 37

1.3.3 Электромагнитное поле в слоистой среде..............39

1.4 Единственность решения задач в ограниченных и неограниченных областях..............................................40

1.4.1 Краевые условия и единственность решения задач в пространстве (кх,у,г) ..................................41

1.4.2 Задачи в пространстве (кх,ку,г) ......................46

1.4.3 Единственность решения задачи в бесконечной полосе 50

2 ЧИСЛЕННОЕ РЕШЕНИЕ МЕТОДОМ КОНЕЧНЫХ

ЭЛЕМЕНТОВ 55

2.1 Численное решение 2.5-D задач в ограниченной области пространства......................................................55

2.1.1 Дискретное преобразование Фурье......................55

2.1.2 Обобщенное решение в ограниченной области .... 56

2.1.3 Приближенное решение задачи методом Бубнова-Галеркина..................................................60

2.2 Алгоритмы приближенного решения 2-D задач методом конечных элементов при кх = 0 ..................................61

2.2.1 Приближенное решение в ограниченной области . . 62

2.2.2 Алгоритмы для прямоугольных конечных элементов 63

2.2.3 Алгоритмы для треугольных конечных элементов . 66

2.2.4 Сочетание прямоугольных и треугольных конечных элементов..................................................69

2.3 Результаты численного решения 2.5-D задач..................70

2.3.1 Сходимость на последовательности сеток...... . 71

2.3.2 Сходимость метода SOR в пространстве (кх, у, z) . 72

2.3.3 Модель прямоугольной неоднородности в трехслойной среде..................................................73

2.3.4 Сравнение результатов численных расчетов с данными физического моделирования......................75

2.4 Результаты моделирования магнитотеллурических полей . 80

2.4.1 Модель трех сегментов на проводящем основании . . 80

2.4.2 Модель низменности цилиндрической формы .... 82

2.4.3 Модель возвышенности в форме трапеции............84

3 ПРИМЕНЕНИЕ МЕТОДА ДЕКОМПОЗИЦИИ

ОБЛАСТИ 86

3.1 Алгоритм численного решения в пространстве (kx,y,z) . . 86

3.2 Сходимость альтернирующего метода Шварца..............88

3.2.1 Слоистая среда............................................88

3.2.2 Однородная среда........................................92

3.2.3 Сочетание краевых условий различного типа .... 94

3.3 Сходимость метода декомпозиции без перекрытия............97

3.4 Результаты численных экспериментов........................98

3.4.1 Численное решение 2.5-Б задачи моделирования поля для прямоугольной неоднородности в полупространстве ............................................98

3.4.2 Моделирование электромагнитного поля для двух

2-Б неоднородностей в пятислойной среде......102

Заключение 123

Список литературы..........................125

/

Введение

В работе рассматриваются вопросы математического моделирования пространственных и плоских электромагнитных полей гармонических источников в двумерно- неоднородных (2-Б) средах.

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

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

ческое тело, расположенное, в общем случае, в плоскослоистой среде. 2-Б неоднородность может быть образована изменением глубины раздела слоев (складки и разрывные нарушения со смещением: горст, грабен).

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

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

Задача для диполя в 2-Б среде называется в литературе 2.5-Б задачей, так как дипольный источник создает З-Б поле, а среда задана 2-Б моделью. Задача моделирования магнитотеллурического поля с моделью источника в виде плоской волны в 2-Б средах является 2-Б задачей.

Электромагнитное поле в плоскослоистой среде с 2-Б цилиндрической неоднородностью (или 3-Б телом) можно представить в виде суммы первичного и вторичного полей. Первичное поле является полем заданных источников в плоскослоистой среде. Вторичное поле (иногда используется другой термин - аномальное поле) обусловлено влиянием неоднородности. Раздельное нахождение первичного и вторичного полей позволяет строить универсальные и достаточно точные алгоритмы расчета полного поля. Раздельное нахождение полей широко применяется в

моделировании (см. [88], [105], [106] и др.).

В настоящей работе рассматриваются вопросы моделирования только вторичного поля. Задача расчета первичного поля в слоистой среде является самостоятельной задачей. Применение метода интегральных преобразований позволяет находить компоненты электромагнитного поля дипольного источника в плоскослоистой среде с высокой точностью. При проведении численных экспериментов для моделей 2-Б неоднородности в полупространстве и трехслойной среде была использована программа расчета электромагнитного поля дипольного источника произвольной ориентации в плоскослоистой среде, авторство которой принадлежит Смагину С.И. и Мазалову В.Н. Алгоритмы этой программы развиты в работах [20], [21], [51], [52] и основаны на переносе пути интегрирования в комплексную плоскость. При этом у подинтегральных функций появляются экспоненциальные множители с отрицательными действительными частями, которые обеспечивают отсутствие осциляций и быстрое их убывание на бесконечности, независимо от количества слоев и положения источника.

Для моделирования электромагнитных полей в 2-Б и 3-Б средах используются, в основном, следующие методы: конечных элементов (МКЭ), конечных разностей (МКР), интегральных уравнений (МИУ), аналитический и гибридный метод. Гибридный метод основан на совместном использовании двух методов: МКЭ или МКР в области содержащей неоднородность и МИУ, аналитического или спектрального для решении внешней краевой задачи. Исключительно аналитический способ решения рассматриваемых задач (см. [56] и др.) возможен только для некоторых простейших случаев на плоскости. Они используются, как правило, для контроля точности численных расчетов. Метод МКР начал использоваться раньше других. Формулировка задачи математического моделирования, постановка граничных условий, вопросы решения

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

В последнее время в математическом моделировании физических полей в сложных средах нашел широкое распространение МКЭ. МКЭ по сравнению с МКР хорошо аппроксимирует области со сложными границами, автоматически обеспечивает выполнение условий сопряжения для компонент поля на границах разрыва среды, позволяет выбирать вид базисных функций для наиболее точной аппроксимации искомого решения на отдельных участках области моделирования. Теория МКЭ излагается в [15], [22], [27], [53], [54] и др. Начиная с работы [76] МКЭ начал использоваться для решения задач геоэлектрики.

Моделирование электромагнитных полей в З-Б средах МКЭ представлено в работах [5], [83], [95], [98]. Подход, основанный на минимизации функционала определяющего энергию электромагнитного поля, используется в [5], [95]. В этих работах исходная задача сводится к поиску стационарной точки функционала, который в случае электрического поля записывается в виде:

/ {2ц)~1[(гоЬ Е)2 + к2Е2 - ЯшщГЦЩ (В.1)

V

где /л - магнитная проницаемость, ш - круговая частота, к - волновое число, ,|ст- плотность сторонних токов. Вычисление приближенного значения функционала в прямоугольной подобласти трехмерной сетки с помощью кубатурных формул дает возможность получить коэффициенты матрицы СЛАУ. В [5] приводятся вид коэффициентов и результаты сходимости решения для 13 и 21- точечных шаблонов, полученных при различных видах кубатурных формул. В работе [98] используется другой

подход: уравнение

rot rot Е + fc2E = О

решается методом Бубнова- Галеркина. Для приближенного решения системы трех дифференциальных уравнений второго порядка используется МКЭ. Представлены результаты расчетов для прямоугольной призмы.

В последнее время появились методы моделирования рассчитанные на многопроцессорные ЭВМ. В [113] матрица для 3-Б магнитного поля разлагается на две субматрицы, которые решаются параллельно с использованием итерационого метода.

В работах [5], [85], [97], [99], [105], [106] МКЭ применяется для моделирования электромагнитного поля в 2-Б средах. Используется как вариационная формулировка задачи, так и метод Бубнова- Галеркина. Решение для 2-Б задач в вариационной формулировке [108] находится из экстремума функционала, записываемого аналогично (В.1), где Е = 1ХЕХ, 1Х— единичный вектор. В [106] метод Бубнова- Галеркина применяется для решения уравнения для вторичного поля:

где Е3Х - х-компоненты первичного и вторичного поля, а' = а + iwe , а - проводимость, Дет'- аномальная комплексная проводимость, е- диэлектрическая проницаемость. Приведено уравнение только для случая Е- поляризации. Моделирование вторичного поля позволяет повысить точность расчетов на низких частотах (особенно важно для магнитной моды при расчетах со стандартной точностью представления вещественных чисел). В этой работе применяется прямоугольная неравномерная сетка, каждая ячейка которой делится, в свою очередь, двумя диагоналями на четыре треугольника. Данный способ триангуляции области моделирования позволил в [105] получить хорошие результаты для маг-

д_ ду

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

Задачи геоэлектрики необходимо решать в неограниченной области. Электромагнитное поле на бесконечности должно удовлетворять условиям излучения [18], [29]. Это требование достаточно для получения однозначного решения уравнений Максвелла. В случаях, когда граница внешней области имеет бесконечно удаленные точки, необходимо задавать парциальные условия излучения [50]. При решении многих задач используется принцип предельного поглощения [16], [50], [55], [59].

Проблема краевых условий для МКР и МКЭ является одной из наиболее важных. Во многих работах граница области моделирования берется досточно удаленной, что позволяет использовать в качестве краевых условий значения первичного поля, если моделируется полное поле [13], или нулевые значения, если моделируется вторичное поле [106]. Однако, это приводит к увеличению числа узлов сетки и, следовательно, к усложнению решения задачи. В однородном окружающем пространстве могут быть использованы граничные условия дифференциального типа, полученные в [46], [47]. Особенно остро проблема краевых условий стоит в непроводящей атмосфере, так как вторичное поле здесь затухает в наименьшей степени. Асимптотические краевые условия Дирихле для 2-Б среды в случае Е- поляризации получены в [71]. В более поздней работе [109] теми же авторами получены краевые условия дифференциального типа, содержащие производные первого порядка. В работах [13], [14], [116] формулируются обобщенные асимптотические краевые условия дифференциального типа, содержащие производные произвольного порядка, которые позволяют сократить размеры области моделирования до десятых долей эффективной длины волны. В то же время, в рамках итерационного подхода к решению СЛАУ, скорость сходимости при использовании таких краевых условий невысока. Предложенные в [71]

краевые условия используются в настоящей работе для моделирования магнитотеллурического поля (см. раздел 1.3.2).

При моделировании динамических и акустических полей в океане также требуется решать задачи в неограниченной области. В [2], [62], [80] задача в неограниченной области сводится к эквивалентной задаче в ограниченной области введением искуственной границы, на которой ставится точное нелокальное граничное условие, моделирующее условие излучения Зоммерфельда. Полученная краевая задача решается МКЭ или МКР, в результате необходимо решать плохо обусловленную систему уравнений. Эффективный прямой метод решения такой системы для области типа протяженного волновода развивается в [2], [62].

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

Большинство авторов отдает предпочтение прямым методам решения систем в случае 2-D задач геоэлектрики как более эффективным, надежным и технологичным [6]. Итерационные методы решения систем, активно применявшиеся в ранних алгоритмах, имеют ряд достоинств: меньшие требования к объему оперативной памяти, возможность снижения вычислительных затрат путем улучшения начального приближения, возможность корректирования системы уравнений в процессе решения. В работе [6] при моделировании МКР сравнивались три метода: итерационный метод верхней релаксации (SOR), метод исключения Гаусса и метод блочной прогонки. Установлено, что прямые методы не уступают, а в задачах моделирования сложных геоэлектрических разрезов с числом узлов 3 — 5 • 103 превосходят метод SOR по быстродействию и позво-

ляют получить точное решение с большей надежностью. Для данного класса задач применялись также и другие итерационные методы: переменных направлений [10], с использованием алгоритмов сверхрелаксации [7], циклической редукции [69].

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

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

В последнее время активн�