Радиотехника и электроника, 2020, T. 65, № 1, стр. 31-48

1. ПОСТАНОВКА ЗАДАЧИ

Интерес к фазированным антенным решеткам (ФАР) с большим коэффициентом усиления (КУ) обусловлен их широким использованием в различных системах локации, связи, диагностики, научных исследованиях и т.д. Реализация высокого КУ неизбежно связана с увеличением числа элементарных излучателей (ЭИ) и, следовательно, электрических размеров ФАР. Несмотря на бурное развитие вычислительной техники и программного обеспечения, электродинамический анализ таких структур представляет сложную задачу, требующую больших затрат компьютерных ресурсов.

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

При всех достоинствах данного подхода нельзя не отметить его недостатки. К их числу относятся:

— невозможность расчета диаграммы направленности (ДН) ФАР, которая для бесконечной структуры вырождается в дельта функцию;

— невозможность оценки влияния краевых эффектов на параметры решетки;

— невозможность анализа чувствительности характеристик ФАР по отношению к отклонениям ее параметров от номинальных значений;

— невозможность статистического анализа характеристик ФАР.

Говоря о недостатках, мы имеем в виду возможность выполнения разных видов анализа на электродинамическом уровне строгости в рамках модели бесконечной решетки. Безусловно, используя разные приближенные методы, можно оценить ДН структуры и, например, влияние краевых эффектов. Однако, подобные оценки всегда будут приближенными и их использование для оценки таких “тонких” эффектов как уровень дальних боковых лепестков вызывает оправданные опасения.

Поэтому, не ставя под сомнение полезность модели ФАР в виде канала Флоке, мы должны отметить ограниченность области ее применения. В связи с этим актуальной задачей остается электродинамический анализ решеток с большим КУ. Как отмечалось выше, непосредственное применение таких методов электродинамического моделирования, как метод конечных элементов [2], метод интегральных уравнений [3], метод конечных элементов во временной области, не является оптимальным подходом, так как получаемые при их помощи вычислительные алгоритмы требуют чрезмерно больших затрат компьютерных ресурсов. Увеличение времени решения граничной задачи делает проблематичным или, как минимум, весьма затратным выполнение этапов инженерного проектирования, связанных с перебором большого числа вариантов.

В связи с ситуацией, описанной выше, в начале ХХ столетия был предложен ряд новых методов моделирования многоэлементных структур больших электрических размеров. Можно отметить ряд подходов, основанных на методе моментов [3]. Основным из них является метод быстрых мультиполей [4] и его развитие – многоуровневый быстрый мультипольный алгоритм [5, 6]. Этот эффективный метод расчета матричных форм используется во многих других методах, таких как метод сложных мультиполей [6], метод локализации матрицы импеданса [7], адаптивный метод [8], многоуровневый алгоритм матричной декомпозиции [8] и метод интегральных уравнений в сочетании с быстрым преобразованием Фурье [10].

Существуют методы, направленные на сокращение числа неизвестных без значительной потери точности. Концепция сокращения размерности матричных уравнений преобладает в недавно разработанных итерационных алгоритмах, например, в методе характеристических базисных функций [11], методе синтетических базисных функций [12], методе базисных функций подобластей [13], методе собственных токов [14] и многоуровневом методе подобластей [15]. Эти методы имеют общую основу в виде применения обобщенных базисных функций для представления поверхностных электрических и/или магнитных токов. Разница между ними прежде всего в способах определения базисных функций.

В ряде методов электродинамического анализа ФАР больших электрических размеров используется в качестве основы решение для бесконечной ФАР, которое затем разными приближенными способами адаптируется таким образом, чтобы новое решение учитывало конечные размеры структуры. Так, в работе [16] коэффициенты связи между элементами ФАР рассчитываются для бесконечной решетки, а затем их используют для анализа конечной структуры. Более строгий подход аналогичного типа предложен в работе [17]. Он может быть эффективно использован для определения ДН элемента ФАР ограниченных размеров. Однако возможность определения коэффициентов связи в рамках данного метода остается неопределенной.

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

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

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

Назовем группу методов, реализующих описанную выше схему решения граничной задачи методами ОМР (МОМР). Наиболее близким подходом к данной работе является вариант МОМР, в котором в качестве возбуждающих и рассеянных волн используются векторные сферические гармоники (ВСГ). Идея разложения по ВСГ (РВСГ) антенных полей была использована в работах [20–23]. При этом необходимо отметить, что развитые в них методы использовались не только для анализа ФАР, но и для решения задач рассеяния на апериодически расположенных одинаковых объектах. Подходы, использованные в [23], наиболее близки к данной работе. Отметим эффективное применение дополнительной теоремы сложения для сферических функций при решении задачи представления полей в смещенных друг относительно друга сферических системах координат.

Нашей задачей является развитие МОМР – РВСГ применительно к важному классу объектов – щелевым и ленточным антенным решеткам с большими электрическими размерами. В первую очередь, с практической точки зрения нас интересуют щелевые решетки. Ленточные структуры представляют методический интерес, так как их удобно использовать в качестве тестовых объектов для сравнения результатов, полученных разными методами.

Отметим, что щелевые решетки активно развиваются, в том числе и в настоящее время [24, 25]. Интерес к ним обусловлен простотой конструкции и технологичностью. В последнее десятилетие дополнительный интерес к щелевым решеткам связан с созданием антенн, сфокусированных в ближней зоне, которые находят применение в разных областях науки и техники. Сфокусированные волноводно-щелевые решетки [26, 27] относятся к числу апериодических решеток. Их электродинамический анализ особенно сложен, поскольку использование модели бесконечной решетки в данном случае невозможно. Единственным средством моделирования в такой ситуации оказывается численное решение граничной задачи для всей решетки со всеми отмеченными выше недостатками. Поэтому разработка новых эффективных методов анализа щелевых ФАР является актуальной задачей.

2. ВЕКТОРНЫЕ СФЕРИЧЕСКИЕ ГАРМОНИКИ

Определение ВСГ следует из общего решения уравнения Максвелла в сферических координатах $r,\theta ,\varphi $ [28]:

где $k$ – волновое число свободного пространства, $\rho $ – безразмерный электрический радиус, $P_{n}^{m}$ – присоединенная функция Лежандра, ${{\psi }_{n}}(\rho ),$ $\xi _{n}^{{(2)}}(\rho )$ – сферические функции Бесселя и Ганкеля второго рода (индекс $n$ меняется от нуля до бесконечности, а индекс $m$ от $ - n$ до $n$). Верхние индексы $e,m$ соответствуют электрическим и магнитным волнам, нижние индексы задают тип волны.

Назовем ВСФ в верхней формуле (2) рассеянными волнами, а в нижней – возбуждающими. Рассеянные волны описывают поле при радиусах бóльших радиальных координат создающих их источников, а возбуждающие волны соответственно при меньших радиусах. Для обозначения возбуждающих и рассеянных волн мы будем использовать нижние индексы $i,s$ от английских терминов incident и scattered.

Соотношения (1) записаны для нормированных компонент поля ${\rm E},{\rm H}$. Они связаны с компонентами векторов $E,H$ следующим образом:

где ${{W}_{0}}$ – волновое сопротивление свободного пространства.

Отметим, что нормированные ВСГ не имеют размерности и зависят только от угловых координат и электрического радиуса. Кроме того, важно указать, что у электрических волн равна нулю компонента поля ${{H}_{r}}$, а у магнитных – ${{E}_{r}}$.

Сферические гармоники удовлетворяют соотношениям ортогональности:

где $N_{{rn,m}}^{{e,m}}$ – нормы радиальных компонент, ${{\delta }_{{n,p}}}$ – символ Кронекера.

При анализе щелевых и ленточных решеток при $\theta = {\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$ нам будет удобно разместить плоскость симметрии, которая для щелевых структур является идеальной электрической стенкой, а для ленточных идеальной магнитной стенкой. Поэтому целесообразно разделить ВСФ на четные и нечетные относительно указанной плоскости:

Верхние индексы $e{\text{,}}o$ соответствуют четным и нечетным волнам (от английских слов even, odd – четный, нечетный). Четные волны удовлетворяют в плоскости $\theta = {\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}$ граничным условиям на магнитной стенке:${{H}_{\tau }} = 0,$ а нечетные на электрической: ${{E}_{\tau }} = 0,$ где индекс $\tau $ соответствует тангенциальным к плоскости компонентам.

Четным электрическим и магнитным волнам соответствуют нечетные и четные значения $n + \left| m \right|\,:$ $n + \left| m \right| = 2l + 1,$ $n + \left| m \right| = 2l,$ а нечетным электрическим и магнитным волнам соответствуют четные и нечетные значения $n + \left| m \right|$: $n + \left| m \right| = 2l,$ $n + \left| m \right| = 2l + 1,$ $l = 0,1,....$ Данные соотношения следуют из свойств присоединенных функций Лежандра [29].

3. СХЕМА МОМР–РВСГ

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

Положение ЭИ в решетке обычно описывают при помощи двух индексов, например, $n,m.$ Такая нумерация удобна в случае периодической структуры. Мы будем использовать сквозную нумерацию ЭИ при помощи индекса $\nu $, который меняется от единицы до ${{M}_{R}}$, где ${{M}_{R}}$ – общее число ЭИ в решетке. Также используем один индекс $N$для задания типа ВСГ вместо двух индексов $n,m,$ которые применялись в формулах (1)–(5).

Введем в пространстве ВСГ следующие векторы:

Здесь $\vec {E}_{{i,s\nu ,N}}^{{e,m}}$ – вектор электрического поля сферической волны электрического или магнитного типа (верхние индексы $e,m$), возбуждающей или рассеянной (нижние индексы $i,s$) заданной в локальной сферической системе координат, связанной с центром ЭИ с номером $\nu $. Индекс $N$задает номер волны. Аналогично задаются векторы магнитного поля $\vec {H}_{{i,s\nu ,N}}^{{e,m}}.$ Величины $A_{{i,s\nu ,N}}^{{e,m}}$ – это амплитуды возбуждающих и рассеянных электрических и магнитных волн ЭИ. Нам необходимо отличать векторы в пространстве ВСГ и вектора отраженных/падающих волн от векторов, заданных в физическом пространстве. Для обозначения векторов первых типов используем символы вида ${\mathbf{A}}$, а для обозначения векторов физического пространства используем стрелки $ \to $. В этом случае символ ${\mathbf{\vec {E}}}$ обозначает вектор в пространстве ВСГ, элементами которого являются векторы из физического пространства, символ ${\mathbf{A}}$ описывает вектор из пространства ВСГ, элементами которого являются скалярные величины. Произведение вида ${\mathbf{A\vec {E}}}$ – это скалярное произведение в пространстве ВСГ, результатом которого является вектор из физического пространства.

Каждый ЭИ окружен сферой с электрическим радиусом $R$, которая не пересекает ни выделенный ЭИ, ни соседние с ним ЭИ. В случае щелевых излучателей сферы заменяются полусферами, что не вносит принципиальных изменений в схему метода. Источником возбуждающих волн для $\mu $-го ЭИ являются волны, рассеянные всеми элементами решетки, за исключением μ-го. Источником рассеянных волн этого излучателя является дифракция возбуждающих волн на нем, а также падающие волны на его входы. Мы предполагаем, что таких входов может быть не один, а несколько. Объединяя амплитуды падающих волн ${{U}_{{i\mu ,j}}}$ ($j$ – номер входа), мы получаем вектор падающих волн μ-го ЭИ: ${{{\mathbf{U}}}_{{i\mu }}}$. Аналогично можно получить вектор отраженных волн ${\mathbf{U}}_{{r\mu }}^{{}}$.

Отметим, что при отсутствии пересечений сфер, окружающих ЭИ, они оказываются расположенными вне областей, занятых источниками, что позволяет нам представить поля на них в виде разложений по ВСГ:

Следующий этап – определение связи между амплитудами волн, рассеянных ν-м ЭИ, и амплитудами волн, падающих на μ-й ЭИ. Для определения указанной связи нам необходимо осуществить разложение поля рассеянной волны ν-го ЭИ в системе координат μ-го ЭИ. Воспользуемся следующими соотношениями, связывающими координаты двух систем:

Здесь под $\Delta x,{{y}_{{\mu ,\nu }}}$ мы понимаем смещение ν-го ЭИ относительно μ-го соответственно вдоль осей 0x и 0у. Их можно также выразить в полярных координатах:

Также нам необходимо учесть поворот векторов при переходе от одной системы к другой. Матрица ${{\hat {M}}_{\nu }}$ связывает компоненты вектора в декартовой системе координат с компонентами в сферической системе ν-го ЭИ:

Поскольку переход от одной декартовой системы координат к другой состоит в ее параллельном переносе, который не меняет ориентацию вектора, то при помощи соотношения (10) можно получить матрицу ${{\hat {M}}_{{\mu ,\nu }}},$ связывающую компоненты векторов в двух сферических системах координат:

Полученные выше соотношения позволяют нам записать падающее на μ-й ЭИ поле:

Суммирование в (12) ведется по всем элементам кроме μ-го. Соотношения (12) описывают радиальные компоненты электрического и магнитного полей. Они являются функциями сферических координат ${{r}_{\mu }},{{\theta }_{\mu }},{{\varphi }_{\mu }}$ согласно выражениям (8). Вместе с тем радиальные компоненты можно представить в виде разложений по возбуждающим волнам μ-го ЭИ:

Приравнивая (12) и (13) и используя ортогональность ВСГ, находим векторы ${\mathbf{A}}_{{i\mu }}^{{e,m}}$:

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

С учетом соотношений (15) запишем выражения (14) в компактном виде

Матрицы ${\mathbf{K}}$ описывают искомую связь между амплитудами рассеянных на ν-м и возбуждающих на μ-м элементах волн. Назовем их матрицами связи.

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

где ${\mathbf{L}}$ – искомые операторы рассеяния.

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

которые позволят нам записать уравнения (17) в более простом виде:

Уравнения (19) представляют собой систему линейных алгебраических уравнений (СЛАУ), которую можно решать численно любым из известных методов. Если векторы ${{{\mathbf{U}}}_{{i\mu }}}$ известны, то достаточно решить первое из уравнений (19). Решение будет иметь вид линейной комбинации векторов ${{{\mathbf{U}}}_{{i\mu }}}$:

Подставим соотношение (20) во второе уравнение (19) и найдем векторы отраженных волн ${\mathbf{U}}_{{r\mu }}^{{}}$:

Из соотношения (21) несложно найти матрицу рассеяния всей решетки:

Матрица рассеяния ${{{\mathbf{S}}}_{a}}$ имеет блочную структуру. Ее размерность равна ${{M}_{R}} \times {{M}_{R}},$ а элементами являются матрицы размерностью $J \times J,$ где $J$ – число входов ЭИ.

Система (19) может быть решена итерационно по следующей схеме:

где $q$ – номер итерации. Расчеты показали, что данный процесс сходится весьма быстро. Отметим, что такое решение не требует вычисления обратной матрицы.

Для определения поля в дальней зоне при $\rho \to \infty $ воспользуемся асимптотическим представлением сферической функции Ганкеля:

Также отметим, что в дальней зоне нас интересуют главные члены асимптотических представлений поля. В рамках такого приближения производная по координате $\rho $ в выражениях (1) эквивалентна умножению на $ - i$. С учетом этих замечаний получаем асимптотические выражения для полей ВСГ в дальней зоне:

В выражениях (25) опущен общий множитель ${{\exp ( - i\rho )} \mathord{\left/ {\vphantom {{\exp ( - i\rho )} \rho }} \right. \kern-0em} \rho }.$

Далее аналогично формуле (6) формируем векторы ${\mathbf{\vec {E}}}_{{s\nu }}^{a},$ ${\mathbf{\vec {H}}}_{{s\nu }}^{a},$ являющиеся функциями угловых координат $\theta ,\varphi $:

Используя выражения (26), запишем поле излучения решетки в дальней зоне:

Соотношения (27) будут использованы нами для определения ДН решетки.

Результаты данного раздела инвариантны по отношению к виду ФАР. Смена типа излучателя приводит только к изменению вида оператора рассеяния ${\mathbf{L}}$, не меняя общей схемы решения. Видно, что МОМР-РВСГ включает два основных этапа: определение матрицы связи и определение оператора рассеяния. Рассмотрим их подробнее.

4. ВЫЧИСЛЕНИЕ МАТРИЦЫ СВЯЗИ

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

Во-первых, необходимо отметить, что на этапе вычисления матриц ${\mathbf{K}}_{{\mu ,\nu }}^{{}}$ удобно использовать нормированные ВСГ (1), так как они не зависят от частоты, что позволяет уменьшить число свободных параметров. Все расчеты можно проводить для электрических расстояний. Переход к ненормированным гармоникам с помощью соотношений (3) не составляет труда.

Во-вторых, следует отметить, что элементы матриц ${\mathbf{K}}_{{\mu ,\nu }}^{{}}$ зависят от угла смещения $\Delta {{\varphi }_{{\mu ,\nu }}}$ (10) тривиальным образом. Допустим мы нашли при $\Delta {{\varphi }_{{\mu ,\nu }}} = 0$ и некотором радиусе смещения $\Delta {{r}_{{\mu ,\nu }}}$ элемент ${{{\mathbf{K}}}_{{\mu ,\nu ,M,N}}}(0),$ где $N$ – номер рассеянной гармоники, а $M$ – номер падающей гармоники. Тогда тот же элемент при произвольном угле и том же значении $\Delta {{r}_{{\mu ,\nu }}}$ можно найти следующим образом:

где ${{m}_{{N,M}}}$ – азимутальные индексы ВСГ. Доказательство соотношения (28) можно выполнить, используя известные свойства экспоненциальных функций, описывающих зависимость ВСГ от координаты $\varphi $.

Таким образом, матрицы ${{{\mathbf{K}}}_{{\mu ,\nu }}},$ по существу, являются функциями одной переменной – электрического расстояния. Такое положение служит основанием для использования следующего алгоритма вычисления матрицы связи. Отметим, что расчет матриц ${{{\mathbf{K}}}_{{\mu ,\nu }}}$ привязан к конкретной решетке, геометрия которой задает значения расстояний между элементами. При большом их числе приходится многократно проводить вычисление матриц связи. Важно также иметь в виду, что выполнение расчетов по формулам (15) является достаточно трудоемким процессом, требующим значительных затрат компьютерных ресурсов. Поэтому необходимо максимально его оптимизировать.

Для этого вместо множества матриц ${\mathbf{K}}_{{\mu ,\nu }}^{{}}$ мы рассматривали матричную функцию ${\mathbf{K}}(\Delta r),$ которая находится при $\Delta \varphi = 0.$ Ее рассчитывали при фиксированных значениях электрического расстояния $\Delta {{r}_{s}}$ с шагом, достаточным для построения с приемлемой точностью аппроксимации функции ${\mathbf{K}}(\Delta r)$:

При малых расстояниях между ЭИ элементы матричной функции быстро меняются при изменении аргумента, поэтому на участке $\Delta {{r}_{s}} \leqslant 2\pi $ шаг изменения расстояния уменьшен: при $\Delta {{r}_{s}} > 2\pi $ он равен ${\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}.$ Максимальное значение электрического расстояния равно $20\pi .$ Для вычисления матричной функции ${\mathbf{K}}(\Delta r)$ при больших расстояниях удобно воспользоваться ее асимптотическим представлением при $\Delta r \gg 1$:

Коэффициенты ${{\left( {a,b} \right)}_{{M,N}}}$ можно найти аналитически, применяя асимптотические представления сферических функций Ганкеля. На практике удобнее оказался численный вариант, при котором указанные постоянные определяли из условий совпадения асимптотических ${\mathbf{K}}_{{M,N}}^{a}(\Delta r)$ и точных функций ${\mathbf{K}}_{{M,N}}^{{}}(\Delta r)$ при $\Delta r = 10,20\pi .$

Если решетка составлена из одинаковых ЭИ с одинаковыми операторами рассеяния ${{{\mathbf{L}}}_{\mu }},$ то лучшие результаты с точки зрения скорости вычислительного процесса дает аппроксимация не матрицы связи, а матриц взаимодействий ${\mathbf{W}}$, являющихся произведением матрицы связи и операторов рассеяния:

Матрицы взаимодействий напрямую входят в СЛАУ (19), описывающую решетку. При анализе решеток щелевых и ленточных излучателей использовалась сплайн-аппроксимация элементов матриц (31).

Независимо от способа аппроксимации переход к матричной функции ${\mathbf{K}}(\Delta r)$ дает существенное повышение эффективности метода. Массивы значений матричной функции рассчитываются однократно и сохраняются в виде файлов данных. Как отмечалось выше, эти данные инвариантны к типу излучателя и конфигурации решетки. При расчетах частных структур файлы считываются, преобразуются, и на их основе строятся те или иные аппроксимации, по которым определяются матрицы ${\mathbf{W}}_{{\mu ,\nu }}^{{}}$, входящие в СЛАУ (19). Таким образом, матричная функция K(Δr) играет роль специальной функции метода МОМР-РВСГ, которая однократно вычисляется и далее без существенных затрат ресурсов многократно используется при анализе ФАР.

Важным вопросом является выбор числа ВСГ. При этом следует иметь в виду, что число рассеянных и возбуждающих гармоник может быть различным. Число возбуждающих волн выбираем так, чтобы обеспечить необходимую точность представления поля на сфере, окружающей ЭИ. Наибольшее число гармоник необходимо для минимального расстояния $\Delta r$, которое в нашем случае равно π.

Число рассеянных волн существенно зависит от вида излучателя. Оно определяется из условия достижения приемлемой точности описания поля на сфере ЭИ суммой ВСГ при его возбуждении по входу. Числом учитываемых угломестных гармоник ${{N}_{\theta }}$ определяется общее число ВСГ ${{N}_{h}}\,:$

Были приняты следующие количественные критерии точности представления поля:

Поля в соотношениях (33) брались на сфере, окружающей ЭИ с радиусом ${\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}.$ Индекс a соответствует приближенным аппроксимациям. Было принято условие, что погрешности ${{\delta }^{{e,m}}}$ не должны превышать 10^–3 ни при каких условиях, т.е. при минимальном расстоянии $\Delta r$ для любой из рассеянных волн. Эти условия удовлетворялись при ${{N}_{{\theta s}}} = 5,$ ${{N}_{{\theta i}}} = 12$ (${{N}_{{\theta i,s}}}$ – число учитываемых угломестных гармоник возбуждающих и рассеянных волн). Таким образом, мы видим, что параметры ${{N}_{{\theta i,s}}}$ неодинаковые и, следовательно, матрицы ${\mathbf{K}},{\mathbf{L}}$ прямоугольные. При этом, однако, матрица ${\mathbf{W}}$ оказывается квадратной.

5. ОПЕРАТОР РАССЕЯНИЯ

Определение оператора рассеяния – следующий этап МОМР-РВСГ. На этом этапе учитывается конструкция ЭИ. Рассматривается одиночный ЭИ, который многократно возбуждается со стороны своих входов, а также со стороны свободного пространства ВСГ. В работе были исследованы три типа излучателей: щель в волноводе, щель в резонаторе, ленточный вибратор. Они показаны на рис. 2а–2в. Щель и лента имеют длину $L$ и ширину w. Ширина волновода a, а высота b. Волновод и волноводный резонатор (см. рис. 2а, 2б) заполнены диэлектриком с относительной проницаемостью ε. Щелевые излучатели имеют бесконечный бесконечно тонкий экран, в котором выполнена щель. Все металлические элементы конструкций идеально проводящие.

Волноводно-щелевой ЭИ предназначен для использования в составе волноводно-щелевой решетки бегущей волны, поэтому он имеет два входа 1, 2, роль которых выполняют полубесконечные волноводы, расположенные слева и справа от щели соответственно. Резонаторно-щелевой ЭИ возбуждается при помощи коаксиальной линии передачи, которая моделируется цилиндром, радиус которого ${{R}_{{\text{ц}}}}$. Цилиндр имеет контакт с верхней стенкой резонатора. В зазор между нижней стенкой и цилиндром включен источник напряжения с внутренним сопротивлением ${{R}_{{\text{г}}}}$, которое равно характеристическому сопротивлению коаксиальной линии. Данный источник выполняет функцию единственного входа 1 излучателя.

Рассмотрим процедуру определения оператора рассеяния в простейшем случае ленточного ЭИ. Решение для всех рассмотренных структур было получено с помощью стандартной методики, которая неоднократно описывалась в литературе. Она основана на записи интегральных уравнений (ИУ) относительно электрических токов на ленте или магнитных на щели и последующем решении ИУ методом Галеркина [30]. Поэтому подробно излагать ее не имеет смысла, остановимся только на некоторых существенных для нас моментах решения.

Излучатель имеет один вход 1 в виде сосредоточенного источника с сопротивлением ${{R}_{{\text{г}}}}$. Мы решаем граничную задачу в приближении узкой ленты, полагая, что

где $\lambda $ – длина волны в свободном пространстве.

Условия (34) позволяют нам учитывать только продольные компоненты электрических токов, текущих по ленте $I_{y}^{e}$ и пренебречь поперечными токами. Для записи ИУ воспользуемся представлением функции Грина свободного пространства в виде интеграла Фурье [28]:

Здесь ${{e}_{0}}$ – электродвижущая сила (ЭДС) источника напряжения, а $J$ – текущий через него полный ток:

В соответствии с методом Галеркина представляем ток в виде ряда по базисным функциям ${{f}_{{bn}}}(x,y)$:

где ${{X}_{p}}$ – неизвестные коэффициенты.

Алгебраизация ИУ (35) проводится в соответствии с методом Галеркина:

Соотношение (38) представляет собой СЛАУ относительно коэффициентов ${{X}_{p}}.$ ЭДС источника ${{e}_{0}}$ может быть выражена через параметры рассеяния:

где ${{U}_{i}}$ – амплитуда падающей на вход волны. Считаем, что она равна единице, ${{e}_{0}} = 2$. Тогда амплитуда отраженной волны ${{U}_{r}}$ представляется следующим образом:

Отметим, что система (39) и соотношение (40) получены для случая возбуждения вибратора со стороны его входа. В выражении (40) мы использовали индекс $w$, чтобы подчеркнуть, что данный параметр соответствует возбуждению ЭИ со стороны его входа.

Полученный результат позволяет нам найти оператор ${{{\mathbf{L}}}_{{ww}}}$ (см. разд. 3), который в данном случае является скалярной величиной из-за наличия у ЭИ одного входа. Поскольку амплитуда падающей волны равна единице, то ${{{\mathbf{L}}}_{{ww}}} = {{U}_{r}}.$

Также мы теперь можем определить оператор ${{{\mathbf{L}}}_{{sw}}}.$ Для этого необходимо найти амплитуды ВСГ, возбуждаемых электрическим током на вибраторе. Для этого воспользуемся результатами [28]. При решении данной задачи был принят ряд упрощающих допущений.

Воспользуемся еще раз условиями (34). Если лента узкая, то мы можем приближенно заменить ток, распределенный по координате $x$ нитью тока:

Представление тока (41) нельзя использовать при решении интегрального уравнения, однако его применение для вычисления полей пространственных гармоник допустимо, так как поля ВСГ достаточно низкого порядка сравнительно слабо меняются по ширине ленты. Поэтому при интегрировании по координате $x$ допустимо поля ВСГ положить равными постоянной, например, значению при x = 0. Такое приближение соответствует представлению (41). Отметим, что при решении интегрального уравнения замена листка тока нитью недопустима, так как поле, возбужденное листком тока может существенно меняться в области, занятой током, т.е. на поверхности ленты.

Выражение (41) можно без труда записать в сферических координатах:

При записи формулы (42) учтены соотношения

Теперь для определения амплитуд рассеянных волн можем воспользоваться соотношением из [28, ф-ла (2.108)]:

где $\theta {\kern 1pt} '$ равно ${\pi \mathord{\left/ {\vphantom {\pi 2}} \right. \kern-0em} 2}.$ Амплитуды магнитных волн равны нулю, так как электрический ток имеет только радиальную компоненту.

Для того чтобы найти оператор ${{{\mathbf{L}}}_{{sw}}}$, достаточно разделить выражение (44) на амплитуду падающей волны ${{U}_{i}}$ (но в силу того, что ток пропорционален ${{U}_{i}}$, она не войдет в соотношение для оператора).

Для того чтобы найти операторы ${{{\mathbf{L}}}_{{ws}}},{{{\mathbf{L}}}_{{ss}}},$ нам необходимо возбудить вибраторы полями ВСГ. Интегральное уравнение в этом случае приобретает вид

Здесь $E_{{iyn,m}}^{e}(x,y)$ – поле возбуждающей ВСГ электрического типа на ленте. Уравнение (45) записано для случая единичной амплитуды ВСГ.

Алгебраизация уравнения (45) приводит к необходимости вычисления интегралов ${{i}_{{n,m,q}}},$ в которые входит функция $E_{{iyn,m}}^{e}(x,y)\,:$

Применим для вычисления параметров ${{i}_{{n,m,q}}}$ прием, аналогичный приему, использованному выше (см. формулу (41)):

Поскольку компонента поля $E_{{irn,m}}^{e}$ выражается с помощью соотношений (1), (2), то для определения параметров ${{i}_{{n,m,q}}}$ необходимо взять интегралы, аналогичные интегралам, стоящим в формуле (44). Их определяли численно. После этого СЛАУ для коэффициентов ${{X}_{p}}$ при возбуждении ВСГ полностью определена и может быть решена. В результате решения находим параметры ${{X}_{{sp}}}$, где индекс $s$ показывает, что они получены при возбуждении вибратора ВСГ.

Для определения оператора ${{{\mathbf{L}}}_{{ss}}}$ мы должны воспользоваться соотношением (44), заменив в нем ${{X}_{{wp}}}$ на ${{X}_{{sp}}}$. Для оператора ${{{\mathbf{L}}}_{{ws}}}$ имеется следующее соотношение:

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

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

На рис. 3 представлена зависимость модуля элемента матрицы взаимодействия $\left| {{{{\mathbf{W}}}_{{6,6}}}} \right|$ от электрического расстояния ρ. Матрица взаимодействия получена для щелевой решетки. Она не зависит от типа щелевого излучателя (см. рис. 2а, 2б). Как отмечено выше, для обозначения ВСГ была использована сквозная нумерация с индексом N. Номеру $N = 6$ соответствуют индексы ВСГ в двойной нумерации: $n = 1,$ $m = - 1.$

Отметим, что учет конкретного типа излучателя позволяет дополнительно уменьшить размерность матрицы взаимодействия, так как многие из ВСГ возбуждаются щелью или лентой с нулевыми или пренебрежимо малыми амплитудами. Выше было отмечено, что при описании рассеянных полей мы ограничились пятью угломестными гармониками ${{N}_{{\theta s}}} = 5.$ При этом в соответствии с формулой (32) общее число гармоник равно 35. Из них существенные амплитуды имеют волны с номерами N в сквозной нумерации, представленными в таблице 1.

Как следует из таблицы, для достижения приемлемой точности расчетов нам достаточно учитывать шесть ВСГ.

Представленные на рис. 3а–3в зависимости получены при $L = 10,$ $w = 1,$ $f = 9.15$ ГГц. Все размеры на рис. 3 и далее приведены в миллиметрах. Сплошная кривая рассчитана с использованием сплайн-аппроксимации матричной функции ${\mathbf{W}}(\rho )$. Точки соответствуют расчету по формуле (30). Как видно из рис. 3а, аппроксимация при больших ρ описывает функцию $\left| {{{{\mathbf{W}}}_{{6,6}}}} \right|$ с достаточно высокой точностью уже при $\rho > 2\pi $. На рис. 3б показана зависимость фазы функции ${{{\mathbf{W}}}_{{6,6}}}$ от электрического расстояния $\rho $. Видно, что для фазы совпадение между ними еще лучше, чем для амплитуды. На рис. 3в показана зависимость модуля функции ${{{\mathbf{W}}}_{{2,2}}}$ также от параметра $\rho $. В данном случае асимптотическая формула хуже описывает исходную функцию. Причина в том, что номеру $N = 2$ соответствует ВСГ с большим угломестным индексом $n = 3$, для которого условия дальней зоны начинают выполнять при больших ρ. Очевидно, что расхождение асимптотики и точной функции будет увеличиваться с ростом индекса $n$.

Из рис. 3а–3в может показаться, что вклады ВСГ с номерами $N = 2,6$ в матрицу взаимодействия сильно отличаются друг от друга. Однако это не совсем так, поскольку при изменении индекса существенно меняются сами пространственные гармоники (см. соотношение (2)). Поэтому изменения в амплитудах не дают полного представления об изменении вклада ВСГ в рассеянное поле.

Рассмотрим далее результаты численного решения граничной задачи для ленточного ЭИ. На рис. 4 показана частотная зависимость модуля коэффициента отражения от входа вибратора ${{S}_{{11}}}$ при $L = 13.85$, $w = 1$, ${{R}_{{\text{г}}}} = 70$ Ом. Для представления тока на ленте мы использовали три базисные функции (37). Число базисных функций выбиралось из условия сходимости решения. При этом использовался следующий критерий сходимости:

где $n$ – число учитываемых базисных функций.

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

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

где

– мощность, падающая на вход вибратора.

В формуле (50) ${{A}_{{sn,m}}}$ – амплитуды четных рассеянных ВСГ электрического типа. При выводе уравнения (51) принято во внимание асимптотическое поведение сферических функций Ганкеля (24).

На рис. 5 показана частотная зависимость параметра $\Delta $, показывающего точность выполнения закона сохранения энергии ($\Delta $ – модуль разницы между выражением в левой части (51) и единицей). Видно, что энергетический дисбаланс весьма мал, что является свидетельством корректности использованного метода решения.

Точность решения проверялась также путем сравнения с результатами, полученными при помощи стандартной системы. Расчеты в ней, как и в случае одиночного вибратора, проводили методом ИУ. На рис. 6 показана решетка из пяти вибраторов, рассчитанная в режиме синфазного возбуждения, расстояние между вибраторами равно 17.

На рис. 7а–7в показаны частотные зависимости модулей коэффициентов отражения от входов вибраторов с номерами 1–3 (см. рис. 6). Под коэффициентом отражения вибратора в составе решетки ${{R}_{j}}$ ($j = 1,...5$) мы понимаем отношение амплитуды отраженной волны к амплитуде падающей волны при условии возбуждения всех элементов решетки. Параметры вибраторов: $L = 13.85,$ $w = 1,$ ${{R}_{{\text{г}}}} = 70$ Ом. Сплошные кривые получены методом МОМР-РВСГ, а точки – методом ИУ. Расхождение в параметрах рассеяния не превышает 0.01.

На рис. 8 показана частотная зависимость погрешности $\Delta $ выполнения закона сохранения энергии для пятиэлементной решетки. На рис. 9 показана нормированная ДН решетки по вектору Пойнтинга в плоскости $ZOY$.

В целом следует отметить достаточно хорошее совпадение результатов, полученных методом МОМР-РВСГ, с расчетами с помощью стандартной системы моделирования. Исключение составляет только уровень дальних боковых лепестков ДН, который в случае МОМР-РВСГ на 2.5 дБ выше, чем при использовании метода ИУ. При этом наблюдается хорошее совпадение ДН в области главного и первого бокового лепестков.

Рассмотрим далее применение МОМР-РВСГ для анализа решетки щелевых ЭИ, изображенных на рис. 2б. Мы рассматривали ЭИ со следующими параметрами: $a = 14,$ $b = 5,$ $L = 10,$ $w = 1.2,$ ${{R}_{{\text{ц}}}} = 0.62,$ ${{L}_{c}} = 4,$ ${{L}_{1}} = 7.06,$ ${{L}_{2}} = 1,$ $\varepsilon = 2.5.$ Решетка состояла из 26 излучателей (${{N}_{x}} = 26$) расположенных вдоль оси 0х с периодом ${{P}_{x}} = 17.$

Параметры ЭИ выбраны из условия его согласования на частоте 10 ГГц. Частотная характеристика модуля коэффициента отражения ${{S}_{{11}}}$ излучателя приведена на рис. 10.

Контроль точности решения для щелевого излучателя также осуществлялся путем расчета погрешности выполнения закона сохранения энергии $\Delta $. Ее зависимость от частоты показана на рис. 11. Из кривых на рис. 10, 11 видно, что решение по МОМР-РВСГ имеет достаточно высокую точность и в случае щелевого ЭИ.

На рис. 12 представлена частотная зависимость модуля коэффициента отражения. Кривая 2 соответствует одиночному ЭИ (см. также рис. 10), а кривая 1 соответствует коэффициенту отражения от решетки. Отметим, что решетка работает в режиме синфазного возбуждения. При этом каждый из ее элементов может быть описан при помощи коэффициента отражения излучателя в составе решетки ${{R}_{j}},$ $j = 1,...,{{N}_{x}}.$ Под коэффициентом отражения всей решетки в целом мы понимаем корень из отношения суммарной мощности всех волн, отраженных от входов ЭИ, к мощности всех падающих волн. Из рис. 12 видно, что резонансная частота, соответствующая минимуму ${{S}_{{11}}}$ сместилась вверх, а минимальное значение коэффициента отражения увеличилось, что является следствием взаимодействия излучателей в решетке.

На рис. 13 представлено распределение коэффициентов ${{R}_{j}}$ вдоль решетки на частоте 10 ГГц. Видно, что на границах структуры наблюдаются достаточно резкие изменения параметров ${{R}_{j}},$ являющиеся следствием краевых эффектов. При смещении к центру вариации коэффициента отражения уменьшаются, а его значение стремится к коэффициенту отражения в составе бесконечной решетки.

На рис. 14 представлена нормированная ДН структуры на частоте 10 ГГц в плоскости $XOZ$. Она соответствует ДН линейной антенны с равномерным амплитудно-фазовым распределением источников. Таким образом, мы видим, что влияние краевых эффектов на ДН незначительно.

В заключение рассмотрим применение МОМР-РВСГ для анализа волноводно-щелевой линейной решетки. Ее излучатель показан на рис. 2в. Особенностью данной структуры является схема возбуждения ЭИ. Если раньше каждый из них возбуждался индивидуально, то теперь все они возбуждаются последовательно общим волноводом. Таким образом, в составе решетки появляется дополнительный блок, выполняющий функцию связи ЭИ.

Для анализа решетки целесообразно рассмотреть схему двух многополюсников (рис. 15). Первый из них с матрицей рассеяния ${{S}_{a}}$ является решеткой, состоящей из ${{N}_{x}}$ щелевых ЭИ. Входы ЭИ не связаны между собой. Матрицу ${{S}_{a}}$ мы можем найти в соответствии с методикой, изложенной в разд. 3, последовательно возбуждая решетку со стороны всех ее входов. Второй многополюсник с матрицей рассеяния ${{S}_{c}}$ состоит из отрезков волноводов длиной ${{P}_{x}}$. Определение матрицы ${{S}_{c}}$ не составляет труда, и мы не будем останавливаться на этом вопросе. Из рис. 15 видно, что два входа решетки с номерами $1,2{{N}_{x}}$ остаются свободными. Они являются входами интересующей нас волноводно-щелевой решетки. Присвоим им соответственно номера $1,2$.

Найти результирующую матрицу рассеяния волноводно-щелевой решетки по заданным матрицам ${{S}_{{a,c}}}$ можно любым из известных способов. Мы использовали методику анализа сложных составных многополюсников описанную в книге [31].

На рис. 16 показана частотная зависимость параметров рассеяния волноводно-щелевой решетки от частоты. Здесь $R$ – модуль коэффициента отражения от входа 1, а $T$ – модуль коэффициента передачи со входа 1 на вход $2{{N}_{x}}.$ Кривые получены при следующих параметрах структуры: ${{N}_{x}} = 30,$ ${{P}_{x}} = 17,$ $a = 16,$ $b = 5,$ $L = 9,$ $w = 1,$ $\varepsilon = 2.5.$

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

Эволюцию ДН анализируемой структуры можно понять из рис. 17а–17в, на котором показаны нормированные ДН в плоскости XOZ, рассчитанные соответственно на частотах 9, 10 и 10.5 ГГц. Положительные углы $\theta $ соответствуют полуплоскости $x > 0,$ а отрицательные $x < 0.$ Видно, что максимум ДН находится в области отрицательных углов. Если учесть, что возбуждение решетки ведется со входа 1, т.е. в положительном направлении по оси 0х, то можно сделать вывод, что излучение происходит на минус первой пространственной гармонике, что возможно благодаря диэлектрическому заполнению волновода.

Главный луч ДН смещается при повышении частоты в направлении угла $\theta = 0,$ т.е. к нормали к плоскости, в которой расположена структура. При этом изменяется амплитудное распределение источников в апертуре решетки. На рис. 18 показано распределение нормированного модуля амплитуды доминирующей рассеянной пространственной гармоники ${{A}_{{s6,j}}}$ ($j$ – номер щели). Номер гармоники в сквозной нумерации соответствует индексам двойной нумерации $n = 1,$ $m = - 1$ (см. таблицу). Нормировка осуществлялась на максимальную амплитуду гармоники.

Видно, что распределение имеет экспоненциальный характер, что соответствует элементарной теории решеток бегущей волны [32] с одинаковыми ЭИ. Скорость затухания волны в решетке сильно зависит от частоты. Кривые 1–3 рассчитаны при $f = 9.5,10,10.5$ ГГц. Амплитудное распределение определяет форму ДН. В частности, на частоте 10.5 ГГц скорость затухания столь велика, что практически работает лишь небольшой участок структуры, расположенный в ее левой части. При этом ДН сильно искажается (см. рис. 17в).

Программа расчета щелевых и ленточных решеток была реализована в среде Маткад. Полный расчет характеристик решетки из 30 элементов на 10 частотных точках занимает 20…30 с при использовании компьютера со средней производительностью. При этом существенная часть времени расходуется на вычисление ДН. Если анализ не связан с изменением параметров ЭИ, то оператор рассеяния определяется однократно.

ЗАКЛЮЧЕНИЕ

Из представленных результатов видно, что алгоритм на основе МОМР-РВСГ обеспечивает при анализе решеток щелевых ЭИ существенно более высокую вычислительную эффективность, чем стандартные системы электродинамического моделирования и при этом дает близкие к ним результаты, что является подтверждением его точности.

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

Современные системы электродинамического моделирования предусматривают решение задач возбуждения сложных структур заданными источниками. Особенностью МОМР-РВСГ является то, что источником служит поле ВСГ. Напрямую современные системы моделирования таких источников не имеют, что заставляет искать те или иные обходные пути решения данной задачи.

В заключение имеет смысл остановиться на ограничениях на конструкцию ЭИ и параметры решетки, которые накладывает МОМР-РВСГ. Главным из них является то, что излучатели должны быть вписаны в сферы некоторого радиуса и эти сферы не должны пересекаться. Несложно увидеть, что выполнить данное условие не всегда возможно даже для таких структур как щели. Например, если они расположены достаточно близко.

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