Full Text

Введение Уровень выбросов вредных веществ из камеры сгорания авиадвигателя должен соответствовать международным стандартам ICAO, что определяет актуальность задачи управления эмиссией для системы автоматического управления двигателем (САУ), надежность которой во многом зависит от валидности поступающей в САУ информации о состоянии двигателя [1]. Во многих случаях прямые измерения физических величин во время полета недоступны. Возможный способ решения проблемы - создание моделей подсистем или процессов двигателя на основе искусственного интеллекта [2; 3]. Камера сгорания авиационного газотурбинного двигателя (ГТД) представляет собой существенно стохастический нелинейный объект, идентификация параметров которого является очень сложной задачей. Эти свойства обусловлены, с одной стороны, требованием компактности конструкции, приводящим к интенсификации процессов перемешивания паров топлива и воздуха посредством наведения турбулентности с поперечным переносом импульса и концентрации топливовоздушной смеси в потоке и последующего химического взаимодействия (горения), скорость которого также в значительной степени зависит от турбулентности. С другой стороны, экспоненциальная зависимость скорости химических реакций от температуры смеси определяет существенную нелинейность горения. С учетом жестких норм на эмиссию вредных веществ все это приводит к узкому диапазону устойчивых режимов работы, ограниченных, с одной стороны, нормируемым уровнем эмиссии окислов азота, а с другой - неприемлемыми в эксплуатации режимами срыва пламени или виброгорения (термоакустическими колебаниями) [4; 5]. При традиционном подходе решение задачи требует большого объема экспериментальных исследований. В связи с недоступностью ряда измерений возникает проблема интерпретации экспериментальных данных, полученных в различных условиях, и экстраполяции эмиссионных характеристик на условия реальной эксплуатации [6]. В настоящее время отсутствует теория подобия режимов малоэмиссионного горения, которая позволяла бы надежно приводить экспериментальные данные, полученные в различных условиях, к неким стандартным условиям, аналогично технологии приведения к МСА параметров лопаточных машин. В связи с этим предлагается создание в составе САУ модели процессов генерации вредных выбросов в камере сгорания ГТД на основе технологии нейронных сетей, обладающих свойством самообучения на базе экспериментальных данных [7]. Кроме того, представляется целесообразной разработка гибридной полуэмпирической модели генерации оксидов азота, состоящей из двух частей. Первая подмодель сочетает в себе физическое уравнение Зельдовича и замыкающий эмпирический коэффициент, являющийся функцией нескольких определяющих переменных. Эта подмодель позволяет генерировать обучающие выборки для второй подмодели в случаях, когда экспериментальных данных недостаточно. Вторая подмодель основана на нейронной сети [8; 9]. Статья посвящена разработке полуэмпирической модели генерации оксидов азота на основе механизма Зельдовича. Теория В полуэмпирической модели генерации оксидов азота моделирование взаимодействия диффузионного и гомогенного факелов в камере сгорания ведется на основе принятой гипотезы возможности применения принципа суперпозиции [10-12]. Для описания процесса генерации окислов азота используется уравнение Зельдовича, позволяющее определить зависимость скорости реакции SNO от состава топливовоздушной смеси с помощью функции распределения плотности вероятности пульсаций концентрации топливовоздушной смеси в зависимости от величины интегрального эквивалентного отношения , характеризующего концентрацию топлива в смеси. Камера сгорания газотурбинного авиационного двигателя представляет собой жаростойкую трубу, куда подают воздух и горючее, которые потом смешиваются. При идеальном смешении образуется топливовоздушная (горючая) смесь, характеризующаяся соотношением «топливо/воздух» - величиной интегрального эквивалентного отношения . Величина интегрального эквивалентного отношения обратна коэффициенту избытка воздуха, т.е. и характеризует количество топлива в смеси. В зависимости от режима работы меняется. Рабочим диапазоном считается диапазон устойчивого горения смеси 0,15 ≤ ≤ 2. При выходе за границы этого диапазона происходит экспоненциальный рост эмиссии вредных веществ в атмосферу (при ≤ 0,15 - «бедная смесь» - и происходит «бедный» срыв пламени, а при ≥ 2 - «богатая смесь» - топлива слишком много и начинается виброгорение - термоакустические колебания) [13-15]. В камере сгорания предусмотрены два факела - диффузионный и гомогенный [16]. Состав горящей смеси (соотношение «топливо/воздух») в факелах имеет стохастический характер (случайный) и изменяется по некоторому случайному закону - то есть случайным образом происходит изменение доли (концентрации) топлива в составе топливовоздушной смеси. Считается, что оба факела - независимы друг от друга и несовместимы в пространстве [17; 18]. Для получения функции распределения плотности вероятности изменения (пульсаций) концентрации (состава топливовоздушной) смеси для диффузионного и гомогенного контуров используются вероятностные характеристики для различных факелов на одном аргументе . При этом учитывается, что функция распределения вероятности [19] в рабочем диапазоне состава смеси изменяется от 0 до 1. Возрастание аргумента соответствует увеличению доли топлива в воздушно-топливной смеси, подаваемой в камеру сгорания (чем больше , тем меньше воздуха в смеси и тем больше топлива, т.е. при имеем чистый воздух без топлива). На диапазоне устойчивого горения смеси 0,15≤ определяется величина математического ожидания (задается дозаторами топлива в камере сгорания). В связи с тем, что технически на длине жаростойкой трубы, ограниченной заданными габаритами камеры сгорания, идеально перемешать горючее и воздух невозможно, имеем суммарную (пространственную и пульсационную) дисперсию нормального распределения (среднеквадратическое отклонение = 20…30 %). Кроме того, в трубе наличествуют продольные акустические волны (давление и скорость), возбуждаемые выделением тепла при горении с собственной частотой газовоздушного столба, приводящие к колебаниям скорости (а следовательно, расхода воздуха и ) на входе в камеру сгорания. То есть продольные акустические волны приводят к гармоническим колебаниям математического ожидания состава смеси с амплитудой 2,5 %. Поскольку эмиссия оксидов азота экспоненциально зависит от температуры пламени и времени пребывания при этой температуре (как всякая химическая реакция), а температура пламени, в свою очередь, зависит от состава смеси , то в математической модели генерации оксидов азота необходимо учитывать эти факторы. Для этого необходимо иметь функции плотности вероятности и вероятности с учетом акустики камеры сгорания. В предлагаемом исследовании решена задача определении функции плотности распределения вероятности изменения концентрации топлива в топливовоздушной (горючей) смеси для диффузионного и гомогенного факелов с учетом всех вышеперечисленных факторов. Данные и методы определения плотности распределения вероятности изменения концентрации топлива в топливовоздушной смеси В качестве исходных данных на заданном режиме работы используются функции распределения плотности вероятности пульсаций концентраций, полученные для диффузионного и гомогенного контуров. Дисперсия концентраций формируется из нескольких составляющих пульсаций параметров потока с соответствующими функциями плотности распределения вероятности [20]: - начальной неоднородности распределения концентраций в факеле по сечению, проявляющейся в турбулентном потоке и в осевом направлении благодаря корреляции поперечного и продольного переноса импульса и концентрации; - регулярной пульсации потока, наводимой завихрителем фронтового устройства за счет прецессирующего вихревого ядра (ПВЯ) при больших числах Рейнольдса, характерных для камер сгорания ГТД; - термоакустических пульсаций потока различных мод колебаний газового столба, фоновой изотропной турбулентности (Колмогоровский масштаб). Для построения адекватной модели как генерации окислов азота, так и определения границы «бедного» срыва необходимо знать перечисленные выше параметры пульсирующего потока. Частично эти данные, например начальная неоднородность распределения концентрации, получены на построенной 3D-метамодели. Для описания взаимодействия диффузионного и гомогенного факелов принята гипотеза возможности применения принципа суперпозиции. Операция суперпозиции гомогенного и диффузионного факелов производится, исходя из допущения об их независимости и несовместности расположения в пространстве. При этом используется теорема сложения вероятностей на диапазоне устойчивого горения. Для диффузионного факела, в силу его природы, принимается гауссова функция плотности распределения вероятности с экспериментально полученными математическим ожиданием и дисперсией. Для гомогенного факела, кроме нормально распределенной вероятности измерения состава смеси, учитывается пространственная неоднородность концентрации топливовоздушной смеси и пульсационная составляющая давления, обусловленная термоакустическим взаимодействием. Средние значения (математические ожидания) состава смеси обоих факелов определяются, исходя из известных (управляемых) расходов топлива по контурам, и расходов воздуха, пропорциональным проходным сечениям контуров. Пример премиксера (предварительного смесителя) приведен на рис. 1, где представлены типичные значения неоднородности концентрации на выходе из премиксера гомогенного контура, полученные на 4D (с учетом нестационарности) математической модели процесса перемешивания, - продольная и окружная неоднородность концентрации топливовоздушной смеси на выходе из премиксера (осредненная по времени и среднеквадратичное отклонение). В качестве исходных данных для гомогенного факела использовались результаты численного моделирования распределения неоднородности и пульсации концентрации топливовоздушной смеси по площади гомогенного факела, полученные с помощью коммерческого пакета программ CFX. Исследование распределения фракции смеси после премиксера (вход в зону горения) для гомогенного факела позволило получить среднеквадратичное отклонение от математического ожидания, которое составило 15 %. Рис. 1. Распределение неоднородности концентрации топливовоздушной смеси по длине премиксера (осредненное по времени и среднеквадратичное отклонение) с помощью коммерческого пакета программ CFX Задача нахождения собственных частот объекта моделирования, формы стоячих волн, фазовых сдвигов и акустического импеданса фронта пламени, решенная для 3D-модели численным методом Навье - Стокса, решалась также с помощью коммерческого пакета программ CFX для неподвижного потока (например, с типичными параметрами Tв = 800К, Tг = 1800 К) при подаче как гармонического возмущения теплового потока (источник) в домене интерфейса между холодным и горячим газом, так и сигнала белого шума. Следует отметить, что решить численным методом задачу определения акустических характеристик для реальной конфигурации камеры сгорания, необходимых для получения данных по термоакустическим пульсациям давления и скорости газового потока, в настоящее время достаточно затруднительно. Поэтому в ходе эксперимента использовался модельный объект, имеющий размеры, соответствующие реальной камере сгорания, но плоской формы, и было принято допущение о плоском фронте пламени и плоском фронте акустической волны. Для нестационарного давления, полученного при расчете на рассматриваемой поверхности, было выполнено преобразования Фурье по времени. На рис. 2 показаны максимальные по рассматриваемой поверхности амплитуды гармоник. В частности, наблюдается пик на частоте вынужденных колебаний f =1000 Гц (q = 100) и на резонансной частоте f = 1240 Гц (q = 124). Рис. 2. Максимальная амплитуда давления Диффузионный факел имеет гауссову функцию плотности распределения вероятности и большую дисперсию состава смеси, чем гомогенный. Диапазон изменения интегрального эквивалентного отношения φ (значения, соответствующие зоне устойчивого горения) для диффузионного факела φ = [0,15; 2,0]. Рассматривая этот диапазон, можно найти среднее значение или математическое ожидание Как показали данные экспериментов, для диффузионного факела справедливо = 2,0 - 1,075 = 0,925. Таким образом, среднеквадратичное отклонение = 0,30833 (30 %). Дисперсия, соответственно, равна D = = 0,095. Гомогенный факел также имеет функцию распределения вероятности, близкую к нормальному закону, и малую дисперсию состава смеси. Диапазон изменения коэффициента избытка топлива (значения, соответствующие зоне устойчивого горения) для гомогенного факела = [0,35; 0,4]. Рассматривая этот диапазон, можно найти среднее значение или математическое ожидание Среднеквадратичное отклонение указывают равным 5 %, а осредненное по времени равным 15 %. Для оценки взято 10 % от математического ожидания. Тогда среднеквадратичное отклонение . Дисперсия, соответственно, = 0,0014063 ≈ ≈ 0,0014. Кривые плотности вероятности изменения состава топливовоздушной смеси для диффузионного и гомогенного факелов без учета акустики камеры сгорания при = 0,4, = 0,04 показаны на рис. 3. Рис. 3. Кривые плотности вероятности изменения состава топливовоздушной смеси для диффузионного и гомогенного факелов без учета акустики камеры сгорания при = 0,4, = 0,04 Модель процессов в гомогенном факеле Как уже отмечалось, при моделировании процессов в гомогенном факеле следует учитывать пульсационную составляющую давления, обусловленную термоакустическим взаимодействием. Как показал модельный эксперимент (см. рис. 2), наличие продольных акустических волн приводит к гармоническим колебаниям (пульсациям) с частотой f = 1 кГц. Концентрация топлива в составе топливовоздушной (горючей) смеси в гомогенном факеле описывается гармоническими колебаниями с шумом : (1) Таким образом, математически задача сводится к расчету плотности вероятности двух совместных (одновременно существующих) независимых процессов: гармонических колебаний и нормально распределенного шума. Результирующая плотность вероятности суммарного процесса вычисляется на основе свойств плотности вероятности [18; 19]: (2) Процесс представляет собой цветной шум(его математическое ожидание не равно нулю, как у белого шума, ), имеющий нормальный закон распределения, плотность вероятности которого описывается известной кривой с максимумом при = 0,4, среднеквадратичной ошибкой = 0,04 и дисперсией (см. рис. 3). Процесс представляет собой регулярные (гармонические) пульсации (колебания) состава смеси (условимся считать, что они изменяются по синусоидальному закону) частотой f = 1÷1,1 кГц и амплитудой Asin = 0,01. Поскольку о величине начальной фазы колебаний концентрации топлива нам достоверно ничего неизвестно, рассматриваемый процесс не может считаться детерминированным. Поэтому условимся считать его случайным, описывающимся гармоническими колебаниями со случайной фазой [21]. Ансамбль (совокупность) реализаций, задающий процесс гармонических колебаний со случайной фазой, показан на рис. 4. Рассмотрим одну из реализаций случайного гармонического процесса (см. рис. 4) и условимся считать ее синусоидальной функцией времени, описывающей пульсации концентрации топлива в горючей смеси: (3) Начальная фаза синусоиды есть случайная величина, которая с одинаковой вероятностью может принимать любое значение за период гармонического колебания - то есть на интервале . Отсюда следует, что плотность вероятности начальной фазы равна (4) Рис. 4. Ансамбль (совокупность) реализаций, задающий процесс гармонических колебаний со случайной фазой Полная фаза колебаний также является случайной величиной, равновероятной на интервале следовательно, ее плотность вероятности равна (5) Для определения одномерной плотности вероятности случайного процесса изменения состава топливовоздушной смеси по гармоническому закону со случайной фазой выделим интервал на одной из его реализаций (рис. 5) из ансамбля (см. рис. 4). Рис. 5. Реализация гармонического колебания со случайной фазой Далее определим вероятность того, что при измерении, проведенном на промежутке времени мгновенное значение сигнала окажется на интервале Эта вероятность соответствует элементу вероятности где - искомая плотность вероятности. Очевидно, что вероятность совпадает с вероятностью попадания фазы гармонических колебаний на один из двух заштрихованных на рис. 5 фазовых интервалов. Поскольку попаданий два, эта вероятность в два раза больше вероятности . Тогда из полученной ранее формулы (5) получим (6) Отсюда искомая функция плотности вероятности (7) Поскольку по условию известно, что (8) Из формулы (6) следует, что где A - амплитуда гармонических колебаний. Тогда (9) Общий вид графика функции плотности вероятности приведен на рис. 6. Рис. 6. Плотность вероятности гармонического колебания амплитуды со случайной фазой Для амплитуды 0,01 пульсаций концентрации топлива в горючей смеси в гомогенном факеле плотность вероятности (10) Кривая плотности вероятности, соответствующая полученным в результате эксперимента параметрам пульсации, изображена на рис. 7. Рис. 7. Плотность вероятности для исследуемого гармонического колебания амплитуды 0,01 со случайной фазой Плотность вероятности изменения состава смеси в гомогенном факеле с учетом акустики камеры сгорания (наличия продольных акустических волны, возбуждаемых выделением тепла при горении с собственной частотой газовоздушного столба, приводящих к расхода воздуха и, соответственно, на входе в камеру сгорания) приведена на рис. 8. Рис. 8. Плотность вероятности изменения состава смеси в гомогенном факеле с учетом гармонических колебаний амплитуды 0,01 со случайной фазой и наличия цветного шума с математическим ожиданием = 0,4, среднеквадратичной ошибкой 0,04 Полученные результаты Полученная плотность вероятности изменения состава топливовоздушной смеси в диффузионном факеле изменяется по нормальному закону распределения с математическим ожиданием = 1,075, среднеквадратической ошибкой = 0,30833: (11) Полученная плотность вероятности изменения состава топливовоздушной смеси в гомогенном факеле с учетом акустики камеры сгорания определяется для гармонических колебаний амплитуды 0,01 со случайной фазой, с наложением цветного шума, имеющего нормальный закон распределения с математическим ожиданием = 0,4, среднеквадратической ошибкой 0,04: (12) Результирующая плотность вероятности состава смеси в диффузионном и гомогенном факеле в соответствии с принципом суперпозиции (рис. 9): (13) Рис. 9. Кривые плотности вероятности изменения состава топливовоздушной смеси для диффузионного и гомогенного факелов с учетом акустики камеры сгорания при = 1,075, = 0,30833 и = 0,4, = 0,04 Далее производится осреднение скорости генерации окислов азота по объему камеры сгорания с использованием зависимости полученной на 3D-метамодели камеры сгорания. Таким образом, осредненная скорость по объему: (16) Индекс эмиссии определяется: (17) где - расход топлива. Эмиссия окислов азота определяется: (18) В соответствии со стандартом Международной организации гражданской авиации (ICAO - International Civil Aviation Organization), параметр эмиссии оксидов азота за стандартный взлётно-посадочный цикл (СВПЦ), определяющий массовый выброс оксидов азота за цикл к установленной взлётной тяге двигателя , определяется: (19) где - время работы двигателя на характерном режиме СВПЦ: руление - 26 мин; заход на посадку - 4 мин; набор высоты - 2,2 мин; взлёт - 0,7 мин, - расход топлива на i-м режиме СВПЦ, - индекс эмиссии оксидов азота на i-м режиме СВПЦ. В рассматриваемой эксперт - модели принято, что массовый выброс эмиссии окислов азота за цикл «взлет - посадка» определяется: (20) Общий вид математической модели оценки генерации окислов азота в среде программирования Matlab Simulink представлен на рис. 10. Рис. 10. Общий вид математической модели оценки генерации окислов азота в среде программирования Matlab Simulink При моделировании в среде Simulink работы малоэмиссионной камеры сгорания выполняется расчет скорости генерации окислов азота в точке по механизму Зельдовича (рис. 11). Рис. 11. Скорость генерации окислов азота в точке по механизму Зельдовича. Моделирование в среде Simulink Далее выполняется расчет эмиссии окислов азота путем интегрирования скорости генерации в точке по объему камеры сгорания с учетом функции плотности распределения вероятностей пульсаций концентрации в камере сгорания (рис. 12). Рис. 12. Результат расчета в среде Simulink индекса эмиссии на заданном режиме работы камеры сгорания Рис. 13. Анализ влияния дисперсии шума на вид кривых плотности вероятности изменения состава топливовоздушной смеси в гомогенном факеле при учете акустики камеры сгорания Кроме того, в рамках данного исследования произведен анализ влияния дисперсии шума на вид кривых плотности вероятности изменения состава топливовоздушной смеси в диффузионном и гомогенном факелах при учете акустики камеры сгорания. В частности, была рассмотрена плотность вероятности состава смеси в гомогенном факеле случаев = 0,04 и = 0,08 (см. рис. 12). Влияние увеличения среднеквадратичной ошибки (дисперсии) шума на вид кривых плотности вероятности изменения состава смеси для гомогенного факела с учетом акустики камеры сгорания иллюстрируется рис. 13. Заключение Расширение пространства состояний двигателя влечет за собой требование увеличения количества прямых измерений параметров, включая такие трудноизмеряемые, как эмиссия вредных веществ, запасы устойчивости к бедному срыву и др. Однако в ряде случаев, где технологии измерения на сегодня еще не достигли требуемого для внедрения уровня готовности, предлагается использование прогностических математических моделей подсистем двигателя. Такие модели, являющиеся функциями измеряемых параметров, имея различный уровень сложности описания объекта идентификации, являются неотъемлемой частью программного обеспечения современных систем автоматического управления (САУ) ГТД. Следует отметить, что математические модели высокого уровня для объектов рассматриваемого класса сложности требуют запретительно больших ресурсов. Наиболее перспективным видится создание моделей МЭКС на основе нейронной сети, эффективно позволяющей увеличить качество решения задачи идентификации объектов высокой размерности на основе синтеза гибких и простых с математической точки зрения алгоритмов [22]. Однако, возвращаясь к проблеме несовершенства технологий измерения, следует отметить неполноту получаемых обучающих выборок для нейронных сетей. Поскольку в условиях полета имеются ограничения на объем вычислительных мощностей (порядка 100 КБ) при сохранении достаточно высоких требований к точности отображения в статике и динамике в большом диапазоне режимов работы, полетных условий и состояний двигателя, была разработана специальная полуэмпирическая эксперт-модель, пригодная как для обучения нейронной сети, так и для встраивания ее в структуру САУ ГТД. Как показал анализ полученных результатов экспериментов, применение разработанных алгоритмов моделирования процессов генерации эмиссии в МЭКС по механизму Зельдовича позволяет достичь необходимых показателей точности прогноза эмиссии и надежности алгоритмов управления малоэмиссионной камерой сгорания ГТД.

About the authors

T. A Kuznetsova

Perm National Research Polytechnic University

V. G Avgustinovich

Perm National Research Polytechnic University

References

Statistics

Views

Abstract - 145

PDF (Russian) - 110

Refbacks

• There are currently no refbacks.