УДК 519.63
БОТ: 10.18698/2309-3684-2024-3-6580
Математическое моделирование распространения пульсирующей волны газовой детонации в водородно-воз-душной смеси с использованием детальной кинетики химических реакций
© А.И. Лопато
ИАП РАН, Москва, 123056, Россия
Работа посвящена численному исследованию распространения пульсирующей волны газовой детонации. Математическая модель основана на системе уравнений Эйлера, записанной для многокомпонентного газа и дополненной моделью детальной кинетики химических реакций Petersen-Hanson. Данная модель кинетики является эффективной и работоспособной при описании процессов в водородно-воздушной и водо-родно-кислородной смеси. Вычислительный алгоритм основан на применении метода конечных объемов, ENO-реконструкции величин газодинамических параметров в расчетных ячейках, расчете потоков через грани ячеек с использованием метода Л USM, а также методов Рунге-Кутты для интегрирования по времени. Рассматривается случай прямого инициирования детонации у закрытого конца канала, заполненного стехиометрической водородно-воздушной смесью. Проведено математическое моделирование распространения пульсирующей детонационной волны. Исследуются механизмы, отвечающие за формирование высокочастотных и высокоамплитудных режимов пульсаций параметров лидирующей волны.
Ключевые слова: математическое моделирование, прямое инициирование детонации, модель Petersen-Hanson, высокочастотный режим пульсаций
Введение. Под детонацией понимают гидродинамический волновой процесс распространения по веществу экзотермической реакции со сверхзвуковой скоростью. Детонационная волна (ДВ) представляет собой самоподдерживающийся ударный разрыв, за фронтом которого непрерывно инициируется химическая реакция вследствие нагрева при адиабатическом сжатии. Структурно ДВ состоит из лидирующей волны (ЛВ) и зоны химических реакций.
Среди множества работ, посвященных исследованиям процессов детонации в газах, выделяют несколько направлений. Одно из направлений связано с работами по исследованию распространения детонации с точки зрения пожаро- и взрывобезопасности в промышленных комплексах, шахтах и других объектов, состоящих из связанных участков и туннелей. В объектах такого типа возникновение процессов детонации представляет угрозу как для конструкций, обеспечивающих корректную работу персонала, так и для здоровья персонала. С этой точки зрения, проектирование конструкций данных объектов и наличие вспомогательного оборудования может обеспечить больший
уровень безопасности и свести последствия от возникновения процессов детонации к минимуму. Другое направление включает работы по изучению инициирования детонации и взаимодействию ДВ между собой в каналах, звездах, сверхновых и других объектах с научной точки зрения. Работы из третьего направления направлены на исследования проблем детонации, и ее применении в промышленных отраслях, включая газовые горелки с импульсной детонацией и двигатели следующего поколения, такие как импульсные детонационные двигатели (ИДД) [1] и двигатели на непрерывной детонации. Стремительное развитие работ по созданию и проектированию детонационных двигателей обусловлено тем фактом, что детонационное горение является термодинамически выгодным способом сжигания топлива и превращения химической энергии топлива в полезную механическую. Из данного утверждения и отмеченных выше особенностей, связанных с детонацией, вытекает актуальность исследований, направленных на исследование и применение процессов газовой детонации, как в фундаментальной, так и прикладной науке.
Проведение натурных экспериментов с исследованием особенностей ДВ и механизмов их инициирования, распространения часто сопряжено с определенными проблемами. Отметим, что одной из таких проблем является измерение параметров течения газовой смеси в рассматриваемой области. Присутствие датчиков в области может приводить к возмущениям течения, которые не возникли бы в случае отсутствия датчиков. Другой проблемой является сама процедура измерения параметров течения. На протяжении выполнения экспериментов параметры, газа, такие как давление и температура, могут меняться в широких диапазонах. Например, минимальное и максимальное значение измеряемых величин могут различаться на несколько порядков. Поэтому прибор измерения должен обладать достаточно большим диапазоном измерения. Кроме того, если прибор находится внутри рассматриваемой области, то важным является также вопрос устойчивости его конструкции к широкому диапазону и перепаду величин таких параметров, как давление и температура.
Численный эксперимент, как метод изучения процессов с помощью математического моделирования, как правило, лишен большей части проблем, свойственных натурному эксперименту (см., например, [2]). Численные расчеты дают возможность исследовать широкий ряд проблем, связанных с течением волн детонации. Расчеты могут проводиться в широком наборе областей разного типа, включая неод-носвязные области, области с периодическими граничными условиями и т.д. Степень точности полученного решения, определяется рядом факторов, включая адекватность математической модели и численного метода, порядков аппроксимации численной схемы, и устойчивость численного метода к ошибкам машинного счета.
Хорошо известно, что ДВ, вообще говоря, является многомерным объектом. Поэтому математическое моделирование ДВ, как правило, предполагает учет многомерных эффектов и подробное описание химических процессов, происходящих в смеси. С другой стороны, математическая модель, соответствующая одномерной структуре ДВ, хоть и является более простой с точки зрения математической модели, но обеспечивает относительно широкий спектр эффектов и особенностей, исследование которых является важным с точки зрения фундаментальной и прикладной науки. Кроме того, ряд эффектов, наблюдаемых в одномерной детонации, имеет отношение к особенностям многомерной, таким, как ячейстая и спиновая детонация.
Как известно из многих численных и экспериментальных исследований, распространение ДВ связано с образованием сложного нелинейного колебательного процесса, включая пульсации параметров за фронтом ДВ. Эффекты, связанные с пульсациями, являются предметом исследования в работах многих научных коллективов. Так, в работах [3,4] представлены различные режимы пульсаций параметров одномерной ДВ в зависимости от величины энергии активации рассмотренной модельной смеси. Математическая модель включает одностадийную кинетику с необратимой реакцией и значениями параметров, предложенных, по-видимому, впервые в [5]. Проведен спектральный Фурье анализ пульсаций пикового давления в зависимости от времени. Выделены значения, соответствующие главным частотам пульсаций. Проведено сравнение с результатами аналитических и численных исследований по разным характеристикам, включая размеры предельного цикла для случая слабо неустойчивого режима пульсаций.
В работах [6,7] рассматривалась детонация в модельной водо-родно-воздушной смеси [7]. В [6] получены высокочастотные и низкочастотные моды пульсаций параметров при распространении ДВ. Показано, что с увеличением порядка аппроксимации численного метода фронт ДВ меньше размазывается по пространству, и зона химических реакций разрешается лучше. Это приводит к более точному расчету эффектов, связанных с неустойчивостями течения. Подобная тенденция подтверждается, в том числе, и результатами из [3]. В работе [7] отмечен механизм, сопровождающий процесс пульсаций параметров за фронтом ДВ. Механизм основан на формировании градиента температуры и спонтанных волн в несгоревшей смеси, прилегающей к фронту ДВ. Влияние химической кинетики на моделирование инициирования детонации посредством температурного градиента в водородно-воздушной смеси рассмотрено в [9]. Авторами предложено считать модель химической реакции пригодной для объяснения механизмов, соответствующих результатам, наблюдаемым в натурных экспериментах, если число компонент смеси и число реакций не меньше
67
некоторых пороговых величин. В противном случае, механизмы, сопровождающих химические превращения в расчетах, могут сильно отличаться от механизмов, соответствующих результатам натурных экспериментов. Так, в [9] используется детальная кинетика реакций, которая содержит 9 компонент и 19 реакций. Модель корректно описывает некоторые результаты натурных экспериментов, включая время задержки воспламенения и характеристики ламинарного течения для широкого диапазона начальных параметров. Помимо детальной кинетики, в [9] анализируются также результаты применения глобальной кинетики. Глобальная кинетика Аррениуса [8] горения модельной во-дородно-воздушной смеси представляет собой одностадийную модель химической кинетики. Модель воспроизводит часть характеристик, включая скорость пламени и ширину ламинарного пламени. Показано, что критический размер «горячих точек», определяющих временные и пространственные координаты инициирования детонации, оказывается больше в случае использования детальной кинетики. Различия объясняются тем фактом, что одностадийная модель кинетики является экзотермической для любых температур. Детальная модель кинетики имеет, как эндотермические, так и экзотермические участки в цепи реакций. Качественные различия в формировании горячих точек приводят к тому, что времена индукции, полученные с использованием одностадийной кинетики, имеют значения на несколько порядков ниже, чем результаты натурных экспериментов. Таким образом, можно отметить, что подробная кинетика обеспечивает достаточно широкий диапазон параметров, на котором полученные результаты хорошо соотносятся с экспериментом. Кроме того, детальная кинетика принимает во внимание ряд факторов, которые не учитывает совсем или учитывает не в полной мере одностадийная. Однако, нельзя не отметить, что расчетное время, выделенное на учет химических превращений по механизмам детальной кинетики, оказывается значительно выше, чем время, соответствующее учету одностадийной кинетики. Что касается результатов расчетов, полученных в ходе применения одностадийной кинетики, то здесь ситуация также неоднозначная. Параметры одностадийной кинетики часто выбираются в соответствии с рабочим диапазоном давлений и температур исходя из специальных таблиц баз данных. При калибровке параметров кинетики и проведении численных исследований с кинетикой в окрестности выбранных параметров имеется возможность получить ряд адекватных результатов и полезных рекомендаций, пригодных для изучения проблем детонации, таких как пульсирующие режимы газовой детонации.
В работе [10] проведен подробный анализ нелинейной динамики детонации в водородно-воздушной смеси с использованием детальной модели кинетики. Математическая модель включает систему уравнений Эйлера, записанной для случая многокомпонентной смеси. Модели кинетики включает 9 компонент и 38 элементарных реакций.
Численный метод повышенного порядка точности включает монотонную схему пятого порядка для повышения порядка по пространству, TVD-схему Рунге-Кутты третьего порядка для интегрирования по времени, поток Роу, метод Gaussian elimination решения системы линейных уравнений для разрешения нелинейной системы, полученной при использовании неявного метода для учета химической кинетики. Авторы рассматривали случай прямого инициирования одномерной ДВ в канале, с использованием расчетных сеток, размеры ячеек которых составляли 2.5 мкм и 12.5 мкм. Получен переход от пересжатой детонации к самоподдерживающейся с формированием двух пульсирующих мод. Для обоих размеров сетки наблюдается картина, когда сначала формируется режим высокочастотных пульсаций, а затем происходит переход к режиму высокоамплитудных пульсаций. Отмечены особенности каждого режима, включая значения частоты. Механизмы, соответствующие процессам в зоне индукции за фронтом ДВ описаны с точки зрения акустических и энтропийных волн, подобно анализу из работы Mcvey и Toong [11], и представляется адекватным описанием механизма пульсирующей детонации. В работе также показана чувствительность результатов к значениям параметров начальных условий, разрешения сетки и свойств численного метода.
Целью работы является математическое моделирование распространения пульсирующей детонационной волны в водородно-воздуш-ной смеси с использованием численного метода второго порядка аппроксимации и детальной кинетики химичеких реакций.
Математическая модель. Математическая модель основана на одномерной системе уравнений Эйлера, записанных в лабораторной системе для случая многокомпонентной смеси и дополненных детальной моделью кинетики химических реакций горения водородно-воз-душной смеси. Определяющая система уравнений имеет вид:
cU SF
dt dx
= S,
U =
p pu " 0 "
pu e , F = pu2 + p (e + p)u , S = 0 0
_pYs _ _ pYsu _ _p®s _
" _ . Pu 2 _ ^ PYS
e = IPP4 {T)~P
s=l Hs 2
С NS
-s=m z z
auci
j=i V I=1
(K -r's)
P=z'
RT,
1 M
NS NS
KjIK' -K,Це;
i =1
i =1
PY
M
p
Здесь р — плотность, и — скорость, р — давление, е — полная энергия на единицу объема, Я — универсальная газовая постоянная, у — массовая доля компоненты 5 смеси, (оз — скорость изменения
массовой доли 5, К — молярная энтальпия, ¡л,, — молярная масса, сг. — молярная концентрация, ац - коэффициент учета взаимодействия с третьим телом в ходе реакции, у'р и у"з — стехиометрические коэффициенты, К^ — коэффициент прямой скорости реакции, Кь. —
коэффициент обратной скорости реакции, N3 — полное число компонент смеси, NR — полное число реакций. Молярная энтальпия рассчитывается, как
К (Т ) = ЯТ
Г а + ^ Т + ^Т2 + ^Т3 + ^Т4 + у 2 3 4 5 Т у
где коэффициенты аь,а2^,...,а6я имеют значения, указанные в [12]. Показатель адиабаты многокомпонентной смеси в такой модели многокомпонентной смеси определяется, как
N3
К/лs
у(Т ) = 1 + Я-
5=1
Е к (Т)-я )л
5
5=1
где С — молярная теплоемкость при постоянном давлении компо-
р5
ненты 5
Ср5 (Т) = Я (а15 + а2Т + азТ2 + а^Т3 + а5^4). Скорость звука определяется для многокомпонентной смеси, как
с =
N3 у С (Т) N3 у
^ 15Ср5 (Т ) ^ У5
^^ 5=1 ¡5 5=1 ¡5
N3 у
ЕЛ(с„ (Т)-Я)
I 5=1 ¡5
В качестве детальной кинетики, описывающей горение водо-родно-воздушной смеси, рассматривается модель Рйегееп-Напвоп (РН) [13]. Модель включает N3 = 9 компонент смеси (Н2,02,Н,0,
ОН, НО, Н202, НО и N ) и NR = 18 элементарных реакций. Значения стехиометрических коэффициентов, скоростей прямых и обратных реакций, коэффициентов взаимодействий с третьими телами приведены в [13]. Адекватность и применимость данной кинетики для
численного моделирования химических реакции в водородно-воздуш-ных и водородно-кислородных смесеИ подтверждается рядом работ (см., например, ссылки в [14]).
Вычислительный алгоритм. Вычислительный алгоритм основан на методе расщепления по физическим процессам Strang [15]. При переходе с одного временного слоя на другой сначала проводится учет конвективной части уравнений без учета химической кинетики, а затем рассчитывается вклад химической кинетики без учета движения газа.
Дискретизация пространственной части определяющей системы уравнений проводится с использованием конечно-объемного метода
dU
д t
F+1/2 F-1/2
Ах
- L, (Q)
Здесь I — индекс ячейки расчетной сетки. Индексы I +1/2 и I -1/2 соответствуют правой и левой грани ячейки I соответственно, Ах — размер ячейки, Q — сеточная функция, значения элементов которой должны быть определены, Е — численный поток. численный поток Е-+1/2 рассчитывается с использованием схемы АИБМ [16], расширенной для случая многокомпонентной смеси:
F,+V2 = \М1/2
где
W (Q+)+w (Q-1 )]+i M[W (Q;)- W (Q-+1)]+P1/2,
ml/2 = m;+ m- 1,
w (Qi ) =
(ре > (0 >
puc Р1/2
рНе P = , P1/2 0 , Р1/2
PYC J г V0 J
Рг -
Здесь Н = е + р — это энтальпия смеси на единицу объема. Верхний индекс "+" (при нижнем индексе I) относится к параметрам на правой грани ячейки I, а индекс "-" (при нижнем индексе I +1) соответствует левой грани ячейки. Элементы Q+, Q-+1 определяются с использованием БКО-реконструкции второго порядка аппроксимации консервативных переменных [17]. Значения параметров М +, МТ+Х, р+, рм рассчитываются в соответствии с [16]. Заметим, что использование схемы АИБМ не является обязательным условием для пригодности численной схемы к исследованиям. Так, в работах [18, 19] численные исследования течений с ДВ и горением смеси проводятся с использованием схемы Куранта-Изаксона-Рис расчета потоков.
Интегрирование системы уравнений по времени осуществляется с использованием метода Рунге-Кутты второго порядка точности [20]:
-
Q« = Q" + At" • L (Q" ), Q;+1 = I Q; +1 Qf ) +1 д f- ц (Q« ),
2
2
2
где At" — шаг по времени, который выбирается динамически из условия устойчивости численной схемы. Верхний индекс "~" указывает на то, что численное решение является результатом первого этапа расщепления по физическим процессам.
Эффект от химических реакций без учета конвекции принимается во внимание на втором этапе расщепления по процессам. На данном этапе в каждой ячейке расчетной сетки происходит решение системы обыкновенных дифференциальных уравнений (ОДУ), описывающих изменения значений молярных концентраций компонент и температур за счет химических превращений. Система уравнений имеет вид
dCs. dt
NR f NS
i II
j=1 V 1=1
aljcl
Г -r; )
NS
K; П С - K; П c
i=1 NS f
dT_ dt
NS J NS f J
rt 11-i [ * (t ) t
NS , .
I ( С (C„ (T )-R ))
s=1
, = 1,2,..., NS,
Система интегрируется на промежутке времени At". В соответствии с методом расщепления, в качестве начальных условий системы ОДУ используются результаты, полученные на предыдущем газодинамическом этапе. Система решается с использованием неявного метода Эйлера и линеаризации по методу Ньютона.
Описанный вычислительный алгоритм основан на алгоритме, построенного для случая двухкомпонентной смеси (реагент и продукт) и постоянного значения показателя адиабаты [21]. Приведенная выше методика представляет обобщение, учитывающее особенности детальной кинетики. Вычислительный алгоритм распараллелен методом декомпозиции расчетной области с использованием библиотеки MPI. Для верификации вычислительного алгоритма и методик были проведены расчеты для ряда задач (см., например [21-23]), включая определение задержки воспламенения водородно-воздушной смеси в нульмерной постановке.
Результаты вычислительных экспериментов. Рассмотрим задачу о формировании и распространении пульсирующей ДВ, которая решается с применением описанной выше математической модели и вычислительного алгоритма. Рассматривается канал длины L = 3 м, заполненный покоящейся стехиометрической модельной водородно-
воздушной смесью. В области размера l — 10 см, примыкающей к левому концу канала, задаются повышенные по сравнению с остальной частью канала значения давления pl — 10 атм и температуры T - 3000 К. В другой части канала заданные параметры имеют значения p0 — 0.1 атм и T — 300 К. Приведенные величины представляют собой начальные условия, необходимые для численного решения системы уравнений. Область канала покрывается расчетной сеткой с размером ячейки Дх = 25 мкм.
12
10
Е
S
Transilition from the overdriven DW to the
self-sustaining one
/
/ HA Mode 1
f |
V
HF Mode \
/ Г
1 1
1 i 1
i A 1 1
4 \ 1 1
1
-L
200
400
600
800
1000
Рис. 1. Схема турбинной лопатки
Формирование профилей, качественно соответствующих профилям ДВ, происходит в момент времени около 50 мкс, который может рассматриваться как конец первого этапа расчета, соответствующего процессу инициирования ДВ. Профиль полученной ДВ содержит такие фрагменты, как лидирующий скачок, зона химических реакций и волна разрежения Тейлора. ДВ является пересжатой и при ее распространении происходит ее переход из одного состояния в другое. ДВ становится самоподдерживающейся, с параметрами, которые коррелирует с теоретическими значениями параметров Чепмена-Жуге (ЧЖ), полученными для кинетики [7]. Так, скорость ДВ достигает значения 1900 м/с, тогда как соответствующее значение скорости из [7] составляет 1993 м/с. По мере приближения детонации к состоянию ЧЖ, что соответствует моменту времени около 240 мкс, начинают
8
6
4
2
0
формироваться высокочастотные пульсации пиковых параметров, как видно из рис. 1. Другими словами, наблюдается переход к этапу высокочастотных пульсаций. Пульсации параметров обусловлены взаимодействием между волнами горения, образующимися в зоне реакции, и ЛВ. Наблюдаемый этап высокочастотных пульсаций характеризуется относительно небольшой продолжительностью по времени и сравнительно небольшими расстояниями от ЛВ, на которых формируются волны горения. Дальнейшее развитие детонационной волны приводит к переходу из высокочастотного режима к высокоамплитудному в момент времени 520 мкс. Наблюдаемые пульсации отличаются периодичностью сигнала, с периодом около 100 мкс. Пространственные распределения плотности и температуры через каждые 10 мкс на временном промежутке 620 - 660 мкс, соответствующем росту пикового давления, приведены на рис. 2, где Т^а1е = 2000 К.
1.2
з? 0.8
0.6
0.4
0.2
1.2
0.8
0.6
0.4
0.2
I I I——I I ——г—;— I—1—|——I I I—|—I—I—г
0.99 0.992 0.994 0.996 0.998 1 1.002
Рис. 2. Профили приведенной температуры (сверху) и плотности (снизу) в системе отсчета ЛВ в последовательные моменты времени, интервал времени 620 - 660 мкс (цикл акустической волны), стрелками отмечено направление движения волн с увеличением времени
1
1
Увеличение скоростей химических реакций за счет подвода теплоты от взаимодействия с внутренними волнами горения приводит к увеличению скорости выделения тепла. Происходит формирование новой волны горения, распространяющейся в сторону зоны индукции и ЛВ. Ускоряющийся фронт пламени проходит по смеси в зоне индукции, сжигая смесь в зоне индукции, до столкновения с ЛВ. Таким образом, на данном этапе наблюдается процесс формирования волны горения, распространяющейся в сторону ЛВ и, что приводит к уменьшению протяженности зоны индукции. Рассмотренная для временного 74
промежутка 620 - 660 мкс стадия распространения ДВ соответствует области, границы которой отмечены зеленой пунктирной линией на рис. 1, и называется циклом акустических волн. Затем начинается следующая стадия, которая соответствует временному интервалу 660 -700 мкс, отмеченному синей пунктирной линией на рис. 2, и называется циклом энтропийной волны.
На рис. 3 приведены пространственные распределения плотности и температуры через каждые 10 мкс на временном промежутке 660 - 700 мкс, соответствующем уменьшению пикового давления.
1.2
wu0.8 . £
0.6
0.4
0.2
—- time
___
—. Ч,
ч.
\
\ \
\ \
\ л
\
\
_ -
J\
/
/ /
¿С — ti me -
—1—1— 1 - 1 1 -
0.99 0.992
0.994 0.996 0.998 х, m
1
1.2
■0.8
■0.6
0.4
0.2
1.002
Рис. 3. Профили приведенной температуры (сверху) и плотности (снизу) в системе отсчета ЛВ в последовательные моменты времени, интервал времени 660 - 700 мкс (цикл энтропийной волны), стрелками отмечено направление движения волн с увеличением времени
1
На данной стадии возмущения распространяются от ЛВ в сторону зоны индукции, что сопровождается ослаблением ЛВ. Цикл характеризуется также увеличением размера зоны индукции и подготовкой условий для наступления нового цикла - акустического. Отметим, что полученные профили плотности и температуры качественно хорошо коррелируют с результатами работы [9].
Полученный в ходе расчетов высокоамплитудный режим включает эффект удвоения периода. Другими словами, в окрестности одного периода истории пикового давления наблюдаются двойные локальные максимумы. Например, временной интервал 700 - 800 мкс содержит два локальных максимума пикового давления. Первый максимум соответствует моменту времени около 750 мкс и значению давления 3.3 атм и отвечает механизму, рассмотренному выше. Рассмотрим
временной интервал 760 - 767 мкс, который является окрестностью второго локального максимума. На рис. 4, представлены распределения давления и массовой доли компонента И2 через каждую микросекунду из временного интервала. Масштаб давления р5са1в = 50 атм. После взаимодействия волны горения с ЛВ горение смеси происходит непосредственно за ЛВ с образованием так называемого "кармана" -области с несгоревшей до конца смесью. Горение смеси в кармане приводит к формированию новой волны горения и повышению температуры и давления газа. В результате такого горения в момент времени около 766 мкс образуется второй локальный максимум давления. Двойные колебания в других временных интервалах на этапе высокоамплитудных пульсаций происходят в соответствии с механизмом, описанном выше.
х, т
Рис. 4. Профили приведенного давления (сверху) и массовой доли компоненты Н2 (снизу) в последовательные моменты времени, интервал времени 760 - 767 мкс, красные линии отвечают моменту времени 760 мкс
Заключение. Разработан вычислительный алгоритм повышенного порядка точности для численного исследования детонационных волн в рамках детальной кинетики химических реакций. Описаны основные этапы вычислительного алгоритма, включая процедуру расщепления по физическим процессам, расчет потоков через грани ячеек с использованием схемы ЛиБМ, учет вклада химических реакций в смеси. Для верификации вычислительного алгоритма были проведены расчеты для ряда задач, включая определение задержки воспламенения водородно-воздушной смеси в нульмерной постановке. 76
Проведено численное исследование распространения пульсирующей детонационной волны с использованием детальной кинетики Petersen-Hanson горения стехиометрической водородно-воздушной смеси. Выбран размер ячейки расчетной сетки, достаточный для формирования пульсаций различных режимов. Прямое инициирование детонации осуществлялось у левого конца канала. Продемонстрирована картина детонации с четко выраженными стадиями и переходами между ними. В ходе расчетов получены, как классические стадии, связанные с формированием и распространением детонационной волны, так и переходы к пульсирующим режимам детонации, включая высокочастотный и высокоамплитудный. Рассмотрены особенности обоих режимов, а также процессы, сопровождающие каждый из режимов. Размер зоны индукции в высокочастотном режиме колеблется незначительно, в отличие от высокоамплитудного режима. Наблюдаемый сигнал пикового давления в высокоамплитудном режиме оказывается близким к периодическому, с наличием двойных локальных максимумов. Описаны механизмы образования локальных максимумов высокоамплитудного режима.
При проведении исследований использовались суперкомпьютеры МСЦ РАН МВС-10П и МВС-10П ОП.
Работа выполнена в рамках госзадания ИАП РАН.
ЛИТЕРАТУРА
[ 1 ] Kailasanath K. Review of propulsion applications of detonation waves. American Institute of Aeronautics and Astronautics, 2000, 38(9), pp. 1698-1708.
[2] Димитриенко Ю.И., Коряков М.Н., Юрин Ю.В., Захаров А.А., Сборщиков С.В., Богданов И.О. Сопряженное моделирование высокоскоростной аэротермодинамики и внутреннего тепломассопереноса в композитных аэрокосмических конструкциях. Математическое моделирование и численные методы, 2021, № 3, с. 42-61.
[3] Sharpe G.J., Falle S.A. Numerical simulations of pulsating detonations: I. Nonlinear stability of steady detonations. Combustion Theory and Modelling, 2000, 4, pp. 557-574.
[4] Lopato A.I., Utkin P.S. Toward Second-Order Algorithm for the Pulsating Detonation Wave Modeling in the Shock-Attached Frame. Combustion Science and Technology, 2016, 188, pp. 1844-1856.
[5] Lee H.I., Stewart D.S. Calculation of linear detonation instability: one-dimensional instability of plane detonation. Journal of Fluid Mechanics, 1990, 216, pp. 103-132.
[6] Lopato A.I., Utkin P.S. Mathematical modeling of pulsating detonation wave propagation using monotone numerical methods of different approximation orders. In: Transient Combustion and Detonation Phenomena: Fundamentals and Applications. Moscow: Torus Press, 2014, pp. 261-268.
[7] Xiao H., Oran E.S. Shock focusing and detonation initiation at a flame front. Combustion and Flame, 2009, vol. 203, pp. 397-406.
[8] Gamezo V., Ogawa T., Oran E. Flame acceleration and DDT in channels with obstacles: Effect of obstacle spacing. Combustion and Flame, 2008, vol. 155, pp. 302-315.
[9] Liberman M., Wang C., Qian C., Liu J. Influence of chemical kinetics on spontaneous waves and detonation initiation in highly reactive and low reactive mixtures. Combustion Theory and Modeling, 2019, vol. 23, iss. 3, pp. 467-495.
[10] Cole L.K., Karagozian A.R., Cambier J.-L. Stability of flame-shock coupling in detonation waves: 1D dynamics. Combustion Science and Technology, 2012, vol. 184, pp. 1502-1525.
[11] McVey J.B., Toong T.Y. Mechanism of instabilities of exothermic hypersonic blunt-body flows. Combustion Science and Technology, 1971, vol. 3, pp. 63-76.
[12] Visit. https://burcat.technion.ac.il/. Last accessed 10 February 2024. Last accessed 10 February 2024.
[13] Petersen E.L., Hanson R.K. Reduced kinetics mechanisms for RAM accelerator combustion. Journal of Propulsion and Power, 1999, vol. 15, iss. 4, pp. 591-600.
[14] Togashi F., Lohner R., Tsuboi N. Numerical simulation of H2/air detonation using unstructured mesh. Shock Waves, 2009, vol. 19, pp. 151-162.
[15] Toro E.F. Riemann Solvers and Numerical Methods for Fluid Dynamics: a Practical Introduction. 3rd edition. Springer, Berlin, 2009, 724 p.
[16] Liou M.-S., Steffen C.J.Jr. A new flux splitting scheme. Journal of Computational Physics, 1993, vol. 107, pp. 23-39.
[17] Shu C.-W. Essentially non-oscillatory and weighted essentially non-oscillatory schemes for hyperbolic conservation laws. NASA/CR-97-206253, ICASE Report No. 97-65, 1997, 78 p.
[18] Lopato A.I., Eremenko A.G., Utkin P.S., Gavrilov D.A. Numerical simulation of detonation initiation: the quest of grid resolution. Advances in Theory and Practice of Computational Mechanics. Smart Innovation, Systems and Technologies, 2020, 173, pp. 79 - 89.
[19] Lopato A.I., Utkin P.S. The usage of grid-characteristic method for the simulation of flows with detonation waves. Smart modeling foe engineering systems (GCM50), 2018, vol. 133, pp. 281-290.
[20] Shu C.W, Osher S. Efficient implementation of essentially non-oscillatory shock-capturing schemes. Journal of Computational Physics, 1988, vol. 77, pp. 439471.
[21] Lopato A.I., Utkin P.S. Numerical study of detonation wave propagation in the variable cross-section channel using unstructured computational grids. Journal of Combustion, 2018, act. 3635797.
[22] Lopato A.I., Favorskaya M.N., Favorskaya A.V., Petrov I.B., Lakhmi C.J. Numerical simulation of shock-to-detonation transitions using one-stage and detailed chemical kinetics mechanism. Smart Modeling for Engineering Systems, 2021, pp. 79-88.
[23] Лопато А.И. Математическое моделирование инициирования детонации в канале с профилированным торцом с использованием параллельных вычислений. Математическое моделирование и численные методы, 2023, № 4, c. 15-26.
Статья поступила в редакцию 22.03.2024
Ссылку на эту статью просим оформлять следующим образом:
Лопато А.И. Математическое моделирование распространения пульсирующей волны газовой детонации в водородно-воздушной смеси с использованием детальной кинетики химических реакций. Математическое моделирование и численные методы, 2024, № 3, с. 65-80.
Лопато Александр Игоревич — канд. физ.-мат. наук, с.н.с. ИАП РАН. e-mail: [email protected].
Mathematical Modeling of the Propagation of a Pulsating Detonation Wave in a Hydrogen-Air Mixture Using Detailed Kinetics of Chemical Reactions
© A.I. Lopato ICAD RAS, Moscow, 123056, Russia
The work is dedicated to the numerical study of pulsating gaseous detonation wave propagation. The mathematical model is based on the Euler equations written for the multicom-ponent gas and supplemented by the detailed chemical reactions model of Petersen and Hanson. to describe the combustion of the hydrogen-air mixture. This kinetics model is effective and efficient in describing processes in hydrogen-air and hydrogen-oxygen mixtures. The numerical algorithm is based on the finite volume approach, essentially non-oscillatory scheme, A USM numerical flux and the Runge-Kutta method for time integration. Direct initiation of detonation at the closed end of a channel filled with a stoichio-metric hydrogen-air mixture is considered. Mathematical modeling of the propagation of a pulsating detonation wave was carried out. The peculiarities of high-frequency and highamplitude pulsations modes are discussed.
Keywords: mathematical modeling, direct detonation initiation, Petersen-Hanson kinetics, high-frequency pulsating mode
REFERENCES
[ 1 ] Kailasanath K. Review of propulsion applications of detonation waves. American Institute of Aeronautics and Astronautics, 2000, 38(9), pp. 1698-1708.
[2] Dimitrienko Yu.I., Koryakov M.N., Yurin Yu.V., Zakharov A.A., Sborschikov S.V., Bogdanov I.O. Coupled modeling of high-speed aerothermodynamics and internal heat and mass transfer in composite aerospace structures. Mathematical Modeling and Computational Methods, 2021, no. 3, pp. 42-61.
[3] Sharpe G.J., Falle S.A. Numerical simulations of pulsating detonations: I. Nonlinear stability of steady detonations. Combustion Theory and Modelling, 2000, 4, pp. 557-574.
[4] Lopato A.I., Utkin P.S. Toward Second-Order Algorithm for the Pulsating Detonation Wave Modeling in the Shock-Attached Frame. Combustion Science and Technology, 2016, 188, pp. 1844-1856.
[5] Lee H.I., Stewart D.S. Calculation of linear detonation instability: one-dimensional instability of plane detonation. Journal of Fluid Mechanics, 1990, 216, pp. 103-132.
[6] Lopato A.I., Utkin P.S. Mathematical modeling of pulsating detonation wave propagation using monotone numerical methods of different approximation orders. In: Transient Combustion and Detonation Phenomena: Fundamentals and Applications. Moscow: Torus Press, 2014, pp. 261-268.
[7] Xiao H., Oran E.S. Shock focusing and detonation initiation at a flame front. Combustion and Flame, 2009, vol. 203, pp. 397-406.
[8] Gamezo V., Ogawa T., Oran E. Flame acceleration and DDT in channels with
obstacles: Effect of obstacle spacing. Combustion and Flame, 2008, vol. 155, pp. 302-315.
[9] Liberman M., Wang C., Qian C., Liu J. Influence of chemical kinetics on spontaneous waves and detonation initiation in highly reactive and low reactive mixtures. Combustion Theory and Modeling, 2019, vol. 23, iss. 3, pp. 467-495.
[10] Cole L.K., Karagozian A.R., Cambier J.-L. Stability of flame-shock coupling in detonation waves: 1D dynamics. Combustion Science and Technology, 2012, vol. 184, pp. 1502-1525.
[11] McVey J.B., Toong T.Y. Mechanism of instabilities of exothermic hypersonic blunt-body flows. Combustion Science and Technology, 1971, vol. 3, pp. 63-76.
[12] Visit. https://burcat.technion.ac.il/. Last accessed 10 February 2024. Last accessed 10 February 2024.
[13] Petersen E.L., Hanson R.K. Reduced kinetics mechanisms for RAM accelerator combustion. Journal of Propulsion and Power, 1999, vol. 15, iss. 4, pp. 591-600.
[14] Togashi F., Lohner R., Tsuboi N. Numerical simulation of H2/air detonation using unstructured mesh. Shock Waves, 2009, vol. 19, pp. 151-162.
[15] Toro E.F. Riemann Solvers and Numerical Methods for Fluid Dynamics: a Practical Introduction. 3rd edition. Springer, Berlin, 2009, 724 p.
[16] Liou M.-S., Steffen C.J.Jr. A new flux splitting scheme. Journal of Computational Physics, 1993, vol. 107, pp. 23-39.
[17] Shu C.-W. Essentially non-oscillatory and weighted essentially non-oscillatory schemes for hyperbolic conservation laws. NASA/CR-97-206253, ICASE Report No. 97-65, 1997, 78 p.
[18] Lopato A.I., Eremenko A.G., Utkin P.S., Gavrilov D.A. Numerical simulation of detonation initiation: the quest of grid resolution. Advances in Theory and Practice of Computational Mechanics. Smart Innovation, Systems and Technologies, 2020, 173, pp. 79 - 89.
[19] Lopato A.I., Utkin P.S. The usage of grid-characteristic method for the simulation of flows with detonation waves. Smart modeling foe engineering systems (GCM50), 2018, vol. 133, pp. 281-290.
[20] Shu C.W, Osher S. Efficient implementation of essentially non-oscillatory shock-capturing schemes. Journal of Computational Physics, 1988, vol. 77, pp. 439471.
[21] Lopato A.I., Utkin P.S. Numerical study of detonation wave propagation in the variable cross-section channel using unstructured computational grids. Journal of Combustion, 2018, act. 3635797.
[22] Lopato A.I., Favorskaya M.N., Favorskaya A.V., Petrov I.B., Lakhmi C.J. Numerical simulation of shock-to-detonation transitions using one-stage and detailed chemical kinetics mechanism. Smart Modeling for Engineering Systems, 2021, pp. 79-88.
[23] Lopato A.I. Mathematical modelling of detonation initiation in the channel with the profiled end using parallel computations. Mathematical Modeling and Computational Methods, 2023, no. 4, pp. 15-26.
Lopato A.I., Cand. Sc. (Phys.-Math.), Senior Researher ICAD RAS. e-mail: [email protected].