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

УЧЕНЫЕ ЗАПИСКИ И А Г И Т о м VII 197 6

№ 6

УДК 533.6.011.55

РАСЧЕТ ТЕЧЕНИЯ В СОПЛЕ ЛАВАЛЯ МЕТОДОМ УСТАНОВЛЕНИЯ ВТОРОГО ПОРЯДКА ТОЧНОСТИ

А. П. Мазуров

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

Применение численных методов для решения задачи о смешанном течении идеального газа в двумерных плоских и осесимметричных соплах позволило получить ряд важных результатов в этой области. Наибольший интерес для практических целей представляют расчеты течений в осесимметричных соплах Лаваля, имеющих большие углы наклона образующей дозвуковой части и малые радиусы кривизны горла. Для построения контуров подобных сопл в работе 11] использовался численный метод, основанный на решении обратной задачи, в которой форма сопла определяется по заданному распределению скорости вдоль оси симметрии. Для сквозных расчетов полей течения в соплах произвольной заданной формы применяются методы установления, основанные на численном интегрировании нестационарных уравнений газовой динамики. Результаты ранее выполненных работ показывают, что методы установления позволяют с достаточной точностью рассчитывать поля течений [2, 3] и коэффициенты расхода осесимметричных сопл [4].

В настоящей статье для расчета полей течения и коэффициентов расхода осесимметричных сопл применяется метод установления второго порядка точности. Основные уравнения газовой динамики записываются в преобразованных пространственных переменных и решаются численно с использованием конечно-разностной схемы Маккормака [5]. В применяемой схеме пространственные производные аппроксимируются разностями вперед в предикторе и разностями назад в корректоре. Этим достигается значительное

упрощение алгоритма расчета по сравнению со схемами, использующими центральные разности, такими, как, например, схема Лакса — Вендроффа [4], что в свою очередь сокращает требуемое время счета до выхода на стационарный режим. Для примера приведены результаты расчета осесимметричного сопла Лаваля с радиусом кривизны горла (отнесенным к радиусу горла) /?с = 0,625, углом сужения дозвуковой части 0* = 45° и углом коничности сверхзвуковой части 6Д=15°. Значения коэффициентов расхода приводятся для сопл с относительными радиусами кривизны горла, изменяющимися в диапазоне /?с = 0,35ч-1,0.

Численные расчеты выполнены на ЭВМ БЭСМ-6 по программе, составленной на языке ФОРТРАН.

1. Газодинамические уравнения осесимметричного течения невязкого нетеплопроводного газа приведем к безразмерному виду и запишем в цилиндрической системе координат х, г в дивергентной форме:

dW . dF(W) . dG(W) rr dt + dx + dr — ’

(1)

где W, F(W), G(W), H — четырехмерные векторы-столбцы, определяемые формулами

W =

Р Г

prU

Р rV

er

, F{W) =

prU (р£Л + p) r prUV (e + p) Ur

, G(W) =

Р rV 0

prUV , и= 0

(Р V2 + p)r р

(е + р) Vr 0

Векторное уравнение (1) соответствует уравнению неразрывности, двум уравнениям сохранения количества движения и уравнению сохранения энергии, в которых компоненты вектора скорости и по оси х и V по оси г отнесены к критической скорости а*, плотность р отнесена к плотности торможения р0, давление р и полная энергия единицы объема е отнесены к произведению р0а2. Координаты х, г отнесены к некоторому характерному размеру сопла Ь, время Ь отнесено к Ь\а.^.

Система (1) замыкается уравнением состояния совершенного

газа

Р = (*— ще

р pV2'

где х — отношение удельных теплоемкостей.

Для того чтобы упростить расчеты параметров потока в ячейках, прилегающих к стенкам сопла, запишем векторное уравнение (1) в новой системе координат £■»), в которой уравнение стенки сопла имеет простой вид ■») = const. С этой целью рассмотрим в общем случае кольцевое сопло, которое показано на фиг. 1, а. Пусть г = г+(х) — уравнение внешней стенки, а г = г_(х) — уравнение внутренней стенки сопла. Тогда преобразование

Ах)

г.у (х) — (X)

переводит контур сопла в линии t\—=0 и 1 (фиг. 1,6).

Соотношение между частными производными в старой и новой системах координат определяется формулами

_д____д_

дх — д;

в 4.

(3)

А = —

1

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

В системе координат Ъ\ уравнение (1) принимает вид

6W

dt

, AdF{W)

a? + di1

+ вдат=н

1 dri

(4)

2. В плоскости l-q введем прямоугольную сетку с размерами ячеек Д£ по оси \ и Ат] по оси т], шаг по времени обозначим М. Пусть W" . обозначает функцию W (tn, £г, т)у.) на временном слое tn = пМ в узловой точке с координатами = y];. =

=уДг](г = 0, 1,...,/А; j = О, 1,..., JW). Аналогично обозначим F?f = F(W?/), G" j = 0{ W’l j), HI i = H (Wl i), Ai,= = Л ft, -/,,), 5, = В (У.

Для численного интегрирования уравнения (4) используем конечно-разностную схему второго порядка точности:

Фиг. 1

№?+/=* Г",

у) + ц- [А,., Пн 1 - К. /) +

+ ИР?.}1 -

■n+l

А5?-!1. /) +

+Д/-Я-+1

(5а)

(56)

где

ПУ = Г( ГП1), О"у = 0( Г?,У), Я"у - Я (и^Г/ )■

В приведенной разностной схеме переход с временного слоя /!" = «• М на слой Ьп+1 — (п+ 1)-Д£ осуществляется в два этапа: сначала по формуле (5 а) определяются величины №"-1,/ и

Затем, подставляя эти величины в формулу (56), находим значение функции = И7(^п+1, ?г, тг)у-).

Для расчета ноля течения в сопле начальные значения параметров и, р и р задаются из одномерного приближения, а компонент скорости V определяется с помощью линейной интерполяции:

У=и[г1 + (г+ — г1)тч] (0<7]<1).

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

№„ = ?>(№!— XV,) + 1^з с учетом условия непротекания.

В выходном сечении сопла (г = /Л) параметры течения определяются по линейной экстраполяции:

ИР/л = 2 \Р/А-і - Гм_2.

Схема расположения узловых точек, завязанных в определении параметров на границах, показана на фиг. 1,6.

Во входном сечении (г = 0) сопла параметры течения определяются из того условия, что расход газа т"+х в этом сечении

в момент времени ?п+1 равен расходу газа от” в горле сопла в момент Ґ. При этом, так как мы предполагаем, что сопло имеет цилиндрический входной участок, и, следовательно, на достаточном удалении от горла сопла статическое давление постоянно по сечению, все параметры потока на входе находятся из уравнения Бернулли и условия то+1 = /га".

Стационарное течение получается интегрированием нестационарного уравнения при неограниченном увеличении времени £->оо.

3. Изложенный метод был применен для расчетов течений идеального газа с отношением удельных теплоемкостей х=1,4 в осесимметричных соплах Лаваля. За характерный размер сопла принимался радиус горла г*.

На фиг. 2 и 3 представлены результаты расчета сопла с относительным радиусом кривизны горла /?с = 0,625, углом сужения дозвуковой части 6*= 45° и углом коничности сверхзвуковой части = 15°. Область течения разбивалась на 61X21 узловых точек. Считалось, что стационарное состояние достигалось тогда, когда изменение массового расхода в горле сопла на каждом шаге по времени становилось меньше Д^2^3- 10~г>. При этом требовалось выполнение приблизительно 500 шагов по времени, что занимало 1 ч 20 мин машинного времени.

----часленныи расчет

----эксперимент [б\

0 100 200 300 400 500 л

Фиг. 3

На фиг. 2 построено распределение чисел М внутри сопла. Здесь же черными кружками нанесены результаты экспериментального исследования [6]. Можно отметить, что проведенный расчет достаточно хорошо согласуется с экспериментом по определению картины течения в сопле. Пунктирными линиями нанесены кривые М = const, расчет которых выполнен по методу [3]. На фиг. 3 показано изменение по времени чисел М на стенке и на оси сопла в плоскости горла. Стационарные значения чисел М составляют М= 1,3925 на стенке и М = 0,8116 на оси. Соответствующие значения, полученные в работе [6], составляют М —1,4 и 0,8. Приблизительно такие же значения чисел М были получены в расчете [4], а именно: М= 1,3935 на стенке и М = 0,8024 на оси. Некоторое расхождение с экспериментальным значением числа М на оси, получаемого в приведенном расчете, по-видимому, связано с тем, что ось симметрии заменяется некоторым тонким цилиндром.

Одной из основных интегральных характеристик реактивных сопл, по которой оценивается качество дозвуковой части сопла, является коэффициент расхода, определяемый отношением

а 99

0 96

0,97

0.96

• метод авп \ м і [7] [S] * Ы юра

\ ... s' s'

т / У / / т •-

~ / - / / /

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

Асимптотическое значение коэффициента расхода сопла с радиусом кривизны горла ^с = 0,625, полученное в данном расчете, составляет ц = = 0,9832. Это значение практически совпадает с экспериментальным определением [а = 0,983+0,008 [6]. Изменение коэффициента расхода по времени в процессе установления проиллюстрировано фиг. 3.

С целью сопоставления результатов, полученных по данному методу и по другим методам, были проведены расчеты коэффициента расхода осесимметричных сопл, имеющих угол сужающейся части 9* = 30° и угол расходящегося участка ва—15°. Радиус кривизны горла имел значения /?с = 0,35; 0,55; 0,75 и 1,00. Во всех вариантах отношение площади входного сечения сопла к площади горла равнялось Гф/Т:,: = 4. Полученные результаты, обозначенные черными кружками, и результаты других авторов [6—9] представлены на фиг. 4. Сопоставляя эти данные, можно отметить, что численный расчет достаточно хорошо согласуется с расчетом методом асимптотических разложений [8] и вполне удовлетворительно с экспериментальными данными [7]. Расчеты показывают, что метод установления второго порядка точности с использованием конечноразностной схемы Маккормака хорошо подходит для расчетов как полей течения, так и коэффициентов расхода осесимметричных сопл.

0,2 0,¥ О,С 0,8 1,0 /гс

Фиг. 4

ЛИТЕРАТУРА

1. Пирумов У. Г. Расчет течения в сопле Лаваля. Изв. АН ГХСР, МЖГ, 1967, № 2.

2. И в а и о в М. Я., К р а й к о А. И. Численное решение прямой задачи о смешанном течении в соплах. Изв. АН СССР, МЖГ, I969, №5.

3. Киреев В. И„ Лифшиц Ю. Б., Михайлов Ю. А. О решении прямой задачи сопла Лаваля. „Ученые записки ЦАГИ*, т. I, № I, 1970.

4. Л а в а л ь П. Нестационарный метод расчета трансзвуковых течений в соплах. В сб. „Численные методы в механике жидкостей*. М., „Мир*, 1973.

5. Mac Cormack R. W. The effect of viscocity in hypervelocity impact cratering. AIAA Paper N 69—354.

6. Back L. H., Cuff el R. E. Flow coefficients for supersonic nozzles with comparatevely small radius of curvature throat. Journal of Spacecraft and Rockets, vol. 8, N 2, 1971.

7. Norton D. J., Shelton S. Performance of rocket nozzles with low radius of curvature radius. JPL Space Programs Summary, vol. Ill, 1969, Pasadena, Calif.

8. К liege l T. D., Levine J. N. Transonic flow in small throat radius of curvature nozzles. AIAA J., vol. 7, N 7, 1969.

9. Соколов В. Д., Яг удин С. В. Коэффициент расхода осесимметричных сужающихся сопл с произвольным контуром. „Ученые записки ЦАГИ*, т. 6, № I, 1975.

Рукопись поступила 23jl 1976 г.