CC BY
73
УЧЕНЫЕ ЗАПИСКИ ЦА Г И Т о м II 197 1
М 5
УДК 533.6.011.5
ОБТЕКАНИЕ СФЕРИЧЕСКОЙ ПОВЕРХНОСТИ СВЕРХЗВУКОВЫМ ПОТОКОМ РАВНОВЕСНО-ДИССОЦИИРУЮЩЕГО ВОЗДУХА
А. П. Косых, А. Н. Минайлос
Получены решения задачи обтекания сферической поверхности гиперзвуковым потоком воздуха в большом диапазоне значений скорости и высоты Н. Задача решалась численно методом установления с использованием явной конечноразностной схемы.
Проведены сравнения и анализ результатов, полученных различными численными методами.
1. Расчету обтекания тел сверхзвуковым потоком газа, нахо^ дящегося в состоянии термодинамического равновесия, посвящено большое число работ. Однако результаты систематических расчетов, соответствующие полету в атмосфере Земли при различных скоростях 1^оо и на различных высотах Н, можно найти только в работах [1]—[3]. Отсутствует сопоставление результатов, полученных различными численными методами, почти нет также анализа полученных результатов. Настоящая работа частично восполняет эти пробелы и расширяет область имеющихся численных решений.
2. Расчеты проведены по явной конечноразностной схеме методом установления, близким к описанному в работах [4] и [5]. Параметры во внутреннем поле рассчитывались по двухшаговой схеме, представляющей собой вариант схемы Лакса—Вендрова. При расчете параметров на поверхностях тела и ударной волны использовались характеристические соотношения. В поле течения с учетом граничных точек располагалось (23X11) узлов сетки. Рассчитывалось обтекание сферической поверхности до центрального угла X, равного 110° включительно (лучи располагались равномерно через 5°). Для учета реальных свойств воздуха применялась методика, предложенная Ю. Н. Дьяконовым [6]. По этой методике скорость звука определяется интерполяцией по таблицам. Зная скорость звука, можно проводить расчеты во внутренних узлах поля. Чтобы избежать потерь времени при итерациях на ударной волне, для каждого варианта использовалась табличная зависимость плотности за ударной волной от нормального к волне
компонента скорости набегающего потока. Эта зависимость рассчитывалась один раз до счета варианта с помощью стандартизированной программы, использующей уточненные аппроксимации термодинамических функций, предложенные В. В. Михайловым [7].
Точность получаемых численных решений оценивалось по значениям интеграла Бернулли в поле течения, по расходу в контрольном сечении, по степени установления варианта. Для большинства вариантов отклонения интеграла Бернулли от точного значения не превосходят 2%. Для. трех наиболее трудных для счета вариантов, соответствующих большим скоростям и высотам полета, отклонения достигают 4%.
Различия в расходе массы через контрольное сечение в сверхзвуковой области (через луч с центральным углом X, равным 110°) и через участок ударной волны, расположенный до этого сечения, не превосходят 2,5%.
3. Основным критерием достоверности результатов численных методов при отсутствии доказательств существования и единственности решения служит сравнение этих результатов с данными эксперимента, с накопленным опытом. Методы расчета обтекания затупленных тел совершенным газом сравнивались неоднократно как с экспериментальными данными, так и между собой (см., например, работы [3], [5], [8]). Эти сравнения показывают, что для совершенного газа имеются достаточно точные результаты (параметры поля течения, полученные различными методами, отличаются не более, чем на 2—4%). При расчете течений равновесно-диссоциирующего газа к ошибкам метода добавляются ошибки расчета термодинамических функций газа (в случае течения совершенного газа их практически нет). В тех диапазонах изменения параметров набегающего потока, где результаты наиболее важны (случай высоких температур) и где точность расчета термодинамических функций чаще всего наиболее низкая, практически нет возможности детально сравнить результаты с экспериментальными данными (экспериментальные данные отсутствуют, либо их точность очень низка), поэтому чрезвычайно существенно сопоставление результатов, полученных различными численными методами.
Таблица 1
Моо Роо [атм] Метод Обо- зна- чения Ссылки
10 0,01 300 МИС, 1 схема, N = 2 1 [9]
20 0,054546 216,65 МИС, I схема, N = 2 2 [10}
20 0,01 250 МИС, 1 схема, N=1 3 [И]
40 0,01 250 МИС, I схема, N = 1 4 [И]
20 0,01 250 Метод прямых О [12]
20 0,05 300 МИС, 11 схема ■ [13]
20 0,00025 250 МИС, 11 схема □ [14]
20 0,0563 (//=20 км) 216,65 (//=20 км) Метод сеток, новление уста- А [3]
20,06 0,01092 231,89 Обратный метод 0 [1]
К сожалению, в литературе отсутствуют расчеты различными методами одного и того же варианта течения. Поэтому помимо систематических расчетов, соответствующих движению в атмосфере Земли, были проведены методические расчеты ряда вариантов, для которых условия в набегающем потоке соответствовали рассмотренным в работах [9]—[14]. Кроме этого, проводится сравнение результатов расчета течений при полете в атмосфере Земли.
4. Рассмотрим результаты методических расчетов. В табл. 1 приведены данные об условиях в набегающем потоке для тех вариантов, с которыми сравниваются полученные результаты. Там же указаны принятые для этих вариантов на фиг. 1 обозначения и соответствующая литература. Буквой N обозначен номер приближения в методе интегральных соотношений (МИС).
На фиг. 1 проведено сравнение результатов для распре- 3 деления давления на поверхно- Ро сти тела. Величина давления отнесена к его величине в критической точке. Это одна из самых консервативных величин ' поля течения при изменении условий набегающего потока и свойств газа. Действительно, исходя из теории гиперзвуко-вых течений, свойства газа ска- о зываются на распределении давления в основном при учете так называемой поправки Буземана (влияние центробежных сил). Однако влияние центробежных сил, уменьшающее давление на теле, оказывается существенным (примерно 4—5%) только для значений х<1,1 [15], а при полете в атмосфере Земли "эффективное значение * изменяется в пределах 1,1-7-1,4 [1]. Все результаты методических расчетов по распределению давления при числе Мсо = 20 легли на одну кривую (сплошная линия на фиг. 1).
Кривая 4, соответствующая Моо = 40, расположена ниже. В окрестности звуковой точки на теле отличие составляет 3%. Результаты работ [1], [3], [12] и [14] (четыре счетных луча) совпадают с результатами наших методических расчетов. Распределения давлений на теле, полученные в работе [13] (три счетных луча) и в данной статье, отличаются на 4% в районе звуковой точки.
Во всех этих методах термодинамические функции рассчитывались на основе различных аппроксимаций. Существенное отличие дает сравнение с первой схемой МИС. Первое приближение МИС (расчеты проведены по методике работы [11]) дает в окрестности звуковой точки отличие в 8,5% при Моо = 20 и в 14% при М» —40. Второе приближение [10] отличается в том же районе на 19% (при Моо = 20). Поскольку результаты остальных методов близки между собой, следует считать, что первая схема МИС определяет давление на теле с большой ошибкой. Сопоставляя результаты при числе Моо=Ю и 20 во втором приближении и числах Моо = 20 и 40 в первом, можно сделать вывод о том, что ошибка увеличивается с ростом числа М» и не зависит от метода определения термодинамических функций воздуха (в работах [10] и [И] эти методы различны). Следовательно, она вызвана не различным за-
Фиг. 1
данием термодинамики, а методикой основного расчета. Возможно, это объясняется большими трудностями прохода особых точек (см. пункт 10 в работе [11]). Сравнения результатов, полученных в работе [11], с результатами других методов проводились до этого либо для совершенного газа, либо при меньших сверхзвуковых скоростях [8], что не позволило выявить рост ошибки. Отметим, что для тел с угловыми точками ошибки метода должны быть меньше, так как скорость в угловой точке задается априори. В связи с анализом I схемы МИС напомним также один известный факт: первое приближение этого метода в случае равновесно-диссодиирующего газа неверно определяет градиент скорости в критической точке.
За исключением первой схемы МИС, результаты всех остальных методов вполне удовлетворительно согласуются между собой. В табл. 2 дано сравнение полученных параметров поля течения (Н =20 км, Моо = 20) между поверхностью сферы (я = 0) и поверхностью ударной волны (п= 1) в сверхзвуковой области течения на луче X = те/2 с данными из работы [3]. Компоненты скорости отнесены к предельной скорости Ушах, плотность —к плотности в набегающем потоке роо, а давление—к р», У^ах- Максимальное отличие результатов, полученных вблизи поверхности тела, не превосходит 3,6%.
Таблица 2
п =0 п — 0,25 п = = 0,5 п— 0,75 п = 1
а 0,6461 0,6455 0,6350 0,6308 0,6576 0,6548 0,7143 0,7123 0,7825 0,7833
V 0 0,0008 0,1224 0,1250 0,2037 0,2054 0,2841 0,2851 0,3714 0,3717
Р 0,03608 0,03686 0,06241 0,6467 0,08943 0,09139 0,1292 0,1299 0,2110 0,2111
Р 0.7290 0,7539 1,246 1,291 1,967 1,964 3,483 3,427 7,832 7,819
Примечание. В правых колонках приведены результаты работы [3].
5. Рассмотрим теперь обтекание сферы сверхзвуковым потоком на различных высотах в атмосфере Земли. На фиг. 2 сплошными линиями показаны границы областей влияния разреженности,
неравновесности и излучения воз-"I духа на течение у сферы радиу-\ сом 2 м [16]. Кривая 1 отделяет область, в которой весь ударный ^ слой является вязким. Для кривой
2 толщина ударной волны имеет порядок О,1е0(е0 — расстояние от волны до критической точки на теле). Можно считать, что ниже этой кривой расположена область,
------где течение в основном невяз-
„[км/сех] кое Кривые 3 и 4 характеризуют зону влияния неравновесности.
Зона релаксации на оси течения для кривой 3 охватывает весь ударный слой е0, а для кривой 4— €,1в0. Несколько выше кривой 3 течение можно считать замороженным, ниже кривой 4—равновесным. Кривая 5 соответствует течению, для которого конвективный поток тепла в критическую
Фиг. 2
точку равен лучистому. В области ниже этой кривой учет излучения может существенным образом изменить ряд характеристик течения, таких как температура, плотность, потоки тепла, градиент скорости, толщина ударного слоя. Штриховыми линиями показаны характерные траектории гиперзвуковых аппаратов: баллистической головной части ракеты (кривая 6), гиперзвукового самолета (кривая 7), спутника (кривая 8), космического корабля типа „Аполлон11, входящего в атмосферу со второй космической скоростью (кривая 9) и, наконец, пилотируемого аппарата, входящего в атмосферу со скоростью, превышающей вторую космическую (кривая 10). Значительная часть траекторий расположена в области, где течение неравновесно, однако аппараты с траекториями 6 и 10 проходят в области с равновесно-диссоциирующим газом, а траектории 8 и 9 имеют большие участки, где течение „околоравновесное“.
Кроме того, в области течений равновесно-диссоциирующего воздуха расположена значительная часть траекторий антиракет, обладающих уже на малой высоте большой скоростью полета.
Наконец, треугольниками показаны рассчитанные варианты. Ряд вариантов рассчитан для течений, в которых необходим учет неравновесности или излучения. Наиболее существенную поправку в области проведения расчетов дает учет излучения. Это необходимо иметь в виду при практическом использовании результатов. Влияние излучения рассмотрено ниже на основе работ [16] и [17].
6. Сводные результаты расчетов представлены на фиг. 3—5.
Сплошными линиями изображе- ' ны результаты настоящей работы £а (Я=0; 5; 10; 15; 30; 45; 57; 60 км); штриховыми с черными кружками — работы [1] (Я=^30,48; 45,72;
60,96 км)-, белыми кружками—
работы [3] (Я=10; 20; 30; 50 км)-, крестами—работы [2] (Я^Ю; 30;
60 км); белыми треугольниками — данные, полученные по методике работы [14] (//=60 км). Наконец, № влияние излучения показано ромбами (работа [16], 67 км,
/? = 2 м) и квадратами (работа [17],
Н=57 км, /?=3 м). При этом белые ромбы и квадраты соответствуют течению без учета излучения^ черные — с учетом излучения.
На фиг. 3 приведено расстояние от сферы до ударной волны е„, отнесенное к радиусу сферы /?. д Отличия в значениях е0, полученных в различных работах, не превосходят 3%.
Для всех кривых фиг. 3 характерна известная ранее немонотонная зависимость значения е0 от скорости (см. также фиг. 4 и 5). Она вызвана процессами диссоциации кислорода и азота в ударном слое. На фиг. 3 внизу изображены молярные концентрации ц основных компонентов воздуха в окрестности критической точки на высоте 30 км в зависимости от скорости полета. Сопоставле-
Фиг. з
ние изменения р и е0 показывает, что для областей с сильной диссоциацией характерно уменьшение значений е0 по сравнению со значениями для течения совершенного газа, достигающее 10% при диссоциации кислорода (1/оо = 3 4,5 км/сек) и 18% при диссо-
циации азота (У» — 6 -5- 9 км/сек). С ростом высоты области сильной диссоциации сдвигаются в сторону меньших значений Ут.
Влияние излучения на величину е0 начинает сказываться со скорости набегающего потока 10 км/сек. Уменьшение значений е0, вызванное высвечиванием ударного слоя, достигает 23% при 1^00= 15 км/сек (при //=45 км), 21% (при Н—Ы км) и 18% (при Н = 67 км).
На фиг. 4 приведены значения градиента скорости ди/дх в критической точке. Зависимость ди/дх от скорости Ут качественно
напоминает поведение функции е0. Проведено сравнение с результатами работ [1], [3] и [14]. Таблицы работы [2] не позволяют определить градиент с приемлемой точностью.Максимальное отличие,
031
6 10 м/сех]
Фиг. 4
037.
Г^/У= 4-5 мм
^ 20 км
, 10км \ г*
!(] 0±!0нм
\ / 15 нм ЗОхм 00 км
У . А
\ ^1 А V
V X \ /, '45хм 57 нм 60 нм
\ \
X
10
[хм/сех]
Фиг. 5
полученное при сравнении с результатами работы [3], составляет 3,2 % {Н == 10 км, Уоэ = 4,5 км/сек).
Учет излучения уменьшает градиент на 19,5% на высоте 57 км и на 27% на высоте 67 км (данные взяты при значении скорости Уоа = 14 км/сек). На высоте 45 км соответствующее уменьшение градиента составляет 15,2% [17] (для этой высоты точки на графике не приведены, чтобы излишне его не перегружать). Следовательно, влияние излучения на градиент скорости возрастает с ростом высоты полета (в противоположность влиянию на значения е0).
Для значений сх полусферы (Л' = тс/2) тоже характерна немонотонная зависимость от скорости полета (см. фиг. 5). Однако она выражена значительно слабее: все изменение в диапазоне скоростей 3 — 15 км/сек не превосходит 4%. Значения сх для сравнения имеются только в работах [2] и [3], причем в работе [2] некоторые величины сх даны с недостаточной для наших целей точностью (примерно 1%) и поэтому на фиг. 5 эти значения не наносились. Точность значений сх определяется точностью расчета давления на поверхности тела и точностью интегрирования. Интеграл вычислялся методом трапеций и методом Симпсона. Оценки показы-
вают, что погрешность интегрирования не превосходит 0,1% величины сх. Погрешность при расчете давления на поверхности тела значительно ниже максимальной погрешности решения (см. разд.
2 и 4). Так, наши результаты отличаются от данных А. Н. Любимова и В. В. Русанова [3] не более чем на 0,001 р'о {Н = 20 км, Моо = 20). По-видимому, это отличие и характеризует ошибку в давлении на теле. Порядок величины этой погрешности такой же, как и погрешности при расчете интеграла. Суммарная погрешность, полученная на основании такого анализа, не превосходит 0,2% величины сх.
Колебательный характер зависимости сх(Уоо) объясняется немонотонностью значений р/р'0, а в конечном итоге диссоциацией: при диссоциации падает температура, растет плотность, уменьшается эффективное значение х и значение р/ро в сечении Л^сопб^ Два минимума в зависимости сх(Уао) на фиг. 5 соответствуют диссоциации кислорода и азота. Они несколько сдинуты в сторону больших значений Уаа по сравнению с минимальными значениями
3 и ди\дх. Так, второй минимум расположен при значении Ут = = \\ ым/сек, а не 8 км\сек. Этот сдвиг объясняется тем, что сх не локальная, а интегральная характеристика, и экстремум достигается при интегральной оценке диссоциации у поверхности тела.
Для скоростей, меньших 6 км!сек, значения сх увеличиваются с ростом высоты, а для скоростей, бблыних 6,5 км/сек, уменьшаются. Рост величины сх с увеличением высоты при малых скоростях объясняется в основном преобладающим влиянием множителя Ро/Роо Ум В формуле ДЛЯ Сх.
В окрестности значений = 5; 6,5; 15 км/сек зависимость
сх — сх{Н) ослабевает. При числе Моо= 20, соответствующем скорости 1/00 = 6,1 — 6,5 км/сек, до больших значений X величина р/р'0не зависит от давления и температуры в набегающем потоке (см. фиг. 1). В связи с этим кривые сх(Усо) в этом диапазоне скоростей совпадают с точностью до 0,1%.
Были рассмотрены также зависимости гзв. т и гзв. в от скорости набегающего потока. Так обозначены расстояния от оси симметрии течения до звуковых точек на поверхностях тела и ударной волны. Эти расстояния характеризуют размеры области влияния. Зависимость гзв. в немонотонна и качественно напоминает графики на фиг. 3 и 4. Величины гзв. в лежат в пределах от 0,26 до 0,44 при изменении значений Уоо в диапазоне от 3 до 15 км/сек. Функции Гзв. т не имеют перегибов. Эти функции значительно консервативней и изменяются в диапазоне от 0,595 до 0,632.
7. В работе [15] предложены аппроксимационные зависимости для определения в потоке совершенного газа значений е0, да/дх, Сх(Х = к/2), Гзв. Т и Гзв. В как функций от отношения плотностей в прямой ударной волне к = Р00/Р1, числа Моо, давления за прямой волной р1/роо Ушах. Вид зависимостей предполагал возможность применения их и в случае равновесно-диссоциирующего газа.
Проверка этих зависимостей, проведенная на основе материалов настоящей статьи, показала, что точность аппроксимаций в этом случае в диапазоне изменения скоростей набегающего потока от 3 до 15 км\сек составляет для функций е0 — 3,5%, ди/дх — 2,7%0, гзв. т —0,17%, гзв. в — 3,6%, с, — 1,8%. Поскольку диапазон изменения параметров подобия, введенных в работе [15], значительно уже для равновесно-диссоциирующего газа, чем для
совершенного, точность аппроксимации для случая равновесно-дис-социирующего газа можно несколько повысить, улучшив аппроксимацию без изменения параметров подобия.
Необходимо отметить, что рассматриваемые аппроксимации не учитывают влияние излучения.
В заключение авторы благодарят В. К. Душина, Ю. Г. Елькина и В. К. Молодцова за предоставленные для сравнения материалы расчетов.
ЛИТЕРАТУРА
1. Lomax Н., Jnouye М. Numerical analisis of flow properties about blunt bodies moving at supersonic speeds in equilibrium gas. NASA TR-NR-204. 1964.
2. Лунев В. В., Магомедов К. М., Павлов В. Г. Ги-перзвуковое обтекание притупленных конусов с учетом равновесных физико-химических превращений. М., Изд. ВЦ АН СССР, 1968.
3. Любимов А. Н., Русанов В. В. Течения газа около тупых тел., М., „Наука", 1970.
4. Moretti G., Abbett М. A time dependent computational method for blunt body flows. A1AA J., 1966, v. 4, No 12.
5. Косых А. П., Минай л ос А. Н. О явных схемах метода установления в задаче сверхзвукового обтекания затупленного тела. Ж. вычисл. матем. и матем. физ., 1970, т. 10, № 2.
6. Д ь я к о н о в Ю. Н. Пространственное обтекание затупленных тел с учетом равновесных физико-химических реакций. ДАН СССР, 1964, т. 157, № 4.
7. Михайлов В. В. Приближенное аналитическое представление термодинамических функций воздуха. ,Инж. сборник", т. 31,
М., Изд. АН СССР, 1961.
8. Обтекание затупленных тел сверхзвуковым потоком газа.
Сб. под ред. О. М. Белоцерковского. М., ВЦ АН СССР, 1966.
9. Белоцерковский О. М., Голомазов М. М., Ш у л и ш-н и н а Н. П. Расчёт обтекания затупленных тел с отошедшей ударной волной потоком равновесно-диссоциирующего газа. Ж. вычисл. матем. и матем. физ., 1964, т. 4, № 2.
10. Белоцерковский О. М. Симметричное обтекание затупленных тел сверхзвуковым потоком совершенного и реального газа. Ж. вычисл. матем. и матем. физ., 1962, т. 2, № 6.
11. Козлова И. Г., Минайлос А. Н. Несимметричное обтекание лобовой части тела вращения сверхзвуковым потоком совершенного или реального газа. Ж. вычисл. матем. и матем. физ., 1967, т. 7, № 3.
12. Те ленин Г. Ф., Тиняков Г. П. Исследование сверхзвукового обтекания сферы воздухом и углекислым газом при термохимическом равновесии. ДАН СССР, 1964, т. 159, № 1.
13. Молодцов В. К. К расчету равновесных течений газа.
Сб. .Вычисл. матем. и матем. физ.“, М., „Наука", 1966.
14. Белоцерковский О. М., Душин В. К. Сверхзвуковое обтекание затупленных тел неравновесным потоком газа. Ж. вычисл. матем. и матем. физ., 1964, т. 4, № 1.
15. Базжин А. П., Минайлос А. Н., Благоскло-нов В. И., Пирогова С. В. Обтекание сферы сверхзвуковым потоком совершенного газа. .Ученые записки ЦАГИ“, 1971, т. 2, № 3.
16. СтуловВ. П., Шапиро Е. Г. Излучение ударного слоя при гиперзвуковом обтекании затупленных тел воздухом. Известия АН СССР, Мех. жидкости и газа, 1970, № 1.
17. Елькин Ю. Г. Гиперзвуковые течения селективно-излу-чающего газа. Труды ЦАГИ, вып. 1258, 1970.
Рукопись поступила 12/1 1971 г.
