CC BY
17
УДК 517.946
А.С. Ильинский, Т.Н. Галишникова
МЕТОД ИНТЕГРАЛЬНЫХ УРАВНЕНИЙ В ЗАДАЧАХ ДИФРАКЦИИ НА КОНЕЧНОМ ИМПЕДАНСНОМ УЧАСТКЕ ГРАНИЦЫ РАЗДЕЛА СРЕД
(лаборатория вычислительной электродинамики факультета ВМиК,
e-mail: tgalish@cs.msu.su)
Рассматривается трехмерная задача отражения поля плоской электромагнитной волны от локального импедансного участка взволнованной поверхности. Краевая задача для системы уравнений Максвелла в области с нерегулярной границей сведена к решению систем гиперсингулярных интегральных уравнений. Предложен численный алгоритм для их решения. Приведены результаты численных расчетов.
При электромагнитном зондировании поверхности Земли в радиочастотном диапазоне для решения широкого круга прикладных задач в области экологии, диагностики состояния окружающей среды большое значение имеют вопросы математического моделирования и численного решения задач отражения электромагнитных волн от нерегулярной границы раздела двух сред. Эти исследования развиваются по двум направлениям: граница раздела рассматривается как строго периодическая волнистая поверхность, так и как поверхность, нерегулярная на конечном зондируемом участке. Выбор периодической модели отражающей поверхности оправдан в том случае, когда протяженность освещенного участка велика по сравнению с длиной волны падающего поля.
В работах [1, 2] дан обзор существующих методов решения таких задач отражения, при этом чаще всего рассматриваются двумерные математические модели. Исследованию двумерных задач отражения и разработке на основе метода интегральных уравнений численных алгоритмов посвящены и работы авторов настоящей статьи [3-5], в которых для Е- и Н-поляризации приведены результаты вычислительных экспериментов для случаев, когда облучаемая среда является хорошо проводящей или прозрачной, а ее граница есть волнистая периодическая поверхность. В работах [3, 6] в качестве отражающей поверхности рассматривается идеально проводящая плоскость, имеющая достаточно протяженный конечный участок импедансной взволнованной поверхности, и на основе анализа распределения полей на взволнованном участке границы установлено, когда достаточно протяженную, но локальную периодическую неоднородность можно рассматривать как бесконечную.
В настоящей работе рассмотрена трехмерная задача отражения поля плоской волны от конечного участка взволнованной поверхности, проводимость которого достаточно велика. Отражающая поверхность не зависит от координаты z, которая совпадает с направляющей цилиндрической границы раздела. Для построения численного алгоритма, как и для двумерных задач, использован метод интегральных уравнений. Но интегральные уравнения, описывающие волновой процесс в такой трехмерной среде, являются сингулярными. Для их решения в данной работе разработан вычислительный метод, использующий явный алгоритм выделения особенности. В отличие от работы [7] рассмотрен случай не только /•.'-. по и Н-поляризации.
Сформулируем математическую постановку задачи отражения поля плоской трехмерной волны от нерегулярной границы раздела двух сред, которая состоит из идеально проводящей плоскости S, имеющей локальный цилиндрический вдоль оси z импедансный участок S, полностью расположенный в полупространстве у ^ 0. Через D\ обозначим область, лежащую над границей раздела, в которой распространяется падающее поле. Пусть £i, ¿¿i, k\ являются соответственно диэлектрической и магнитной проницаемостью, а также волновым числом среды D\. Здесь к\ = Шу/еГДГ, ш — круговая частота. Зависимость от времени: exp(—iut).
Математически исследование рассматриваемой задачи сводится к решению систем уравнений Максвелла
rot Н (ж, у, z) + гшеiE (ж, у, z) = 0,
rot Е (ж, у, z) — ш(j,iH (ж, у, z) = 0
в полупространстве Di с неоднородной границей раздела сред S |J S, на локальном взволнованном участке поверхности S которой выполняются условия Леонтовича
[п х Е (ж, у, z)} = -IV j [n х [n х Н (ж, у, z)}}, (2)
а на плоском участке Е границы раздела сред условия идеальной проводимости
[п х Е (ж, у, г)] = 0. (3)
Здесь И^ = у//¿г/^г! с2, М2 — характеристики участка 5; п — внешняя нормаль к области с нулевой координатой по оси г. В качестве источника возбуждения берется поле плоской волны
Е0 (ж, у, г) = Е0 ехр (га0ж - г/30у + ¿70г),
(4)
Н0 (ж, у, г) = Н0 ехр (ш0ж - г/30у + г70г),
где «о = вт^о йш^о, А) = кхсоз^озтсро, 70 = ¿^соб^о- Здесь (тг/2 — сро) — угол между осью г и направлением распространения падающей из области 1>1 плоской волны. Если = 7г/2, то соб^о = О, и получаем двумерную задачу отражения. Когда $о = 0, в плоской задаче имеем случай нормального падения. Отраженное поле удовлетворяет условиям излучения.
Учитывая характер падающего поля (4) и тот факт, что геометрия границы раздела сред не зависит от координаты г, решение сформулированной краевой задачи (1)-(4) для системы уравнений Максвелла будем искать в предположении, что полное отраженное поле имеет такую же зависимость от координаты г, что и падающее, т. е.
Е(ж,у,.г) =Е (ж, у) ехр (¿7о г), Н (ж, у, ¿) = Н (ж, у) ехр . (5)
Обозначим через и(х,у) = Ег(х,у), и(ж, у) = Нг(х,у) проекции искомого электромагнитного поля на ось г. Тогда рассматриваемая краевая задача для системы уравнений Максвелла (1)-(5) сводится к решению однородных уравнений Гельмгольца в области 1>1
(<Эж2 + ду2 (г>(ж, у)) ^
Граничные условия Леонтовича (2) на поверхности Б примут вид [1]
(7)
, ч Щ ( . ду(х,у) . ди(х, у)
и(х,у) = -щ. 1Ш£\ —
1 / . д'о(х, у) . ди(х,у)
где 01 = к\ — 7д, а условия (3) на поверхности I! запишутся следующим образом:
Цж,у) = 0, = (8)
В предыдущих формулах использованы обозначения: д/дп = пхд/дж + пуд/ду — производная по нормали; д/дт = пхд/ду — пуд/дж — производная по касательной; пх, пу — направляющие косинусы нормали п.
Для получения систем интегральных уравнений, эквивалентных краевым задачам (6)—(8), введем функции Грина дЕ,н (М, Р), удовлетворяющие в области у ^ 0 неоднородному уравнению Гельмгольца
Э2 д2 2\ /дЕ(М,РУ
а* ■ * Р') \(1н{м,р)) =-2*г>(м<р)< V»
граничным условиям идеального металла на плоскости у = О
Лм,Р) = о, = 0 (10)
и имеющие вид
„В/
Здесь и в дальнейшем в приводимых формулах верхний знак будет соответствовать ^-поляризованному полю, а нижний — Деполяризованному. В формуле (11) использованы обозначения: Нд1-5 (ж) — функция Ханкеля нулевого порядка первого рода; гм,р = у/ Аж2 + Ау2; гм,р = л/Аж2 + А«2; Аж = ж м — жр;
Ду = Ум — Ур; $У = Ум + Ур; (хм,р',Ум,р) — координаты точек М и Р соответственно. Далее применим в области формулы Грина и учтем свойства поверхностных потенциалов, в результате получим две независимые системы из двух интегральных уравнений [3, 7] относительно неизвестных на конечном участке Б функций 9, в случае ^-поляризации:
Иц\У2 ди(М) ^оРУз ду(М) 2ШФ1 дпм + дтм + 2тт
ддЕ(М,Р) Эу(Р)
д (М, Р)
Иц\У2 дд (М, Р)\ Эи(Р)
И'1 ?
дп}
дп}
(3\
1 ди(М) _
2 дпм 2ж
ёяр = (М,), М € Б,
дпр дтр \
ЭдЕ(М, Р) ък^УУ,2 д2дЕ(М, Р)\ &и(Р)
<Эп
А/
И7!/?! дпмдпр ) дпр
¿70РУ2 д2дЕ(М,Р) ду(Р) дпмдпр дтр
й.3р =
диЕ(М)
дп
м
и функций у(Р), в случае поляризации:
дт
кгШг
■ди(М,Р)
кгШг
ди(Р)
дп}
М е 5,
г>(Р) +
9гр 7о д«(М)
йвР = Уо(М), Ме5,
___
7о ддн(М,Р)ди(Р)
ЬШ1
<Эп
А/
(Этг
1 г \iiW2Pl
2ж ] я
й.3р = дуЛ(м)
дпм
/11. (М, Р) д2дн (М, Р)
дпм Ме5,
дпмдпр
у(Р) +
где = — импеданс свободного пространства и
Щ {Х1 у)
= ехр (га0ж — «Аш) Т ехр (га^х + г/Зоу).
(12)
(13)
(14)
(15)
(16)
Если падающее поле не зависит от координаты г, то 70 = 0, и вместо двух систем интегральных уравнений (12), (13) и (14), (15) получаем два независимых интегральных уравнения Фредгольма второго рода по конечному импедансному взволнованному участку границы раздела сред, а именно
ддЕ(М,Р)'
1 ди(М)
1
Ъх
дь(М,Р) + аЕ-
дпт
ди(Р) дпг
¿8Р = -и%(М), М € 5,
(17)
в случае ^-поляризации, схе =
ИУ2 к 1^1 '
И
1
у(М)
1
2ж
дд (М, Р) дпр
андн(М,Р)
у(Р)йзР = у^(М), МеБ,
(18)
в случае Деполяризации, «я = Результаты вычислительных экспериментов для двумерных
математических моделей, описываемых соотношениями (16)—(18), приведены в работах [2, 3, 6].
Численные алгоритмы для решения систем интегральных уравнений (12), (13) и (14), (15) основаны на сведении их к системам линейных алгебраических уравнений. В связи с этим возникает вопрос о вычислении интегралов, входящих в расчетные формулы для получения элементов матриц.
Пусть у = /(ж) — уравнение цилиндрической неоднородной поверхности Б, 0 ^ ж ^ а. Разобьем отрезок [0, а] точками хц на N равных частей: жо = 0, ждг = а, Хг = га/Ы, г = О,..., N. Тогда точками разбиения кривой Б будут точки с координатами (ж*, / (ж*)). Интегралы, входящие в системы уравнений (12), (13) и (14), (15), заменим суммами интегралов по отрезкам разбиения и положим
М = i = 0,...,-/V — 1. В результате получим системы линейных алгебраических уравнений. При вычислении тех интегралов на отрезке разбиения [$р ; , //' 0...../V ^ I. подынтегральная
функция в которых зависит либо от дЕ,н (М,Р), либо от 9д , переходим от интегрирования
по я к интегрированию по ж, заменив йвр = + (/'(ж))2с^жр, а в тех интегралах, которые зависят
д2дЕ'н(М,Р)
от —дпмдпр ' сохраним интегрирование по дуге 5.
Заметим, что слагаемые в ядрах интегральных уравнений (12)—(15), зависящие от дЕ,н (М,Р) (см. (9)—(11)), имеют логарифмическую особенность при совпадении аргументов, поэтому при построении численного алгоритма она выделяется в явном виде, и интегралы, содержащие особенность,
а Е.Н /^ р\
вычисляются аналитически. В слагаемых, зависящих от Р ' ' п°ДьштегРальная функция при
совпадении аргументов доопределяется через значение кривизны контура в этих точках. Эти вопросы рассмотрены в работе [2] при решении двумерных задач.
При математическом моделировании трехмерных задач отражения системы интегральных уравнений (12), (13) и (14), (15) имеют ядра, зависящие от —|гамэПр ; поэтому полученные интегральные уравнения являются сингулярными. Следовательно, для вычисления элементов матриц, являющихся интегралами от функций с гиперсингулярной особенностью, необходимо построить численный алгоритм. Для этого внешнюю нормаль к области 1>1 запишем в виде п = —у'а1 + ха}, где — орты осей ж и у, а для поверхности 5 используем параметрическую форму ее задания: ж = ж(«), у = у (я). В отличие от результатов работы [8] точки М и Р могут совпадать. Заметим, что в случае задания поверхности Б в виде у = /(ж) связь между направляющими косинусами пх, пу нормали п есть
у'а = + (/'(ж))2, х'а = пу^1 + и>{х))\ где пх = /'(ж)/^1 + (Г(х))2, пу = -1/^1 + (Г(х))2.
Далее можно показать [7, 8], что
/3? I нЫ^гм^р^у'^р + и I 6(М,Р)у'арй8р. (19)
3рз-рз+1
В формуле (19) интеграл, содержащий полный дифференциал, вычисляется в явном виде, а в интеграле, содержащем функцию Нд1-5 {¡3\Гм- х ,р), логарифмическая особенность выделяется в явном
виде.
Аналогично можно показать, что
+ Р\ I Н™ (Рггм^^х'^йзр-Ы I 6(М,Р)х'арй8Р. (20) Следовательно, учитывая (19), (20) и тот факт,
что при вычислении —дпмдпр —~ точки М и Р не совпадают, получаем численный алгоритм для вычисления вторых нормальных производных от
функций Грина д ' (М, Р), а именно
Эргрз +1
н!(1) (&ГМ р, + 1)
, ,2/А/ 1 - УР] + 1 ) У'ам 1 + (хм I - ХР]+1 ) х'ам г(1) о „ \ 4 2_' ^ 4 2_' ^
3 + 1
-Я
(ум1+, ~ Ур) y'sMi+h + К, - xPj) x'SMi i [PirM Pj)-
rMiMl,Pj
H
(1)
i rM.+ i,pj+1
(:Умг+h +ypj+1)y'SM + (*Mt+i-xPj+1)x'
rMi+i,PJ + 1
-H
i(1)
{yMi+1 + УР,) y'sM + (xMi+i - XPj) X> -
+
Ш n2
Sp3'p3 +1
я,
1 Гм.
SM i sP i+b
ТН{01] (PirM.
<+i'P) (ysMi+iysp+<
+ 2ТГ I S (M, P) (21)
5pJ'pJ+l
Заметим, что наличие последнего слагаемого в (21), содержащего дельта-функцию, при совпадении точек М и Р дает выброс при счете соответствующих диагональных элементов матрицы.
Разработанный алгоритм был реализован для решения системы интегральных уравнений (12), (13) в случае Я-поляризации. Пусть уравнение контура S в плоскости XOY задается уравнением f(x) = 1 — cos ж, 0 ^ х ^ 2тт. Расчеты проводились при N = 20 для следующих электродинамических параметров среды: к± = 0,5; £2 = 50 +4г. На рис. 1 при фиксированном значении угла = 0 показано распределение функции
д u(xj(x))
1{х) =
дп
на импедансной границе 5 при изменении угла (ро. Кривые 1-5 соответствуют (ро = 85; 65; 45; 25; 5°. Кривая 1, полученная при (ро = 85°, является тестовой и близка к результатам расчета двумерной задачи дифракции поля плоской волны в случае ее нормального падения, которые приведены в работах [2, 7]. Как и следовало ожидать, при (ро = 5° (кривая 5), когда угол между образующей цилиндрической поверхности 5 и падающей волной мал, величина наводимого тока 1(х) на импедансном участке 5 значительно меньше по сравнению с результатами, соответствующими кривым 1~4-
Рис. 1
Рис. 2
1(х)
__I_,_I_,_I_,_I_,_V.......:........г1_,
0 1 2 3 4 5 6 х
Рис. 3
На рис. 2, в отличие от результатов, представленных на рис. 1, показано распределение функции 1{х) на границе S при изменении угла <ро5 но значение угла = 45°. Кривые 1—5 также соответствуют <Ро = 85; 65; 45; 25; 5°. Токи, как и следовало ожидать, уже не являются симметричными.
Результаты расчетов, представленные на рис. 3, получены при фиксированном значении сро = 45°, а значение угла менялось. Кривые 1—5 соответствуют = 85; 65; 45; 25; 5°. Приведенные на рис. 13 результаты вычислительных экспериментов показывают возможности разработанного численного алгоритма для решения задач отражения, когда падающее поле имеет трехмерный характер.
СПИСОК ЛИТЕРАТУРЫ
1. Ильинский A.C., Галишникова Т. Н. Численные методы в задачах дифракции. М.: Изд-во МГУ, 1987.
2. Ильинский A.C., Галишникова Т.Н. Рассеяние плоской волны на волнистой поверхности // Математические модели естествознания. М.: Изд-во МГУ, 1999. С. 86-111.
3. Ильинский A.C., Галишникова Т.Н. Математическое моделирование процесса отражения плоской электромагнитной волны от волнистой поверхности // Радиотехника и электроника. 1999. 44. № 7. С. 773786.
4. Ильинский A.C., Галишникова Т.Н., Бережная И. В. Характеристики отражения от волнистой границы раздела прозрачных сред // Прикладная математика и информатика. № 3. М.: Изд. отдел ф-та ВМиК МГУ, 1999. С. 33-42.
Б.Ильинский A.C., Галишникова Т.Н., Бережная И.В. Исследование характеристик отражения электромагнитного поля от волнистой периодической границы раздела прозрачных сред / / Прикладная математика и информатика. № 6. М.: Изд. отдел ф-та ВМиК МГУ, 2000. С. 5-11.
6. Ильинский А. С., Галишникова Т.Н., Бережная И. В. Сравнение двух математических моделей в задаче дифракции ii-поляризованной волны на нерегулярной границе раздела сред // Вестн. Моск. ун-та. Сер. 15. Вычисл. матем. и киберн. 2001. № 1. С. 8-13.
7. Ильинский A.C., Галишникова Т.Н. Рассеяние электромагнитного поля на конечном импедансном участке границы раздела сред // Прикладная математика и информатика. № 25. М.: Изд. отдел ф-та ВМиК МГУ, 2007. С. 56-69.
8. Вайникко Г.М., Лифанов И. К., Полтавский Л.Н. Численные методы в гиперсингулярных интегральных уравнениях и их приложения. М.: Янус-К, 2001.
Поступила в редакцию 21.01.08
