ТЕПЛОФИЗИКА ВЫСОКИХ ТЕМПЕРАТУР, 2015, том 53, № 4, с. 556-560
УДК 533.6
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ТОНКОЙ СТРУКТУРЫ ЦИЛИНДРИЧЕСКОЙ ДЕТОНАЦИОННОЙ ВОЛНЫ В ВОДОРОДНО-ВОЗДУШНОЙ ГОРЮЧЕЙ СМЕСИ © 2015 г. В. Ю. Гидаспов, Н. С. Северина
Московский авиационный институт (национальный исследовательский университет)
E-mail: severina@mai.ru Поступила в редакцию 29.01.2015 г.
В работе исследуются режимы распространения в водородно-воздушной горючей смеси цилиндрической детонационной волны, сформированной в результате распада разрыва на границе камер высокого и низкого давления. Моделирование проводится в квазиодномерной нестационарной постановке с помощью оригинального сеточно-характеристического метода с явным выделением сильных и слабых разрывов. Показано, что при приближении к скорости детонации Чепмена—Жуге детонационная волна распространяется в пульсирующем режиме; что частота пульсаций на больших расстояниях слабо зависит от параметров в камере высокого давления; что шаг колебаний температуры и плотности по пространственной координате за головной ударной волной близок к измеренным в экспериментах размерам детонационной ячейки.
Ключевые слова: ячеистая структура детонационной волны, сеточно-характеристический метод, выделение разрывов.
DOI: 10.7868/S0040364415040122
ВВЕДЕНИЕ
Результаты расчетно-теоретических и экспериментальных исследований режимов распространения и структуры детонационных волн (ДВ) приведены в большом числе отечественных и зарубежных работ, в частности в [1—7]. Исследованы ячеистая структура, энергетические и асимптотические свойства ДВ. В последнее время особое внимание уделяется созданию двигательных установок, использующих детонационное сжигание топлива. Как правило, ячеистую структуру ДВ связывают с многомерными эффектами, сложной ударно-волновой картиной течения и наличием вязкости и турбулентности. Однако возможности квазиодномерного нестационарного моделирования далеко не исчерпаны. В [5] при рассмотрении течений со сферической и цилиндрической симметрией в зависимости от подведенной энергии описаны три режима распространения ДВ на больших расстояниях от места инициирования: пульсирующий, галопирующий и срыв детонации. В настоящее время в связи с развитием вычислительной техники возникла возможность моделировать описанные классиками режимы течения с использованием детальных термодинамически согласованных кинетических механизмов на разностных сетках, пространственный шаг которых существенно меньше длин самовоспламенения в рассматриваемых горючих смесях. Особенностью
используемой авторами методики моделирования [8—10] является, во-первых, то, что расчетными линиями являются траектории частиц газа, вдоль которых уравнения химической кинетики записываются как обыкновенные дифференциальные, во-вторых, ударные волны (УВ), контактные разрывы (КР) и характеристики, образующие веера разрежения, явно выделяются и рассчитываются по специальным алгоритмам. Такой подход позволяет отслеживать распространение возмущений, обеспечить точное выполнение законов сохранения массы и элементного состава, свести к нулю сеточные диффузию и теплопроводность.
ПОСТАНОВКА ЗАДАЧИ
Осесимметричное нестационарное течение реагирующего газа в областях непрерывности параметров можно описать системой дифференциальных уравнений в частных производных, записанной в ее характеристическом виде:
— = u ± a: du ± — dp ± dt
N
pa
рyfr - eyppt
ua +
x ,=1 p2a(erpp - epPr)
W
dt = 0,
У
— = u dh -1 dp = 0, dji -1W(dt = 0, i = 1,...,N. dt p p
Здесь х — радиус-вектор точки, р, и, р, е, к, а — плотность, скорость, давление, удельные внутренняя энергия и энтальпия газа, скорость звука соответственно, у! — мольно-массовые концентрации, N — число рассматриваемых компонентов, Щ — скорость образования 1-го компонента в единице объема в результате химических реакций. Нижние индексы "Т", "р", "у(" обозначают частное дифференцирование по соответствующему параметру. Влиянием эффектов вязкости, теплопроводности и диффузии пренебрегаем.
При переходе через УВ и КР выполняются соотношения типа Ренкина—Гюгонио, которые являются следствием интегральных законов сохранения и могут быть записаны в виде
Рл (Р - ил) = Рп (Р - ып), Рл +Рл (Р - ил) = Рп +Рп (Р - "и), Рл (Р - ил)( + 0.5(Р - "л)2) = = Рп (Р - "п)(( + 0.5(Р - "п)2).
Индексом "л" помечены параметры слева от разрыва, индексом "п" — справа. В случае УВ концентрации непрерывны: уш- = уп ¡, I = 1, ..., N Б — скорость УВ. В случае КР Б = ил = ип, рл = рп; изменение остальных величин произвольно.
Термодинамические свойства реагирующего газа описываются с помощью модели многокомпонентного совершенного газа в рамках допущения о равновесной заселенности энергетических уровней, отвечающих всем внутренним степеням свободы молекул и атомов [11]. Удельный термодинамический потенциал Гиббса здесь имеет вид
выражения для Щ имеют вид
G (p,T, у ) = X Y i
i=i
RT ln
PY t Po X Y J
j=i
+ g0(T )
гдер0 = 101325 Па — стандартное давление, 0° (Г) — температурная часть стандартных молярных потенциалов Гиббса отдельных компонентов [11], Я — универсальная газовая постоянная. Другие термодинамические величины, используемые при математическом моделировании, выражаются через потенциал Гиббса и его частные производные.
В случае произвольного механизма из N обратимых химических реакций
N
N
X V (r M о X V TMi, r = 1,..., Nr
Nr
w = X ((r) -v (r) )x
i=1
i=1
r =1
*(rП (PYj)Г) - K(r)I! (PYj)"
V j=i j=i
Здесь Mi — символ /-ого вещества, V(r), V(r) — стехиометрические коэффициенты, константы
скоростей прямых K(r) и обратных K(r) реакций связаны через константу равновесия [11].
МЕТОД РАСЧЕТА
Численное моделирование проводилось оригинальным сеточно-характеристическим методом, позволяющим рассчитывать квазиодномерные нестационарные течения многокомпонентного реагирующего газа с явным выделением произвольного числа взаимодействующих разрывов [8—10]. Кратко опишем особенности используемой методики моделирования.
— Расчет ведется маршевым сеточно-характе-ристическим методом второго порядка точности по слоям t = const.
— Разностная сетка может состоять из подвижных и неподвижных узлов. Каждому расчетному узлу ставится в соответствие набор параметров течения (скорость, давление, температура и вектор состава газовой смеси и др.), а также пространственная координата и скорость движения сеточного узла. Сильному разрыву (УВ, КР) ставятся в соответствие два расчетных узла, имеющих одинаковые координаты и скорости движения, и отличающиеся наборы характеристик течения, которые соответствуют левому и правому пределам гладких решений, существующих по обе стороны от разрыва, и удовлетворяют условиям динамической совместности.
— При пересечении траекторий двух сильных разрывов пределы справа и слева перестают удовлетворять условиям динамической совместности. В этом случае происходит распад разрыва с образованием новых сильных и слабых разрывов, параметры которых находятся из решения известной задачи Римана.
— Шаг интегрирования выбирается стандартным образом из условия Куранта. Однако если пересечение какой-либо пары сеточных линий происходит раньше, то шаг интегрирования выбирается из условия точного выхода на время этого пересечения.
— Учитывается возможность зарождения УВ в области, где решение было изначально гладким. Координата и время точки зарождения определяются по пересечению звуковых характеристик одного семейства.
t, мс
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3
х, м
Рис. 1. Временная развертка течения: 1 — УВ, 2 — КР, 3 — характеристики, 4 — траектории, 5 — головная УВ.
— Параметры потока в каждом узле расчетной сетки определяются из характеристических соотношений, записанных в разностной форме. Опорные точки для соответствующих характеристик могут быть найдены как на предыдущем слое, так и на выделяемых разрывах. При восстановлении параметров течения между узлами разностной сетки используются интерполяционные полиномы от первого до третьего порядка, при этом не допускается существование экстремумов между узлами [8].
Соответствующие системы конечно-разностных уравнений, дополненные в случае необходимости условиями динамической совместности, решаются итерационным методом до достижения некоторой наперед заданной точности. Если за максимально допустимое число итераций сходимость во всех точках разностной сетки не достигается, шаг по времени уменьшается.
ЗАДАНИЕ НАЧАЛЬНЫХ ДАННЫХ.
РЕЗУЛЬТАТЫ МАТЕМАТИЧЕСКОГО МОДЕЛИРОВАНИЯ
Рассматривалась цилиндрическая камера высокого давления (КВД), заполненная гелием (7Квд = 3000 К), и низкого давления (КНД), заполненная стехиометрической водородно-воз-душной смесью при нормальных условиях (^КНд = = 101325 Па, ТКвд = 298.15 К). Для описания химических превращений в горючей смеси ис-
пользовались восемь обратимых химических реакций [12], в которых участвовали семь компонентов: Н2, ОН, Н20, 02, Н, О, теплофизиче-ские свойства которых заимствованы из [11].
На рис. 1 приводится временная развертка течения при давлении в КВД, равном 6 МПа. В начальный момент времени происходит распад разрыва в х = 0.1 м. Влево распространяется веер волн разрежения, вправо — УВ. Горючая смесь воспламеняется справа от КР, образуется волна горения, которая догоняет головную УВ и взаимодействует с ней. В результате образуется пересжатая ДВ, распространяющаяся по горючей смеси.
В КНД образуется висячая УВ, распространяющаяся в направлении левой границы расчетной области, в качестве которой используется жесткая стенка, отстоящая на 1 см от оси симметрии (разрешение особенностей течения, возникающих при взаимодействии УВ с осью симметрии, в данном случае не имеет принципиального значения). Происходит серия взаимодействий УВ со стенкой и КР, в результате несколько ударных волн распространяются по горючей смеси вслед за головной УВ. Необходимо отметить, что возмущения от стенки на момент окончания численного моделирования не достигли головной УВ.
При численном моделировании границами расчетной области являются левая стенка и головная УВ, используется подвижная разностная сетка, со
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.