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

УДК 533.6.011.8

В. В. Кузенов

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

Рассмотрена задача, связанная с численным моделированием взаимодействия потока, возникающего вблизи движущегося в газовой среде реактивного снаряда на твердом топливе, и продуктов горения твердого топлива. Особенностью предлагаемого подхода является описание газодинамических процессов внутри реактивного снаряда на твердом топливе и в спутной струе на основе единой расчетной методики. Численное исследование газодинамических течений выполнено с применением нестационарных уравнений Навье-Стокса.

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

В настоящей работе приведена расчетно-теоретическая модель моноблока, позволяющая рассчитывать характеристики газодинамических процессов в РСТТ, а также проводить расчеты взаимодействия струи продуктов горения с окружающей газовой средой на основе единой численной методики. Одной из целей работы является расчетно-теоретическое изучение особенностей структуры и пространственных распределений газодинамических параметров потока внутри сопла и при его взаимодействии со спутным потоком газа вне соплового аппарата. При этом расчетные исследования проводятся в рамках вязкого течения газа [1-4].

Постановка задачи. Численное исследование газодинамических течений с учетом различного рода особенностей, которые возникают в газовых потоках как внутри РСТТ, так и при истечении газовых струй из них в спутный газовый поток, может быть выполнено на основе решения нестационарных трехмерных уравнений Навье-Стокса.

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

Введем преобразование координат вида

г = г(Ь п С); г = п С); р = п, С).

Нестационарные уравнения Навье-Стокса в векторной полудивергентной форме в произвольной криволинейной системе координат £, п, С имеют следующий вид:

дЩ д£ dG dQ dt дп д(

dFn да

дОг ' дп д(

+

+ Sv; Fv — Fi + Sv.

При этом вектор зависимых переменных задачи и, векторы потоков Б, С, (( и вектор Б в криволинейной системе координат £, п, С выражаются следующим образом:

и = 3 [р, ри, рю, рь, ре]т ; Г = 3 [рУ^, риУ^, рюУ%,рьЩ, реУ^]т ;

С = 3 [рУп, риУп, рюУп, рьУп, реУп]т ; (( = 3 [рУ^, риУ^, рюУ^, рьУ^, реУс];

S — J

а-

ри r

р(и2 — w2) др др др

а-—-+ Щ + п дп + Z дС

puv др др др

а— ++ п дП + с Щ а 2puw + Cfi др + др + i др.

r ra ra дп ra дС ареи + р [д (Щ) + д (JVn) + д (JV()\ +

r J\ д6 дп д( )

ри pw2 puw

+ а--+ au--aw-

r r r

Здесь £ — время, р — плотность газа, р — давление, е — внутренняя энергия, 7 = д (г, г, <^)/д (£,п, С)- якобиан преобразования координат,

£ п С

V = и + ^V + -тШ, Кп = Пги + ПгV + -7-Ш, КС = (ги + С*V + Ш —

^ г^ га га

контравариантные компоненты вектора скорости в новой криволинейной системе координат, Уг = и, Кг = V, = ш — проекции вектора скорости на оси г, г и ^ цилиндрической системы координат г, г, а = 0 — соответствует плоскому, а а = 1 — осесимметричному случаю течения газовой смеси.

Компоненты тензора вязких напряжений а^ и вектора в системе координат г, г, для сжимаемого газа имеют вид:

TTrr — P

T p p р

Tr p — р

du 2 2 — - - div V ör 3

^ 1 öw ^ u 2 2—— + 2a— - - divV ra ö^ r 3

öw 1 du w

—I--ö--a—

ör ra ö^ r

Tzz — P

Trz — P

Tzp — P

öz - 2 div V 3

du dz öv

dr ;

dw dz 1 öv

r» ö^ ;

^ du dv 1 öw u

dr dz ra ö^ r'

0,a

Trr - T

7 rz

T

a—, a r

Tr p + Tpr a а г dT+n dTöT r ' r \ r ö? r ön r öC

где р — молекулярная вязкость.

Преобразуем тензор вязких напряжений а^- с учетом криволиней-ности системы координат £, п, С:

Trr ^^

; (^+^+с^ I - 2 i% dw +% dw + % dw'i -

4 ( du du 3 ö? + nr ön

du

dC

3 V ra ra dn ra öC

2 / öv öv Л 2u

- 316 de +dn +Сг dcj- a 37.

т

Tzz — P

4 / öv öv öv\ 2 / du du öu\ 3 V de + nz ön +Cz de; - 3 V de + nr dn +Cr öcy -

ёр öw + n p dw

dn

+ Cpd^

ra öC

a

2u 3r

Tp p p

ёр öw + npdw + Cpdww^ + a 4u_

ra ö? ra ön ra öC / a 3r

r

r

2 / ди ди ди\ 2 / дv дv dv^

—з д6 + пг дП + Crдс)—з д6 + nzдП + Ж.

Orz — V Orv р

ди ди ди\ / дv дv дv

д6 дп ди V д4 дп дС

dvj dvj dvj\

+ пг dn + с öCJ +

/6v^ +^^ + Cvd^ _ aw \ ra д6 ra дп ra д^ r

Ozv — V

^+п^+с. + С % öv+* ?+% dv

dw dvj

д4 + п-~öq ' d^ ' Vra dQ ' ra дп ' ra dC

Запишем компоненты векторов уу ,Су,((у, связанные с симметричным тензором вязких напряжений а^, в криволинейной системе координат £, п, С:

Fy —

G у —

Q v —

Fy (1 Fv (2

Fy (3 Fy (4

Gy (1

Gv (2

Gy (3 Gy (4

Qy (1

Qy(2 Qy (3 Qy (4

—J

0

6 i 6 |

6r Ozr + 6z Ozz + ~~ O

6 i 6 | ^^

Sr Ovr + 4z Ovz + a Ovv

—J

пv

пг Orr + ^ Orz + ra O Пг Ozr + nz O zz + ra O.

rv

+ + ^

nrO vr + nzO vz + yOtO VV

—J

CrO rr + CzO rz + aOrV

Cr Ozr + CzO zz + ~ZO

C + C + ^

CrO vr + CzO vz + aOvv

Укажем компоненты векторов Г, (С^, ((^ и уу = ^^ + б^, обусловленные вектором теплового потока I = —Лgrаd (Т) + аУ (здесь Л — коэффициент молекулярной теплопроводности; Т — температура газовой смеси), в криволинейной системе координат £, п, С:

0

0

Fv (5) — AJ

?r2 + Ja + + ( nr + + nz )Tn +

p

+ UrCr + ^ + Cz )TZ

Gv (5) — AJ

(nr + npl + nz б) T + (n2 + + nz2) Tn+

+ (nr Cr + npCp + nz Cz )Tf

Qv (5) — AJ

Cr + %0p + Cz IT, + ( Cr nr + + Cz nz )Tn +

Cpn p

+ (Cr2 + + C2 )TC

S2(1)—0, S2(2) —0, S2(3) —0, S2 (4)— 0;

52(5) = -7{и$}(2) + ^(3) + ^5,(4)} +

^ , ч ди ^ . . д-у ^ . . дш ^ . . ди ^ . . д-у + Я. (2) - + я,(3) - + я, (4) - + О.(2) щ + «V (3) ^+

+ О (4) ^ + «V (2) ди + «V (3) | + (4)

При известных в физическом пространстве координатах узлов сетки в расчетной области £, п, С метрические коэффициенты в общем случае могут быть найдены путем численного дифференцирования по формулам:

т . х/о/А ör öz öw örözöw ör öz öw

J — ö (r,z, w)/ö (?,n,C) — -önöC + öCö?ön + öndCö?

ör öz öw ör öz öw ör öz öw

ö? öC ön ön ö? öC öC ön ö? '

/öz öw öz öw

r \ön öC öC ön

— J

-1

öw ör

ön öC

ör öw

ön öC

T, . öw öz öz öw Пг — J 77 77 - 7777

ö? öC ö? öC T-1 f ör öw ör öw

nz — V öC - öC ö^

-1 ör öz öz ör -1 öz ör ör öz

— J VdnöC - ö^do' ; np — J Vö^dC - ö^ÖÖ';

Cr — J

1

öz öw öw öz ö? ön ö? ön

Cz — J

1

ör öw ör öw

ön ö? ö? ön

ör öz öz ör \ \ö?ön ö? ön /

Cp — J

Для численного решения, как правило, используется безразмерная форма уравнений Навье-Стокса. Данная форма уравнений позволяет выделить систему безразмерных комплексов: число Маха — М, число Рейнольдса — Яе, число Прандтля — Рг. При этом все переменные, присутствующие в уравнениях Навье-Стокса, нормируются на характерные величины и изменяются в пределах от 0 до 1.

При переходе к безразмерным уравнениям в качестве масштабов времени, длины, скорости, температуры, плотности, давления и внутренней энергии выбираются некоторые характерные значения этих величин:

— р

V' ^

К,

V2

T* — VT,

с*

р*р* — p*v ,

— V2.

В качестве масштабов коэффициентов переноса и физических свойств газа используются соответствующие значения, определенные при характерных величинах температуры Т* и давления р*:

P*

A*,

где с*, ср — характерные значения теплоемкости при постоянном объеме и давлении; р* — характерное значение коэффициента динамической вязкости; Л* — характерное значение коэффициента теплопроводности. Числа Рейнольдса и Прандтля определяются соотношениями:

Re —

P*

Рг — с*—

i i ср .

A

Далее переходят к безразмерным переменным, используя следующие связи размерных и безразмерных переменных:

ш

ш —

ш*

где с — безразмерная величина, соответствующая размерной форме исходной переменной с.

t

e

с

с

Тогда в безразмерных переменных нестационарные уравнения Навье-Стокса в произвольной криволинейной системе координат п, С примут вид

dU dF dG dQ * 1 дЬ д* дп ß( Re

ÖF dGv dQ v *

~ I T7Z I ~ Г Oy

д* дп d(

где безразмерные компоненты ]?»5,<О5,(?»,£1;5 безразмерных векторов

]?»,С?»,((» и Б» = + Б1, обусловленные вектором теплового потока, находятся с помощью соотношений:

5(< 5(< 5(< 5(<

= уг(5), О» = -Р;О»(5), ((5 = ^(V(5), ?1;5 = ^

с*

здесь 7* = — характерное значение показателя адиабаты; компо-

сч_

ненты ]?»(5), С?»(5), (?»(5), ¿'»(5) имеют вид исходных размерных компонентов ]» (5), О» (5), (»(5), (5), но записаны относительно безразмерных переменных с?.

Граничные условия, определяющие характеристики продуктов сгорания, поступающих в камеру двигателя, следующие:

Р КГа

иа = ит (рт/рз) n, ра = е8 = _а1).

Здесь и8, р8, Т8, е8, Я — соответственно скорость, плотность, температура, внутренняя энергия и газовая постоянная поступающих газообразных продуктов сгорания с поверхности твердого топлива; рт — плотность твердого топлива; ит — скорость перемещения поверхности топлива при его выгорании. Предполагается, что вдув газа в камеру двигателя происходит по нормали к поверхности топлива. Величина скорости выгорания ит известна из эксперимента и задается в виде зависимости [5]

ит = иТо р»,

где иТо — величина скорости перемещения границы выгорающего топлива при давлении р = 1, V — константа, зависящая от типа используемого топлива (0 < V < 1).

На границе расчетной области, совпадающей с твердой поверхностью, для динамических переменных, присутствующих в уравнении Навье-Стокса, ставятся условия прилипания, а также условие постоянства температуры на обтекаемой газом поверхности. На границе расчетной области, через которую вдувается спутная струя, задаются следующие газодинамические параметры: 7 = 1,4, рю, рж, W(Xl. Здесь рю, рж — соответственно скорость, плотность, давление поступающего через граничную поверхность газа. На внешних (находящихся

на достаточно больших расстояниях от среза сопла) поверхностях рас-

д/

четной области ставятся условия свободного истечения - = 0, где

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

Приведенная система уравнений дополняется начальными условиями

и (г, 0) = 0, р (г, 0) = ро, е (г, 0) = е0.

Численная реализация разработанной модели основана на многоблочной многосеточной технологии расчетов на неортогональных структурированных сетках с использованием схем расщепления по физическим процессам и направлениям. Система уравнений Навье-Стокса решалась на основе явного метода Годунова [6, 7]. При этом "вязкая" часть системы уравнений Навье-Стокса определялась также явным образом. При аппроксимации конвективной составляющей векторов потоков на границах расчетной ячейки применялась процедура расчета распада разрыва, разработанная С.К. Годуновым, с использованием реконструкции сеточной функции (внутри расчетной ячейки), которая была взята из работ [8, 9]. Более подробная информация об используемой в работе численной методики приведена в работе [10].

Результаты численного моделирования. При численном моделировании газодинамических процессов в рассматриваемой задаче расчетная область представлялась на основе двух сопряженных между собой областей:

первая область состоит из внутренней части конструкции РСТТ и прямоугольной зоны, являющейся продолжением предыдущей зоны, располагающейся за соплом (вниз по течению) и имеющей поперечные размеры, равные диаметру среза сопла;

вторая область представляет собой многоугольник, нижняя граница которого сначала проходит по оси симметрии моноблока, далее по его наружной поверхности и оканчивается границей сопряжения с первой областью. Этот многоугольник слева ограничен вертикальной плоскостью, через которую внешний поток газа с параметрами Тте, 7 =1,4 втекает в расчетную область.

В затопленном пространстве в качестве начальных условий в спут-ном потоке и краевых условий на границе расчетной области, в зоне вдува спутного потока, использовались параметры, соответствующие следующим вариантам расчетов:

— = 72000 —, рте = 0,036 атм, Тте = 270 К, 7 =1,4

с

(высота Н = 25 км);

— woo = 147300 —, рте = 0,00065 атм, Тте = 270 К, 7 =1,4

с

(высота Н = 55 км).

Рис. 1. Схематичное изображение реактивного снаряда на твердом топливе

Процессы горения твердого топлива описывались следующим образом: с боковой поверхности камеры сгорания с момента времени £ = 0 и до окончания расчета в камеру сгорания поступает газ со скоростью 500 см/с, температурой Т = 800 К и давлением р = 5 атм.

Схематичное изображение и основные размеры РСТТ приведены на рис. 1. При этом моноблок состоит из цилиндра длиной 1060 см, диаметром, равным диаметру среза сопла Дмоноб = Дсопл = 215 см, и конической головной части с углом раскрытия 54°.

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

Первая стадия — формирование стационарного газового потока внутри РСТТ (время выхода на режим £ « 3 • 10-2 с). Вторая стадия — это рабочий режим двигателя.

Момент времени (£ « 0,05 с) выхода газодинамических параметров внутри РСТТ на установившийся режим характеризуется следующими величинами основных газодинамических параметров: в камере сгорания РСТТ число Маха изменяется в диапазоне А М € [2 • 10-3; 2 • 10-1], продольная компонента скорости находится в диапазоне значений Уг € [10 ... 250] м/с, температура и давление принимают характерные значения Т « 800 К, р « 5 атм.

Рассмотрим процессы взаимодействия внешнего газового потока со сверхзвуковой выхлопной струей продуктов горения твердого топлива (ТТ), вытекающей из соплового аппарата РСТТ с конической формой расширяющейся части.

Проведенное в работе моделирование взаимодействия струи продуктов горения ТТ со спутным потоком газа вне соплового аппарата на основе разработанной численной методики показало, что вариант расчета, соответствующий высоте к = 25 км (степень нерасчетности , Ра

— п ~ 1, п = —, где ра — давление на срезе сопла, рх — давление во

р^

внешнем потоке газа), является режимом течения с перерасширением, а вариант расчета, соответствующий к = 55 км (степень нерасчетно-сти — п ~ 80), относится к течению с недорасширением.

Из результатов работ [10, 11] известно, что в зависимости от степени нерасчетности п можно выделить некоторые характерные режимы течения на начальном участке смешивающихся со спутным потоком струй.

При сильном недорасширении (вариант к = 55 км, п ~ 80) реализуется картина течения, подробно описанная в работах [1-3, 10-12] и представленная на рис.2,а и 3,а. Данные рисунки демонстрируют волновую структуру сильно недорасширенной струи, истекающей в спутный поток воздуха. Из приведенных распределений видно, что при достаточно больших значениях степени нерасчетности п вблизи выходной кромки сопла и в спутной струе окружающего воздуха из-за столкновения расширяющейся струи выхлопных газов РСТТ и спут-ного потока возникают косая ударная волна (УВ) и падающая на ось струи висячая УВ, которая характеризуется регулярным отражением

31 8.27Е+02

29 7 53Е+02

27 6 Е9Е+02

25 5.21Е»02

23 5 52Е»02

21 4 83Е+С2

1S 4 I4E+C2

17 3 45E+Ü2

15 i 34E+Q2

13 liZL'02

11 2 0&Е+02

9 I Jlc»02

7 I 55Е+02

5 1 2SE+C2

г 1 05Е+02

1 3 60Е+01

о юоо 2000 г

б

Рис. 2. Пространственное распределение температуры Т в продуктах горения и спутной струе воздуха (К = 55 км, Ма « 4, = 270 ф:

а — п = 80; = 4,4; б — п = 1; = 2

от оси струи (осевая координата области падения х ~ 1700 см — коническая форма сопла). Вследствие того, что размеры первой бочки растут, как л/п, наблюдается только первая бочка, размеры которой существенно превышают характерный поперечный размер РСТТ. Напомним, что термин "бочка" используется в работах [1-3, 11, 12] для описания периодической волновой структуры начального участка смешения струй из-за того, что на начальном участке граница смешения струй имеет характерную бочкообразную форму с выпуклостью, обращенной в сторону внешнего потока.

На рис. 2, б и 3, б приведены пространственные распределения температуры и числа Маха, возникающих вблизи реактивного снаряда на твердом топливе и соответствующих высоте полета Н = 25 км. В этом случае реализуется картина течения, соответствующая степени нерас-четности, близкой к единице (Н = 25 км, п ~ 1). Здесь в силу достаточно большой величины давления в спутной струе (рж ~ ра) радиальное расширение центральной струи выхлопных газов сильно ограничено. Это ограничение приводит к падению на ось струи висячего скачка уплотнения (осевая координата области падения х ~ 1700 см — коническая форма сопла) с последующим регулярным отражением от нее, характерный поперечный размер центральной струи близок к диаметру среза сопла, распределение температуры вдоль оси струи изменяется в следующем диапазоне значений: Т « 130 ... 200 К. В этом

Рис. 3. Пространственное распределение числа Маха М в продуктах горения и спутной струе воздуха (Ма « 4, = 270 Ю:

h = 55 км; n и 80; M„

4,4; б — h = 25 км; n и 1; M„

2

а

случае при взаимодействии спутного потока и центральной струи количество бочек увеличивается (если сравнивать с результатами варианта Н = 55 км) и становится больше 1.

Работа выполнена в рамках Программы фундаментальных исследований РАН № 1 (координаторы программы академик Г.Г. Черный и академик Д.М. Климов) и при поддержке гранта РФФИ № 05-01-00780.

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

¡.Сверхзвуковые струи идеального газа. Ч. 1, 2. Г.И.Аверенкова, Э.А. Ашратов, Т.Г. Волковская и др. - М.: Изд-во МГУ им. М.В. Ломоносова, 1970.

2. К о в а л е в Б. Д., Мышенков В. И. Расчет вязкой сверхзвуковой струи, истекающей в спутный поток // Уч. зап. ЦАГИ. - 1978. - Т. 9. - № 3. - С. 125-130.

3. Мышенков В. И. Расчет течения вязкой ламинарной сверхзвуковой струи в спутный поток // ЖВМ и МФ. - 1979. - Т. 19. - № 2. - С. 474-485.

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

5. С о р к и н Р. Е. Теория внутрикамерных процессов в ракетных системах на твердом топливе. - М.: Наука, 1983.

6. Численное решение многомерных задач газовой динамики / Под ред. С.К.Годунова. - М.: Наука, 1976.

7. К о л г а н В. П. Применение принципа минимальных значений производной к построению конечно-разностных схем для расчета разрывных решений газовой динамики // Уч. зап. ЦАГИ. - 1972. - T. 3. - № 6. - 376 с.

8. Куликовский А. Г., Погорелов Н. В., Семенов А. Ю. Математические вопросы численного решения гиперболических систем уравнений. - М.: Физматлит, 2001.

9. И в а н о в М. Я., Крупа В. Г., Н и г м а т у л л и н Р. Э. Неявная схема С.К. Годунова повышенной точности для интегрирования уравнений Навье-Стокса // ВМиМФ. - 1989. - Т. 29. - № 6. - С. 888.

10. Кузенов В. В., Суржиков С. Т. Численное моделирование внутренней баллистики РДТТ и структуры спутных струй. - М.: ИПМех РАН. Препринт № 780, 2005.

11. А в д у е в с к и й В. С., Ашратов Э. А., Иванов А. В., П и р у м о в У. Г. Сверхзвуковые неизобарические струи газа. - М.: Машиностроение, 1985.

12. П и р у м о в У Г., Росляков Г. С. Газовая динамика сопел. - М.: Высшая школа, 1987.

Статья поступила в редакцию 27.06.2006

Виктор Витальевич Кузенов окончил в 1983 году Московский государственный университет им. М.В. Ломоносова, старший научный сотрудник лаборатории "Радиационная газовая динамика" ИПМех РАН, канд. техн. наук, доцент кафедры "Теплофизика" МГТУ им. Н.Э. Баумана. Автор более 50 научных работ в области теплофизики, радиационной газо- и плазмодинамики.

V.V. Kuzenov graduated from the Lomonosov Moscow State University in 1983. Ph. D. (Eng.), assoc. professor of "Thermal Physics" department of the Bauman Moscow State Technical University. Senior researcher of "Radiation Gas-dynamics" laboratory of the Institute for Applied Mechanics of the Russian Academy of Sciences. Author of more than 50 publications in the field of thermal physics, radiation gas and plasma dynamics.