Математическая модель рабочего процесса пульсирующего детонационного двигателя Текст научной статьи по специальности «Математика»

УДК 621.452

МАТЕМАТИЧЕСКАЯ МОДЕЛЬ РАБОЧЕГО ПРОЦЕССА ПУЛЬСИРУЮЩЕГО ДЕТОНАЦИОННОГО ДВИГАТЕЛЯ

© 2007 Р.Б. Сейфетдинов

Самарский государственный аэрокосмический университет

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

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

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

смесь заполняет определенный объем камеры, инициируется детонационная волна, которая распространяется по смеси, сжигая ее.

Временной интервал между началом зарядки камеры горючей смесью и началом инициирования детонационной волны рассчитывается с учетом скоростей и1 и и2 фронтов смеси и детонационной волны соответственно таким образом, чтобы фронты достигли открытого конца камеры одновременно. Далее начинается истечение продуктов сгорания в атмосферу и понижение давления в камере до атмосферного, при этом в камеру распространяется волна разряжения [1].

За последние десятилетия проведено множество теоретических и экспериментальных исследований в области детонационного го-

Рис. 1. Схематичное представление рабочего цикла простейшего ПуДД

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

Выражение законов сохранения массы, импульса и энергии в дифференциальной форме в представлении Эйлера для одномерного случая можно записать [2,3]

д 1п р д 1п р ды д 1п Г

-— + ы-— + — + ы-= 0

дг дх дх дх

ды ды 1 др

— + ы — + —— = 0 дг дх р дх '

Т

дS дS

--+ ы —

дх дх

дО- - ы Ъ = 0

д1 дх '

(1)

(2) (3)

д ы д ы 2

--ъ ы--1--1

д г д х k - 1

д 1п а а2 д 8

кЯ дх

дх

(4)

(5)

д( 2 ^ д ( 2 1 д 1п Г

— I-а + ы 1+ (ы + а) — I-а + ы 1 + ыа--

дг {к-1 ) дх {к -1 ) дх

а2 д8 а (д8 д8

-----1 — + ы —

kR дх R {дг дх

, (6)

д ( 2 д ( 2 1 д 1п Г

— I-а - ы 1 + (ы - а)—I-а - ы 1+ ыа--+

дг {к -1 ) дх {к -1 ) дх

а2 д8 а (д8 д8Л „

+----1 — + ы— I = 0

kR дх Я { дг дх )

, (7)

Здесь

2

а + ы = и ^ к -1 +

2

к-1

а - ы = J

- инварианты Римана.

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

д/, , . ди

■ + (ы + а)

дг дх

д 1п Г а2 д8 а ( д8 д8 Л, (8)

= ыа-+--I— + ы— I

дх кЯ дх Я \ дг дх )

дJ , лди — +(ы - а)

дг

= -ыа-

д 1п Г

дх

дх

а2 д8 а (д8 д8

--^ —I — + ы —

кЯ дх Я {дг дх

. (9)

Инварианты Римана J+ и J являются функциями двух переменных х и г. Тогда для инвариантов Римана можно записать выражение

СЛ д^ д^ дх

где г - время, х - координата, р - плотность, ы - скорость, р - давление, Г - площадь, Т-температура, 8 - энтропия, Q - количество теплоты.

Уравнения (1) и (2) можно переписать, введя в качестве независимых переменных а и 8

2 д 1п а 2 д 1п а ды д 1п Г 1 (д8 д8 Л Л

---+ ы---1---+ ы---1--+ ы— I = 0

к-1 дг к-1 дх дх дх Я {дг дх) '

dг

си_ dг

дг

и

дг

- + -

дх дг

+ ди - дх

дх дг

Сравнивая эти выражения с (8) и (9) можно заметить, что

сх

сг

= ы + а

(10)

Используя преобразования предложенные Риманом, получаем систему из двух уравнений

Выражения (10) называются характеристическими уравнениями. Данные уравнения описывают траектории движения фронтов элементарных волн. Зная форму волны в начальный момент времени можно описать ее поведение в последующие моменты времени. Для этого можно использовать плоскость х-г координат, отмечая на ней траектории проходимые различными точками волн и + . Направления этих траекторий можно определить по заданным начальным условиям исходя из того,

0

что котангенс угла наклона равен (10).

Анализ движения волн конечной амплитуды показал, что задача может быть решена, если уравнения в частных производных (8) и (9) заменить двумя парами уравнений в полных производных

йх

= и + а

йг

дШ + . . дJ+ + - + (и + а)-

дг

дх

д 1пГ а2 дS а (дS дS^

= иа--1----1—I--+ и —

дх kR дх R \ дг дх

йх

йг

(11)

= и - а

дШ , дШ — +(и - а)

дг

= -иа-

д 1п Г

дх

дх

а2 дS а (дS дS

--^ —I — + и —

kR дх R { дг дх

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

В этой системе уравнений три неизвестные зависимые переменные и, а и S. Систему необходимо замкнуть уравнением энтропии (3), записав его в следующем виде

йх дS дS

— = и--+ и —

йг ' дг дх

Выражения типа

д

1 д_а + и д_а

Т дг дх

(12)

-+сА

г х

представляют собой производную функции двух переменных по^хаправлению, задаваемому выражением — = с , характеризующую изменение какой либо величины, связанной с данной движущейся частицей вещества. Она так же называется субстанционной или индивидуальной производной.

Система уравнений (11) и (12) в общем случае не имеет аналитического решения и должна решаться численно.

Для проведения исследований рабочего процесса пульсирующего детонационного двигателя был выбран одномерный нестаци-

онарный метод характеристик (МХ). Одномерный подход к задаче был предпочтен двумерному подходу по причине обеспечения высокой скорости расчета и как следствие снижение затрат временных ресурсов при исследовании. Например, решение задачи об ударной трубе. Численная схема типа Годунова может решать задачу Римана сотни или даже тысячи раз. МХ получит решение, решая задачу Римана всего один раз.

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

Основным недостатком МХ является его "жесткость", т.е. неспособность трансформироваться для решения пространственных задач или задач с измененным набором уравнений. Другими словами составленная система уравнений и алгоритм, написанный для решения какой-то частной задачи трудно адаптировать для решения другой частной задачи.

Таким образом, МХ может быть удачно применен для исследования рабочего процесса ПуДД. Математически это означает, что могут быть определены кривые - характери-

стики, поперек которых производные физических параметров могут иметь разрывы, в то время как сами параметры остаются непрерывными. Вдоль характеристик дифференциальные уравнения в частных производных приводят к полным дифференциальным уравнениям. Физическая суть алгоритма МХ - это решение волн и их взаимодействий.

МХ может быть применен для решения системы дифференциальных уравнений при выполнении следующих условий:

1) Система уравнений должна быть линейной или квазилинейной системой дифференциальных уравнений в частных производных (квазилинейное дифференциальное уравнение в частных производных линейно в производных наивысшего порядка зависимых переменных).

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

3) Зависимые переменные должны быть непрерывны во всей области решения. Другими словами МХ неприменим для переходов через ударные волны и поверхности разрыва. Однако отсутствуют ограничения на непрерывность зависимых переменных.

Как говорилось ранее, используя МХ газодинамические разрывы необходимо отделять и решать отдельно в явном виде /5/. В рабочем процессе ПуДД можно выделить три основных типа газодинамических разрывов: ударные и детонационные волны, контактные поверхности.

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

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

выражению [2, 4, 6]

а,

и 2 =

К

I Р

V Л

2 -1

2К1

к +1

1

Л 2

В2_ + К -1

^ 1

К 2 + 1 у

± Щ

(13)

где нижний индекс 1 идентифицирует параметры газа перед волной, а нижний индекс 2 идентифицирует параметры газа за волной.

В случае, когда давление за волной меньше давления перед волной, то есть имеет место волна разряжения, распространяющаяся вправо (влево), то (13) примет вид

Г

U 2 =

2а1

К -1

1 -

К -1 л ( pl Л щ

2 у

± и

(14)

V у

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

(

( К + 1 + Р2

р2

Р1

1

М 2

К -1

Рз

+ 1 Р2

1 +

- 1 Р

(15)

1 уу

а2 = а1

К2 -1

( р3 Л ^ — • (16)

V Р1 у

Основными видами газодинамических разрывов в ПуДД являются детонационные волны. Согласно рис. 1, можно сделать допущение, что детонационная волна всегда распространяется в однородном потоке. Тогда можно использовать модели детонации Чэп-мена-Жуге или модель Зельдовича-Ноймана-Доринга (ЗНД). Единственным взаимодействием детонационной волны является взаимодействие с контактной поверхностью разделяющей горючую смесь и воздух (продукты сгорания) на выходе из камеры. Однако данная задача может быть решена как задача Римана. Подробное рассмотрение этих моделей детонации изложено во множестве

Волна разряженги (Тэйлора) Р3^0.4РЧЖ ^У^ р =р 2 гчж ичж

Pi

1 1 I

) *=L/2 л—L координата по оси

и=о и>0 -— ичж и=0

Рис. 2. Представление детонационной волны в р-х координатах по модели Чэпмена-Жуге

книг и публикаций, в том числе в [2, 4].

Отметим основные особенности данных моделей. Модель Чэпмена-Жуге рассматривает детонационную волну как разрыв, в котором происходи мгновенное сгорание (рис. 2). То есть, непосредственно за детонационной волной находятся продукты сгорания. Модель Чэпмена-Жуге дает в процессе решения информацию о скорости детонации, давления за детонационной волной и т.д., не требуя при этом информации о скоростях протекания химических реакций в детонационной волне.

ЗНД модель детонации является усовершенствованной моделью Чэпмена-Жуге. Она принимает во внимание скорость протекания химических реакций. ЗНД модель описывает детонационную волну как ударную волну и следующую за ней зону реакций (рис. 3). Зона реакций определяется скоростью протекания химических процессов.

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

Рис. 3. Представление детонационной волны в p-х координатах по модели Зельдовича-Ноймана-Доринга

дут определяться характеристиками.

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

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

СПИСОК ЛИТЕРАТУРЫ

1. Bussing, T. R., Hinkey, J. B., and Kaye, L. "Pulse Detonation Engine Preliminary Design Considerations". AIAA 94-3200, 30th AIAA/ ASME/SAE/ASEE Joint Propulsion Conference, June 1994.

2. Курант Р., Фридрихе К. Сверхзвуковое течение и ударные волны. Пер. с англ. Ком-панейца А.С. М.: Ин. лит., 1950.

3. Рудой Б.П. Основы теории газообмена

ДВС. Уч. пос. Уфа, УАИ, 1977.

4. ЗельдовичЯ.Б., Райзер Ю.П., Физика удар-

ных волн и высокотемпературных гидродинамических явлений. М.: Наука, 1966.

5. Магомедов К.М., Холодов А.С. Сеточно-

характеристические численные методы. М.: Наука, 1988. 6. Anderson, J. D. Jr. Modern Compressible Flow with Historical Perspective. McGraw-Hill, Inc., 2nd edition, 1990.

THE THEORETICAL MODEL OF PULSE DETONATION ENGINE

© 2007 R.B. Seyfetdinov

Samara State Aerospace University

The analysis and consideration of the modeling issues of the simple pulse detonation engine was carried out in this paper. Also, the main hyperbolic differential equations, that describe the wave-like nature of working process of the pulse detonation engine was considered. The selected method for a numerical simulation was justified from the point of view of accuracy and accessibility of this examination method for application purposes.