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

Вычислительные технологии

Том 18, № 4, 2013

Численная схема высокого порядка для

моделирования динамики в смеси реагирующих газов

>к

и инертных частиц*

Д. А. Тропин, А. В. Фёдоров Институт теоретической и прикладной механики СО РАН, Новосибирск, Россия e-mail: d.a.tropin@yandex.ru, fedorov@itam.nsc.ru

Для расчёта ударно-волновых и детонационных процессов, протекающих в одномерных нестационарных течениях инертных/реагирующих газовых смесей и твёрдых инертных частиц, предлагается схема высокого порядка точности, основанная на подходе TVD и использовании решателя жёстких систем RADAU5. Верификация результатов расчётов по данной схеме проводится путём решения тестовых задач распространения ударных и детонационных волн в двухфазной смеси (задачи о структуре ударной волны в газовзвеси воздуха и инертных частиц и о детонационной динамике в смеси водорода — кислорода с мелкими инертными частицами песка). Сопоставление результатов расчётов, полученных с использованием предложенного метода, и стационарного решения показало вполне удовлетворительную точность метода. В качестве приложения рассмотрена задача о подавлении детонации водород-кислородной смеси путём вброса инертных частиц.

Ключевые слова: газовзвеси, ударные и детонационные волны, численное моделирование, схемы TVD.

Введение

Проблема математического моделирования течений инертных/реагирующих газовзвесей представляет большой интерес при исследовании ослабления и подавления детонации газовых смесей. Это в свою очередь актуально при изучении вопросов взрыво- и пожаробезопасности в шахтах, а также на производствах, использующих реагирующие газовые смеси в качестве рабочих тел. Один из способов ослабления и подавления нежелательной детонации — введение инертных компонентов в реагирующую газовую смесь, например, вброс инертных частиц. Добавление химически инертных твёрдых частиц является эффективным способом контролирования и модификации процессов горения и детонации в газовых системах. Этот способ может быть использован и для того, чтобы уменьшить скорость детонации взрывчатых веществ в некоторых технологических приложениях. Ниже будут рассмотрены численный метод и проблемы, которые возникают при представлении данных процессов в рамках модели механики гетерогенных реагирующих сред.

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

*Работа выполнена при частичной финансовой поддержке РФФИ (коды проектов 11-08-00144-а, 11-01-00634).

Не останавливаясь на ранних работах, ссылки на которые можно найти в [1-6], упомянем лишь несколько работ последних лет, близких к решению интересующих нас задач. Так, в [7] использовалась схема Годунова первого порядка, а в [8] её модификация (схема второго порядка) для решения задач механики гетерогенных сред в двухскоростном и двухтемпературном приближении, учитывающая давление и объёмную концентрацию фаз. Решение подобных начально-краевых проблем механики гетерогенных сред с учётом сил межфазного взаимодействия рассматривалось в работе [9]. Здесь наряду со схемой Годунова применялась схема МИБСЬ второго порядка точности. Для разрешения производных от объёмной концентрации использовалось предложение [10] об однородности потока. В качестве тестовой в [9] была решена задача о распаде произвольного разрыва, разделяющего ударную трубу на две части, заполненные жидким и газообразным (паром) додеканом соответственно.

Заметим, что в работах [3, 4] решалась начальная задача для системы уравнений механики гетерогенных сред в двухскоростном, двухтемпературном приближении без учёта объёмной доли частиц и с простой кинетикой химических превращений в твёрдой фазе. В этом случае использовались ТУВ-схема для аппроксимации уравнений газовой фазы и схема Джентри — Мартина — Дэйли для аппроксимации уравнений твёрдой фазы. На основе предложенного метода авторами данных работ были решены различные задачи о распространении волн гетерогенной детонации в алюминиево-кислородной смеси.

Во всех упомянутых работах уравнения математических моделей механики гетерогенных сред (МГС) либо не учитывают объёмные концентрации фаз, либо — при учёте последних — решаются с использованием некоторых предположений для определения вида аппроксимации производных от объёмных концентраций.

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

1. Математическая модель детонационных процессов в смесях реагирующих газов и твёрдых инертных частиц

Рассмотрим ударную трубу, заполненную стехиометрической смесью водорода — кислорода и облаком частиц песка, длиной Ь, расположенную на расстоянии I = 0.3 м от диафрагмы, разделяющей камеру высокого давления (КВД, примыкает к торцу ударной трубы и имеет длину 0.2 м), заполненную в начальный момент времени кислородом с высокими параметрами (р = 20 атм, Т = 3000 К), и камеру низкого давления (КНД), заполненную стехиометрической водородокислородной смесью с начальными параметрами р0 = 1 атм, Т0 = 300 К. Передняя кромка облака находится на расстоянии 0.5 м от торца ударной трубы. После разрыва диафрагмы в КНД происходит распространение ударной волны (УВ) и инициирование детонационной волны (ДВ), которая до попадания в облако частиц выходит на режим Чепмена — Жуге (Ч — Ж).

Система уравнений, описывающая математическую модель механики реагирующей смеси газов (водород — кислород, продукты химических превращений) и инертных частиц в рамках неравновесной динамики смеси газа и твёрдых частиц, имеет вид

др1 + д (р1«1) = 0

дг

дх

д (Р1М1) д (р!«? + т,1р) дш,1

+---= Р^Г —

дг

дх

дх

д ( р1«1 ( Е1 + т1 —

д (рЕ) + \ \ р1

дг

дх

др2 + д (Р2М2)

дт1

дг

дх

дг

дх

0,

д (р2«2) + д (р2«2 + т2р)

дт2

+1

д ( Р2«2 ( Е2 + т,2 — д (Р2Е2) + \ \ р2

дг

дх

дт2 , , {

+я +

п~

Здесь индекс 1 относится к смеси газов, индекс 2 — к частицам; Ег = ег + —--полная

энергия фазы или компонента; рг, иг, ег — средняя плотность, скорость и внутренняя энергия г-й фазы; р — давление смеси; тг — объёмная концентрация г-й фазы; рг = тгргг, где ргг — истинная плотность г-й фазы, ш = в«1 + (1 — в) и2, в € [0, 1] — параметр, показывающий, какая часть энергии за счёт сил трения передается газу, какая — частицам. Внутренняя энергия реакционноспособной газовой смеси определяется из соотношения

е1 = 111

Т1 + ^ £,аЬоа 1T00,

а=1

где с€, ср — удельные теплоёмкости смеси при постоянном объёме и давлении соответ-

8

ственно, £а — относительная массовая концентрация компонента а, с€ = ^ с0а£а, к0а —

а=1

энтальпия образования компонента. Величина а меняется от 1 до 8, Т00 = 298.15 К. Внутренняя энергия частиц е2 = сю,2Т2.

Законы силового взаимодействия и теплообмена между газом и частицами выражаются известными формулами

I

3т2ри 4в

Сп \щ — «2| («1 — «2) , Я

6т2Л1 в2

Ш (Т2 — Т1).

При этом зависимость числа Нуссельта от чисел Рейнольдса и Прандтля принимается в виде Ыи = 2 + 0.6Ке1/2Рг1/3. Для коэффициента сопротивления с учётом его зависимости от чисел Рейнольдса и Маха относительного движения частиц используется формула

Сп =

Р11в \и1 — «2 \ где Ке =-, М12

1+е

0.43 "м4.67 м12

24 4 0.38+ — + —=

\«1 — «2\ —р11 —ЪР '

Система (1), дополненная уравнением состояния для газовой смеси в целом

Со

р = РиТ 1

а=1

а

Шп

и уравнениями детальной кинетики [11, 12]

¿г

Р

Ма

Мт

Г=1

в=1

Мв

V вт

кы Д

в=1

Мв

'вт

где Мг — порядок г-й реакции, а Е [1, 8], позволяет после постановки соответствующей начально-краевой задачи инициирования/гашения детонации рассчитать процесс распространения ударной/детонационной волны в канале, заполненном смесью реагирующих газов и движущимся облаком инертных частиц.

1

"аг "аг

2. Численный метод решения уравнений, учитывающих объёмную долю частиц

При решении системы одномерных нестационарных уравнений динамики смеси газов для аппроксимации по времени использован конечно-разностный метод типа универсального алгоритма [13]. Для пространственной аппроксимации уравнений газовой смеси применена ТУВ-схема третьего прядка точности с расщеплением вектора потоков по Ван Лиру. Детали алгоритма подробно изложены в [14]. Отметим, что учёт объёмной концентрации частиц приводит к появлению этой величины в векторах потоков. В настоящей работе проведён переход к новой искомой переменной Р = т1р, при этом вид векторов потоков становится таким же, как в случае чистой смеси газов.

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

др1 д (р1П1)

дг

дх

W1

1 ,Р,

д (р1и1) д (р1п\ + Р) Р дт1

+

дг дх т1 дх

д (рЕ) д [(Р1Е1 + Р) щ] Р дт1

- I = + Р1,и + Я1>и,

дг

дх

т1 дг

- I (П1 - Щ) - q = Щ1,е + Р1,е + #1,,

др2 дг

д (р2П2) дх

Щ

2 ,р,

д (Р2Щ2) д (Р2Щ2 + (т2/т1) Р) Р дт2 _ утг , „

+ —---Г I = Щ2,и + Р2,и + #2,'

дг дх (1 - т2) дх

д (Р2Е2) д [(Р2Е2 + (т2/т1) Р) Щ2] Р дт2

дг

дх

(1 - т2) дг

+ I (щ - Щ2) + q =

Щое + Р2,е +

Здесь

в

1,и

-в

2,и =

Р дт1

т1 дх Я1,

Р

дт2

(1 — т 2 ) дх

, Р1

Р2

Р дт 1

Р дт 2

е = —р2,е

т1 дг (1 — т2) дг

-Я

2,и

— ¡, Я1,е = —Я2,е = —¡Ш — д.

Как сказано выше, для построения разностной схемы аппроксимации по времени уравнений (12) используется метод универсального алгоритма [13]:

др1

дг

р™+1

рП

т

де1

~дг

д«1 ~дг ="

е?+1 - еП

и^1 — и?

-1 (Щ^ — и?^) + 1 Р1ии + р1- Я?,и, р1 р1 р1

т

р~ — (е? + (и?)2 /2) — и? (Ш"и — и?Щд] + р1

1

1

+ — Р"е + — Я?

р

др2 дг

р?+1

р р2

т

1,е

ди2 ~дг

и?+1 — и?

- Щи — и?Щр) + р-Р?и + Р1-Я?,и, р2 , ,Р р2 р2

2

де

~дг

е"+1 _ е?

е2 е2 т

р- [Щ,е — (е? + (и?)2 /2) — и? (Ш2,и — и?Ш2,р)] + -Р™е + -Я^. р2 р2 р2

(5)

Временная аппроксимация (5) имеет первый порядок точности. Максимально допустимое значение шага по времени определяется условием Куранта:

тах ^ + аг) х ^ ^Х^ ^ 1

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

В векторной форме уравнения механики гетерогенных сред для газовой смеси имеют вид

дЕ дЦ

дг + дх

и,

где

Е

Ц (

и

Я1,Р

Я1,и Я1,е

0

Р дтл

р1и1

р1и\ + Р

(рЕ + р) и1 \

т1 дх Р дт 1

— I

т1 дг

+ — Я

Введём полную энтальпию И1 = е1 + Р/р1 + п\/2. Тогда Q можно представить как

д

Р1П1

Р1Щ + Р

Р1И1П1

Для построения пространственной аппроксимации первых трёх уравнений системы (4) с применением ТУВ-подхода необходимо расщепить вектор потоков д (6). В данном случае для получения устойчивой противопотоковой аппроксимации правых частей разностной схемы (5) использован метод расщепления вектора потоков по Ван Лиру [14]. Разобьём вектор потоков д на положительную и отрицательную составляющие: д = д+ + д-.

Пусть

Я+ = р1а (М1 + 1)2 к_

4 '

Р1а (М1 - 1)2

~4

где М1 = п1/а — число Маха, а = л/грТр! — замороженная скорость звука, 7 Ср, 1 /су, 1 показатель адиабаты газа. Тогда

р1п1' п1 > а' Я+' —а < п1 < а, Q О, п1 < —а'

О,

п1 > а'

1,р

Я ' —а <п1 < а' р1п1' п1 < —а'

Q+

1,«

Р1П1 + Р'

Я+ О'

(7 — 1) п1 + 2а

1

п1 > а' — а < п1 < а' п < а'

Q-

Q+e

Q

1,е

О'

Я

(7 — 1) п1 — 2а

1

Р1 п21 + Р' (Р1Е1 + Р) Щ'

Я+ О'

[(7 — 1) щ + 2а]

2(72 — 1)

О'

[(7 — 1) т — 2а]2 Я 2 (72 — 1) (Р1Е1 + Р) Щ'

п1 > а' — а < п1 < а' п1 < — а'

2 п1 > а' ' — а < п1 < а' п1 < — а'

п1 > а' — а < п1 < а' п1 < —а.

Далее для построения ТУВ-схемы введём в направлении х равномерную дискретную сетку {хг' г = 1'... N} с шагом Ах = хг — хг-1 и аппроксимируем производные потоков по пространственным переменным:

д_ дг+1/2 — дг-1/2 (8)

дх) г ~ Ах ' ()

^±1/2 = д-±1/2 + д+±1/2'

а

0+1/2 = - 4 [(1 - к) 6+ + (1 + к) г] (д-+1),

0+1/2 = + а [(1 - к) б- + (1 + к) 6+] (0+),

где 6+ и 6- определяются по формулам

6+

0,

ш1п(|Д+| , 0 |Д-|)

signД+signД < 0, sigпД+sigпД- > 0,

0, sigпД+sigпД < 0,

шт(|Д-| , 0 |Д+|), signД+signД- > 0.

Параметр 0 изменяется в пределах 1 < 0 <

3 — к 1 - к

и определяет порядок схемы.

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

Для пространственной аппроксимации уравнений инертных частиц применена схема Джентри, Мартина и Дэйли [15], которая обладает свойствами консервативности и обеспечивает второй порядок точности для члена, описывающего конвекцию. Как уже было сказано, эта схема использовалась в [3, 4] для аппроксимации уравнений МГС без учёта объёмной доли частиц. В настоящей работе данная схема будет применяться для расчёта уравнений МГС, учитывающих объёмную концентрацию частиц. В этой схеме по каждую сторону от узловой точки пространственной сетки находятся средние значения скоростей на границах ячейки. Знак скоростей определяет, из какого именно узла ячейки надо взять значение я для написания разностей против потока. В случае одномерных уравнений будем иметь

/1 / ^ _ /1, П

щи _ _ иЯЯЯ иЬЯЬ

где

Щ? _ щгаи + Я _ 1/2 «+1 + и?) ,

Дх

д (т2/т1) Р

'12)

дх

и.

1/2 (и? + иП-1)

я _ (р2, Р2и2, Р2Е2 + (т2/т1) Р). Значения я определяются следующим образом:

Яя

Яг, ия > 0, Г Яг—1, иЬ > 0,

Яг+1, ия < 0, ЯЬ \ Яг, иЬ < 0.

Для пространственной аппроксимации производных, содержащих объёмную концентрацию частиц, используем выражение

д ( т2 /т1 ) Р дх

Р дт2 т2 дР

(1 — т2) дх т1 дх

Производные от давления будем рассчитывать разностями против потока:

т2 дР ^ (т2)г Рг - Р—

т1 дх (т1)г Дх

6

Производные от объёмной концентрации по времени выражаются через производные по пространству через уравнение неразрывности континуума частиц:

дш.2 дт2п2 дЬ дх

Производные по пространству в свою очередь рассчитываются по формуле (12), в которой с = т2.

Рассмотрим один из вариантов аппроксимации конвективной производной по данной схеме (пд > 0, п< 0). Аппроксимация производной в этом случае имеет вид

С«+1 + О — С« + <-1)

2Ах .

Первое дифференциальное приближение будет следующим:

с (п + п/Ах + п" Ах2/2 + п///Ах3/6 — п + п/Ах — п" Ах2/2 + п///Ах3/6) =

2А =

= сп/ + сп///Ах2/6.

Таким образом, видно, что представленная схема даёт второй порядок аппроксимации уравнений МГС.

После определения динамических и термодинамических параметров газовой смеси решалась задача по вычислению относительных массовых концентраций компонентов смеси. Для этого на каждом временном шаге решалась задача Коши для обыкновенных дифференциальных уравнений (ОДУ) химической кинетики, записанных вдоль траекторий газовых частиц. Реализация осуществлялась решателем жёстких систем ИАБАШ. Для интерполяции значений относительной массовой концентрации в узлы сетки использовался кубический сплайн. Оказалось, что применение последнего для интерполяции разрывных решений, возникающих на контактных разрывах, разделяющих различные газы, приводит к возникновению осцилляций, так как кубический сплайн выстраивает только гладкую функцию. Эти осцилляции удалось подавить использованием сглаживающего фильтра, который работает по принципу шт-шо^

На равномерных сетках данный метод приводит к большим временным затратам, поскольку для задач распространения и подавления детонации необходимо разрешать структуру ДВ, т.е. нужно, чтобы на структуру ДВ приходилось достаточное количество ячеек разностной сетки. В результате для равномерной сетки следует задавать очень большое количество ячеек на всю расчётную область. Поэтому в данном случае вполне уместно применить метод адаптивных сеток с привлечением методологии, описанной в [16] для одного уравнения переноса (уравнения Хопфа). Обобщение данного метода на уравнения неравновесной газовой динамики на примере задачи о распространении ДВ предложено в [17]. Эффективность метода адаптивных сеток путём сравнения с результатами решения задачи о распространении ДВ на равномерной сетке показана в работах [17, 18].

3. Расчёты распространения одномерных ударных и детонационных волн в двухфазной смеси

Проведём тестовые расчёты распространения ударных и детонационных волн в смеси реагирующих газов и твёрдых инертных частиц.

3.1. Структура ударных волн в смеси газа и частиц

Применим предложенную численную технологию для определения структуры УВ в смеси воздуха и инертных частиц песка диаметром 10 мкм на основе решения нестационарной задачи о распаде разрыва. В качестве начального будем задавать разрыв между двумя равновесными состояниями. Первое (справа от разрыва) — начальное состояние, соответствующее атмосферным условиям: р = р0 = 1 атм, Т1 = Т2 = Т0 = 300 К, п1 = п2 = 0, второе (слева от разрыва) — равновесное состояние за возникающей УВ, движущейся со скоростью О, когда температуры и скорости газа и частиц далеко вниз по потоку от УВ приближаются к равновесным значениям п1 = п2 = пе, Т1 = Т2 = Те. Очевидно, что при распаде такого разрыва в одну сторону распространяется замороженная УВ, за фронтом которой реализуется неравновесное течение с непрерывным переходом к равновесному в области релаксации параметров. В качестве точного решения данной задачи используем расчётные данные, полученные из решения краевой задачи для ОДУ [6] в стационарной модели МГС.

Структуры УВ, полученные по предложенным в работе моделям, приведены на рис. 1 в виде распределения температур и скоростей фаз. Очевидно, что такие структуры состоят из замороженной УВ, за которой следует зона неравновесного течения, в котором происходит релаксация к конечному равновесному по скоростям и температурам состоянию. Расчёты показали, что структура УВ, полученная в нестационарных расчётах, практически повторяет таковую, полученную в стационарной задаче (см. сплошные и штрихпунктирные линии на рис. 1): значения замороженных и равновесных параметров совпадают, длины зон релаксации по скоростям практически одинаковы, по температурам — различаются незначительно. Тем самым показано, что в нестационарном случае реализуется и сохраняется структура УВ, полученная в расчётах по стационарной модели.

Кроме того, были проведены расчёты по предложенной схеме на последовательности сгущающихся сеток (с различными шагами разностной сетки — Н Н/2). Полученные

Г, К

450

400

350

300

3.2 3.3 3.4 3.5 3.6 X, М 3.2 3.3 3.4 3.5 3.6 х. М

Рис. 1. Распределения термодинамических параметров в смеси газа (воздуха) и частиц песка (а — скорости газа и частиц, б — температуры газа и частиц), рассчитанные по нестационарной модели методом ТУБ и стационарной модели из решения ОДУ

; I

N

N -

- \ -

> :

г -«1, Дх 1 ^ = 1 ММ -

-----щ, Дх = 2

: -------н2, Дх = 1 :

- ........................... «2, Дх = 2 -

- :

т, к

400

350

300

250

200

150

100

50

0

450

400

350

300

Т\, Дх = 1 мм Ти Ах = 2 Т2, Ах = 1 Т2, Ах = 2

2.7

2.8

3.1

3.2

— 300

_I_-

3.2 X. М

Рис. 2. Распределения термодинамических параметров в смеси газа (воздуха) и частиц песка (а — скорости газа и частиц, б — температуры газа и частиц), рассчитанные для различных шагов разностной сетки

а

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

3.2. Распространение детонационных волн в смеси реагирующих газов и твёрдых инертных частиц

Рассмотрим приложение представленного расчётного метода к решению задачи о подавлении детонации реагирующей газовой смеси путём вброса инертных частиц [19, 20]. На примере задачи о распространении ДВ в смеси водорода — кислорода и инертных частиц песка в трубе с геометрией, описанной выше. На рис. 3, а представлены распределения скоростей газа и частиц диаметром 100 мкм, т2,0 _ 10-4. В такой газовзвеси детонационная волна движется стационарно [20]. Как известно, для подавления ДВ необходимо, чтобы химически инертные частицы "сработали" в зоне химической реакции. После выхода частиц за плоскость Чепмена — Жуге значительного воздействия на процесс подавления детонации они не оказывают (см., например, [21]). Из рис. 3, а видно, что в зоне реакции частицы остаются практически неподвижными (их скорость в структуре ДВ достигает и2 ~ 50 м/с при скорости газа и1 _ 2100 м/с), а разгоняются лишь за этой зоной, где подавляющего влияния уже не оказывают.

На рис. 3, б приведены распределения скоростей газа и частиц, а также объёмная концентрация частиц диаметром 10 мкм в области течения. Очевидно, что более мелкие частицы в зоне реакции разгоняются сильнее (и2 ~ 300 м/с) по сравнению с частицами 100 мкм. Кроме того, оказалось, что уменьшением размера частиц при сохранении их объёмной концентрации распространение ДВ по газовзвеси замедляется (это обусловлено возрастающим сопротивлением частиц), т. е. за то же время £ _ 5 • 10-4 с ДВ при уменьшении диаметра частиц проходит меньшее расстояние по газовзвеси. Распределение объёмной концентрации частиц (рис. 3, а) показывает, что за фронтом УВ

и, м/с

----частицы

-------объёмная концентрация частиц

2000

1500

1000

500

I ... I

0.00014

II, м/с

2500

0.00012

00001 о

а

0 8

1.2

1.4

1 .0

Рис. 3. Распределения скоростей газа (водород —кислород) и частиц и объёмной концентрации частиц в ударной трубе; начальная объёмная концентрация частиц Ш2,о = 10-4, момент времени £ = 5 ■ 10-4 с; частицы диаметром 100 мкм (а), 10 мкм (б)

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

Представляется интересным привести интегральные зависимости, характеризующие процесс подавления детонации. Так, в работе [22] были получены зависимости скорости детонации от диаметра частиц для различных значений объёмной концентрации, на которых с ростом диаметра частиц при их неизменной доле в объёме скорость детонации увеличивалась и все расчётные данные стремились к определённому асимптотическому значению, соответствующему замороженному течению. Замороженное течение реализуется при больших диаметрах частиц, когда частицы не успевают разогнаться и их скорость в детонационной волне намного меньше скорости газа, как было видно выше для частиц диаметром 100 мкм. Полученные кривые по мере роста диаметра частиц сближаются, но не пересекаются, и чем больше объёмная концентрация частиц, тем ниже лежит соответствующая кривая. Следует отметить, что при уменьшении объёмной концентрации частиц этот асимптотический предел приближается к скорости детонации в чистой водород-кислородной смеси.

Кроме того, были рассчитаны концентрационные пределы детонации. Получено, что с увеличением диаметра частиц их начальная объёмная концентрация в облаке, необходимая для подавления ДВ, возрастает. Эти результаты качественно согласуются с данными работы [23].

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

Подобные зависимости скорости детонационной волны от диаметра частиц для различных объёмных концентраций были описаны ещё в [1], а позднее в [24]. Представленные в настоящей работе результаты уточняют эти данные, полученные в рамках более ранних математических моделей с приведённой кинетикой. В [24] с использованием модельной кинетики было показано, что скорость ДВ монотонно возрастает с увеличением диаметра частиц независимо от начальной объёмной концентрации частиц, однако этот рост очень мал. Например, при m2,0 = Ю-4 и увеличении диаметра частиц на порядок (от 5 до 5G мкм) рост скорости детонации составил только Б.Б %. Изложенные выше результаты, полученные с учётом детальной химической кинетики в данном диапазоне диаметров частиц, показывают более значительное увеличение скорости ДВ (на 25%).

Выводы

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

2. Численный метод верифицирован на основе показанного качественного и количественного соответствия расчётных данных по структуре УВ в двухфазной смеси, полученных в нестационарном и стационарном течениях.

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

Список литературы

[1] Казаков Ю.В., Фёдоров А.В., Фомин В.М. Режимы нормальной детонации в ре-лаксирующих средах // Физика горения и взрыва. 1989. № 1. C. 119-127.

[2] Казаков Ю.В., Миронов Ю.В., Фёдоров А.В. Расчёт детонации газовой смеси при наличии инертных твёрдых частиц // Моделирование в механике. 1991. Т. 5(22), № 3. C. 152-162.

[3] Хмель Т.А. Численное моделирование двумерных детонационных течений в газовзвеси реагирующих твёрдых частиц // Математическое моделирование. 2004. Т. 16, № 6. C. 73-77.

[4] Хмель Т.А., Фёдоров А.В. Численные технологии исследования гетерогенной детонации газовзвесей // Там же. 2006. Т. 18, № 8. C. 49-63.

[5] Рычков А.Д. Математическое моделирование газодинамических процессов в каналах и соплах. Новосибирск: Наука, 1988. 222 с.

[6] Фёдоров А.В., Федорченко И.А. Численное моделирование распространения ударной волны в смеси газа и твёрдых частиц // Физика горения и взрыва. 2010. Т. 46, № 5. C. 97-107.

[7] Deledicque V., Papalexandris M.V. An exact Riemann solver for compressible two-phase flow models containing non-conservative products // J. of Comput. Phys. 2007. No. 222. P. 217-245.

[8] ScHWENDEMAN D.W., Wahle C.W., Kapila A.K. The Riemann problem and a highresolution Godunov method for a model of compressible two-phase flow // Ibid. 2006. No. 212. P. 490-526.

[9] Zein A., Hantke M., Warnecke G. Modelling phase transition for compressible two-phase flow applied to metastable liquids // Ibid. 2010. No. 229. P. 2964-2998.

[10] Abgrall R. How to prevent pressure oscillations in multicomponent flow calculations: A quasi conservative approach // Ibid. 1996. No. 125. P. 150-160.

[11] Tien J.H., Stalker R.J. Release of chemical energy by combustion in a supersonic mixing layer of hydrogen and air // Combustion and Flame. 2002. No. 130. P. 329-348.

[12] Бедарев И.А., Фёдоров А.В. Сравнительный анализ трёх математических моделей воспламенения водорода // Физика горения и взрыва. 2006. Т. 42, № 1. C. 26-33.

[13] Ковеня В.М., Яненко Н.Н. Метод расщепления в задачах газовой динамики // Новосибирск: Наука, 1981. 384 с.

[14] Van Leer B. Flux-vector splitting for the Euler equations // Lecture Notes in Phys. 1982. Vol. 170. P. 507-512.

[15] Gentry R.A., Martin R.E., Daly B.J. An Eulerian differencing method for unsteady compressible flow problems // J. of Comput. Phys. 1966. Vol. 1. P. 87-118.

[16] Ceniceros H.D., Hou T.Y. An efficient dynamically adaptive mesh for potentially singular solutions // Ibid. 2001. No. 172. P. 609-639.

[17] Бедарев И.А., Фёдоров А.В. Тестирование метода адаптивных сеток на расчётах одномерных детонационных волн // Вычисл. технологии. 2009. Т. 14, № 3. C. 14-24.

[18] Тропин Д.А. Физико-математическое моделирование ослабления и подавления детонации в реагирующих газовых смесях инертными частицами: Дис. ... канд. физ.-мат. наук. Новосибирск: ИТПМ СО РАН, 2012. 118 с.

[19] Фёдоров А.В., Фомин П.А., Фомин В.М. и др. Физико-математическое моделирование подавления детонации облаками мелких частиц. Новосибирск: НГАСУ (Сибстрин), 2011. 156 с.

[20] Фёдоров А.В., Тропин Д.А., Бедарев И.А. Математическое моделирование подавления детонации водород-кислородной смеси инертными частицами // Физика горения и взрыва. 2010. Т. 46, № 3. C. 103-115.

[21] Фомин П.А., Чен Дж.-Р. Влияние химически инертных частиц на параметры и подавление детонации в газах // Там же. 2009. Т. 45, № 3. C. 77-88.

[22] Фёдоров А.В., Тропин Д.А. Моделирование прохождения детонационной волны через облако частиц в двухскоростной, двухтемпературной постановке // Физика горения и взрыва. 2013. Т. 49, № 2. С. 61-70.

[23] Wolanski P., Liu J.C., Kaufman C.W. et al. The effect of inert particles on methan-air detonations // Archivum Combust. 1988. Vol. 8, No. 1. P. 15-32.

[24] Papalexandris M.V. Numerical simulation of detonations in mixtures of gases and solid particles // J. of Fluid Mech. 2004. Vol. 507. P. 95-142.

Поступила в 'редакцию 1 апреля 2013 г., с доработки — 20 мая 2013 г.