МОДЕЛИРОВАНИЕ ОТРАЖАТЕЛЬНЫХ ХАРАКТЕРИСТИК ОСЕСИММЕТРИЧНЫХ РАДИОЛОКАЦИОННЫХ ОБЪЕКТОВ МЕТОДОМ ИНТЕГРАЛЬНЫХ УРАВНЕНИЙ Текст научной статьи по специальности «Физика»

-

https://doi.org/10.38013/2542-0542-2021-4-76-93 УДК 621.391.812.7

Моделирование отражательных характеристик осесимметричных радиолокационных объектов методом интегральных уравнений

А. А. Рунов

Публичное акционерное общество «Межгосударственная акционерная корпорация «Вымпел», Москва, Российская Федерация

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

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

Для цитирования: Рунов А. А. Моделирование отражательных характеристик осесимметричных радиолокационных объектов методом интегральных уравнений // Вестник Концерна ВКО «Алмаз - Антей». 2021. № 4. С. 76-93. https://doi.org/10.38013/2542-0542-2021-4-76-93

For citation: Runov А. А. Modelling of the reflectance profile of axisymmetrical radiolocation objects using the integral equation method // Vestnik Koncerna VKO "Almaz - Antey". 2021. No. 4. P. 76-93. https://doi.org/10.38013/2542-0542-2021-4-76-93

Поступила 04.08.2021 Отрецензирована 23.08.2021 Одобрена 01.12.2021 Опубликована 27.12.2021

см о см

тг

^ Оценка эффективной поверхности рас-

£ сеяния (ЭПР) радиолокационных объектов является одной из актуальных задач в радиоло-

$ кации. Одним из широко применяемых теоре-

^ тических методов расчета ЭПР является метод

^ интегральных уравнений (ИУ) [1]. Метод ИУ обладает принципиально более высокой точно-

Щ стью по сравнению с геометрической теорией

« дифракции (ГТД) [2] и методом краевых волн

о (МКВ) [3] во всем диапазоне своего примене-

| ния и монопольно используется в области ре-

ь зонансного и рэлеевского рассеяния, где ГТД

т и МКВ неприменимы.

™ Формулировка задачи дифракции ме-

ср тодом интегральных уравнений возможна

™ в виде ИУ I или в виде ИУ II рода. В общем

2 _

(П -

- © Рунов А. А., 2021

случае она является некорректно поставленной задачей по Тихонову, то есть неустойчивой к погрешностям в исходных данных. Таким образом, возникает проблема точности получаемых решений. Решения электродинамических задач на основе ИУ II рода в большинстве практических случаев (например, для идеально проводящих тел) являются более устойчивыми по отношению к решениям на основе ИУ I рода, поэтому исторически при решении практических задач давалось предпочтение ИУ II рода.

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

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

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

Достаточно представительной выборкой радиолокационных объектов являются металлические тела вращения, для которых учет осевой симметрии при математической формулировке задач дифракции позволяет дополнительно снизить размерность до одномерной [1].

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

Задача дифракции электромагнитных волн (ЭМВ) на поверхности идеально проводящего тела может быть сформулирована в виде ИУ I рода для электрического поля. Граничные условия для электрического поля на идеально проводящей поверхности определяются как [4]:

[п, Е]=0,

где: п - единичный вектор нормали к идеально

проводящей поверхности;

Е - вектор электрического поля на идеально

проводящей поверхности;

[п, Е] - векторное произведение векторов п и Е.

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

В более наглядном для рассматриваемой задачи виде граничные условия можно представить как:

Ех = 0,

где Ет - тангенциальная составляющая электрического поля на идеально проводящей поверхности.

С учетом данных граничных условий для суммарного поля Е соотношение для падающего (первичного) и рассеянного (вторичного) полей на поверхности идеально проводящего тела определяется как:

ЕТ + ЕТ = 0, (1)

где: ЕТ - тангенциальная составляющая падающего (первичного) электрического поля; ЕТ - тангенциальная составляющая рассеянного (вторичного) электрического поля.

В формуле (1) и далее по тексту верхний индекс г для векторов электрического поля Е означает падающее поле, верхний индекс ^ -рассеянное поле.

Для комплексных гармонических величин рассеянное электрическое поле можно определить как функцию от векторного А и скалярного Ф потенциалов на поверхности тела:

Е'= -/юЛ - VФ, (2)

Г Г е~1кя где: А(5) = ц]] 1(5') —(3)

5'

1 Г Г е~>ш

(4)

5' _

та

В формулах (3), (4) используются следу- | ющие обозначения:

О

|(5") - распределение плотности электрическо- £ го тока на поверхности проводника;

CM

о

CM

< I

со го

s

о со

о.

ф

о

о

V

со

см ■ci-io

с?

см ■ci-io см

(П (П

o(S') - распределение плотности электрических зарядов на поверхности проводника; S - поверхность проводника, соответствующая точкам наблюдения (измерения) потенциалов Л^) и Ф^);

S' - поверхность проводника, соответствующая точкам интегрирования плотности тока и плотности заряда о(^), в описываемой модели S и S' геометрически совпадают; s - точка наблюдения (измерения) потенциалов Л^) и Ф(s);

R - расстояние между точками наблюдения и интегрирования;

ц, 8 - абсолютные магнитная и электрические проницаемости вакуума; к = 2л/Х - волновое число; X - длина волны в свободном пространстве; V - векторный дифференциальный оператор набла, операция VФ(s) определяет градиент скалярной функции Ф^) в точке наблюдения s; / - мнимая единица.

С учетом (1) и (2) задача дифракции в виде ИУ для электрического поля формулируется как:

Ег(?) = ,/шЛ^)+ УФ^). (5)

Уравнение (5) с учетом (3) и (4) содержит две неизвестные величины - |(У) и о(^), которые в каждой точке s' поверхности S' связаны уравнением непрерывности:

У'т = -/шо^'), (6)

где V - скалярный дифференциальный оператор набла, операция определяет градиент векторной функции в точке интегрирования s'.

С учетом (6) уравнение (5) может быть представлено в виде, содержащем одну неизвестную функцию

. Г Г

!(*') —И'-

(7)

1 ГГ е~]кК

/СОЕ ]] } 4пЯ

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

Для сведения интегрального уравнения к дискретному обычно используются проекционные методы, например метод моментов (МОМ) [5] или его разновидности. В МОМ для ИУ типа (7) искомая плотность электрического тока J(S') (неизвестная функция) и первичное поле E'(S (известная функция) аппроксимируются наборами базисных и весовых функций. Наборы функций должны быть полными и по возможности содержать ортогональные составляющие.

Искомая плотность электрического тока J(S') на поверхности S' аппроксимируется суммой базисных функций с весовыми коэффициентами:

к,

7=1

где: 1/^) -/-я базисная функция; I/ - весовой коэффициент/-й базисной функции;

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

Используемые для аппроксимации базисные функции являются известными (выбранными) функциями, а значения весовых коэффициентов I/ неизвестны.

Первичное поле Ег^), соответственно, аппроксимируется суммой весовых функций с соответствующими коэффициентами:

Ки,

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

i=i

где: - /-я весовая функция;

V - весовой коэффициент /-й весовой функции;

г - при записи как нижний индекс в Wi(S) или нижний индекс в V означает номер весовой функции или ее весового коэффициента; Кцг количество используемых весовых функций.

В аппроксимации Ег(5) известными являются как весовые функции Wi(S), так и их весовые коэффициенты V.

Верхний индекс г в Ег(5), как отмечалось ранее, означает падающее (первичное) поле.

Для вывода дальнейших соотношений введем скалярное произведение весовых и базисных функций по поверхности S:

5

Далее по тексту выражение вида <А(5), В(5)> означает скалярное произведение векторных функций А(5) и В(5) на поверхности

Тогда:

1] = <|(5'), Щ')>, (8)

V = (Е'(^, Wi(S)>. (9)

Для компактности и наглядности дальнейшего изложения представим (7) в операторном виде:

Е^ = Ь{|^}.

Далее, подставив в данное операторное уравнение дискретную аппроксимацию |(5")

набором базисных функций 1]-(5"), получим: (к] у

Е = L

С учетом линейности оператора Ь данное уравнение можно представить как:

к,

С учетом (9) и свойства линейности операции скалярного произведения полученное уравнение можно представить как:

К1

^ (10) У=1

На основании (10) получается система линейных алгебраических уравнений (СЛАУ), являющаяся дискретной аппроксимацией ИУ (7):

2Уг = XII(Щ(Ю) 1] +

+ ХЦ Ь 1Д5') ) ) М.

I 5

\]=1

Второе слагаемое системы уравнений (11) содержит операцию производной УФ(|(5")) по координатам точки наблюдения, что достаточно неудобно для последующего численного решения. Для приведения (11) к модификации более приемлемой в вычислительном отношении воспользуемся свойством интегрирования по частям и теоремой о дивергенции векторной величины на замкнутой поверхности. В условно краткой форме второе слагаемое (11) можно записать как:

JJw(V 4>)dS = JJ V(W<D)dS-

- JJ(VW)0 dS.

s

(12)

7=1

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

с весовыми функциями Wi(S), получаем:

К]

<Ш4(5),Е'(5)> = Ъ{]Д5')}).

У=1

Подынтегральное выражение первого слагаемого (12) содержит дивергенцию векторной величины WФ. Поскольку уравнение (11) соответствует дифракционной задаче на идеально проводящем теле, в соответствии с используемыми граничными условиями все входящие в него векторные величины являются тангенциальными к поверхностям 5 и 5'. По этой причине второе слагаемое (12) равно нулю, так как поверхность 5 замкнута. В этих _ условиях СЛАУ (11) может быть преобразова- | но к виду более удобному в вычислительном й

отношении:

см о см

< I

со та

г

о со

о.

<и

о

о <и со

см ■ч-ю

с?

см ■ч-ю см

(П (П

I 5

В отличие от уравнения (11) во втором слагаемом уравнения (13) вместо трудно вычисляемой операции градиента скалярного потенциала УФ(|^')) появилась операция дивергенции весовой векторной функции V(Wi(S)). Как следствие, векторная функция Wi(S) должна быть дифференцируемой, а скалярная функция V(Wi(S)) - относительно легко вычисляемой.

В матричной форме СЛАУ (13) имеет

вид:

= |^|| ||1||,

где: ||У|| = ||<Ei(S), Wi(S)>|| - матрица-столбец коэффициентов разложения первичного (падающего) поля в ряд по весовым функциям; |ОД| = ||ЬЩ^)}^ - прямоугольная матрица взаимных сопротивлений Wi (на поверхности S) по отношению к ]/ (на поверхности S'); 1|1|| = 11<- матрица-столбец коэффициентов разложения искомого тока в ряд по базисным функциям

В приведенном выражении, расшифровках и далее по тексту обозначение вида ||А|| (символы, исполненные прямым жирным шрифтом, в двойных прямых скобках) означает матрицу.

Р-»|| -

11.2

В дальнейшем будет рассматриваться случай, когда размерности ||У|| и ||1|| совпадают (матрица |^|| - квадратная, т.е. К/ = К№).

С учетом вышеизложенного выражение для элементов матрицы взаимных сопротивлений тела |^|| с произвольной геометрией в развернутом виде имеет вид:

гц И \ViCS) +

1 ГГ ГГ е~'ш (14)

где нижние индексы / и / обозначают номера строк и столбцов матрицы |^||.

Приведенный вывод ИУ (7) или СЛАУ

(13) не претендует на оригинальность. Варианты получения подобных ИУ, а также выражений для коэффициентов СЛАУ типа

(14) известны и широко приведены в литературе, например в [1], стр. 23-26; [5], стр. 25-30; [6] и др.

Для тел вращения СЛАУ (13) можно модифицировать. Если в (13) использовать Wi(S) и |у(^), которые являются функциями азимутальных координат Ф и Ф' поверхности тела вращения (рис. 1), то тангенциальные составляющие Е;Х$), Е£(£) (1) и |(^) (3) можно разложить в ряд Фурье по азимутальной координате на азимутальные гармоники. В этих условиях, учитывая свойство ортогональности тригонометрических рядов, матрицу |^|| можно преобразовать к блочно-диагональному виду относительно номеров гармоник п:

Цг-Л

Р„11 II о || игл

РН =

- Ш\

где: п - максимальный («верхний») номер используемой азимутальной гармоники, общее количество гармоник = 2п + 1;

|^п|| - матрица взаимных сопротивлений, соответствующая азимутальной гармонике с номером п; ||0|| - нулевая матрица;

нижний индекс в |^п|| означает номер азимутальной гармоники.

Соответственно, матрицы иметь вид:

IIV-J

l|V_2||

IlV-ill

IlVoll

IlVill

HV2||

и ||I|| будут

l|Vn||

IIUII

i-i II

IUI = llloll IHill

inj

III', II

Учитывая ортогональность тригонометрических функций, вычисление матриц ||Z||, ||V|| и ||I||, а также решение СЛАУ (13) можно производить отдельно для каждой гармоники.

При введении системы цилиндрических координат на поверхности тела вращения согласно рисунку 1, матричное уравнение для одного блока, соответствующего n-й азимутальной гармонике, будет иметь вид:

7t,t' 7t,4>'

71 TL

V |Ф'

п

vn

17Ф

где: n - номер азимутальной гармоники;

г t,t'

•п

Ф ,t'

'п

k,l

,t,4>' ■n

Ф,Ф'

'n

- квадратная матрица ||Zn||,

||Z„ || - матрица взаимного сопротивления

Рис. 1. Система координат на поверхности тела вращения и сферические координаты падающей ЭМВ: и иф; и^; - единичные векторы, касательные к поверхности тела вращения в точках наблюдения и

интегрирования; t, ? - координаты вдоль образующей тела вращения

в точках наблюдения и интегрирования; ф, ф' - координаты, ортогональные образующей тела вращения;

0, ф - сферические координаты направления прихода падающей ЭМВ

поверхностного тока (первый верхний индекс относится к точке наблюдения, второй - к точке интегрирования);

Г

ф'

ТФ

- матрица-столбец ||In||;

- матрица-столбец ||Vn||.

Формулы для элементов матриц ||Znk,l||

к-го компонента касательной составляю- для различных азимутальных гармоник п бу-щей электрического поля к /'-му компоненту дут иметь вид: {=т с'=г

,9п+1 + дп-1>

щ J m J /y(t') (sine' sinß0n+1 * + f}(t') (cose'cosfl'gj dt'dt ' t-o t'=0

t=T .

^ Ü _T

+ 17Z f f'M f f j(.t')gndt' dt, t=0

(15)

та

та 5 <D

t=T t =Т t=T ,_

(zf). . =-соц I m J fj(t') (sinQdn+1 ~ J fwf^fjWgndt'dt,

' t=o t'=o t=o 4

(16)

t=T t -T t=T ,_T

(zf').. = +c0(i J nit) J fjQtT) (srnd3n+1 ~ У"'1) dt'dt + |/,(t)[" f'^gndfdt, (17) ' t=o t'=o t=o f

t=T t'=T t=T ,_

(zf')..=Mi J m \ fM9n+1l9n-^^p f fiVf^fjWandt'dt, (18)

' t=o t'=o t=o f _0

где:

0n

и

-s

,-jkR

R

■ cos пфйф;

(19)

CM

о

CM

<

1

(0 та

2

О CO

.

<D

О

о

<D CO

CM ■Clio 9

CM ■Clio

CM

w w

п - номер азимутальной гармоники; 0 - угол между вектором и и осью Z (рис. 1);

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

иф - единичный орт, касательный к поверхности тела вращения и направленный в направлении вращения образующей (перпендикулярно и);

^ - координата на поверхности тела вращения вдоль его образующей; ф - координата на поверхности тела вращения в направлении вращения образующей; £ (¿) - функция, имеющая первую производную//(¿) по ^ (координате образующей).

В формулах (15)—(18) нижние индексы / (рис. 1), то матрицу коэффициентов возбужде-

и обозначают номера строк и столбцов, к которым относится данный элемент матрицы

I |Z | | .

В используемой модели в качестве f (t) иf (t') применялись кусочно-треугольные функции с весовыми коэффициентами, являющимися функциями от р и р' в целой степени (рис. 2).

В качестве весовых функций в (15)—(18) выбраны:

Щп = Ut f (t) j (20)

Wf,„ = иф f (t) (21)

В качестве базисных функций выбраны:

Jj,n = U f (t') e

M'

(22)

ния весовых функций как:

можно представить

Ег >|

n,n ^ф' I

Ee> Ег >

Г,п = иф'£ (О ^пф, (23)

Если обозначить ориентацию вектора поляризации падающей на тело вращения ЭМВ в плоскости, проходящей через оси X и I, как 9-поляризацию, а поляризацию в перпендикулярной плоскости, как ф-поляризацию

| | УП,0|| = № ||УП,ф|| = ll<wn ||Уф,0|| = ll<wфn ||Уф,ф|| = ||^ф„

где: Е0 - 0-й компонент электрического поля на поверхности тела вращения (рис. 1); Еф - ф-й компонент электрического поля на поверхности тела вращения (рис. 1).

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

Используя представление функции Бесселя I рода целого индекса действительного аргумента интегральным представлением [7]:

tl.t'l t2,t'2 t3,t'3 t4,t'4

t,t'

t'2

б

Рис. 2. Выбор узлов численного интегрирования при вычислении элементов главной диагонали матриц а - координаты узлов численного интегрирования весовых и базисных функций совпадают. Интегрирование по t и V по 4-м точкам; б - координаты узлов численного интегрирования весовых и базисных функций не совпадают. Интегрирование по t по 4-м точкам, интегрирование по ^ по 3-м точкам.

а

Jn = Jn (к р sin 0),

а также учитывая ориентацию векторных составляющих весовых функций W в формулах (20), (21), полные выражения для элементов матрицы ||V|| можно представить в виде:

= 2njn+í J pfi(t)eJkzcose (cosGsinfl ]n+1 + p/¿(t)e^zcos0O'sine cosfl/J dt, (24)

t=о

t=T

= -2njn+1 J pfi (t)e;fcz cos 0 (cos 0 /n+1^/n~1) dt, (25)

t=о

t=T

V^. = 2njn+1 J pfi(t)eJkzcose (sinfl /n+1^/n~1) dt, (26)

t=о

t=T

V*». = +2izjn+1 J p/¿(t)^fczcos9 (/n+1 y11'1) dt, (27)

t=о

где: i - номер коэффициента возбуждения; 0 - угол между вектором uf и осью Z.

В представленных формулах VP'qi = <up, Wqn > является коэффициентом возбуждения q-го компонента n-й азимутальной гармоники i-й весовой функции касательной составляющей электрического поля p-й поляризацией первичной (падающей) ЭМВ. Верхний индекс i в up означает падающее (первичное)

введя обозначение:

t=T

электрическое поле, нижний индекс г в Wqn г означает номер весовой функции. Замена в УрАг скалярного произведения <Ер, Wn г> на (и'р, Wqn > означает использование в формулах для Рр,(1г первичного электрического поля единичной —

та

амплитуды. |

С учетом полученных ранее результатов § рассеянное осесимметричным объектом поле ^ определяется как: —

см о см

< I

со та 5

о со

.

<и

о

о <и со

см ■ч-ю

с?

см ■ч-ю см

(П (П

Е5 =

ц

4пг

п=к

II

п=—к

и5

и

е I

>ф.ф

ф|

,-}кг

(28)

где: г - расстояние от точки наблюдения до объекта, г >> R;

и0, иф - орты векторов поляризации рассеянного поля Е*, перпендикулярные направлению распространения рассеянной электромагнитной волны;

||К£,9|| - матрица-строка коэффициентов излучения р-й поляризации электрического поля q-м компонентом п-й азимутальной гармоники, элементы матрицы определяются как: Rp,qj = <ир, }%) - коэффициент излучения р-й поляризации электрического поля q-м компонентом п-й азимутальной гармоники /-й базисной функции плотности тока.

Выражение для элементов п-й гармоники поляризационной матрицы рассеяния будет иметь вид:

СР.Ч _

4тс

>р-Ф|

м

Г<М|||^<М>1

\КЧ

где: ЦУ^'Ц ||УП'ф|| ||Упфд|| и ||Уф,ф|| - блоки матрицы ||Уп||, обратной к матрице взаимных сопротивлений |^п||. Вывод формул (15)—(18) и выражений для коэффициентов возбуждения и излучения V и R аналогичен приведенному в [8].

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

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

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

Аргументом функций Бесселя в (24)-(27) является kpsin 0, таким образом необходимое количество азимутальных гармоник определяется поперечными волновыми размерами тела р и сектором углов 0, отсчитываемых от оси тела вращения, в котором необходимо вычислять ЭПР. В [1], стр. 169, предложена соответствующая формула, которую приведем в обозначениях, принятых в настоящей статье:

п = 2Ртах 0

(29)

0

2

4

6

8 10 12 14 16 18

Рис. 3. Функции Бесселя первого рода целого индекса действительного аргумента

--ЛМ;--J\(x);--J2(x)

где: п - требуемое число азимутальных гармоник;

ртах - максимальный нормированный радиус тела вращения;

0 - угол падения электромагнитной волны к оси тела вращения.

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

х

Для практической оценки точности аппроксимации (28) предлагается следующий способ.

ЭПР металлической сферы имеет одинаковое значение независимо от углов падения плоской волны. Однако если представить сферу как тело вращения, то в зависимости от угла прихода падающей волны 0 для определения рассеянного поля в соответствии с (29) потребуется учитывать разное количество азимутальных гармоник. Эталонным значением рассеянного поля можно считать случай, когда облучение происходит вдоль оси симметрии тела вращения, при этом коэффициенты возбуждения всех гармоник, за исключением п = 1, -1 равны нулю (рис. 4). Предложенный способ также применим для оценки точности аппроксимации азимутальными гармониками решений на основе ИУ II рода.

Сравнивая результаты расчета ЭПР сферы при облучении с направления 0 = 0 град. (учитываются гармоники с п = +1, -1) и при облучении с направления 0 = 90 град. (используется различное количество гармоник), можно определить их необходимое количество для фиксированных X и Ятах, обеспечивающее требуемую точность вычисления ЭПР во всем секторе углов облучения. Ниже,

на рисунках 5, 6, приведены результаты определения необходимого количества азимутальных гармоник для различных 0, Ятах, обеспечивающих точность определения ЭПР 5 и 10 %. В качестве показателя точности использовалось максимальное отклонение значения ЭПР от эталонного (при 0 = 90 град.). Также из рисунков 5, 6 можно определить размеры секторов по углам 0, в которых обеспечивается гарантированная точность вычисления ЭПР 5 и 10 % при фиксированном п. На рисунке 7 приведены значения максимальных ошибок ЭПР при ее расчете вне секторов гарантированной точности.

Как видно из рисунков 5, 6, для тела вращения с максимальным радиусом, равным одной длине волны, точность расчета ЭПР во всем диапазоне углов 0 не хуже 5 % обеспечивается при использовании 8 азимутальных гармоник падающего поля. Для обеспечения такой же точности расчета ЭПР для тела с максимальным радиусом 2,5 длины волны потребуется 17 гармоник. На рисунке 7 приведены зависимости максимальных ошибок определения ЭПР при фиксированном числе азимутальных гармоник. Данные ошибки определялись за границами гарантированных секторов расчета ЭПР, приведенных на рисунках 5 и 6,

Рис. 4. Аппроксимация падающего поля азимутальными гармониками при различных углах падения те

первичной ЭМВ Е.

100

=5 Э

Ёв т

и 3 ° й

& Й £ §

X [Й

о. м

й О

V Е «

й а

80

60

40

20

100

=5 Э я Л

и 3 ° й & й [3 ^

X г-Й

о. и

й О

V % «

й а

80

60

40

20

3 6 9 12 15 18 21 - требуемое количество азимутальных гармоник

3 6 9 12 15 18 21 - требуемое количество азимутальных гармоник

Рис. 5. Зависимость размера секторов гарантированной точности вычисления ЭПР 5% от количества азимутальных гармоник и радиуса тела вращения (сферы) для X = 1 м

Радиусы сферы, м:--Я = 0,1;--Я = 0,25;

--Я = 0,5;--Я = 1,0;--Я = 2,5

Рис. 6. Зависимость размера секторов гарантированной точности вычисления ЭПР 10% от количества азимутальных гармоник и радиуса тела вращения (сферы) для X = 1 м

Радиусы сферы, м:--Я = 0,1;--Я = 0,25;

--Я = 0,5;--Я = 1,0;--Я = 2,5

см о см

< I

со те

г

I

о ^

со

о.

<и

о

о <и со

см ■ч-ю

с?

см ■ч-ю см

(П (П

м

о

я у

£

12 10 8

6

—ЛЛ

1 у— л

V, л

3 6 9 12 15 18 21 п - требуемое количество азимутальных гармоник

Рис. 7. Зависимость максимальной погрешности расчета ЭПР дБ/м2 от количества учитываемых азимутальных гармоник и максимального радиуса тела вращения (сферы)

Радиусы сферы, м:--Я = 0,1;--Я = 0,25;

--Я = 0,5;--Я = 1,0;--Я = 2,5

то есть для случаев, когда значения ошибок принимают максимальные значения.

Как следует из рисунка 7, для достижения погрешности расчета ЭПР не хуже 3 дБ для тела с максимальным радиусом 2,5 длины волны достаточно использование 6 азимутальных гармоник.

Вычисление элементов матрицы взаимных сопротивлений по формулам (15)—(18) имеет некоторые сложности. Поскольку поверхности интегрирования 5 и S' в уравнении

(7) совпадают, у ядра ИУ (7) имеется логарифмическая особенность вида 1/Я. Как следствие, при вычислении элементов главной диагонали матриц \\Ъкп1'\\ (то есть для I = ]) при ис-

пользовании одинаковых способов численного интегрирования по 5 и координаты узлов численных интегралов по t и V в формулах (15)—(18) будут совпадать. В этих случаях Я в формуле (19) равно нулю и возникает машинное переполнение.

На практике данную проблему решают использованием формул численного интегрирования с различным количеством узлов для t и t'. Иллюстрация данной проблемы и способа ее решения приведена на рисунке 2.

В разработанной программе для вычисления (15)—(18) использовалась формула численного интеграла Гаусса в вариантах с 2, 3, 4, 5, 6, 7, 8, 9 и 10 узлами.

Для вычисления (19) были опробованы варианты численного интегрирования по формулам прямоугольников, трапеций, парабол (Симпсона), Чебышева и Гаусса. Во всех случаях необходимая точность вычисления обеспечивается подбором необходимого количества сегментов численного интеграла. Наиболее быстрая сходимость интеграла обеспечивалась при использовании формулы Гаусса.

Следует отметить, что факт наилучшей (по скорости) сходимости при использовании

0

0

п

п

4

2

0

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

Рассмотрим источники и оценки погрешностей вычисления элементов матриц ||Z||, ||V||

и ||R||.

Расчет на ЭВМ коэффициентов возбуждения ||V|| и рассеяния ||R|| осуществляется по формулам, содержащим стабильно вычисляемые элементы с требуемой точностью. Точность обеспечивается размерностью разрядной сетки ЭВМ в среде программирования (7 или 14 значащих цифр при ординарной или двойной точности, например для FORTRAN), а также точностью выполнения арифметических операций и вычисления встроенных (sin x, cos x, e~mt, x2) и внешних специальных функций (в рассматриваемом случае - Jn(x)). При численном интегрировании вдоль образующей тела вращения при ее корректном разбиении на сегменты (по переменной t) текущие значения р, sin cos 0 и z, входящие в выражения для ||R|| и ||V||, в пределах координат одной весовой или базисной функции изменяются незначительно. По этой причине значения слагаемых, входящих в численные интегралы, имеют сходные порядки. Это обеспечивает устойчивость процесса численного интегрирования, так как при суммировании величин одного порядка не происходит потерь большого числа значащих цифр.

Аналогичный вывод можно сделать и в отношении расчета элементов матрицы взаимных сопротивлений ||Z||, за исключением одного отличия - в формулах ее элементов ||Z|| используется функция Грина кольцевого источника gn(p, р', z, k), которая не вычисляется аналитически, не табулирована и не входит в перечень встроенных или внешних функций прикладных математических библиотек алгоритмических языков. Вычисление данной функции отдано на «откуп» пользователей,

при этом, как правило, используется численное интегрирование.

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

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

При вычислении функции gn(р, р', г, к) (формула (19)) подобный алгоритм не всегда себя оправдывает. Так, не удается обеспечить сходимость численного интеграла при малых р и р' и больших г. Также сходимость нарушается при вычислении gn(р, р', г, к) для гармоник с большим п. Причиной оказывается явление, известное в вычислительной математике как «катастрофическая потеря верных знаков» [9]. Проблема проявляется для любых формул численного интегрирования, в том числе при использовании вычислений с двойной точностью (14 значащих цифр).

Решением данной проблемы является:

- использование асимптотических приближений gn(р, р', к) для малых р и р' и больших г;

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

В качестве асимптотического прибли- — жения функции gn(р, р', г, к) для малых р и р' | и больших г в разработанной модели использо- | валась функция Бесселя I рода целого индекса £

см о см

< I

со та

г

о со

о.

<и

о

о <и со

см ■ч-ю

с?

см ■ч-ю см

(П (П

действительного аргумента, в которую при данных условиях вырождается gn(p, р', 2, к).

Для обеспечения необходимой чувствительности интегральной суммы при вычислении gn(р, р', к) к входящим в нее малым слагаемым был использован следующий способ. Вычисляемый интеграл представлялся суммой отдельно вычисляемых интегралов (слагаемых) с интервалами интегрирования по ф, кратными п/2п. В этих условиях значения cosnф в пределах одного слагаемого содержит величины одного знака; как следствие, частные интегральные суммы и вычисляемое значение gn(р, р', 2, к) оказываются более устойчивыми к катастрофической потере верных знаков.

Использованные в разработанной модели вышеописанные способы улучшения сходимости вычислений gn(р, р', 2, к) обеспечили интегрирование с требуемой точностью при использовании не более 3-х итераций.

Ниже приводятся результаты расчетов тестовых задач с помощью разработанной модели. На рисунках 8 и 9 приведены графики моностатической ЭПР кругового цилиндра радиусом 7,5 см и длиной боковой поверхности 25 см на длине волны = 15 см. Рассматриваемый пример взят из [10]. Приведенное на этих же

о -10 -20 -30 -40

эпр ее, дБ/м2

50

100

150 200 Угол е, град.

250

300 350

Рис. 8. Результаты расчета ЭПР кругового цилиндра, дБ/м2, радиусом 7,5 см и длиной боковой поверхности 25 см на длине волны 15 см, сравнение с данными из [10]. Наложение результатов расчета на оригинальный график из [10]. В расчетах использовалось 13 азимутальных гармоник, 99 - поляризация

--эксперимент [10];--расчет методом ИУ

[9];--расчет методом ИУ по методике модели,

приведенной в настоящей статье

графиках сравнение с [10] показывает удовлетворительное совпадение по угловым направлениям максимумов и провалов индикатрис, а также приемлемое совпадение по их уровням. Предположительно, минимальные значения приведенных в [10] экспериментальных данных ограничены динамическим диапазоном или фоном измерительного оборудования на уровне -35 дБ.

На рисунке 10 приведены результаты расчета ЭПР методом ИУ модели кругового конуса радиусом р = 3,5 и углом при вершине 20 = 23 градуса. Края и носик конуса скруглены с г = 0,01, то есть 0,00286р. Также приведен тестовый пример - расчет ЭПР данного конуса методом ИУ, взятый из [1], стр. 178, носик и края конуса скруглены радиусом гкрая = гнос = 0,0125р, при этом автор в [1] приводит данные, что дальнейшее уменьшение радиуса закругления носика и края не приводят к изменению рассеянного поля и, следовательно, к изменению ЭПР.

На рисунке 11 приведено сравнение результатов расчета ЭПР кругового конуса с углом при вершине 15 град. и ка = 3,08 (левая половина графика) с результатами расчета и измерений ЭПР такого же конуса, опубликованных

0

ЭПР фф, дБ/м2

-40

50 100 150 200 250 Угол град.

300 350

Рис. 9. Результаты расчета ЭПР кругового цилиндра, дБ/м2, радиусом 7,5 см и длиной боковой поверхности 25 см на длине волны 15 см, сравнение с данными из [10]. Наложение результатов расчета на оригинальный график из [10]. В расчетах использовалось 13 азимутальных гармоник, фф - поляризация.

--эксперимент [10];--расчет методом ИУ

[9];--расчет методом ИУ по методике модели,

приведенной в настоящей статье

0

0

10

£

С -10

-20

-30

/-

\\ \\ \\ \\ t \ / \ /' \ i \

Y 'íÚ / MV '' |l v^ 1 1 Kv ' \ ' у

í i w ч u " и » V

и i i

10

£

&

is -10

-20

-30

\ \ \ \ f- 0^ iL .................... a

\ \ f > \ /У> / 1 \ ¿j Ii V 1 1 V'ч - ^..... / \ T

\l \ J

0 30 60 90 120 150 180 угол град.

а

0 30 60 90 120 150 180 0, град. б

20 10 0 -10 -20 -30

Рис. 10. Расчет ЭПР конуса: а — методом ИУ, дБ/м2, б — сравнение с результатом из [1]

офф, дБ 10 log {-т

180 150 120 90 60 угол град.

а00, дБ

30

10 0 -10 -20 -30

\ yrtttttfix / <1 Vßf >Ül*

V nil N V

20 10 0 -10 -20 -30

0 30 60 90 120 160 ф, град 10 log (J)

180 150 120 90 60 угол град.

а

30

10 0 -10 -20 -30

A

V" SiiííM1 ' \ Äif

V / Ii/ 1

1

0 30 60 90 120 160 ф, град б

Рис. 11. ЭПР конуса, дБ/м2, рассчитанная методом ИУ на разработанной модели (а) в сравнении с приведенной в [11] (б)

в [11] (правая половина графика). Результаты расчета лучше совпадают с экспериментальными данными из [11] по сравнению с расчетными данными из [11]. Основные отличия — в глубине «провалов» индикатрис рассеяния, что, по-видимому, объясняется большим шагом дискретизации расчетов, приведенных в [11] (из-за чего точки «провалов» могут быть пропущены), либо недостаточным количеством использованных азимутальных гармоник. Меньший по сравнению с рассчитанным

уровень «провалов» в экспериментальных данных [11], по-видимому, можно объяснить ограничением динамического диапазона использованной для измерений безэховой камеры на уровне —25—28 дБ либо наличием «фона» камеры на этом уровне.

Метод ИУ может успешно использоваться для решения электродинамических задач в случаях, когда применение МКВ и ГТД невозможно. В качестве примера ниже приводятся результаты расчета ЭПР конуса с различными

<я

s ф

0

0

0

0

радиусами закругления ребра «боковая поверхность - дно» вплоть до геометрии - «капля». На рисунке 12 приведена геометрия конуса с углом при вершине 20 носа = 15 град., радиус закругления носика конуса (г носа) постоянный и равен 3 мм, радиус закругления стыка дугового ребра «боковая поверхность - дно» (г дна) изменяется от 0,04 до 23 см (значение г дна = 23 см соответствует геометрии «капля»).

На рисунках 13.1-13.10 приведены моностатические индикатрисы рассеяния (ЭПР) данного конуса для длины волны X = 0,1 м для различных радиусов закругления дна (г дна). На всех рисунках направление 0 = 180 град. соответствует облучению конуса с «носика», направление 0 = 0 град. соответствует облучению конуса со стороны «дна».

На рисунке 14 приведено полученное на основании данных рисунков 13.1-13.10 изменение моностатической ЭПР рассматриваемого конуса от радиуса ребра «боковая

поверхность - дно» при осевом облучении - со стороны «носика» конуса и со стороны «дна».

Как видно из полученной зависимости, моностатическая ЭПР конуса при осевом облучении со стороны «носика» в существенной степени зависит от радиуса скругления («стыка») боковой поверхности и донной части. Это обстоятельство объясняется известным фактом - отраженная от конуса волна при его облучении со стороны носа в значительной степени формируется за счет интерференции волн, отраженных от «стыка», и волн, огибающих конус вдоль его образующей. Амплитуда данных волн может значительно превышать амплитуду прямых отражений от носовой части конуса. Амплитуды волн, отраженных от «стыка», уменьшаются при увеличении радиусов скругления.

Приведенные в данном разделе результаты расчетов и сделанные на их основании выводы не являются новыми, а иллюстрируют

г носа

см о см

< I

со та

г

I

о ^

со

.

<и

о

о <и со

см ■ч-ю

с?

см ■ч-ю см

(П (П

Рис. 12. Геометрия конуса, г носа = 0,003 м, 20 носа = 15 град.

20 10 0 -10 -20 -30 -40 -50

ЭПР аОО, афф, дБ

г дна = 0,004 м, X = 0,1 м

Л

\f\f-i 1Л

П | /

1 1 1 \ /

0 30 60 90 120 угол О, град.

150 180

Рис. 13.1. ЭПР конуса, дБ/м2, гноса = 0,003 м, гдна = 0,004 м, X = 0,1 м --ЭПР 99;--ЭПР фф

20 10 0 -10 -20 -30 -40 -50

ЭПР аОО, афф, дБ

г дна = 0,006 м, X = 0,1 м

0 30 60 90 120 угол О, град.

150 180

Рис. 13.2. ЭПР конуса, дБ/м2, гноса = 0,003 м, гдна = 0,006 м, X = 0,1 м

дна

--ЭПР 99;--ЭПР фф

BfA BIB

20 10 0 -10 -20 -30 -40 -50

ЭПР а00, афф, дБ

r дна = 0,009 м, X = 0,1 м

0 30 60 90 120 угол 0, град.

150 180

Рис. 13.3. ЭПР конуса, дБ/м2, гноса = 0,003 м, гдня = 0,009 м, 1 = 0,1 м

дна 7

--ЭПР 00;--ЭПР фф

20 10 0 -10 -20 -30 -40 -50

ЭПР а00, афф, дБ

r дна = 0,02025 м, X = 0,1 м

0 30 60 90 120 угол 0, град.

150 180

Рис. 13.5. ЭПР конуса, дБ/м2, гноса = 0,003 м,

„ = 0,1 м - ЭПР фф

гдна = 0,02025 м, 1 = 0,1 м

20 10 0 -10 -20 -30 -40 -50

--ЭПР

ЭПР а00, афф, дБ

r дна = 0,0455625 м, X = 0,1 м

0 30 60 90 120 угол 0, град.

150 180

Рис. 13.7. ЭПР конуса, дБ/м2, гноса = 0,003 м, гдна = 0,0455625 м, 1 = 0,1 м --ЭПР 00;--ЭПР фф

20 10 0 -10 -20 -30 -40 -50

ЭПР а00, афф, дБ

r дна = 0,0135 м, X = 0,1 м

30 60 90 120 угол 0, град.

150 180

Рис. 13.4. ЭПР конуса, дБ/м2, гноса = 0,003 м,

= 0,1 м - ЭПР фф

гдна = 0,0135 м, 1 = 0,1 м

дна 7

20 10 0 -10 -20 -30 -40 -50

- ЭПР

ЭПР а00, афф, дБ

r дна = 0,030375 м, X = 0,1 м

Л

W у Л л \л

п лЛл 0 /7

Ii I г 1 v'

f

0 30 60 90 120 угол 0, град.

150 180

Рис. 13.6. ЭПР конуса, дБ/м2, гноса = 0,003 м, гдна = 0,030375 м, 1 = 0,1 м --ЭПР 00;--ЭПР фф

20 10 0 -10 -20 -30 -40 -50

ЭПР а00, афф, дБ r дна = 0,06834375 м, X = 0,1 м

0 30 60 90 120 угол 0, град.

150 180

Рис. 13.8. ЭПР конуса, дБ/м2, гноса = 0,003 м, гдна = 0,06834375 м, 1 = 0,1 м --ЭПР 00;--ЭПР фф

та

та

s ф

0

20

ЭПР о00, Офф, дБ

r дна = 0,102515625 м, X = 0,1 м

см о см

< I

to те

s

о

CÛ

о.

V

о

о ф

CÛ

см ■clin

с?

см ■clin см

(П (П

20

ЭПР аОО, афф, дБ

r дна = 0,102515625 м, X = 0,1 м

-30 -40

0 30 60 90 120 угол 0, град.

150

180

Рис. 13.9. ЭПР конуса, дБ/м2, гноса = 0,003 м, гдня = 0,102515625 м, X = 0,1 м

дна ?

--ЭПР 00;--ЭПР фф

ЭПР а носа, а дна, дБ

X = 0,1 м

20

10

10

-20

-30

-40

-50

-40

30 60 90 120 угол О, град.

150

180

Рис. 13.10. ЭПР конуса, дБ/м2, гноса = 0,003 м,

м, X = 0,1 м - ЭПР фф

гдна = 0,153773438 м, X = 0,1 м

дна ?

--ЭПР 00;

0,0030 0,0090 0,0203 0,0456 0,1025 г дна - радиур ребра (стыка) «дно - боковая образующая конуса»

Рис. 14. Зависимость моностатической ЭПР конуса

при осевом облучении, дБ/м2, геометрия рис. 12, от радиуса скругления стыка дна с боковой образующей, X = 0,1 м

--ЭПР, отражение со стороны носа;

--ЭПР, отражение со стороны дна

возможности метода ИУ и возможности разработанной модели для решения прикладных задач электродинамики.

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

позволяют получать устойчивые в вычислительном отношении решения для широкого круга практических задач. Предложенный способ оценки точности аппроксимации решений конечным рядом азимутальных гармоник применим как к ИУ I, так и к ИУ II рода. Приведенные примеры решения тестовых задач и сравнение с известными результатами иллюстрируют их хорошее совпадение. Данные примеры рассчитывались на персональной ЭВМ общего пользования, не содержащей спецвычислителей, при этом время расчетов составило от единиц минут до нескольких десятков минут в зависимости от электрического размера объекта и сектора исследуемых углов, что определяется количеством используемых азимутальных гармоник. Модель обеспечивает построение геометрий широкого класса объектов - тел вращения, образующая которых аппроксимируется набором отрезков и дуг различных размеров и кривизны. Использование метода ИУ позволяет получить в резонансном диапазоне точность решения электродинамических задач, недостижимую для МКВ и ГТД.

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

1. Vasilyev Е. N. Vozbuzhdeniye tel vrashch-eniya. М.: Radio i svyaz, 1987. 272 s. (Russian).

2. Borovikov E. A., Kinber V. E. Geometrical theory of diffraction. M.: Svyaz, 1978, 248 p. (Russian).

0

0

3. Ufimtsev P. Ya. Method of edge waves in the physical theory of diffraction. M.: Sovetskoe radio, 1962, 244 p. (Russian).

4. Markov G. T., Vasilyev E. N. Matematicheskiye metody prikladnoi elektrodinamiki. M. Sovetskoe radio, 1970, str. 11. (Russian).

5. Muath Gouda. The Method of Moment for the Electromagnetic Scattering from Bodies of Revolution Final Master's Degree Thesis Subject Category: Electrical Engineering Series Number: 4/ 2008., 74 p.

6. Ra'ed A. Malallah, Dr. Zaki A. Ahmed & Dr. Ahmed H. Abood. Studying the effect of conducting bodies of revolution (bor) shape on computing radar cross section (RCS). Journal of Missan Researchers, Vol (6), No (11), 2009, p. 45.

7. Sveshnikov A. G., Bogolyubov A. N., Kravtsov V. V. Lektsii po matematicheskoi fiz-ike Uchebnoye posobiye. M.: Izdatelstvo MGU, 1993, str. 75. (Russian).

8. J. R. Mautz, R. F. Harrington. Radiation and scattering from bodies of revolution. Appl. Sci. Res. 20, June 1969, p. 405-435.

9. Forsythe G., Malcolm M., Moler C. Computer methods for mathematical computations. Transl. from Eng. by Kh. D. Ikramov. M.: Mir, 1980, 277 p., p. 28.

10. Setukha A. V., Aparinov A.A., Stavtsev S. L., Fetisov S. N. Superkompyuternye tekhnologii pri realizatsii metoda granichnykh integralnykh uravneniy v zadachakh elektromagnitnogo ras-seyaniya. MGU im. M. V. Lomonosova, TsAGI im. N. E. Zhukovskogo, IVM RAN, OKB im. A. Lulki. Prezentatsiya seminara 26 noyabrya 2019. http://agora.guru.ru/ sct/files/2019_10_01_ Setuha_p1.pdf, 51 s. (Russian).

11. Shustikov V Yu. Software package for calculating the reflectance profile of axisymmetrically-shaped targets in resonant and high frequency wavelength ranges. Journal of "Almaz - Antey" Air and Space Defence Corporation, 2018, No. 2, p. 75-81.

Об авторе

Рунов Александр Адольфович - кандидат технических наук, ведущий научный сотрудник публичного акционерного общества «Межгосударственная акционерная корпорация «Вымпел», Москва, Российская Федерация. Область научных интересов: электродинамика, радиолокация.

Modelling of the reflectance profile of axisymmetrical radiolocation objects using the integral equation method

Runov А. А.

Public Joint Stock Company Vympel Interstate Corporation, Moscow, Russian Federation

Based on the first-kind integral equation method for the electric field, the procedure and software for calculating the radar cross-section of axisymmetrical objects, bodies of revolution, are developed. Algorithms are proposed for computation of the matrix of mutual impedances and Green's function of a ring source providing the computation accuracy required for obtaining a stable solution. The method of solution approximation accuracy evaluation by azimuthal harmonics is proposed. Comparison with test examples is carried out and the applicability for solving real-world problems is shown.

Keywords: electromagnetic wave diffraction, radar cross-section, integral equation method

Information about the author |

Runov Aleksandr Adolfovich - Candidate of Engineering Sciences, Leading Researcher, Public Joint Stock Company i

0)

Vympel Interstate Corporation, Moscow, Russian Federation. ¡5

Science research interests: electrodynamics, radiolocation. ^