УДК 517.9
Ю. Г. Смирнов, М. Ю. Медведик, Д. И. Васюнин
МЕТОД КОЛЛОКАЦИИ РЕШЕНИЯ ОБЪЕМНОГО СИНГУЛЯРНОГО ИНТЕГРАЛЬНОГО УРАВНЕНИЯ В ЗАДАЧЕ ОПРЕДЕЛЕНИЯ ДИЭЛЕКТРИЧЕСКОЙ
прониЦаемости материала
Аннотация. Рассмотрена задача дифракции электромагнитного поля на диэлектрическом теле, расположенном в прямоугольном волноводе. Задача сведена к объемному сингулярному интегральному уравнению на теле. Рассмотрен численный метод коллокации для решения этого уравнения. Представлены расчетные формулы для матричных коэффициентов метода коллокации.
Ключевые слова: краевая задача, электромагнитная задача дифракции, интегральное уравнение, численный метод.
Abstract. Electromagnetic diffraction problem on dielectric body located in rectangular waveguide is considered. The problem is reduced to volume singular integral equation on the body. Numerical collocation method for solving the equation is considered. The formulas of matrix coefficients for collocation method are presented. Keywords: boundary value problem, electromagnetic scattering, integral equations, numerical method.
Определение диэлектрических и магнитных параметров нанокомпозит-ных материалов и сложных наноструктур с различной геометрией является актуальной задачей нанотехнологии и наноэлектроники. Однако эти параметры, как правило, недоступны для экспериментального измерения (ввиду композитного характера материалов и малых размеров образцов), что приводит к необходимости применять методы математического моделирования и решать задачи численно с помощью компьютеров. При этом приходится решать трехмерные векторные задачи в полной электродинамической постановке. Решение таких задач является в настоящее время одной из самых актуальных проблем в электродинамике и с приемлемой для практики точностью на электродинамическом уровне строгости математическими методами требует очень большого объема вычислений, что зачастую невозможно даже на самых современных суперкомпьютерах. Особенно остро стоит проблема решения обратных электродинамических задач на сложной системе поверхностей и тел в резонансном диапазоне частот, возникающая при определении параметров нанокомпозитных материалов и наноструктур. Многочисленные дорогостоящие пакеты прикладных программ для решения задач электродинамики (Ansis, Quikwave и т.д.), имеющиеся на рынке программных продуктов, решают задачу традиционными конечно-разностными методами или методами конечных элементов и не дают удовлетворительных по точности результатов [1].
Альтернативным подходом является применение метода объемных сингулярных интегральных уравнений [2-4]. Настоящая статья посвящена разработке численного метода для решения уравнения. Применяется метод коллокации с аналитическим суммированием медленно сходящихся рядов в функциях Грина.
1 Краевая задача дифракции для системы уравнений Максвелла
Рассмотрим следующую задачу дифракции. Пусть в декартовой системе координат Р = {х: 0 < xi < a, 0 < Х2 < b, — го < X3 < го} - волновод с идеально проводящей поверхностью ЭР . В волноводе расположено объемное тело Q (Q с Р - область), характеризующееся постоянной магнитной ^оницаемос-тью Цо и положительной 3 х 3-матрицей-функцией (тензором) диэлектрической проницаемости е(х). Компоненты е(х) являются ограниченными функциями в области Q , e (Q), а также e 1 е L(Q).
Граница ЭQ области Q кусочно-гладкая. Точнее, предположим, что
3
для каждой точки границы хо е ЭQ существует окрестность 0 (в R ) и
2 3
C -диффеоморфизм этой окрестности на R , при котором точка хо переходит в точку 0, а образом множества 0nQ является множество одного из
следующих типов (ниже (х1, х2, хз) - декартовы, (г, 0), r > 0, 0 е S - сфери-
3
ческие координаты в R ). Либо х1 > 0 (х0 - точка гладкости границы), либо
3
х1 > 0, х2 > 0 (х0 - точка на «выходящем» ребре), либо R \ {х1 > 0, х2 > 0} (х0 - точка на «входящем» ребре), либо r > 0, 0 е Q', где Q' с S - односвязная область с кусочно-гладкой границей ЭQ' (х0 - вершина «конуса с ребрами»). В частности, если ЭQ' - гладкая, то х0 - коническая точка; если ЭQ' образована дугами больших окружностей, то х0 - вершина многогранного угла. Пусть Q - ограниченная область и каждая точка х е ЭQ принадлежит одному из этих типов. Тогда будем говорить, что Q - область с кусочногладкой границей. Будем также предполагать, что тело Q не касается стенок волновода, ЭQ пЭР = 0 . В Р \ Q среда изотропна и однородна с постоянными &0 (> 0), М> 0).
Требуется определить электромагнитное поле E, H е L2 loc (Р), возбуждаемое в волноводе сторонним полем с временной зависимостью вида e imt. Источник стороннего поля - электрический ток j0 е L2 loc (Р). В области
3
Р с R стандартные дифференциальные операторы grad, div, rot понимаются в смысле обобщенных функций.
Будем искать «слабые» (обобщенные) решения системы уравнений Максвелла:
rot H = —/юеЕ + jE;
rot Е = /ю^Н. (1)
Эти решения должны удовлетворять условиям на бесконечности: поля Е и Н при |хз| > C для достаточно больших C > 0 имеют представление («+» соответствует +го , «-» соответствует —го )
f E' H
=Z R
Z Q.
p±)eiY
X (і)П peз — iY £4^ —/ЮЄо ( V2Пp )x єз
ГОЦ° ( V2¥p )x єз
^X(pVpeз — iY(p)V2Tp j
(2)
где
Ур" ) = ^ $ “X Р ) , Іт Ур ) > О или Іт у ( ) = О, £ур' ) ^ О и X ^, П р( Ху, Х2) и X?) , *¥р (Ху, Х2) ( &о =ю2ЄоЦо І - полная система собственных значений и ортонормированных в Ь^( П) собственных функций двумерного оператора Лапласа -А в прямоугольнике П := {(х1, Х2): О < Ху < а, О < Х2 < Ь} с условиями Дирихле и Неймана соответственно, и У2 = еу д/дху + Є2 д/дх2 . Для коэф-
фициентов разложений (2) имеют место оценки
Rp
Iq(+) = О(pm )p
(з)
для всех т е N.
С физической точки зрения условия (2) означают, что рассеянное поле является суперпозицией нормальных волн, расходящихся от тела. Условия (3) обеспечивают экспоненциальную сходимость рядов (2), а также возможность их почленного дифференцирования по х) любое число раз.
Для Е, Н должны выполняться краевые условия на стенках волновода:
(4)
Если выполняются уравнения Максвелла, то второе условие в (4) следует из первого, и его можно опустить. Но если рассматривать оператор Максвелла, порождаемый левой частью (1), то надо ставить оба условия.
Для и е H 1(Р) существуют граничные значения из пространства
1/2
H (ЭР) в смысле теории следов. Почти везде на ЭР определен вектор нормали.
Пусть также Е0 и Н0 - решения рассматриваемой краевой задачи в отсутствие неоднородного тела Q , £(х) = £0I, хе Р (I - единичный тензор):
rot Н0 = — 7Ю£0Е° + jE, rot Е0 = /Ю|10Н0 (5)
с краевыми условиями
0 0 (6)
E° ІЭP = О, HV ІЭP = 0.
Эти решения могут быть выражены аналитически через с помощью введенного ниже тензора Грина. Решения не обязаны удовлетворять условиям на бесконечности. Например, ЕО и НО могут быть ТМ- или ТЕ-модой этого волновода.
7з
Имеют место результаты о гладкости решений задач (1)-(4) и (5)-(6) при более гладких данных [3]. Сформулируем один из таких результатов.
Утверждение 1. Пусть ^ е ННос(Р). Тогда Е0, Н0 е НЬос(Р). Пусть,
2 1 — 1
кроме того, дОе С , ее С (О). Тогда сужения Е |д, Н |де Н (0 и Е \р\д,
Н |р\де Н 11ос(Р \ О). Кроме того, справедливы условия сопряжения на дО:
[Ет]|эе = 0, [Нт]|эе = 0,
где [ • ] означает разность следов с разных сторон Э<2.
В предположениях утверждения 1 краевые условия на дР и условия сопряжения на д<2 понимаются в смысле равенства следов элементов из 1/2 1/2
Н 1ос(дР) и Н (д<2). Ясно, что при первоначальных общих предположениях о тензоре е такие условия сопряжения не имеют смысла.
2 Тензорная функция Грина прямоугольного волновода
Построим диагональный тензор Грина О^ , компоненты которого являются фундаментальными решениями уравнения Гельмгольца в Р с коэффициентом к0 =ю2 £оМю и удовлетворяют краевым условиям первого или второго рода на дР , обеспечивающим обращение в нуль тангенциальных составляющих напряженности электрического поля на стенках волновода. Его компоненты имеют вид [3]
_ ^ ^ —у Хо — Уз
^ 2 ^ _ е '«т1 3 -Л31 пп . пт пп . пт
О1 =—Т X X-л , О ч со8~ Х1 вш — х2 со^ У1 вШ —У2 ; (7)
аЬП=0 т=1 У пт(1 + 80п) а Ь а Ь
_ ^ ^ —у Хо — У
„ 2 ^ ^ е пт 3 3 . пп пт . пп пт
О2 =—Т X X------ Т81П х1со^ — х2 81П У1со^—у (8)
аЬп=1 т=0 У пт(1 +80т) а Ь а Ь
~ —у Хо — Уо
„ 2 ^ е пт 3 3 . пп . пт . пп . пт
О3 = — X X-----------------8т—х^ш—Х2Й1П—ляп—У2 . (9)
аЬ^ “ упт а Ь а Ь
п=1 т=1 пгп
В этих выражениях упт =
пп ^2 (пт ^2 2
— I +1---I — к0 , при этом ветвь квад-
а ) I Ь )
ратного корня выбирается так, чтобы 1т упт ^ 0 .
Запишем От с выделенной особенностью при х = у :
1 ек\х—У|
От =~-Т + gm (Х, У), Х, у е Р , (10)
4п | х — у |
где функция gm е С“ (О X Р).
Отсюда и в силу симметрии функций Грина От (х, у) = От (у, х), (т = 1, 2, 3) имеем следующее утверждение.
Утверждение 2. Тензор Грина Ge допускает представление
е 1 eik0 Iх-У1 е
Ge =-л—;------г I + g (х, У), х, у е Р, (11)
4я | х — у |
где матрица-функция (тензор) g е C(Q х Р) и g е C(Р х Q).
Такое представление функции Грина удобно для теоретического исследования задачи дифракции, но непригодно для численных расчетов, т.к. не содержит алгоритма вычисления g . Отметим, что функции Грина имеют
1 eik0| х—У1
единственную особенность вида---------------и не имеют других особенностей
4гс I х — у I
в силу сделанного нами предположения о том, что тело не касается поверхности волновода.
3 Объемное сингулярное интегральное уравнение
Наша ближайшая цель - свести краевую задачу к объемному сингулярному интегральному уравнению и доказать теорему эквивалентности.
Пусть решения краевых задач (1)-(4) и (5), (6) существуют и единственны. Перепишем (1) в эквивалентной форме:
rot Н = — /ю^Е + jE, го!Е = /ЮЦ0Н, (12)
где
jE = jE + jE . (13)
В последнем равенстве jE = -/ю(£(х) — £01)Е - электрический ток поляризации.
Нетрудно проверить, что решение краевой задачи (4), (12) имеет вид
Е = /ЮЦ0 A e ———grad div A e , Н = rot Ae , (14)
/ЮЁ0
где
Ae = {Ge (r)jE (y)dy - (15)
Р
векторный потенциал электрического тока.
Потенциал Ae удовлетворяет уравнению
AAe + k0 AE =— jE . (16)
Таким образом, потенциал Ae есть свертка с тензором Грина прямоугольного волновода для уравнения Гельмгольца, обеспечивающий выполнение требуемых краевых условий для полей.
Однако формулы (14) не дают явного решения задачи (4), (12), т.к. ток jE зависит от Е. Из соотношений (13)-(15) для поля Е следует интегро-дифференциальное уравнение
Е( х) = Е0( х) + к21 Ое (г) <2
ё( У)
-1
ё0
Е(у)йу +
- grad div | Ое (г)
2
ё( У)
ё0
-1
Е(у)<яУ, х є 2.
(17)
Кроме того,
Е(х) = Е0(х) + ко21 Ое (г)
2
ё( У)
ё0
-1
Е(у)йу +
- grad div | Ое (г)
2
ё( У)
ё0
-1
Е(у)оУ, х є Р \ 2.
(18)
Формула (18) дает представление решения Е(х) в области Р \ Q , если Е(у), у е Q - решение уравнения (17). Поле Н выражается через решение (17) в виде
Н( х) = Н 0( х) - /ЮЄ0 тої1 Ое (г)
2
ё( У)
ё0
-1
Е(у)оУ, х є Р.
(19)
Сведем полученное выше интегродифференциальное уравнение к объемному векторному сингулярному интегральному уравнению.
Представим функцию Грина в виде
Ое (г) = О0 (г) + 000 (г) + ^2 (г), г =| х - у |;
(20)
О0(г ) =
ік0г Л
Є 0 -1
. 1, °00(г) = т1--1, ё0(г) = ё^2, &0). (21)
4кг 4кг
Применяя теорему о дифференцировании интеграла с ядром, имеющим слабую особенность, придем к известному представлению:
~~ ----- ип (у)ёу = v.p. Г
дх J дг 4,пг J
дх1 * дхп 4пг
і 2 п
дхідхп 4пг
1 ип (у^у -15^ (х). (22)
Используя полученные соотношения, переходим от интегродифферен-циального уравнения (16) к векторному сингулярному интегральному уравнению:
Е( х) + з
1 - ( 1 Е( х) - -^р.Г Г1(х, у) " ё( у) / Е( У^У
1 ё0 ] •! 2 _ ё0 _
Г( х, у)
ё (У)
ё0
-1
Е(У^у - Г ?2(х, у)
2
ё( У)
ё0
-1
Е( у ^у = Е0( х). (23)
Здесь тензоры Г, Г1, Г2 имеют вид
f( x, y) = klGE (r) + ( - , grad) grad Go(r); fi( x, y) = ( - , grad) grad Goo(r);
(24)
(25)
(2б)
Часто интерес представляют задачи рассеяния в среде, характеризующейся постоянной во всем объеме волновода диэлектрической проницаемостью (ё = ё01) и тензорной магнитной проницаемостью Д в 2 (вне 2 Д = Д01). В этом случае краевая задача сводится к объемному сингулярному интегральному уравнению (такого же типа) для магнитного поля и выражению для электрического поля через решение этого уравнения:
В последних формулах Оц (х, у) - тензорная функция Грина прямоугольного волновода, отвечающая произвольному распределению источников магнитного поля. Как и для рассматривавшейся функции Грина Ое (х, у), имеет место представление в виде суммы сингулярного слагаемого того же вида и гладкой функции.
4 Представление функций Грина с суммированием медленно сходящихся рядов
Рассмотрим выражение
J (GzU + graddiv(GEU))y,
Q
Здесь каждая компонента тензора представляется в следующем виде:
“ “ e~1nm\x3 _Уз\
g = от УУ
ab n=im=0 Ynm (i + 50n )
sin (nXi )cos (mX2 )sin (nYi )cos (»Y);
g3 =— УУ
3 ab
-sin (nXi) sin (mX2) sin (nYi) sin (mY2).
' n=i m=i Y nm
Проинтегрировав компоненты тензора Грина по параллелепипеду
Pji2i3 = j(xi, x2, x3): ii < h < ii + i, i2 < -h2 < i2 + і, із < h. < із + і j, будем иметь
(сохраним прежние обозначения для проинтегрированных компонент тензора Грина)
0 ^ ^ /*0 / \
Gi = — У У пн ( 3) cos (nXi )sin (mX2 )cos (nH“ (( + 0,5)) х
— n =im=i Ynmnm
Xsin
THHljsin(mH2 (2 + 0,5))sin)
mH2
2Hi ^ /o0m (x3) , ( X , , ( ( | 05) , (m2
—2“У—^—^sin(mX2)sin[mH2(.2 + 0,5))smI —^ I;
— m=i YQmm
S ^ ^ j0 (x )
G2=4, уу^нш
sin (nXi) cos (mX2) sin (nHi (ii + 0,5)) х
— n =im=i Ynmnm
х sin
nHi
cos
(mH2 (2 + 0,5)) sin
mH2
^ у /n02 (x3 ) sin (nXi )sin ^nHi (.і + 0,5))^п І;
— n=i Yn0m V 2 J
о ~ ~ J0 (x \
Jnm\x3) •
G3=T УУ
— n=i m=i Ynmnm
sin
(nXi )sin(mX2 )sin(nHi (ii + 0,5))х
х sin
tHl j sin (mH2 (2 + 0,5)) sin ^
mH2
Здесь
-xi
-x2
ЛУі
Xl = —К X2 =-гЧ Yl = -^, Y2 = b
-У2
—h
Hi = ^,H2 =
, > i“------.я2
b a
=—h2 ,
aba
xi = Jihi, x2 = J2h2, Уі = kh, У2 = i2h2 ;
b
Y nm
Jnm(x3) = <
(-exP (-(x3 -(z3 + *) ( Ynm )- exP (-(x3 - i3h3 )Y nm )) x3 > (( +
(exp(-(/'3h3 - (Y nm ) exP(-((( + !)h3 - x3 )Ynm )) x3 < i3h3,
(2 - exp (-(x3 - /3Л3 ) Y nm) - exP (-((z3 + !)h3 - x3 )y nm ))
i3h3 < x3 <(( + !)h3-Введем обозначения:
Лч sin nx cos ny
r (x, y; A) = y- '
=L n(n2 УА2)
• = p(x, y;А)y q(x, y;А) (0 < x, y < —);
p(x, y; А) =
—
4А2 (i _ є"
-2—A
XI e-M2*-x-y) -e-^(x+y)+ sign(x-y)g-AMx-yl) -e-Ax-yl Yl-q (x, y; A) =-^y (71- x - y + sign (x - y )(-| x - y|));
4А
p2( x, y; А) = p( X, y; іА) =
_—
л
4А sin —A
х(іп (A(x У y _—))y sign (x _ y )sin (|x _ y| _—))); pi(x, y; А) = p(x, y; А);
_ sin nx sin ny _
y( X, y; А) = У-
n=1
— e-Al x_y| y Є-А(2—_І x_y) _ e~A(xyy ) _ e-A(2—_x_y)
-2—A
,(x> 0,y >0, x У y <2—);
4А i _ Є
S2 (x, y; А) = s( x, y; іА) = 4а (cos (— _ x _ У))_ cos (A(—_ |x _ y
si( x, y; A) = s( x, y; A);
d (x, y; А) = ds = _ —
х
dx 4 і _ Є 2—A
_e-A(xyy) y e-A(2—_x_y)
,_AIx_yl _e“A(2—_lx_У^sign(x_y)_
(x > 0, y > 0, x У y < 2—);
d2( x, y; A) = d (x, y; iA) =
“ a (sin (A(—_x _ y ))_sign (x _ y) sin (A(—_| x _ y|
n
i
/иш(х3)
М х, у; А) = й (х, у; А);
(-ехР(-(х3 _(73 + 1)к3 )иш )-ехР(-(х3 -/3%3 )иш )) х3 >(( + 1)h3, (ехр(-(/3Й3 - х3 (иш )-ехР(-((73 + 1)И3 - х3 )иш )), х3 < /3%3-’
(-ехР(-(х3 - /3%3 )иш )- ехР(-((73 + 1)И3 - х3 )иш )),
/3%3 < х3 <(/3 + 1)3. В точке коллокации
1
Qj1 }2}Ъ - \ х1 -| Л + Т | h1, х2 -| 72 + Т | h2, х3 -| 73 + Т I %3
1
1
значения компонент тензора Грина после суммирования медленно сходящихся рядов будут иметь вид (снова сохраним прежние обозначения для проинтегрированных компонент тензора Грина в точке коллокации)
го го 1иш | %3 | 73 +
°1 -4II-
К и-1 ш-1 У ишиш
Г 1 I ^ Г 1 I ^ • ин1
008 и I ?1 + — IН1008 и I 71 + — IН1 эт—^ х
хэтш| /'2 +— |Н2э1пш| 72 + — |Н2э1п шН2
2 н /ош ( I 73 + 2 | %3 I А 1 А А 1
+ —211--------^1 япш| /2 + - |ш| 72 + - |
л
ш-1
УОшш
. тН2 Н
2 + 8, ,■ +
/3 73
I 1 сое и I/1 + 2 IН1008 и I 71 + 1 IН1 х
К 1 и
и-1
. иН1 х эт—1
72 + 2 IН2,/2Н2,Аи
72 + ^ |Н 2, (/2 + 1)Н2,
+8,-
26 2 Н1
/3Л \ 4
К
1 кЬ 1 кЬ
Р2(Н2(72 + 2), Н2/2; —) - Р2(Н2(72 + 2), Н2О2 + 1); —)
2 к 2 к
Н2 + 72 - /2 + 2
72 - /2 +1
вщп 72 - /2
./2 - /2
х
X <1
К
р21 Н1 (/1 +l),Н11 71 + ^;— |-р2 | Н1/ЪН1171 + ;—
4к 2 л2
Н1 + зі§п| к - /! - ^
(
л- Нл
к - /і- 2
кі -/і +-
1
(
л- Н
І1 - /1 +■
1
Ні
2л2к2.
./ига І И3 I 73 ^
л п=1 га=1 У птпга
■ (. 1 Л ^ • Г 1 Л ^ • пн1
ЭШ ПI /1 + — IН ЭШ ПI к1 + — IН1 ЭШ ^ X
і ■ 1Л и ( ■ 1 • тН2
X соэ га| /2 + ^ IН2С0Э га I 72 + ^ IН28т—+
+ 2Я2£5ц1 пГ/1 + 2 |н15шпІл +2]н151пПА +
л2 п=1
. тН2 X ЭШ----2
п2 о а
2 <3 ^ 00 £■
л га=1
((
П71 2
2а 2 Н 2
1зк\ 4
л
Р2\Ні \ А+21,Н1/1; л
ка Л
к + | ] Ні, (/1 + і)Ні, А га Н1 ( 71 + 2 ^, Н1 (/1 +1);—
Ні + ві§п | а - /і + 2
71 - /і + 2
-8Щп| кі -/1 -2
(
л- Ні
л
Р2
Н2 (/2 + l),Н2 Г72 + ‘2^;“
Р2
Н2/2, Н2 І 72 +
і Л кь
4к 2 л2
Н2 + ^§п| 72 - /2 - 2
72 - /2 - 2
ві§п І 72 - /2 +■
1
(
72 - /2 + ■
Н2
2л2к 2
О = "г !!■
/ига І й3 І 73 ^
л и=1 га=1 у птпга
■ (■ 1 Л ^ • ( ■ 1 Л ^ • пН1
эт п I /і + ^ I Ні эт п I 7і + -2! Ні эпп 2 X
I . 1Л . mH2
Ц J2 + ^ IН2"'”
X sin ml /'2 + — | H2sin ml j + — | H 2sin—+
+ 8;
8b
f
л4 -1 n
^ n=1
1
2 sinn I ;i + 11 Hisin n 1 ji + 11 Hisin nHHl x
H2 I J2 + ^J,H2/2;An -p H2 |J2 + 2 I’H2 (/2 + i);
1
я a2
+ /з J3 _4 л
H2 + sign| J2 - ;2 + 2
J2 - /2+2
\
sign| J2 - ;2
J2 - ;2
X
P2
ka
\ A
Hi| J1 + “ | , H1/1; -p2 Hi| J1 +— | , H1 (/1 + 1);
ka
л
2
4k 2 a 2
H1 + sign| Ji - ;1 +
1
л- H
J1 - ;1 +■
1
-sign| Ji -;1 -2 Il л-H1
J1 - * - 2
/J
Вычислим последовательно операции дивергенции и градиента: graddiv(G£t7) = graddiv(GiUi, 2, G3U 3) =
= grad
■Ui +—2 U2 +
Эх1 d%2 6x3
U
+
32G. Ui +^^^ U2 +^ u3
\
3xi
ЭХ1ЭХ2
ЭХ1ЭХ3
xi +
ЭХ2ЭХ1
3Х22
ЭХ2ЭХ3
Х2 +
32<G Ui +^2. u2 +33'
ЭХ3ЭХ1
ЭХ3ЭХ2
Х3.
Продифференцируем функции Gi, G2, G3 и вычислим вторые производные
32Gl 32Gl 32Gl 32G2 32G2 32G2 32G3 32G3
3Х2 ’ Эх2Эх^ ЭХ30Х^ Эх10Х^ ЭХ2 ЭХ30Х^ 3x^3 ’ ЗХ2ЗХ3
32G
3
3Х
2
Для вторых производных имеем
3 2G
3x2
1 8a ^ 1 . .. 1 . . . 1
1 = --3 2 — sin m(/2 + -)H sin m(J2 + -H
. HH
л3 “ m m=1
sin
X
X
( (
5/ .- Б
/3 .3
I V
Л ( - Б
.1 + 2 |^,/1H1,Ага
//
1
л3 ^ п/т ^3173 + 2II ( _ 1 Л Н ( . 1 Л г . пН1
+ —з Е------- —2--------соэ п| /і + — | Н^оэ п| 7і + — І Ні эш-—1
а п=1
V пт
э2о2
8Ь 1 . 1 . 1 . пН 1
= —3Е -вт п(/і +2>Ні эш п(71 + -)Н эт——
X
X
( (
5/ 7 Б
/3 .3
I V
л п=1п
1
І2 + 2 ІН2,(/2 + 1)Н2,Ап
- Б
ЛЛ
.2 + 2 ІН2,/2Н2,Ап
/у
л3 - га/пга [ к3 [ 73 + 2
+лг Е—
Ь3 , У2
^ га=1 * пга
і ■ 1 Ли ( ■ 1 Ли • тН2
соэ га| /2 + 2 ІН2С0Э га I 72 + 2І Н2 sin—2^
Здесь функция й(х, у; А) = йу (х, у; А) - для вещественных Ап , и й(х, у; А) = а?2 (х, у; А) - для чисто мнимых Ап ;
Э20
— — у пга
- = -— ЕЕ-
„и ^ ^
/пт І И3 I 73 + ~
у 2 I пга
2 эшп(7і +-п'|HlС0sп(/'і +-2 |НіX
Эх2Эх1 аЬ л л у
z 1 п=1 га=1 I,
■ пН1 . ( 1 Л„ (. 1 Л„ . гаН2
Xsin-----га І /2 +— |Н2 соэ га І 72 + — |Н2Sin
2 И 21 Г 2I 2
гаН2
/373 2т ^
л Ь га=1
-5/373 “^Т Е SІn га I /2 + -2 IН2 соэ га І 72 + -2 IН2 sin
X
X
Ні 171 + 2 Л, Ні (/1 +1); А,
Л
(
Ні І 71 + і I,Ні/і;Аи
V
Далее,
ЭО (. . . . ) Э20і ( .... )
(/Ь/2, 7ъ 72) = ^^~ (7l, 72, /1, /2) •
ЭХ1ЭХ2 ЭХ2ЭХ1
Здесь функция 5(х, у; А) = 5 (х, у; А) - для вещественных Ат , и 5(х, у; А) = ^ (х, у; А) - для чисто мнимых Ат ;
Э20
— —
к=—ЕЕ
V. ТТҐ1
/пга ( Хз)
ЭХ3ЭХ1 ла л л у2 га
3 1 п=1 га=1 іпт'п
- эт (пХі) эт (гаХ2) соэ (пНі (/'і +0,5)) X
X 81П
Щ1 Ц sin (тН 2 (/2 + 0,5)) эт (
гаН
2 1.
Э202 8
ЭХ3ЭХ2 лЬ
=-лЬ ЕЕ
/пга ( Х3)
Х3
п=1 га=1 у пгап
пН
sin (пХі) эт (гаХ2) эт (пНі (/і + 0,5)) X
X 81П
гаН
Э 2G3 = _8_
ЭХ1ЭХ3 ла
ЕЕ
/пга ( Х3)
Х3
008
п=1га=1 у пгага
пН
(пХі )sin(гаХ2 )іп(п( (/'і + 0,5))X
X 81П
гаН
2 Л.
Э 2G3 = _8_ ЭХ2ЭХ3 лЬ
/пга (Х3)
ЕЕ _
п=1 га=1 у пгап
пН
Х3
sin (пХі) соэ (гаХ2) эт (пНі Оі + 0,5)) X
X 81П
^ sin (тН 2 (/2 + 0,5)) эт (
гаН
2 Л.
Э 2о
/пга ( Х3)
______3 = у ^_______________-1 Х3 Х3
-Л 2 2 2-і 2-і 2
3 л п=1 га=1 у пгапга
пН
эт(пХі )sin(гаХ2 )іп(пНі (/'і + 0,5))X
X 81П
гаН2
Полученные формулы позволяют эффективно вычислять приближенные суммы соответствующих рядов, т.к. оставшиеся ряды сходятся быстро (экспоненциально).
5 Метод коллокации
Рассмотрим теперь вопрос о построении схемы для метода коллокации для рассматриваемой задачи дифракции. Будем формулировать метод не для сингулярного интегрального уравнения, а для интегродифференциального уравнения. Этот подход оказывается эффективным в силу более удобного
представления интегралов. Будем предполагать, что матрица
об-
ратима в Q,
ё( х)
е0
N-1
- I
є Ь— (0, I - единичная матрица.
Введем обозначения
£ =
1
е0
и перейдем к следующему уравнению:
е0
Е
АЗ = £/(х) - к010(Х, у)3(у- ега<і <ііу | О(Х, у)3(у)йу = Е0 (х) . (27)
е е
Представим это уравнение в виде системы трех скалярных уравнений:
3 ^
^ и (х) - к01С(х, у)31 (у)йу-—Шух | в(х,у)3(у)йу = Е01 (х), / = 1,2,3 (28)
1=1 б Ц б
Определим компоненты приближенного решения 3 следующим образом:
_ N _ N _ N
31 = ^ Ок/кк х), 3 2 = ^ Ьк/к2( х), 3 3 = ^ ОкЙ(, х),
к=1 к =1 к=1
где /к - базисные функции-«ступеньки», существенно зависящие лишь от переменной хI.
Ниже проводится построение функций /к . Будем считать, что б - параллелепипед: б = {х: о < х1 < О2, Ь>1 < х2 < Ь2, С1 < хз < С2} . Разобьем б параллелепипедами:
Пк/т = {х : х1,к < х1 < х1,к+1, х2,/ < х1 < х2,1+1, х3,т < х3 < х3,т+1};
а2 - Щ 1 1 Ь2 - Ь , с2 - С
х1к = а1 +—-----1 k, х2 / = Ь1 +—-1 /, х3т = С1 +—-1 m,
п п п
где к = 1,..., п -1; /, т = 1,..., п -1.
Обозначив й1 :=| х1 к - х1 к-1 |, получим формулы для /кт :
1, х є ПкІга,
0, х *£ Пк/т.
23
Функции /к[т, /кт , зависящие от переменных х2 и хз соответственно, определяются аналогичными соотношениями. Построенное множество базисных функций удовлетворяет необходимому условию аппроксимации в Ь2 . Перенумеруем базисные функции:
/к1, /к, /к, к = 1,..., N,
где N = !(п3 - п2).
Расширенную матрицу для нахождения неизвестных коэффициентов а^, 0^ удобно представить в блочной форме:
А11 А12 А13 БЛ
А21 А22 А23 B2
А31 А32 А33 Б3 J
элементы столбцов Бк и матриц Aki определяются из соотношений:
4 = 4 (x); 4=ш! (xj) - §k/ko2 j Gk (xj, j) fi (-
Q
-^- j 4^Gl (Xj, y)f (y)dy, (29)
Oxk J OXi
Q 1
где координаты точки коллокации:
xi = (xi1,xi2,xi3 ) Xi1 = (( + 1/2) xi2 = (( + 1/2)h2’ xi3 = (( + 1/2)h3 , к, l = 1,2,3; i, j = 1,..., .
Таким образом, в статье представлены расчетные формулы для матричных коэффициентов метода коллокации для решения объемного сингулярного интегрального уравнения в задаче определения диэлектрической проницаемости материала. Численные результаты расчетов решения интегрального уравнения методом коллокации будут приведены в отдельной статье.
Список литературы
1. Shestopalov, Yu. V. Development of Mathematical Methods for Reconstructing Complex Permittivity of a Scatterer in a Waveguide / Yu. V. Shestopalov, Yu. G. Smirnov, V. V. Yakovlev // Proceedings of 5th International Workshop on Electromagnetic Wave Scattering, October 22-25. - Antalya, Turkey, 2008.
2. Медведик, М. Ю. Применение ГРИД-технологий для решения объемного сингулярного интегрального уравнения для задачи дифракции на диэлектрическом теле субиерархическим методом / М. Ю. Медведик, Ю. Г. Смирнов // Известия высших учебных заведений. Поволжский регион. Физико-математические науки. - 2008. - № 2. - С. 2-14.
3. Смирнов, Ю. Г. О существовании и единственности решений обратной краевой задачи для определения эффективной диэлектрической проницаемости наноматериалов / Ю. Г. Смирнов // Известия высших учебных заведений. Поволжский регион. Физико-математические науки. - 2009. - № 1. - С. 11-24.
4. Смирнов, Ю. Г. Применение ГРИД-технологий для решения нелинейного объемного сингулярного интегрального уравнения для определения эффективной диэлектрической проницаемости наноматериалов / Ю. Г. Смирнов // Известия высших учебных заведений. Поволжский регион. Физико-математические науки. - 2008. - № 3- С. 2-10.
Смирнов Юрий Геннадьевич
доктор физико-математических наук, профессор, заведующий кафедрой математики и суперкомпьютерного моделирования, Пензенский государственный университет
E-mail: mmm@pnzgu.ru
Медведик Михаил Юрьевич
кандидат физико-математических наук, доцент, кафедра математики и суперкомпьютерного моделирования, Пензенский государственный университет
E-mail|: _medv@mail.ru
Васюнин Денис Игоревич аспирант, Пензенский государственный университет
E-mail: mmm@pnzgu.ru
Smirnov Yury Gennadyevich Doctor of physico-mathematical sciences, professor, head of sub-department of mathematics and supercomputer modeling, Penza State University
Medvedik Mikhail Yuryevich Candidate of physico-mathematical sciences, associate professor, sub-department of mathematics and supercomputer modeling, Penza State University
Vasyunin Denis Igorevich
Post graduate student, Penza State University
УДК 517.9 Смирнов, Ю. Г.
Метод коллокации решения объемного сингулярного интегрального уравнения в задаче определения диэлектрической проницаемости материала I Ю. Г. Смирнов, М. Ю. Медведик, Д. И. Васюнин II Известия высших учебных заведений. Поволжский регион. Физико-математические науки. -2009. - № 3 (11). - С. 71-87.