Численный метод решения задач электромагнитного рассеяния на трехмерном магнитодиэлектрическом теле Текст научной статьи по специальности «Математика»

А.Г. Дмитренко, Т.Н. Пастухова

ЧИСЛЕННЫИ МЕТОД РЕШЕНИЯ ЗАДАЧ ЭЛЕКТРОМАГНИТНОГО РАССЕЯНИЯ НА ТРЕХМЕРНОМ МАГНИТОДИЭЛЕКТРИЧЕСКОМ ТЕЛЕ

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

Задачи электромагнитного рассеяния на трехмерных однородных магнитодиэлектрических телах в резонансной области частот уже много лет привлекают внимание исследователей. Интерес к этим задачам обусловлен многообразием применения конкретных результатов их решения в различных областях науки и техники, в частности, в радиолокации, метеорологии, биологии и атмосферной оптике (см., например, [1 - 4]). К настоящему времени разработаны различные методы решения рассматриваемых задач. Эти методы основаны на использовании как уравнений Максвелла в дифференциальной форме, так и интегральных соотношений теории электромагнитного поля. К первой группе методов относятся, например, метод разделения переменных [5], методы конечных разностей и конечных элементов [6 - 7], ко второй группе - например, метод продолженных граничных условий [8], различные схемы метода моментов [9]. Существуют также гибридные методы, объединяющие основные черты конечных методов и методов интегральных уравнений (см., например, [10]). Однако вычислительные алгоритмы, основанные на вышеупомянутых методах, получаются чрезвычайно емкими по затратам ресурсов и времени компьютерных вычислений, особенно в случае рассеивателей, не обладающих симметрией вращения. Для конечных методов это обусловлено необходимостью распространения вычислений на всю рассматриваемую область пространства, для методов интегральных уравнений - необходимостью вычисления большого числа поверхностных или объемных интегралов.

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

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

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

ФОРМУЛИРОВКА ЗАДАЧИ И МЕТОД ЕЕ РЕШЕНИЯ

Геометрия задачи показана на рис. 1. Будем рассматривать стационарную (зависимость от времени выбрана в виде ехр(-,'ю/)) задачу дифракции электромагнитного поля {Е0,Й0} на трехмерном однородном теле Д, ограниченном поверхностью Б, с диэлектрической є, и магнитной д,- проницаемостями. Вне тела (область Бе) среда безгранична, однородна и характеризуется электромагнитными параметрами єе,

де. Требуется найти рассеянное поле {Ее, Йе} в области Бе.

Рис. 1. Геометрия задачи

Математическая постановка задачи имеет следующий вид:

V х Ее = т^еЙе УхЙе = -іюєе Ее

VхEi = гацЙ

б ^хЙ, =-,а£гЕг

;(Е, - Ее )= ПХ Е0

- Йе ) = п х Й0

{Ш ;л/йЄЙе }х к / к+

-{4^еЙ е; -4^еЁе } = О (К-), К

(1)

(2)

где Ее, Не и Е,, - поля в областях Д и Д ; п -

единичный вектор нормали к поверхности Б; Я2 = х2 + у2 + г2; а х Ь - векторное произведение. Будем полагать, что Б е С(2), п х Е0, п х Н0 е Ст(5), 1т е,, 1т ц, > 0, тогда существует единственное решение граничной задачи (1) - (3).

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

Выберем на вспомогательной поверхности Бе конечную совокупность точек {Мп е }}= . В каждой точке Мп е разместим пару независимых вспомогательных элементарных электрических диполей с момен-

хп,е,1 _,п,е,1 ^ п,е,1

п,е,1 _,п,е,1 ^ п,е,1

а также

еТ2 , излучающих в однородную среду

с параметрами єе, це. Единичные направления еП{е,х

^е > г'е

?п,е,2 т2

Ее (М) = £\- Vх (Vх П„,ед) —Vх П„

п=1 [К Єе

Йе (М) = £\-VхП „,ед + V х (V х Пи,е,2 ) ^, (4)

п=1 [М'е ке

П„,е,1 = ¥е (М, МП е )р.

Пп,е,2 = ¥е (М, М„е)р.

п,е,1

п ,е,2

рп е,1 = е,1е1 е,1 + рТ2е лё„2 е,1,

рп^2 = pnl,e,2exn1,e,2 + р^Ч^2, М е Д,

а поле Е,, Hi в Д - в виде суммы полей вспомогательных диполей, расположенных на вспомогательной поверхности Б,-:

п-

N

Е, (М) = £ \ ^ V х (V х Пп,гД) --V х Пп =1 [к Єг

Йг (М) = £ \- Vх ГІпЛІ + ^ Vх (Vх гіп,,2)

п=1 [ь к2

(5)

Пп ,і,1 = ¥ і (М, М п -) р Г1пЛ2 = ¥ (М, мп,) р.

п,і,1

п,і,2

тами Рт1’ = Рт1 ■ ет^ И Рт2’ = Рт2’ е,

пару независимых вспомогательных элементарных магнитных диполей с моментами р”[’е’2 = Р^2 ет”’е’2 и

^ п,е,2 п,е,2 ^ п,е,2

Рт2 , = Рт2 , ет', ,

Р п,і,1 = пЩЛ£ п,і,1 + пп,і,1^ п,і,1 і т -Гт1 ^т1 -г т2 т2 >

Р"’'’2 = р^2?Г1л2 + РтТ2е^2, м є д .

В представлениях (4), (5)

¥е (М,Мп,е ) = ЄХР(ікеКММ„,е ) 7 4пКМ,Мпе ,

¥г (М,Мп,г ) = ЄХР(ікгКММ*, ) 7 4ПКМ,Мп ,■

- фундаментальные решения уравнения Гельмгольца для областей Бе и Д; ке,, = ю^,, це г-, КММ, ^ - расстояние от диполя в точке Мп на Бе до точки М в Бе;

К

М ,М„гі

- расстояние от точки Мп г на Б, до точки М в

п,е,1 ^п,е\ ^п,е,2 ^п,е,2

Рт2 , Рт1 , , Рй ,

РTI^г’1 :

Рxnl’г’2 >

pTI^г’2

и eтn2e, выбраны в плоскости, касательной к Бе в точке Мп е. Аналогично, на вспомогательной поверхности Б, выберем совокупность точек {Мп, 1, в каждой из которых разместим пару неза-

висимых вспомогательных электрических диполей с

моментами р”,гД = pтnl,г,1eтnl,г,1 и Р^1,1 = РтТ1^1 и пару независимых вспомогательных магнитных диполей

с моментами р^’2 = р”!1,2 еТ112 и рп22 = Рхп2 ,2 ет;2’2 , которые расположены в плоскости, касательной Б, в точке Mni, и излучают в однородную среду с параметрами е,, ц,. Представим неизвестное рассеянное поле {Ее, Не} в Де в виде суммы полей введенных вспомогательных диполей, расположенных на Бе:

Д; Ые и Ы, - число точек размещения диполей соответственно на Бе и Б; р^

(п = 1,2,..., Же) и ртп1,и

(п = 1,2,..., Ж,) - неизвестные комплексные постоянные (дипольные моменты).

Представления (4), (5) удовлетворяют уравнениям Максвелла (1) и условиям излучения (3). Для того чтобы удовлетворить граничным условиям (2), необходимо соответствующим образом выбрать значения

постоянных р"^1,рп2еД , р-п^2,Рrn^,2 (п = 1,2,...,Же)

и р^\ рЩ’1 , pтr,2, рЩ’2 (п = 1,2,...,N). Используем для этого метод коллокаций. Пусть Му (у = 1,2,...,Ь)- точки коллокации на поверхности

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

4 Ь х 4( Ые + Ж,):

п1 х (Еі - Ее; ) = п1 х Е0

Я1 х - ЙI) = п1 х Й

(6)

где п - значение единичного вектора нормали в точке Му на поверхности тела; ЕУ, Ну и Е/, Ну - значения компонент внешнего и внутреннего полей в точке М/; Е0, НО - значения компонент возбуждающего поля в этой же точке. Решение системы (6) определяем путем минимизации функционала

Ф = X{| n1 х^Ё1 - Ё^)- n1 х £0| + 1=1

+Ь>|п1 х ( - Hj ) - n1 х H0 |2 } . (7)

Задачу минимизации функционала (7) решаем методом сопряженных градиентов.

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

Ёв,е (M) = (цв/ 8е f He(S) (M) =

= (exp (ikeR )/кеЯ)Ве (е, ф) + O(R-2), (8)

Ёв,ф (M) = -(^/8е )12 Не,е (M) =

= (exp (ikeR )/keR) )ф (е, ф) + O(R-2),

где компоненты диаграммы рассеяния D0 (е, ф) и Бф (е, ф) определяются выражениями

Ne

D0(е,ф) = £Gne(е,ф) {(cosеcosфcosa^1 +

П=1

+ cos е sin ф cos РП,еЛ - sin е cos уП,еЛ)рП^1 +

+(cos е cos ф cos а2,еД + cos е sin ф cos pn,e,1 -

- sinе cos Y2,e,1)p^1 + ((xJSe )1/2 х х[(- sin ф cos a^2 + cos ф cos Pf^2) p^2 +

+(- sin ф cos a 2a2 + cos ф cos Р2а2) p^2]},

Ne

Dp (е, ф) = X G„,e (е, ф){ (- sin ф cos aj^1 +

n=1

+cos ф cos Pn,e,1) p^1 + (- sin ф cos а„,еД + +cos фcos P2,e,1) p^1 + (v-ehe )1/2 х х[(-cosе cosфcos an,e,2 - cosе sinфcos pj^2 + +cos е sin ф cos pn,e,2 + sin е cos Y1„,e,2) p^2 +

+(- cos е cos ф cos a2a2 - cos е sin ф cos P2a2 +

+ sin е cos Y Г2) px„2e,2l}, (9)

Gne (е, ф) = exp{-ike (sin е cos фх„,е +

+ sin е sin фу„е + cos егпе)},

в которых cos a„,e,1,cos p„,e,1,cos Y„,e1 - направляющие косинусы единичного вектора e^1; cos a^,e1,

P„e 1 „ e1

2 , cos Y2 - направляющие косинусы еди-

ничного вектора e^1; cosa„,e,2,cosp„,e,2,cosY1„,e,2 -направляющие косинусы единичного вектора e^^ 2; cos a 2,e,2,cos P2,e,2,cos y „,e2- направляющие косинусы единичного вектора eT„2e 2; x„,e, y„,e, zne - декартовы координаты точки Mn e; 0 и ф - общепринятые угловые сферические координаты точки наблюдения M.

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

Д = (Ф' / Ф0)1/2, (10)

где Ф' - значение функционала (7) при условии, что вместо точек коллокации / выбирается сетка промежуточных точек]' (/' = 1,2,...,Ь'), а Ф0 определяется выражением

Ф0 = !|{|п1 х Ё]/12 + Це\п1 х Н0;12|.

ЧИСЛЕННЫЕ РЕЗУЛЬТАТЫ

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

Входными величинами программы являются геометрические параметры рассеивателей (в длинах волн), значения относительных диэлектрической е = е,/ее и магнитной ц = ц,/це проницаемостей,

возбуждающее поле {Е0, Н0}, параметры подобия Ке

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

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

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

На рис. 2 приведены кривые, иллюстрирующие изменение относительной нормы невязки граничных условий (10) при одновременном увеличении числа точек размещения диполей и числа точек коллокации Ь на поверхности рассеивателя. В качестве последне-

го выбран эллипсоид с полуосями кеа = 5,0, кеЬ = 3,0, кес = 2,0, ориентированными вдоль осей X, У, 2 соответственно. Электродинамические параметры эллипсоида следующие: е = 8, ц = 1. Эллипсоид

возбуждается линейно поляризованной плоской волной, распространяющейся вдоль оси 2; электрический вектор Е0 волны ориентирован вдоль оси X. Положения вспомогательных поверхностей характеризуются следующими значениями параметров подобия:

Ке = 0, 5, К, = 5, 0 . Необходимо отметить, что точность решения задачи зависит от положения вспомогательных поверхностей, поэтому предварительно были проведены соответствующие исследования, которые позволили установить, что наименьшее значение нормы невязки граничных условий достигается при значениях параметров Ке и К,, лежащих в интервалах 0,4 < Ке < 0,6, 3,0 < К, < 6,0, что совпадает с выводами работы [12]. При проведении исследований числа точек размещения диполей на вспомогательных поверхностях Бе и Б, выбраны одинаковыми и равными половине числа точек коллокации: Ые = Ы, = Ь /2. Таким образом, числа точек размещения диполей жестко связаны с числом точек коллокации. Это позволяет исследовать изменение нормы невязки граничных условий при одновременном изменении чисел точек размещения диполей и точек коллокации.

Рис. 2. Зависимость нормы невязки граничных условий от чисел точек размещения диполей и точек коллокации для эллипсоида с полуосями kea = 5.0, keb = 3.0, kec = 2.0 и относительными значениями проницаемостей 8 = 8, ц = 1. Кр. 1 - комбинированная система диполей; кр. 2 - система электрических диполей

Кривая 1 на рис. 2 относится к обсуждаемому методу, в котором в каждой точке размещения диполей на обеих вспомогательных поверхностях располагаются пара электрических и пара магнитных диполей; кривая 2 относится к методу, предложенному ранее в [12], в котором в каждой точке располагается только пара электрических диполей. Как видно из рисунка, использование комбинированных (электрических и магнитных) систем диполей позволяет существенно (в рассмотренном примере на 15%) повысить точность решения задачи, причем выигрыш в точности возрастает с увеличением значения диэлектрической прони-

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

Рис. 3 и 4 иллюстрируют результаты сравнения нормированных амплитудных диаграмм рассеяния в Ё- и H-плоскостях, рассчитанных предполагаемым методом, с такими же результатами, приведенными в работе [13]. Результаты относятся к вытянутому сфероиду с полуосями kea = keb = п/2, kec = п и проницаемостями 8 = 4 + i • 2, ц = 1. Полуось kec ориентированна вдоль оси Z, полуоси keа и keb - вдоль осей X и Y соответственно. Сфероид возбуждается плоской волной, падающей вдоль оси Z. Параметры подобия выбраны следующими: Ke = 0,5, Ki = 5,0.

1Ёе1/ 1Ёе Imax

Рис. 3. Нормированные амплитудные диаграммы рассеяния в Е-плоскости сфероида с полуосями кеа = кеЬ = п/2, кес = п и проницаемостями е = 4 +,• 2, ц = 1. Кр. 1 - данный метод; кр. 2 - результаты работы [13]

KI/I ёф1

Рис. 4. Нормированные амплитудные диаграммы рассеяния в Н-плоскости сфероида с полуосями кеа = кеЬ = п/2, кес = п и проницаемостями е = 4 +,• 2, ц = 1. Кр. 1 - данный метод; кр. 2 - результаты работы [13]

Числа точек размещения диполей на внутренней и внешней вспомогательных поверхностях выбраны одинаковыми: Ne = Ni = 144. Алгоритм размещения этих точек на обеих поверхностях также выбран одинаковым: в каждом из двенадцати полусечений ф = const, отстоящих друг от друга на угловое расстояние Дф = 30°, равномерно по углу 0 расположены двенадцать точек размещения диполей. Число точек коллокации выбрано в два раза большим: L = 288. Алгоритм их расположения по углу 0 выбран таким же, как для точек размещения диполей, но выбирают их как в полусечениях ф = const, определенных для точек размещения диполей, так и посередине между ними.

Кривые 1 на рис. 3 и 4 - это результаты, полученные предлагаемым методом; им соответствует значение невязки Д, равное 7%. Кривые 2 - результаты работы [13], полученные методом расширенных граничных условий. Как видно из рисунков, кривые иллюстрируют очень хорошее совпадение сравниваемых результатов. Небольшие наблюдаемые различия можно объяснить ошибками при графическом съеме информации с рисунка статьи [13], а также погрешностями расчетов.

Рис. 5 и 6 характеризуют влияние отклонения формы рассеивателя от осесимметричной на бистати-ческое сечение рассеяния

ст(е,ф)= lim4nR2 {е (е,ф)2 +|Ёе ф (е,ф)2^^э|2 (11)

R——ад ' '•

для эллипсоидов, характеризующихся одинаковыми значениями проницаемостей: 8 = 8, ц = 1, но различными значениями полуосей. Результаты приведены в Ё-плоскости (ф = 0). По оси абсцисс на рис. 5 и 6 отложен угол 0 (см. рис. 1), по оси ординат - величина ст/х2 , где X - длина волны во внешней среде. Как и ранее, во всех рассмотренных случаях полуоси kea, keb, kec ориентированы вдоль осейX, Y, Z декартовой системы координат; возбуждающая рассеиватели линейно поляризованная плоская волна падает вдоль оси Z, вектор Е0 ориентирован вдоль оси X. Параметры метода для всех рассеивателей выбраны одинаковыми: Ke = 0,5, K = 5,0, L = 2Ne = 2Nt = 392 . Точки размещения диполей на внутренней и внешней вспомогательных поверхностях выбраны следующим образом: в каждом из четырнадцати полусечений ф = const, отстоящих друг от друга на угловое расстояние Дф = 25,7°, равномерно по углу 0 расположены четырнадцать точек размещения диполей. Алгоритм расположения точек коллокации по углу 0 выбран таким же, как для точек размещения диполей, но выбирают их как в полусечениях ф = const, определенных для точек размещения диполей, так и посередине между ними. При данных параметрах метода во всех случаях значение невязки (10) не превышает 11%.

Кривая 1 на рис. 5 относится к осесимметричному рассеивателю - сфероиду с полуосями kea=keb=3,0, kec = 2.0 ; кривые 2 и 3 - к эллипсоидам, которые отличаются от вышеупомянутого сфероида только тем, что полуось kea увеличена на 0,05X и 0,1X соответственно. Аналогично, кривая 1 на рис. 6 относится к

сфероиду с полуосями кеа = кеЬ = 4,0, кес = 2,0, а кривые 2 и 3 - к эллипсоидам, полуось кеа которых увеличена на 0,05Х и 0,1Х соответственно.

ст/х2,дБ

Рис. 5. Бистатические сечения рассеяния в Е-плоскости эллипсоидов с проницаемостями е = 8, ц = 1, характеризующихся различными значениями полуосей. Кр. 1 - сфероид, кеа = кеЬ = 3,0, кес = 2,0; кр. 2 - эллипсоид, кеа = 3,314, кеЬ = 3,0, кес = 2,0; кр. 3 - эллипсоид, кеа = 3,628, кеЬ = 3,0, кес = 2,0

ст/х2,дБ

Рис. 6. Бистатические сечения рассеяния в Е-плоскости эллипсоидов с проницаемостями е = 8, ц = 1, характеризующихся различными значениями полуосей. Кр. 1 - сфероид, кеа = кеЬ = 4,0, кес = 2,0; кр. 2 - эллипсоид, кеа = 4,314,

кеЬ = 4,0, кес = 2,0 ; кр. 3 - эллипсоид, кеа = 4,628, кеЬ = 4,0, кес = 2,0

Сравнение кривых, приведенных на рис. 5 и 6, между собой позволяет сделать вывод, что небольшие отклонения формы рассеивателя от осесимметричной в наименьшей степени влияют на сечение рассеяния в прямом направлении (0 = 0) и в наибольшей степени -на сечение рассеяния в направлениях, прилегающих к направлению бокового рассеяния (0 = 90°). Особенно чувствительными к малому изменению формы рассеивателя являются окрестности резко выраженных минимумов диаграммы рассеяния (рис. 6). Малое изменение формы рассеивателя приводит как к изменению положения минимумов диаграммы рассеяния, так и к изменению значения сечения рассеяния в минимуме.

1. Seker S.S. // IEE Proceed. 1986. V. 133. Pt. H. Ко. 4. P. 305.

2. Апельцин В.Ф., Ильинский А.С., Сабитов Б.Р. // Вестник МГУ. Сер. 15. Вычислит. мат-ка и кибернетика. 1985. № 3. C. 47.

3. Lakhtakia A., IskanderM.F. // IEEE Trans. 1983. V. AP-31. No. 1. P. 111.

4. Havemann S., Baran A.J. // J. Quant. Spectrosc. Radiat. Transfer. 2001. V. 70. P. 139.

5. Rother T. // J. Quant. Spectrosc. Radiat. Transfer. 1998. V. 60. P. 335.

6. Sun W, Fu Q., Chen Z. // Appl. Opt. 1999. V. 38. P. 3141.

7. JensenM.A., Rahmat-Samii Y. // Radio Sci. 1996. V. 31. P. 1823.

8. Mishchenko M.I., Travis L.D., Mackowski D.W. // J. Quant. Spectrosc. Radiat. Transfer. 1996. V. 55. P. 535.

9. Wang J.J. // Generalized Moment Methods in Electromagnetics. New York: Wiley, 1991.

10. ShengXin-Qing, Jin Jian-Ming, Song Jiming et al. // IEEE Trans. 1998. V. AP-46. No. 3. P. 303.

11. Wriedt T. Generalized Multipole Techniques for Electromagnetic and Light Scattering. Amsterdam: Elsevier Science, 1999.

12. Дмитренко А.Г., Мукомолов А.И. // РЭ. 1995. Т. 40. № 6. С. 875.

13. MorganM.A., Chen C.H., Hill S.C., BarberP.W. // Wave Motion. 1984. V. 6. P. 91.

Статья представлена кафедрой исследования операций факультета прикладной математики и кибернетики Томского государственного уни-

верситета, поступила в научную редакцию «Математика» 14 июня 2005 г.