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

МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ Им. М.В. ЛОМОНОСОВА

Факультет вычислительной математики и кибернетики

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

Ивахненко Владимир Игоревич

УДК 537.8

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ДИФРАКЦИИ ЭЛЕКТРОМАГНИТНЫХ ВОЛН НА ТРЕХМЕРНЫХ РАССЕИВАТЕЛЯХ, РАСПОЛОЖЕННЫХ В СЛОИСТОЙ СРЕДЕ.

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

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

Научный руководитель

доктор физико-математических наук

профессор

Ильинский A.C.

Москва, 1998 г.

СОДЕРЖАНИЕ

Стр.

Введение. 3

Глава 1. Моделирование рассеивания электромагнитного

поля на частице, расположенной на границе двух полупространств. 18 П. 1 Вывод объемного интегродифференциального уравнения

для рассеивателя, расположенного в свободном пространстве. 20

П.2 Численный метод решения. 25

П.З Решение системы линейных уравнений. 27 П.4 Вывод объемного интегродифференциального уравнения

для рассеивателя, расположенного на подложке. 29

П.5 Численная реализация метода. 34 П.6. Вычисление электродинамических характеристик в дальней зоне. 36 Глава 2. Результаты математического моделирования

рассеяния электромагнитных волн на диэлектрических телах. 39

П. 1. Рассеиватели, расположенные в свободном пространстве. 40 П.2. Рассеиватели, расположенные на границе двух

полупространств. 59 Глава 3. Моделирование дифракции электромагнитной волны

на локальном идеально проводящем экране. 76

П1. Описание математической модели задачи дифракции. 77

П.2 Численный метод решения системы уравнений. 88

П.З Результаты математического моделирования. 91

Заключение. 100

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

Введение.

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

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

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

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

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

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

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

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

В данной работе проводится математическое моделирование задач дифракции стороннего электромагнитного поля на следующих двух типах рассеивателей:

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

2. тонком идеально проводящем экране, расположенном как на границе раздела двух полупространств, так и в свободном пространстве.

Задача дифракции электромагнитного поля на диэлектрическом теле относится к числу классических в электродинамике. В настоящее время разработаны различные методы решения задач рассеяния электромагнитного поля на трехмерных прозрачных локальных рассеивателях, расположенных как в свободном пространстве, так и в слоистой среде. Существует обширная литература, относящаяся как к теоретическим вопросам существования и единственности решения задачи дифракции на диэлектрическом теле [6,8,14,23, 47,61,71,77], так и к построению и реализации численных алгоритмов [2426,30,33-43,46,49-52,63,65,68,72,74-76,79,81-89,91-94]. Выбор того или иного метода определяется конкретной задачей. Метод дискретных источников (ДИ) [8,46,49] является весьма эффективным для решения задач рассеяния на осесимметричных рассеивателях. Для решения задач дифракции на трехмерных рассеивателях произвольной формы применяются метод Т-матриц (и его модификации) [38,87,91], метод граничных интегральных уравнений [6,65,84],

метод конечных разностей [34,93], метод конечных элементов с использованием различных типов "реберных" функций [33,39,85], комбинированный метод конечных элементов и граничных интегральных уравнений [40,43], а также асимптотические методы [68].

Кроме перечисленных выше методов, существует подход, основанный на методе объемных интегральных уравнений (ОИУ) [6,23,30,32,37,63]. В этом случае в качестве неизвестного рассматривается полное электрическое поле Е внутри рассеивателя, относительно которого выписывается операторное уравнение [32]:

Е - (ТК + к20 ) \С(ег - 1)ЕА> = Е 0 (1)

Здесь в тензор Грина окружающей среды, гг - относительная диэлектрическая проницаемость рассеивателя, Е°- электрическое поле в отсутствии рассеивателя. Это уравнение также может быть записано относительно объемного тока поляризации 1 = (ег- 1)Е внутри рассеивателя:

—---)\вЫу=Е° (1а)

£ -1 и т/ г V

Если внести в (1) или (1а) операторы дифференцирования под знак интеграла, получится сингулярное интегральное уравнение [63], теоретическое исследование разрешимости которого при различных предположениях о свойствах распределения диэлектрической проницаемости было проведено в работах [23,24,77]. Вопросы численного решения сингулярных интегральных уравнений и, в частности, построение кубатурных формул для вычисления многомерных сингулярных интегралов исследованы в книге [19].

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

исходной краевой задаче объемное интегродифференциальное уравнение (ИДУ) относительно тока поляризации и объемного заряда поляризации [71]. Это уравнение может быть записано в виде:

Щ-- l(VQGVPI + kfa)dvP=E°(Q) (16) er V

Для данного уравнения были доказаны теоремы существования и единственности решения, а также теорема эквивалентности решения ИДУ решению исходной задачи дифракции [71].

В работе [79] для случая двумерного рассеивателя в свободном пространстве впервые получены результаты численного моделирования. Вопросы численного решения ОИУ с применением метода моментов и различных координатных и проекционных функций, а также разные схемы сведения к системе линейных алгебраических уравнений рассматривались в работах [25,26,63,74,76,83,94].

Отметим особенности применения метода ОИУ, которые делают данный подход привлекательным для практических приложений. Использование равномерной сетки при дискретизации объема рассеивателя приводит к хотя и плотно заполненным, но хорошо структурированным матрицам, что позволяет эффективно хранить их в памяти ЭВМ и применять быстрое преобразование Фурье (БПФ) для умножения матрицы на вектор в итерационном методе решения системы линейных уравнений [3]. Объем, занимаемый матрицей в памяти ЭВМ, составляет О(N), где N число неизвестных. Применение БПФ позволяет существенно снизить затраты на решение системы уравнений, т.к. необходимое число арифметических операций в этом случае оценивается как NUa log N где Na - общее число итераций и степень а, как правило, меньше 0,4. Аналогичная оценка для метода конечных элементов с абсорбционными

граничными условиями составляет ]\[1,33[85]. Поскольку представления для

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

Следует отметить, что при численном решении ОИУ апроксимационные свойства оператора, а следовательно, и точность получаемого решения существенно зависят от выбранной численной схемы. Большой интерес для приложений представляют методы, позволяющие получить решение с необходимой точностью при сравнительно небольших вычислительных затратах. Последнее означает, что предпочтительны численные схемы, которые позволяют вычислять матричные элементы, используя кубатурные формулы низкого порядка. Простейший подход основан на методе коллокаций с кусочно постоянной аппроксимацией тока поляризации. Зарубежом также широко применяется метод аппроксимации дискретными диполями (ДДА) [41,51,52,72,75,82], который можно рассматривать как модификацию коллокационной схемы для ОИУ с кусочно постоянной аппроксимацией поля внутри рассеивателя. Эквивалентность этих схем показана в работе [62].

Для рассеивателей, относительные диэлектрические проницаемости которых имеют небольшие значения эти методы дают правильные результаты. Однако, как показали результаты вычислительного эксперимента [52,57], для рассеивателей с большой относительной диэлектрической проницаемостью эти методы имеют очень медленную сходимость приближенного решения к точному, при уменьшении шага сетки. Это объясняется следующими причинами:

1. грубая аппроксимация электромагнитного поля на элементах, содержащих границу раздела сред;

2. кусочно-постоянная апрокснмацпя поля Е не позволяет правильно аппроксимировать йп> Е, которая должна быть непрерывной функцией внутри рассеивателя, если относительная диэлектрическая проницаемость рассеивателя является непрерывно-дифференцируемой функцией;

3. высокая сингулярность тензора Грина [8] приводит к неточному вычислению матричных элементов при кусочно постоянной аппроксимации, когда точки интегрирования и наблюдения совпадают или расположены близко друг к другу.

Следует отметить, что учет квадратичных членов разложения поля в ряд Тейлора наряду с нулевым членом разложения [74], а также более точная аппроксимация электромагнитного поля в окрестности границы [75] позволяют улучшить скорость сходимости как коллокационной схемы, так и метода ДДА для рассеивателей с большой относительной диэлектрической проницаемостью.

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

Идея построения слабого решения объемного ИДУ для рассеивателя, расположенного в свободном пространстве, впервые была предложена в работе [94], однако не была реализована полностью. Использованная авторами статьи аппроксимация векторного потенциала объемными функциями-"крышками" приводит к необходимости решать задачу на области существенно большей, чем занимаемая рассеивателем, что не позволяет применить данный подход для рассеивателей, расположенных на границе раздела сред. Для снижения сингулярности ядра можно использовать непрерывно дифференцируемые функции при аппроксимации объемного тока поляризации и перенести оператор дифференцирования с ядра на ток поляризации [30,37,76,83]. В общем случае эта операция приводит к нагруженному интегральному уравнению [76], что

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

Задача дифракции электромагнитного поля на тонком ограниченном идеально проводящем экране Э принадлежит к числу классических в электродинамике. Наиболее естественный подход к решению этой задачи -сведение ее к векторному интегродифференциальному уравнению на экране относительно неизвестных поверхностных токов у (известное в зарубежной литературе как интегральное уравнение смешанных потенциалов - МР1Е), которое впервые было получено А. Мауэ в 1949 году [64] для экрана, расположенного в свободном пространстве. Это уравнение на плоском экране Б имеет следующий вид для временной зависимости ехр(ггуг):

где Grad и D/v- операторы поверхностного градиента и поверхностной дивергенции, А - интегральный оператор:

Grad А (Div j) + к20А j = -io)s0Eincr (2)

Еы- стороннее электромагнитное поле;£й - диэлектрическая проницаемость свободного пространства; со- циклическая частота падающего поля. Индекс т показывает взятие касательных компонент к поверхности экрана. Для случая плоского экрана 8, расположенного на границе раздела слоистой среды, данное уравнение имеет вид [69]:

—вгай(А1 {Шу/)) + 1соА2) = Е°т (2а) со

где

электромагнитное поле, создаваемое в слоистой среде сторонними источниками, при условии, что экран 8 отсутствует; Ах и А2 - интегральный операторы:

Агту№= Цт^ояоьт*»

А2 т= \\т2(а,р)хр)с1*р

Ядра Т1 и Т2, имеют особенность 1/Я и определяются тензором Грина слоистой среды [7,28].

Центральным моментом при теоретическом исследовании уравнения (2) является правильный выбор пространств для решения и для правой части. Пространства решений и правых частей должны быть достаточно широкими, чтобы соде�