Имитационное математическое моделирование излучения выхлопных струй самолетных двигателей Текст научной статьи по специальности «Физика»

УДК 621.384.32

И. С. Григорьев, А. А. Гусева

Имитационное математическое моделирование излучения выхлопных струй самолетных двигателей

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

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

Введение

Более 20 лет специалисты АО «Научно-производственное объединение «Государственный институт прикладной оптики» (ГИПО) выполняют расчеты излучения струй двигательных установок, которые необходимы для решения задач трехмерного имитационного моделирования различных сюжетов на уровне компьютерных имитационных моделей (КИМ), оснащенных реальными физическими характеристиками и необходимых для оптимизации конструкций современных оптико-электронных систем, работающих в сложных погодных о условиях, различных географических широту тах, при различной орографии подстилающих ^ поверхностей и освещенности визируемых сцен с разных ракурсов наблюдения. Моде-ь ли излучения струй (МИС), разработанные в ГИПО в 1971 г., были одобрены Межведомст-те венным советом при ГОИ им. С. И. Вавилова. ^ Однако современные требования к КИМ,

Ф появление новых специализированных иссле-ш дований [1] актуализировали необходимость ?? совершенствования МИС. В статье рассмот-^ рено создание трехмерной имитационной мо-о дели излучения горячих выхлопных струй, | образованных продуктами сгорания в авиаци-ь онных двигателях. Возникает необходимость оз разработки алгоритма, который, во-первых, будет применим по скорости выполнения к задачам трехмерной визуализации на стандартной ™ видеокарте, во-вторых, будет обладать макси-

й мально реализуемой теоретической строгого _

(П

© Григорьев И. С., Гусева А. А., 2018

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

Существуют два подхода к решению указанной задачи. Первый связан со способом визуализации, применяемой в индустрии видеоигр, в которой полупрозрачные среды моделируются текстурами, не несущими в себе реального физического смысла. Второй подход связан с расчетом сначала уравнений Навье - Сток-са (получение распределения температур), а затем уравнения переноса излучения (для вычисления яркости в заданном направлении) [2]. Расчет по этому методу требует большого количества времени. В представленной работе реализован подход, в котором сочетаются высокая скорость расчета и визуализации, с одной стороны, и методологическая корректность - с другой. Используемые модели и приближения Задача расчета струй традиционно включает в себя два этапа: вычисление пространственных распределений газодинамических параметров (температуры, давления, скорости, плотности) и расчет уравнений переноса излучения. Для исключения расчетов уравнений Навье - Сток-са первый этап реализуется с использованием аналитической теории турбулентных струй Г. Н. Абрамовича [3]. Ее преимущества заключаются в том, что она позволяет единым образом описать распределение параметров

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

Граница струи по теории Абрамовича рассчитывается с использованием уравнения, описывающего увеличение толщины слоя смешения (г) от продольной координаты (х):

dr _ с и - и2 йх и '

где С - константа;

щи и2 - скорости потока на границах слоя смеш ения;

и - средняя скорость.

Значение константы определяет величину начального участка струи при фиксированных соотношениях температур и скоростей струи и спутного потока. В расчеты по теории входят поперечная и продольная координаты, отнесенные к радиусу сопла. Таким образом, данная теория включает информацию как о размере сопла, так и о соотношении температур и скоростей на срезе сопла и во внешней среде. Пример расчета распределения температур в струе по данной теории представлен на рис. 1.

Основной участок струи ввиду его низкой температуры берется в расчет частично - конец струи полагается на расстоянии 1,6 длины начального участка.

При вычислении на видеокарте уравнения переноса излучения используется метод трассировки лучей. Решение уравнения переноса излучения осуществляется по групповой модели [4].

Исходное уравнение переноса излучения (УПИ)

дJx

Рис. 1. Распределение температуры, отнесенной к температуре на срезе сопла по модели струи Абрамовича (Т / Т). По осям отложены координаты, отнесенные к радиусу сопла (хД В начальном участке струи располагается изотермическое ядро конусообразной формы

где Jx- спектральная плотность энергетической яркости (СПЭЯ), Вт/ср/мкм/м2;

s,Q - пространственная координата на линии визирования (м) и телесный угол;

Ьх - объемный коэффициент поглощения, м-1;

^аск - СПЭЯ абсолютно черного тела, Вт/ср/мкм/м2.

Данное уравнение имеет решение в виде

Jx (s,Q) = J0exp

ь| kx (s')Jbback (s')exp

-j kx (s')ds'

0 _

s

-j kx (s'')ds

ds

"max

J(s, ü) = J Jx (s, ü)dX;

^min

Jx (s, ü) = J^exp

-J kx (s ')ds'

(1)

+Jbback (1 - exp[-kx s]),

- + kx (s)Jx (s, V) = kx (s )JbacCk (s),

ds

где J0 - СПЭЯ на входной границе (в данном случае в зависимости от конкретной линии визирования это либо излучение сопла, либо ноль).

Уравнение (1) представляет решение в частном случае изотермического участка. Проинтегрируем исходное уравнение в пределах малого спектрального диапазона g = [Я^ Я2]:

е

о р

ст о

оте

к

а р

а

m

о Ч е л с с

к с е

у

и м с о К

ду ^ а)+^ (яу* (я, а)=^ ш*ь1ал (я);

дя

л (я, а) =| у (я, оу4я, у*МаЛ (я) =

Я2 1 Я2 = | уьЬаЛ (8)й я, ^ (8)=—— | ^ Ш я.

Я1 Я 2 Я1 я1

Тогда для СПЭЯ и энергетической яркости (ЭЯ) будем иметь следующие выражения:

У* (я, а) = уехр -[ № (я')йя'

_ о _

+|^(¡гу11^(я^ехр -[^(я''^" Ж',

О _ я' _

у (я, а)=£ у* (я, а). (2)

Спектры коэффициентов поглощения для воды и углекислого газа берутся из современной версии базы данных HITRAN [5]. Пример усреднения спектра в интервале 0,05 мкм пред-

ставлен на рис. 2. Использовались спектры молекул воды и углекислого газа, взятые для температур 300, 550, 1500 К, значения коэффициентов поглощения для промежуточных температур находились линейной интерполяцией.

Расчет СПЭЯ по уравнению (2) осуществляется в дискретной форме для каждой линии визирования, пересекающей расчетный объем, как показано на рис. 3.

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

со о о

о

см

о

тг и

01 §

- я

и

>и ет

ш о

1- и

I о

<1 а

1 0

(1 V

та г Й

§ о

в

о

со

га

I

о.

<и

I о 15

10

I 5

1- О 0

<и

со -5

-10

сч 10

ю

о

сч

ю

см

(П

(П

4,20

4,25 4,30

Длина волны, мкм

4,35

4,40

Рис. 2. Усреднение спектра коэффициента поглощения

-10 0

Рис. 3. Параллельный рачет уравнения переноса излучения вдоль каждой линии визирования

К<1

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

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

Сравнение с экспериментальными данными

Температура на срезе сопла для режимов «мак-симал» и «форсаж» определялась на основании имеющихся экспериментальных данных. В расчетах были приняты значения 1050 К для режима «форсаж» и 450 К для режима «максимал». Значение скорости истечения газов на срезе сопла определялось из характеристики удельной тяги (Н-с/кг) самолетного

двигателя [6]. Значения объемных концентраций углекислого газа и воды составляют, соответственно, 4 и 6 % для режима «форсаж» и 0,4 и 0,6 % для режима «максимал». Расчеты проводились для среднего ИК-диапазона при учете сильного поглощения излучения молекулами углекислого газа. Рассмотрим сначала струю самолета в виде объема излучающей и поглощающей среды с неравномерным распределением температуры и концентраций оптически активных компонент. В этом случае плоскость сопла представляет собой просто слой газа, и можно рассмотреть решение уравнения (2) в отсутствие внешних источников. Угловая зависимость силы излучения и максимальной яркости струи (значение энергетической яркости в наиболее яркой точке струи) в такой постановке представлена на рис. 4. Под силой излучения понимают интеграл от энергетической яркости поверхности по площади этой поверхности, видимой с большого расстояния под заданным углом.

Углы отсчитываются следующим образом: 90° соответствует наблюдению струи вбок, а 180° - наблюдению струи в сопло. Проанализируем графики на рис. 4. Ввиду различия концентраций, поглощающих излучение молекул, струя на режиме «форсаж» обладает большей оптической плотностью, чем в режиме «максимал», что характеризуется различной угловой зависимостью максимальной яркости. При изменении угла от 90 до 180° уменьшается визируемая площадь струи, но вместе с этим ввиду оптической плотности струи снижается и сама энергетическая яр- 5 кость точек струи. Для направления в 180° в § режиме «форсаж» характерно сильное самопо- ^ глощение энергии газового слоя, так что вклад £ в яркость вносят только внешние, относитель- Ц но холодные слои. В то же время для режима * «максимал» снижение на этих углах силы из- | лучения обусловлено в основном фактором § уменьшения видимой площади струи. о

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

160 Угол, град

160 Угол, град

90°

о см

Ом

14 м

<1

I

м га

г

о ^

со

га г о.

<и

о

и <и со

см ■ч-ю о

I

см ■ч-ю см

(П (П

Рис. 4. Угловая зависимость максимальной яркости и силы излучения для струи в относительных единицах (а - «форсаж», б - «максимал») и примеры отрисовки струи под разными углами с близкой дистанции (в);

—>- - сила излучения; —л— максимальная яркость

валась следующая дискретизация: при длине струи около 15 м цилиндрическая область разбивалась вдоль на 15 частей, угловая разбивка по окружности - на 10, выходная покрышка разбивалась по радиусу на 10 частей. Представим такую разбивку в форме (15, 10, 10). Также этому случаю соответствует величина ds разбиения вдоль линии визирования 0,025 м. На рис. 5, а представлен график разницы данной угловой зависимости и зависимости, полученной при разбивке (100, 100, 100), а на рис. 5, б - разницы данной угловой зависимости и зависимости, полученной для значения ds = 0,002 м.

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

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

х10

СЧ

о

•ч

.--V

сч о о

80

100

120 140 Угол, град б

160

180

Рис. 5. Различие решения от уровня дискретизации поверхности струи (а) и величины пространственного шага

расчета УПИ (б)

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

Ввиду экспериментально подтвержденной слабой энергетической яркости струи в режиме «максимал» основной вклад в общее излучение вносит излучение от внутреннего пространства сопла двигателя. Эксперименты

по замеру силы излучения струй проводились на базе АО «Научно-производственное объединение «Государственный институт прикладной оптики» [1]. Для режима «форсаж» график экспериментальных значений определен до угла 120°, для режима «максимал» - до 170°. В обоих случаях графики являются результатом усреднения по нескольким сериям измерений. Заключение

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

I

S cd

1,0

0,8

0,6

3 0,4 л

и

0,2

*ПМУ""Т..................т...................

100

120 140 Угол,град б

160

180

120 140 Угол, град а

Рис. 6. Представление общей силы излучения струи в виде суммы двух компонент и сопоставление с экспериментальными данными:

а - «форсаж», б - «максимал»; —л— - сила излучения струи; —с— - сила излучения сопла;-- суммарная сила

излучения; —<•— - экспериментальные данные

е

о р

ст о

оте

к

а р

а

m

о ч е л с с

к с е

у

и м с о К

ния производился методом трассировки лучей с использованием групповой модели. Модель учитывает распределение концентраций поглощающих излучение молекул и неравномерность температуры в струе. Полученная модель позволяет описать основные закономерности изменения силы излучения при изменении параметров газового потока и углов наблюдения. Она применима к задачам быстрой визуализации модели в графическом приложении. Список литературы

1. Логинов Д. Я., Филиппов В. Л., Яцык В. С. и др. Техническое обеспечение исследований спектроэнергетических сигнатур сценариев применения оптико-электронной аппаратуры, наблюдаемых в атмосфере в ее верхней и нижней полусферах // XXV Международная научно-техническая конф. и школа по фотоэлектронике и приборам ночного видения. Труды конф. Т. 1. М.: Орион, 2018. С. 186-189.

2. Qinglin N., Zhihong H., Shikui D. IR radiation characteristics of rocket exhaust plumes under varying motor operating conditions // Chienese Society of Aeronautics and Astronautics & Beihang University. 2016. Vol. 21. Рр. 1-14.

3. Абрамович Г. Н. Теория турбулентных струй. М.: Наука, 1960. 715 с.

4. Суржиков С. Т. Тепловое излучение газов и плазмы. М.: Изд-во МГТУ им. Н.Э. Баумана, 2004. 543 с.

5. Rothman L. S. The HITRAN 2008 molecular spectroscopic database // Journal of Quantitative Spectroscopy & Radiative Transfer. 2009. Vol. 110. Рр. 533-572.

6. Жернаков В. С., Кривошеев И. А. Прогноз основных характеристик ТРДДФсм для перспективных авиационных комплексов // Вестник УГАТУ 2015. Т. 19. № 2. С. 56-62.

Поступила 12.12.18

о см

Григорьев Иван Сергеевич - инженер-программист 3-й категории Акционерного общества «Научно-производственное объединение «Государственный институт прикладной оптики», г. Казань.

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

Гусева Александра Александровна - старший научный сотрудник Акционерного общества «Научно-производственное объединение «Государственный институт прикладной оптики», г. Казань.

Область научных интересов: имитационное математическое моделирование фоно-целевых обстановок для разно-спектральных оптико-электронных систем.

< Mathematical simulation of aircraft engine

i

£ jet exhausts radiation

^ The paper deals with the problems of mathematical simulation of aircraft engine jet exhausts radiation, the

® simulation being carried out by means of shader subroutines for the concurrent computation of the radiative

§ transfer equation on the video card resources. The combination of an analytical model of an isobaric jet and

™ ray tracing of computation of the radiative transfer equation allows us to develop a flexible model of aircraft jet

£ radiation, the model taking into account the main parameters of streams in the jet and in the co-current flow, the

■j spectral lines of the radiating components, and provides real-time computation. For the graphic implementation

0 of the model, the OpenGL standard is used.

| Keywords: radiative transfer equation, exhaust jets, concurrent computations, shader subroutines, OpenGL.

1 to Grigor'ev Ivan Sergeevich - third category programming engineer, Joint-stock company "Scientific and Production

Association "State Institute of Applied Optics", Kazan.

cm Science research interests: gas dynamics, radiation transfer, concurrent computing on a video card.

w o

cm Guseva Aleksandra Aleksandrovna - Senior Research Fellow, Joint-stock company "Scientific and Production Association

S "State Institute of Applied Optics", Kazan.

jjj Science research interests: mathematical simulation of background and target situations for multi-spectroscopic

— optoelectronic systems.