CC BY
66
УДК 621.396.6(07), 537.8(07)
СУПЕРКОМПЬЮТЕРНОЕ МОДЕЛИРОВАНИЕ ПОЛЕЙ РАССЕЯНИЯ НА ОБЪЕКТАХ СЛОЖНОЙ ФОРМЫ
Р.Р. Салихов, А.Б. Хашимов
SUPERCOMPUTER SIMULATION OF SCATTERED FIELDS ON COMPLEX SHAPE OBJECTS
R.R. Salikhov, A.B. Khashimov
Рассмотрены методы расчета радиолокационных характеристик объектов сложной формы. Предложены математические модели антенных систем, устанавливаемых на таких объектах. Строгие электродинамические соотношения основаны на интегральных уравнениях Фредгольма II рода. Показано, что применение суперкомпьютеров позволяет исследовать различные радиолокационные сценарии и оптимизировать диаграммы направленности бортовых антенных систем.
Ключевые слова: математическая модель, интегральные уравнения, диаграмма рассеяния, эффективная площадь рассеяния.
In this article the methods for determining of radar characteristycs of complex shape objects are considered. Also the mathematical models of antenna systems, which are usually disposed on such objects, are proposed. The rigorous electrodynamics formulations are based on Fredholm’s integral equations of II type. It is shown that the application of supercomputer simulation allows us to investigate different radar scenes and optimize the patterns of antenna systems.
Keywords: mathematical model, integral equations, scattering field, radar cross section.
Развитие современных радиолокационных и радионавигационных систем связано с повышением точности определения радиолокационных характеристик объектов, имеющих сложную форму поверхности. Широкое применение стелс-техноло-гий [1, 2] требует от разработчиков использования строгих электродинамических подходов к построению математических моделей таких объектов, так как анализ полей рассеяния необходимо производить с очень высокой точностью. Кроме того, такие модели должны быть достаточно универсальными с точки зрения управления геометрией поверхности объектов и электрофизических параметров областей размещения композитных радиопоглощающих материалов (РПМ). Очень близки к таким задачам проблемы определения диаграмм направленности (ДН) бортовых антенных систем. Оптимизация ДН с учетом влияния объекта установки резко повышает эффективность работы систем ближней радиолокации [3].
Широко используемый в настоящее время метод физической теории дифракции (ФТД) позволил получить фундаментальные результаты в области оценки рассеивающих свойств таких объектов, как ребра и двугранные вогнутые структуры [4, 5]. Эти оценки позволили при моделировании реальных объектов выработать критерии скругле-ния резких изломов геометрии, при которых полное поле рассеяния, определяемое методом ФТД и методом интегральных уравнений (ИУ), отличаются в пределах контролируемой точности. Общим вопросом для всех методов является оценка точности полученных результатов. Она может быть получена из сравнения расчетных и экспериментальных результатов, но создание экспериментального образца часто затруднено, стоимость натурных испытаний, требующих специальных условий, оказывется очень высокой. Второй способ оценки точности моделирования заключается в решении серии тестовых задач для объектов из-
Салихов Ринат Рафикович - заместитель главного конструктора ООО «НПО РТС», [email protected]
Хашимов Амур Бариевич - канд. физ.-мат. наук, доцент кафедры конструирования и производства радиоаппаратуры, Южно-Уральский государственный университет; [email protected]
Salikhov Rinat Rafikovich - Deputy Chief Designer of JSC «NPO RTS», [email protected]
Khashimov Amur Barievich - Candidate of Science (Physics and Mathematics), Associate Professor of Radio Equipment Design and Production Department, South Ural State University; [email protected]
вестной формы, допускающих аналитическое определение полей рассеяния, например, для идеально проводящей сферы [6]. Здесь основные сложности связаны с реализацией вычислительных процедур, требующих значительных компьютерных ресурсов. Широкое применение суперкомпьютеров позволило значительно расширить круг задач, решаемых строгими методами.
В строгой электродинамической постановке нахождение поля рассеяния от объекта заданной формы и известного электромагнитного поля возбуждения можно сформулировать в виде интегральных уравнений (ИУ) I и II рода. Каждое из этих уравнений имеет свои специфические особенности, которые определяют выбор численных методов решения. В частности, ИУ I рода с гипер-сингулярным ядром используется для задач с идеально проводящими незамкнутыми бесконечно тонкими пластинами [4]:
4-
lim J gradp (js • gradq ф) dsq -
pP0
- k2 J Гф dsq
= -n x Ег ( p),
(1)
где Zc = yj^-a/ea - волновое сопротивление среды распространения; n - вектор единичной внешней нормали к поверхности объекта S ; js - распределение электрического поверхностного тока на S ; Ег - вектор напряженности возбуждающего электрического поля; p, q - координатные точки наблюдения и источника соответственно; k = 2-/X ; X - длина волны электромагнитного поля; ф = exp (-ikrpq ')jrpq . Отметим, что численное решение таких ИУ относится к некорректным задачам математической физики. Используя предельный переход под знаком интеграла и аппроксимацию решения функциями, удовлетворяющими условию Мейкснера на границе пластины, можно получить равномерную сходимость решения в классе функций, удовлетворяющих заданным требованиям гладкости. Это требует предварительного исследования устойчивости численного метода решения ИУ I рода, что ограничивает его применение для плоских конструктивных элементов малой, но конечной толщины (в сравнении с длиной волны). В действительности любые практические конструкции объектов в местах сопряжения имеют конечный радиус скругления, допускающий непрерывность нормали. Для замкнутых металлических поверхностей широко используется ИУ II рода [7]:
eXP (ikrpq )
.1 x i i i q i x giau„
2-
js (p ) = 2- n xJjs (q)x gradq 2% s + 2n x Нг (p ),
dsq +
* nJ^c|p-q|
где Нг - вектор напряженности возбуждающего магнитного поля. В дальнейшем будем предполагать, что геометрические характеристики поверхности объекта £ удовлетворяют условию непрерывности нормали п - условию Гельдера [8]:
(3)
где с, 0 < а < 1 - некоторые постоянные. Важно подчеркнуть, что ИУ II рода очень часто используются для построения итерационных схем вида Гаусса - Зейделя, последовательной верхней релаксации, что особенно эффективно при использовании суперкомпьютерных методов численного решения.
Вместе с тем необходимо отметить, что ИУ II рода сложно использовать для численного исследования поля рассеяния объекта, содержащего в своем составе плоские пластины очень малой электрической толщины (для kd < 0,05, где d -толщина плоского элемента). Это связано с особенностями ядра ИУ (2), так как при малой толщине векторное произведение в ядре ИУ стремится к нулю для достаточно больших расстояний между точками наблюдения и источника. С точки зрения вычислительных процедур такие особенности могут приводить к неконтролируемому поведению численного решения, что требует особо тщательного подхода к вычислению квадратур, которые используются для аппроксимации ядра ИУ. Эти квадратуры следуют из общей схемы численного решения ИУ проекционными методами.
В статье [9] показано, что для решения трехмерных векторных электродинамических задач при выполнении соответствующих условий по поляризации поля излучения и размеров антенн возможен переход к значительно более простым двумерным скалярным задачам. Это означает, что вместо решения ИУ (2) возможен переход к его скалярному аналогу:
2)|
2
-dlq =
= 2H ( p) ,
(4)
(2)
где } - распределение электрического тока на
контуре L , образованном сечением £ плоскостью Ф = 0 ; я{2^ (кгп) - функция Ганкеля первого порядка второго рода от соответствующего аргумента; гм - радиус-вектор, соединяющий точки наблюдения и источника; Н1 - тангенциальная к контуру L составляющая возбуждающего магнитного поля. Такое приближение асимптотических решений двумерных задач позволяет определить оценки возможных искажений поля излучения бортовых антенн и получить примерный уровень эффективной площади рассеяния (ЭПР) для заданных угловых интервалов.
n • r
p pq
r
Численное решение ИУ (2) или его скалярного аналога (4) дает возможность определить ЭПР объекта или ДН антенной системы, установленной на нем. Основные ограничения, возникающие при численном решении, обусловлены компьютерными ресурсами и погрешностями, возникающими при решении задач очень большой размерности. В связи с этим большое распространение получили итерационные методы численного решения ИУ. Используемая платформа MATLAB предоставляет разработчикам широкие возможности для составления программ с использованием эффективного представления матричных форм Function Handle для метода сопряженных градиентов, обобщенного метода взвешенных невязок [10].
Одной из ключевых проблем, возникающих при использовании метода ИУ, является создание моделей объектов и их последующая дискретизация. Необходимо отметить, что универсального метода дискретного представления S не существует, так как практические требования к исследуемым задачам рассеяния весьма разнообразны, например, по степени формализации деполяризирующих свойств объекта, по частотным свойствам ЭПР. Наибольшее применение получили два метода дискретного представления поверхности S , первый из которых состоит в аппроксимации области построения численного решения ИУ тонкопроволочной сеточной моделью, имитирующей
распределение поверхностного тока jx некоторой линейной структурой, с условиями непрерывности в узлах сопряжения сетки. Для численного решения используется система ИУ Поклингтона [7], с помощью такого подхода решен широкий круг задач, имеющих важное практическое значение. Недостаток такого метода обусловлен его ограниченностью в представлении объектов сложной геометрической формы.
Второй метод - метод конечных элементов, получил очень широкое распространение в чис-
ленном решении задач прикладной электродинамики. К достоинствам этого метода следует отнести возможность анализа объектов со сложными геометрическими характеристиками, адаптивное управление количеством и размерами элементов дискретизации исследуемого объекта, а также возможности внутренней манипуляции параметрами самого элемента. Разработаны и детально исследованы различные типы элементов дискретизации исследуемых поверхностей, допускающих соответствующие аппроксимации численного решения ИУ. Общие оценки эффективности применения конечных элементов различной конфигурации основаны на условиях универсальности аппроксимации поверхности £, а также точности определения компонент вектора нормали к элементу декомпозиции, так как находится численное решение векторного ИУ. С точки зрения практического использования наиболее удобным является треугольный плоский конечный элемент, по геометрическим свойствам близкий к равностороннему треугольнику. Отметим, что именно это свойство обеспечивает высокую точность численного решения, так как вырожденные треугольники с одной стороной, значительно меньшей, чем две другие, могут пересекать линии поверхностных токов с большими вариациями, что резко увеличивает общую погрешность решения. На рис. 1 представлен пример декомпозиции поверхности самолета Як-40.
Выбор треугольного элемента определяется простотой нахождения всех геометрических параметров: координаты вершин треугольника; площадь элемента; компоненты нормали к элементу, необходимые для формирования вычислительного ядра задачи. Важным практическим достоинством триангуляции поверхности £ является возможность оперативного управления процессом декомпозиции, при которой определение общего количества элементов задается минимальной площа-
дью элемента, зависящей от требований к точности численного решения задачи и возможностью введения элементов меньшей площади для эффективной аппроксимации резких изломов геометрии поверхности. В ряде работ показано, что скругле-ние таких сингулярных участков цилиндрическими или тороидальными поверхностями не приводит к росту погрешности численного решения [2]. Большое количество вычислительных экспериментов и сравнение их с известными аналитическими решениями, например, для рассеяния плоской волны на идеально проводящей сфере, показывают, что для практических расчетов такой базовый параметр, как сторона треугольного элемента, можно выбирать в интервале X/8...Х/5. Таким образом, перечисленные свойства треугольных элементов декомпозиции позволяют провести серию вычислительных экспериментов и выбрать такую дискретизацию, при которой достигается оптимальное соотношение между временем расчетов и точностью полученных результатов.
Основным методом декомпозиции поверхности объектов треугольными элементами является триангуляция Делоне [11]. Важной отличительной чертой этого метода является то, что триангуляция Делоне строит набор треугольников, которые «стремятся» к равноугольности, то есть этот набор обладает максимальной суммой минимальных углов всех своих треугольников среди всех возможных триангуляций на заданном наборе точек. Кроме того, триангуляция Делоне обладает минимальной суммой радиусов окружностей, описанных около треугольников, среди всех возможных триангуляций на заданном наборе точек. Эти свойства в большинстве случаев используются для построения алгоритмов триангуляции поверхности. Отметим, что эти алгоритмы эффективно реализованы в современных средствах автоматизированного проектирования и трехмерной графики. Это обстоятельство имеет практическое значение для разработчиков, так как позволяет быстро формировать геометрию декомпозиции исследуемого объекта с возможностью адаптивной коррекции для участков резкого изменения геометрии поверхности £ . Следует особо подчеркнуть, что в ряде случаев автоматизированная триангуляция поверхности может приводить к существенной неравномерности декомпозиции, особенно для объектов сложной формы. Это проявляется в появлении «дефектных» элементов с выраженной неравномерностью сторон, углов и площадей треугольников. Поэтому необходима «ручная» коррекция отдельных участков триангуляции.
Процесс получения дискретной аппроксимации поверхности £ с помощью современных программных средств можно разбить на следующие этапы. Вначале модель исследуемого объекта создается в системе трехмерного проектирования. Далее эта модель преобразуется в поверхностную
модель, состоящую из треугольных элементов. Эта процедура является стандартной для программного обеспечения трехмерного моделирования. На выходе этого этапа формируется файл с расширением *. &^1, в котором содержатся данные о координатах вершин треугольников и проекций нормалей к плоскости каждого треугольника. Полученный файл можно корректировать в зависимости от результатов предварительного анализа качества триангуляции, вводить дополнительные элементы для более точной аппроксимации. Программные средства дают возможность локализовать расположение «дефектных» элементов, например, для нахождения номеров элементов с ошибочным определением направления нормали. На рис. 2 представлен пример управления параметрами элементов для повышения точности численного решения ИУ, когда поверхность объекта имеет резкие изменения геометрии в пределах небольшого участка.
Рис. 2. Пример использования элементов с разными параметрами
В настоящее время использование суперкомпьютеров позволило исследовать значительное количество задач, имеющих большое значение для практических задач радиолокации и радионавигации, синтеза бортовых антенных систем. Использование суперкомпьютеров позволяет получить в распоряжение разработчиков несколько десятков вычислительных ядер и десятки терабайт оперативной памяти. Стоит отметить, что это имеет ключевое значение для решения задач рассеяния прямыми численными методами, так как появляется возможность построения и вычисления матричных уравнений сверхбольшой размерности.
Численное решение задач рассеяния проводилось на вычислительном кластере «СКИФ Урал» Южно-Уральского государственного университета с использованием программного комплекса MATLAB. Пакеты расширения MATLAB Distributed Computing Toolbox и MATLAB Distributed Computing Engine позволяют использовать эффективные программные реализации для распределенных вычислений. Вместе с тем, для получения ощутимого эффекта от использования
кластера потребовалось радикально изменить процедуры формирования и распределенного размещения сверхбольших матриц в оперативной памяти, решения распределенного матричного уравнения. В частности, формирование матриц уравнения производится с использованием стандартных программных средств МА^АВ, таких как индексация и комплексное сопряжение матриц. Это позволило создать процедуру блочного формирования для учета векторного характера задачи с последовательным объединением матриц в общую матрицу. Такое блочное формирование существенно сокращает временные затраты на формирование общей матрицы, так как здесь можно использовать эффективный прием работы с вложенными циклами, так называемые неявные циклы.
Распределенное размещение сверхбольших матриц матричного уравнения производится на платформе расширений МА^АВ. При этом процедура размещения матрицы в памяти кластера построена таким образом, что для пользователя
это выражается в использование двух функций -создание распределённой матрицы и сопряжение распределённой матрицы с локальной. Все действия с этой матрицей не отличаются от действий с обычной матрицей. Распределённое формирование и хранение матриц позволило эффективно использовать вычислительные возможности кластера, а также провести ряд численных экспериментов для объектов различной размерности и геометрии.
В качестве примера рассмотрим определение ЭПР самолета Як-40 для частоты зондирующего поля 75 МГц. На рис. 3 представлена зависимость ЭПР объекта для 4096 элементов дискретизации поверхности, при этом на каждом элементе определяются три комплексные скалярные компоненты
вектора поверхностного тока ^ .
Полученные зависимости для различного количества элементов дискретизации хорошо согласуются между собой в рамках равномерной сходимости результатов.
120
150
90
180
1000
800^
600
400
60
30
210
330
240
270
300
Рис. 3. ЭПР самолета Як-40, частота 75 МГц
90
Рис. 4. Геометрия объекта с бортовой антенной
Рис. 5. ДН антенной системы на объекте
0
1
0
Второй пример численного решения задачи рассеяния показывает изменение ДН антенной системы в виде электрического диполя, ось которого параллельна образующей корпуса объекта, приведенного на рис. 4. Размеры объекта: диаметр большого торца 2Х , диаметр малого торца 0,4Х, высота 2,2Х.
На рис. 5 приведена ДН антенны в плоскости, проходящей через оси диполя и объекта. Расчеты показывают существенное влияние поля рассеяния на объекте на ДН диполя, что подчеркивает необходимость проведения математического моделирования таких задач. Проводя серию таких расчетов, можно оптимизировать расположение антенной системы на объекте.
Выводы
1. Применение ИУ II рода позволяет построить эффективные математические модели задач рассеяния электромагнитных полей на объектах сложной формы. Квадратурные формулы повышенной точности и триангуляция поверхности объекта по методу Делоне обеспечивает высокую точность численного решения с возможностью адаптивного управления процессом декомпозиции. Применение принципа соответствия асимптотических решений двумерных и трехмерных электродинамических задач значительно расширяет возможности математического моделирования сложных практических задач радиолокации, радионавигации, проектирования бортовых антенных систем. Сравнение полученных результатов с известными тестовыми решениями показывает высокую эффективность и универсальность предложенных математических моделей.
2. Применение суперкомпьютеров и программного обеспечения платформы МАТЬАВ дает возможность исследовать объекты больших электрических размеров. Программная реализация методов формирования блочных матриц обеспечивает оптимальное использование кластера. Использование распределенных вычислений для итерационных методов численного решения ИУ еще больше расширит круг задач, имеющих важное практическое значение.
Авторы выражают искреннюю признательность сотрудникам лаборатории суперкомпьютер-ного моделирования ЮУрГУ за ценные консультации по использованию прикладного программного обеспечения кластера «СКИФ Урал».
Литература
1. Лагарьков, А.Н. Фундаментальные и прикладные проблемы стелс-технологий / А.Н. Лагарьков, М.А. Погосян // Вестник Российской академии наук. - 2003. - Т. 73, № 9. - С. 848-862.
2. Сотников, А.М. Оценка отражающих свойств наземных и воздушных объектов с пассивной защитой на основе композитных радио-изотопных покрытий / А.М. Сотников, Р.Г. Сидоренко, Г.В. Рыбалка // Системы управления, навигации и связи: сб. /Харьковский университет Воздушных Сил им. И. Кожедуба. - 2009. - Вып. 1 (9). -С. 70-74.
3. Борзов, А.В. Цифровое моделирование входных сигналов систем ближней радиолокации от сложных радиолокационных схем конфигурации / А.В. Борзов, А.В. Соколов, В.Б. Сучков // Журнал радиоэлектроники. - 2004. - № 4. - С. 17-58.
4. Захаров, Е.В. Интегральные уравнения с ядрами типа Адамара в задачах дифракции / Е.В. Захаров, А.Г. Давыдов, И.В. Халеева //Актуальные вопросы прикладной математики / под ред. А.Н. Тихонова, А.А.Самарского - М.: Изд-во МГУ, 1989. - С. 164-179.
5. Борзов, А.В. Анализ радиолокационных характеристик объектов сложной пространственной конфигурации / А.В. Борзов, Р.П. Быстров, А.В. Соколов // Журнал радиоэлектроники. - 1998. -№ 1. - С. 34-63.
6. Хенл, X Теория дифракции /X Хенл, А. Мауэ, К. Вестпфаль. - М.: Мир, 1964. - 333 с.
7. Вычислительные методы в электродинамике / под ред. Р. Митры. - М.: Мир, 1977. - 486 с.
8. Захаров, Е.В. Численный анализ дифракции радиоволн / Е.В. Захаров, Ю.В. Пименов. - М: Радио и связь, 1982. - 184 с.
9. Войтович, Н.И. О соответствии асимптотических решений двумерных и трехмерных задач в антенной технике /Н.И. Войтович, А.Б. Хашимов // Радиотехника и электроника. - 2010. -Т. 55, № 12. - С. 1471-1476.
10. Ануфриев, И.В. MATLAB 7 / И.В. Ануфриев, А.Б. Смирнов, Е.Н. Смирнова. - СПб.: БХВ-Петербург, 2005. - 1090 с.
11. Скворцов, А.В. Эффективные алгоритмы построения триангуляция Делоне / А.В. Скворцов, Ю.Л. Костюк // Геоинформатика. Теория и практика. - Томск: Изд-во Томского ун-та, 1998. -Вып. 1. - С. 22-47.
Поступила в редакцию 5 декабря 2012 г.
