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

УЧЕНЫЕ ЗАПИСКИ Ц А Г И

Т о м IX

197 8

М 3

УДК 629.7.015.3.036:533.697.4

ОБ УВЕЛИЧЕНИИ ТОЧНОСТИ ОПРЕДЕЛЕНИЯ ИНТЕГРАЛЬНЫХ ХАРАКТЕРИСТИК СОПЛ НА ОСНОВАНИИ ЧИСЛЕННЫХ РАСЧЕТОВ ПОЛЯ ТЕЧЕНИЯ

В. Л. Зимонт, С. В. ^гудин

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

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

где С — расход, 1 — проекция потока импульса на ось сопла, р, р, .и — давление, плотность и продольная составляющая скорости, индекс „с“ относится к одномерным параметрам в рассматриваемом сечении с площадью Р.

К точности определения интегральных характеристик потока в соплах реактивных двигателей предъявляются высокие требования (—0,2%). При численном решении прямой задачи о смешанном течении невязкого и нетеплопроводного газа в сверхзвуковых соплах методами установления, изложенными, например, в работах |.6 — 8], ввиду погрешностей аппроксимаций и многоцикличных итераций, параметры потока определяются с некоторой степенью точности, что приводит, как правило, к недостаточной точности определения интегральных характеристик. Более того, анализ конкретных численных расчетов приводит к йарадоксальному на первый взгляд выводу о том, что точность определения интегральных

0)

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

Такое положение связано с тем, что причина погрешностей формул одномерной теории и конечно-разностных расчетов носит различный характер. В первом случае интегральные законы сохранения выполняются точно, но используются схематизированные представления о распределении параметров в сечениях, тогда как во втором случае газодинамическая картина течения определяется относительно правильно, но, во-первых, точность выполнения интегральных законов сохранения массы, количества движения в конце счета зависит от выбранных критериев выхода решения на стационарный режим и часто из-за экономии машинного времени, необходимого для расчета одного варианта, оказывается недостаточно высокой (—1%). Во-вторых, что более важно, нарушается условие постоянства энтропии в непрерывных течениях (без скачков уплотнения). При этом для плавных сопл указанные ошибки оказываются более существенными, чем одномерная схематизация картины течения.

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

По-видимому, наиболее широко используемой работой такого рода является статья [9], в которой решение разлагается в ряд по степеням 1//? (/? — отношение радиуса кривизны контура в горле сопла к радиусу критического сечения сопла /?*). Полученные в [9] результаты применимы при /?>1,5. Несмотря на то, что расширить диапазон применимости аналитических методов для случая ^<1, можно, как показано в работе [5], путем выбора тороидальной системы координат и использования разложений в ряд по степеням 1/(/?+1), такие решения непригодны для сопл с более сложными профилями. Необходимые аналитические результаты для сверхзвуковой области вообще отсутствуют.

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

р=7?- = 1-£*«, 7=-^=1-£л„, (2)

ис I 'с* ;

где К*1, — некоторые величины, много меньшие единицы, чис-

ленные значения которых определяются формой распределения

параметров в критическом сечении (безразмерными профилями, отнесенными к средним по сечению значениям параметров). Здесь и далее * показывает принадлежность к критическому сечению, индекс I относится к группе членов одинакового порядка малости.

Приведенные ниже результаты численных расчетов показывают, что в практически интересных случаях даже при больших неравномерностях распределений параметров достаточно ограничиться лишь главными членами поправок К*\ и Л*!, имеющими вид:

•^*1 = А* + ~ В* + — Е* + ~ С*>

2х! Л* + 4- Ещ + —Ч- С* + -^гВ,

где

(3)

ц2 цЗ р2

Члены вида а, Да, АаАЬ распределений а(г), Ь(г) вычисляются следующим образом:

Я2

R*

f*

а—----- 1 а (г) г dr, ba(r) — a(r)—a, АаАЬ

Rl

Д а (г) Д b (г) г dr.

о

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

Если профили распределения параметров точно удовлетворяют интегральным законам сохранения, значения ц и I*, найденные, согласно выражениям (1) и (2), по параметрам в критическом сечении, совпадают. Однако, если профили параметров лишь приближенно удовлетворяют законам сохранения, соотношения (2) дают значительно меньшую погрешность, поскольку ошибка при определении профилей скажется лишь на величины малых по сравнению с единицей значений К*/ и А.,, в то время как при интегрировании, согласно (1), эта ошибка непосредственно влияет на значения р.

и 7*.

Для определения погрешности, вносимой отбрасыванием в (2) поправок более высокого порядка, использовались аналитические решения для трансзвуковой области, учитывающие три члена разложения по степеням 1/(1+/?) [5]. Сравнивались значения коэффициентов расхода и относительных импульсов, найденных на основе этих решений ц.*!, /*!, согласно выражениям (1), и соглас-

но (2), с учетом лишь членов /С%1 и А*!. При этом рассчитывались также значения ЛГ*г и А^2 (выражения для них приведены в [1]). Отметим, что характеристика неравномерности У приведена в _[1] с опечаткой в первом слагаемом, которое следует читать р'гЛ/ри2. Результаты расчетов для /? = 2 следующие:

[л*! =0,996205, р.^.2 = 0,996190, при этом К%2/К*1 =—0,0068;

Ли =0,997675, /*2 = 0,997677, А*2/Л*, =- 0,0043.

При /? = 0,625 значения р.*! = 0,981225, (а*2 = 0,980921, АГ^2//С*ь =

= — 0,0230, 7*! = 0,987855, /*, = 0,987810, &*г/Ач =—0,0131.

Из соотношений (2) и (3) непосредственно следует, что с точностью до членов более высокого порядка малости для неравномерного течения выполняются неравенства [а<1 и /*<1. Действительно, соотношения (2) и (3) можно переписать в виде:

причем члены В*, С*, Е* и подынтегральное выражение — не отрицательны. (Отметим, что удельный критический импульс для неравномерного изоэнтропического течения г* больше, чем для одномерного г* = /„./0 > = /*с/Ос. Доказательство этого факта

приведено в работе [10]).

Эффективность метода увеличения точности расчета интегральных характеристик в критическом сечении сопла исследовалась с использованием результатов численных расчетов, полученных методами [2-—4]. Выбор этих методов среди многих других не случаен. Метод С. К. Годунова [2], являющийся методом первого порядка точности, развитый и дополненный в [6], применяется при расчете сложных двумерных и пространственных течений [11]. Схема В. П. Колгана [3], представляющая собой модификацию метода С. К. Годунова, ввиду использования кусочно-линейной аппроксимации функций внутри элементарной расчетной ячейки, имеет порядок по пространственным переменным выше первого. В работе [12] отмечается, что оба эти метода обладают чрезвычайно важным свойством монотонности, и приводится подробный анализ погрешностей, возникающих при расчете этими методами. Метод Мак-Кормака [4], являющийся методом второго порядка точности, также применяется для расчета широкого класса задач, причем отмечается хорошее соответствие результатов расчета и экспериментальных данных.

Численно исследовались сопла, входная часть которых состоит из комбинации цилиндра с безразмерным радиусом, равным трем, и конуса с углом наклона образующей к оси, равным 45°, сопряженных дугой окружности радиуса, равного единице. (В качестве единицы длины используется радиус критического сечения /?*). Образующая стенки сопла в критическом сечении является окружностью радиуса Я. Сравнение расчетов, полученных указанными методами для /? = 2 и 0,625, проводится с аналитическим решением [5] и известными экспериментальными данными [13]. При расчетах методами [2] и [3] расчетная сетка состояла из 16X40 ячеек, а при расчетах методом [4] использовалось 21X61 узлов. Показатель адиабаты принимался равным 1,4.

На фиг. 1 приведены профили р, р, и, V в критическом сечении сопл с Я = 2 (сплошные линии) и /? —0,625 (пунктирные линии). Точкой на кривых показаны экспериментальные значения давления

в критическом сечении сопла с /? = 0,625 по данным [13]. Лучшее совпадение с расчетом по формулам аналитического решения [5] {кривые /) для /? = 2 достигается при расчете методом Мак-Корма-ка (кривые 2). Укажем, что максимальная погрешность энтропийной функции при расчете методом Мак-Кормака составляет —0,8%, а при расчете методом Годунова (результаты расчета показаны кривыми 3) составляет 3,5% для сопла с /? = 2 и 5% с Я = 0,625. С. В. Ягудин провел исследования по уменьшению погрешностей, возникающих при расчете сопл методом Годунова. Уточненный расчет на границах, использование сгущающейся к стенкам сопла расчетной сетки, замена в исходной системе уравнений движения уравнения энергии конечным соотношением — интегралом Бернулли, уточнение параметров при расчете распада разрыва и т. д. не дали удовлетворительных результатов. Ввиду ошибок, о чем свидетельствует несохранение, например, энтропийной функции, такие величины, как коэффициент расхода ^ и импульс /, оказываются заниженными.

На фиг. 2 и 3 приведено изменение значений коэффициентов расхода у. (сплошные кривые) и импульса / (пунктирные кривые) сопла с /? = 0,625 и 2 в процессе установления решения. Здесь п — число циклов. На фиг. 2, а и 3, а — результаты интегрирования (1) параметров, рассчитанных методами [2 — 4]. Цифры 1, 2 (см. фиг. 3) и 3 соответствуют использованным методам. На фиг. 2, б и 3, 6 приведены уточненные значения ^ и /*, полученные по соотношениям (2). Сплошной и пунктирной стрелками на фиг. 2, б указаны значения р и /*, полученные, согласно (1), с использованием аналитических решений работы [5]. Экспериментальное значение коэффициента расхода р для сопла с /? = 0,625 равно 0,983 + 0,008 и р = = 0,990 + 0,008 для сопла с Я —2 [13].

Анализ показал, что форма профилей параметров в процессе установления решения достигается довольно быстро, и с ростом числа циклов происходит, в основном, уточнение численных значений параметров. Поскольку при расчете р. и /* по соотношениям (2), согласно выражениям (3), используются лишь безразмерные профили параметров, систематические ошибки, связанные со сдвигом профилей, исключаются. Поэтому расчеты [д. и /* с помощью

ИЛ

1,00

0,98

0,96

0,94

0,92

0,90

иЛ

/ V

'—V' \\ \ • ___ / КУ5

ГХ п

к V т./ !

\ ч 4 N л

) V- V А

* N.

} 1 . 3

2 У

О 200 400 600 800 1000 1200 п а)

’ о 200 т боо воо то то л ________а)

м

0,99

0,98

0,97

0,96

1\ ^ 1\ / =- —

\\)^Л

Л3 —и-—7*

О 200 400 ВОО 800 1000 1200 п

б)

Фиг. 2

■ 0,938 0,996 0,994 0,992 0,990

г

к»-"*’" >1

—-*■ ч

3

~4

О 200 400 600 800 1000 1200 п

б)

Фиг. 3

соотношений (2) оказываются значительно точнее, чем ракеты, согласно формуле (1).

Отметим, ЧТО приведенные В [1] численные значения [1 и для сопла с Я — 0,625 отличаются 'от результатов фиг. 2, что, по-видимому, связано с погрешностями в предварительных расчетах [1].

Для получения выражения для величины импульса в сечении сверхзвукового сопла площади Р представим его, следуя работе [1], в виде:

/=С«(1 + А + В) + рГ. (4)

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

0 = рм/7(1 + А),

уравнения адиабаты (в уравнениях учитываются лишь поправки первого порядка малости)

_Р_

Ро

и уравнения энергии

4[1+Л(х -1)]

1

*(*— 1)

1 + о— Е

иг(1 +2А 4-3В + С)

Р о

где р0, ро — параметры торможения. 58

Поскольку эта система не разрешается аналитически, поступим следующим образом. Обозначим средние значения в рассматриваемом сечении через

Р = Рс(1+Л/), /»=•/»<: 0 + (?). и = ис(1+/?), (5)

где рс, рс, ис — параметры в сечении при одномерном течении,

N, ф, /? —некоторые малые добавки.

Учитывая, что для одномерного течения

Сс = Рс мсГ, Р^-1^

Ро V Ро / 2■*. рс х + 1 (л X — ] 2

1 — Л„

X— 1рс X — 1 \ х+1

где Хс — коэффициент скорости, из перечисленных уравнений получим систему уравнений относительно лУ, С} и Я:

Л^+Я^-О-Л, _ЛГ + -^-(г=^—

х 2

2 X ^ (24 + ЗВ + С)

^Ас /?-<3=— (х— 1)Л с

X— 1

(6)

Значение коэффициента расхода, входящего в первое уравнение, определяется из соотношения (2).

Подставляя полученные из системы уравнений (6) значения <3 и Я в (5) и используя (4), получим:

-----а_| Ч'.+ ч->.) — [,+ 1_Х|(«-1)|Х

1с (х+1)(1 + Х2)1

х(^- + Л+^)-Х2(5 + С)). (7)

Уравнения (6) были ранее получены в работе [1], однако приведенное в ней выражение для удельного импульса, на основании которого был сделан вывод о возможности уточнения величины удельного импульса лишь через безразмерные профили параметров в рассматриваемом сечении сопла, оказалось ошибочным.

Таким образом, для того, чтобы уточнить величину потока импульса сверхзвукового сопла, необходимо знать рассчитанные профили параметров не только в выходном сечении сопла, но также и в критическом сечении [для нахождения уточненной величины [л, входящей в выражение (7)], т. е. нельзя уточнить величину расхода, импульса и удельного импульса на основании лишь рассчитанных профилей в выходном сечении сверхзвукового сопла.

На фиг. 4 и 5 представлены результаты, иллюстрирующие возможность уточнения величины относительного импульса 7 в различных сечениях сверхзвуковой части сопла.

Используемые для анализа контуры сопл, обеспечивающие равномерный и прямолинейный поток на выходе с числами М = 2,197 (фиг. 4) и М = 5 (фиг. 5), рассчитывались методом характеристик и указаны пунктирными линиями (во втором случае рассматривалось укороченное сопло, составляющее 0,11 части от расчетной длины). Кривые 1 соответствуют значениям коэффициента импульса в зависимости от длины сопла I(х), полученные методом характеристик. При расчетах течений в указанных соплах методом

С. К. Годунова расчетная область состояла из 16X46 ячеек (разбиение вертикальных сечений равномерное, шаг в продольном направлении одинаковый).

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

7= '*-+ . ,8)

/с

Кривые 4 показывают уточненные значения 1(х) с помощью соотношения (7). Значения Т(х) при малых длинах сопла, полученные по формуле (8), лежат выше соответствующих значений, полученных методом характеристик, что объясняется погрешностями метода [2] при расчете ускоряющегося потока на начальном участке сопла. С увеличением длины сопла определение / по формуле (8) приводит к заниженным значениям коэффициента импульса по сравнению с методом характеристик. При расчете укороченного сопла (см. фиг. 5) наблюдаются аналогичные, но более ярко выраженные тенденции в поведении кривых 7(х), определенных отмеченными

способами. Для уменьшения ошибок при расчете течения на начальном участке сопла был проведен расчет начального участка укороченного сопла с уменьшенным вдвое шагом вдоль оси х. В этом случае кривые, соответствующие рассматриваемому способу определения I (.х), отмечены цифрами со штрихами. Определение уточненных значений Т {х) с помощью приближенного соотношения (7), как и в случае сопла, рассчитанного на М = 2,2, приводит к неплохим результатам, однако в этом случае, по-видимому, более целесообразно пользоваться точным решением системы уравнений, состоящей из уравнения расхода, уравнения адиабаты и уравнения энергии.

Таким образом, приведенная в [1] и в настоящей статье мето дика, состоящая в использовании обобщения гидравлической теории на случай неравномерных профилей параметров (численные характеристики неравномерностей определяются на основании ко-нечно-разностных решений уравнений газовой динамики), может быть применена как в горле, так и других сечениях сверхзвукового сопла (при этом получаются, конечно, разные соотношения) и оказывается весьма эффективным средством для получения достаточной точности при решении прикладных задач, связанных с расчетом интегральных характеристик сопл двигателей. Этот подход может быть применен также для некоторых типов пространственных сопл, обладающих двумя плоскостями симметрии, например,, эллиптических или прямоугольных, и может быть обобщен на случай совместного истечения через сопло Лаваля двух или нескольких потоков.

В заключение отметим, что для уточнения величины расхода при численных расчетах методом [2] существует другой прием, предложенный для случая расчета трансзвуковых течений в плоских турбинных решетках [14]. Этот прием опирается на допущение о высокой точности расчета давления и угла наклона вектора скорости 0 и заключается в пересчете модуля скорости и плотности на основании уравнения энергии и из условия постоянства энтропии. Мы применяли этот прием для определения потоков массы и импульса в сверхзвуковых соплах. Оказалось, что он дает хорошие результаты в минимальном сечении, хотя, как следует из фиг. 1, указанные в [14] допущения строго не выполняются, что приводит к резкому увеличению ошибки в расходе при незначительных отклонениях от критического сечения. В результате неоднократных обсуждений этого обстоятельства с М. Я. Ивановым, одним из авторов [14], и с А. Н. Крайко было выяснено, что в этом случае причина хорошей точности при определении расхода в критическом сечении связана не с высокой точностью определения давления и угла наклона вектора скорости при расчете методом С. К. Годунова, а с тем, что при точном выполнении уравнения энергии и соотношения для изоэнтропы ошибки при определении расхода и потока импульса малы в силу экстремальных свойств подынтегральных выражений (1) при М = 1 и 9 = 0.

Авторы признательны А. П. Мазурову, предоставившему программу расчета течений в соплах методом Мак-Кормака.

1. 3 и м о н т В. Л. Метод повышения точности определения интегральных характеристик потока в сверхзвуковых соплах. .Ученые записки ЦАГИ“, т. 7, №,2, 1976.

2. Годунов С. К,, Забродин А. В., Прокопов Г. П. Разностная схема для двумерных нестационарных задач газовой динамики и расчет обтекания с отошедшей ударной волной. „Ж. вычисл. матем. и матем. физ.*, т. 1, № 6, 1961.

3. Колган В. П. Конечно-разностная схема для расчета дву-

мерных решений нестационарной газовой динамики. „Ученые записки ЦАГИ“, т. 6, № 1, 1975.

4. М а с-С о г m а с k R. W. Jhe effect of viscosity in hypervelocity

impact cratering, AIAA Paper 69-354, 1969.

5. К liege 1 I. K., Levine 1. N. Transonic flow in small throat

radius of curvative noggles AIAA J., vol. 7, N 7, 1969.

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

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

8. Мазуров А. П. Расчет течения в сопле методом установления второго порядка точности. .Ученые записки ЦАГИ“, № 6, т. 7, 1976.

9. Hall Т. Н. Transonic flow in two-dimensional and axially-symmet-ric nozzles. .Quart. J. Mech. Appl. Math.“, vol. 15, pt. 4, 1962.

10. Крайко A. H., Соколов В. E. Об удельном импульсе потока в минимальном сечении сопла Лаваля. .Изв. АН СССР, МЖГ“, 1976, № 1.

11. Годунов С. К., Забродин А. В., Ив а но в М. Я., Крайко А. Н., П р о к о п о в Г. П. Численное решение многомерных задач газовой динамики. М., „Наука", 1976.

12. Косых А. П., Минайлос А. Н. Исследование методов сквозного счета для задач сверхзвуковой аэродинамики. „Ученые записки ЦАГИ“, т. 7, № 1, 1976.

13. В а с k L. Н., С u f f е 1 R. F, Flow coefficients for supersonic nozzles with comparatively small radius of curvature throats. .Journal of spacecraft and Rockets", vol. 8, N 2, 1971.

14. Богод А. Б., Грановский A. B.p Иванов М. Я. Численное исследование некоторых особенностей трансзвуковых течений в плоских турбинных решетках. „Изв. АН СССР, МЖГ“, 1976, № 2.

Рукопись поступила 25/II 1977