CC BY
52
УЧЕНЫЕ ЗАПИСКИ Ц А Г И Т О м IV 19 7 3
№ 4
УДК 518.5.533.6.011.5
К РАСЧЕТУ ПРОСТРАНСТВЕННОГО ОБТЕКАНИЯ СВЕРХЗВУКОВЫМ ПОТОКОМ ТЕЛ СЛОЖНОЙ ФОРМЫ
М. Я- Иванов, Т. В. Никитина
Проанализирована картина обтекания сложной пространственной конфигурации типа летательного аппарата сверхзвуковым потоком невязкого и нетеплопроводного газа под углом атаки. Для численного интегрирования стационарной системы уравнений, записанной в виде законов сохранения, использована конечноразностная схема сквозного счета, примененная ранее к расчету двумерных и пространственных сверхзвуковых течений в соплах и струях. Проведено сравнение результатов настоящей работы с расчетами по методу характеристик и одному из вариантов схемы Лакса— Вендрова.
Метод сквозного счета [1], [2], использованный в работе [3] для исследования сверхзвукового обтекания различных конических тел, в настоящей статье применяется к анализу картины течения около конфигурации типа летательного аппарата, расположенного под углом атаки в сверхзвуковом потоке идеального газа. Из существующих численных методов исследования пространственных сверхзвуковых течений идеального газа наибольшее распространение получили различные варианты метода характеристик и неявного конечноразностного метода, предложенного К. И. Бабенко и др. [4]. С использованием этих методов решено большое число задач внешнего обтекания гладких тел сверхзвуковым потоком. При наличии в области течения достаточно сложных систем поверхностей разрыва применение указанных методов сопряжено со значительным усложнением численных алгоритмов и программ для ЭЦВМ, что связано с необходимостью выделения разрывов. В последнее время расчет пространственных сверхзвуковых течений все чаще проводится с помощью методов сквозного счета, основанных на дивергентной форме записи уравнений газодинамики. Так, анализ трехмерной картины течения около целого ряда конических тел проведен в работах |3] и [5], причем работа [3] использует метод, предложенный и опробованный в работах [1] и [2], а работа |5] — одну из модификаций схемы Лакса — Вендрова. Все упомянутые работы содержат также обширную библиографию, сравнение и анализ большого числа конечноразностных методов
для расчета двумерных и пространственных стационарных сверхзвуковых течений идеального газа.
Интенсивное развитие работ по численным методам сделало возможным решение весьма сложных и важных для практики задач, например, расчета поля течения около тел типа орбитального летательного аппарата в сверхзвуковом потоке [6] и [7]. В настоящей статье расчет обтекания подобной конфигурации, рассмотренной ранее в работе [б], выполнен по разностной схеме [3]. Проведено сравнение полученных результатов с результатами расчетов, выполненных в работе [6] методом характеристик и методом сквозного счета работы [5].
1. Исходными уравнениями при построении конечноразностной схемы, как и в работе [3], служат интегральные законы сохранения, записанные в сферической системе координат
d
dr
JJ ac?6d<p= C^[(c— c№)db — (b— a£9)d:p]+ JJ fdbdy. (1)
Здесь a, b, с и /— четырехмерные векторы, определяемые как
pv
Ь = Г S І П 9 I 2
Р + pv2
pvw
о
2 р + р (v2 + а»2)
(р -f- рву2) ctg 8 — puv pw (v ctg 0 4- u)
В (1) и далее обозначения для газодинамических параметров общепринятые, — сферическая система координат, Г — некоторый замкнутый контур на сферической поверхности г = const, ограничивающий площадку 5. Контур Г и площадка S считаются функциями г. Деформация контура Г описывается касательным
к сферической поверхности г = const вектором £ = dnjdr, где dti — смещение Г в направлении внешней нормали при переходе от сферической поверхности г --= const к r + afr = const. В каждой точке Г
вектор £, как и dn, перпендикулярен радиус-вектору г и полностью определяется своими проекциями Г?9 И .
Рассматривается стационарное обтекание тела сверхзвуковым равномерным потоком совершенного газа с постоянными теплоемкостями в предположении, что во всей рассчитываемой области поток г—сверхзвуковой, т. е. проекция вектора скорости на ось г больше местной скорости звука (и > а).
Для описания формы рассматриваемой конфигурации (модели) введем также связанную с телом декартову систему координат xyz, причем начало декартовой системы, а также и сферической системы, совместим с вершиной тела. Тогда уголкоторый отсчитывается от отрицательного направления оси г, определяет положение меридиональной плоскости, проходящей через радиус-вектор г и ось х, а 0 — угол между г и осью х. Набегающий по+ок
определяется заданием числа М.ж и угла атаки а, причем предпо-. -*■ лагается, что вектор скорости лежит в плоскости хг.
Введение безразмерных параметров проводится отнесением пространственных координат к характерной длине тела, плотности и компонентов скорости — к своим критическим значениям в набегающем потоке и давления — к произведению критической плотности на квадрат критической скорости.
Система (1) замыкается условием постоянства полной энтальпии, которое в безразмерных переменных имеет вид
2 % р „ *.+ 1 ----г — + и2 + V + ге;2 =------г ,
У. — I р X— 1’
где х — показатель адиабаты.
Метод численного решения системы (1), порядок вычислений и полная система конечноразностных уравнений приведены в работе [3]. Формулы для расчета параметров на границах разностных ячеек даны в работе [1].
2. Как уже указывалось, для выяснения возможностей метода при численном решении задач сверхзвукового обтекания достаточно сложных пространственных тел был проведен расчет течения около тела типа летательного аппарата с дельтавидным крылом, имеющим затупленные передние кромки. Расчет обтекания данной конфигурации с использованием других численных методов проведен в работе [6]. При вычислениях использовалась аналитическая аппроксимация модели летательного аппарата, приведенная в этой работе. Так, поперечные сечения модели получались сопряжением эллипсов с различной степенью эллиптичности для наветренной и подветренной сторон поверхности. Для задания контуров модели в плоскости симметрии хг, т. е. правых частей уравнений г = 2+(х) и г = 2-(х), и формы в плане, т. е. аналогичной функции у+{х), использовались полиномы третьей степени, причем, как и в работе [6], для удобства аналитического описания исследуемое тело разбивалось на шесть сегментов. За характерный размер принималось расстояние от начала координат до концевого сечения пятого сегмента, так что полная длина модели ¿ — 1,275. Аппроксимирующие полиномы и их коэффициенты, взятые из работы [6], приведены в таблице.
*_(•*). 1
г+ (х) I = а0+а, х+а, хг+а3 х'\ х-х-х4
У+ <■*>
се Коэффициенты “і
О- * Ф а> г ^ (подветренная сторона) (наветренная сторона) у+ (боковая сторона)
"Г ^ .X о Оо аі а3 а3 «0 аі | а% а. «о а, аа а3
1 0,05 0,025 0,26795 -4,0678 30,176 0,025 0,26795 0 0 0,025 0,26795 —1,0446 1,7120
2 0,10 0,032 0,08749 -0,4499 0,68725 0,0384 0,26795 0,12770 -0,88725 0,036 0,17633 —0,3383 0,90825
3 0,30 0,037 —0,010 0 0 0,090 0,21256 —0,22505 0,4385 0,065 0,150 0 0
4 0,50 0,035 -0,010 0 0 0,120 0,06993 -0,0336 -0,31525 0,095 0,150 —0,02935 0,77175
5 0,70 0,033 —0,010 0 0 0,130 0,01746 0,02121 —0,07952 0,130 0,23078 1,3697 -1,7605
6 1.00 0,030 —0,010 0 0 0,135 0,00873 0 0 0.275 0,57735 0 0
or
S миф
f -ІИф
g -лиф
о09I oOZl сOS
Z '-іиф
ЧІИф
На фиг. 1—6 представлены результаты расчета обтекания описанного летательного аппарата сверхзвуковым потоком с числом Моо = 5 под углом атаки а = 5°. Носок модели (л:-<0,05), как и в работе [6], был коническим и плавно сопрягался при х = 0,05 с основной частью контура. Задание начальных данных на сферической поверхности, пересекающей ось абсцисс при л = 0,05, для
Фиг. 6
сферы проводится в соответствии с работой [3]. Число расчетных слоев N в направлении координаты 6, т. е. между телом и внешним скачком, равнялось 20, и число слоев К в направлении координаты у между плоскостями симметрии на наветренной (<? = 0) и подветренной (9 = 180°) сторонах также равнялось 20. Результаты работы [6] получены при 11X19 и 29X19 в случае исполь-
зования метода характеристик и метода сквозного счета соответственно.
На фиг. 1 показана * форма контура (линии 1) и скачков уплотнения (линии 2) в двух плоскостях ху и хг. Штриховые линии и светлые кружки здесь и ниже отвечают методу характеристик и методу сквозного счета [6]. При х = 0,7 от крыла отходит вторичный скачок, который в окрестности л:=1,0 достигает носового скачка.
Распределение давления по х вдоль поверхности модели для <р = 0 и 180° приведено на фиг. 2. Изменение давления по поверхности тела в поперечном сечении при л; = 0,7 и 0,765 показано на фиг. 3,а, б. В окрестности корневой хорды крыла наблюдается резкое увеличение давления, вызванное возникновением вторичного скачка. О распределении параметров в ударном слое при <р — 90°, х = 0,7 и 0,765 можно судить по данным фиг. 4, а, б. На этих фигурах по оси ординат отложена нормализованная коорди- . ната С = (у—ут)/(Ув ~Ут)> изменяющаяся от 0 на теле до 1 на ударной волне (_ут и ув — координаты тела и ударной волны соответ-
ственно). Отметим, что в рассматриваемом случае плоскость ®=^90° является наиболее сложной для расчета, что обусловлено формой рассчитываемого тела.
На фиг.5 приведены линии постоянства давления,отнесенного к давлению набегающего потока (р/рт = const), носовой скачок и контур летательного аппарата в сечении х = 0,765. Кривые plpco = const нанесены через интервал 0,2. Сгущение линий отвечает распространяющемуся в потоке вторичному скачку.
Пространственная картина обтекания летательного аппарата изображена на фиг. 6, на которой показано пересечение поверхности модели, носового скачка и поверхностей pjp^ — const с плоскостями у= 0; х—0,2; 0,6 и 1,0.
Расчет обтекания летательного аппарата проведен дол:=1,04. Пересечение втооичного скачка с головной ударной волной (см. фиг. 1) происходит при д: it; 1,02. Метод, использованный в настоящей работе, позволяет рассчитывать интерференцию скачков уплотнения без усложнения вычислительного алгоритма, что было отмечено и проиллюстрировано в работах [1] — [3]. Поскольку при дальнейшем увеличении х толщина ударного слоя, т. е. расстояние между ударной волной и телом в окрестности плоскости ? = 90°, уменьшается, это ведет к уменьшению в указанной области размеров ячеек разностной сетки и, как следствие, шага интегрирования. Данное обстоятельство при фиксированной сетке увеличивает время счета. Поэтому при проведении дальнейших вычислений целесообразно вблизи плоскости 9 = 90° производить объединение ячеек разностной сетки в направлении от поверхности тела к ударной волне, что в случае применения использованного в данной работе метода не вызывает каких-либо трудностей.
Так как расчеты в работе [6] и в данной работе проводились на машинах, существенно различающихся по быстродействию (в работе [6] на IBM 360/67, а в настоящей работе на БЭСМ-ЗМ), то непосредственное сравнение времени счета не имеет смысла. Тем не менее сопоставление эффективности методов все же возможно. При численном интегрировании системы (1), выполненном в настоящей работе, длина шага интегрирования Дг была переменной. Она имеет наименьшее значение Дг = 0,00035 при малых л: = 0,5, затем монотонно увеличивается до значения Дг = 0,0093 при х = 0,7. Далее из-за сгущения разностной сетки в направлении от поверхности тела до головной ударной волны (при 90°) величина Дг уменьшается и при;с=1,0 становится равной 0,0013. Время расчета одной точки разностной сетки на ЭЦВМ БЭСМ-ЗМ равнялось 0,005 мин. Весь расчет потребовал около 16 ч машинного времени, что указывает на целесообразность использования для подобных расчетов машин с большим быстродействием (например, БЭСМ-6). Аналогичные данные, приведенные в работе [6], где, как уже отмечалось, расчеты проводились на ЭЦВМ IBM 360/67, следующие: средняя длина шага интегрирования для метода характеристик 0,0032 и для метода сквозного счета 0,0018; время расчета одной точки 4ХЮ-4 мин (метод характеристик) и 1,2ХЮ~4 мин (метод сквозного счета). Учитывая различие в быстродействии ЭЦВМ БЭСМ-ЗМ и IBM 360/67 можно сделать вывод, что метод работ ,[ 1 ] — [3], использованный в настоящей работе, не менее эффективен, чем методы расчета, использованные в работе [6].
Авторы признательны А. Н. Крайко за полезные советы и внимание к работе.
1. И в а но в М. Я., К р а й к о А. Н., Михайлов Н. В. Метод сквозного счета для двумерных и пространственных сверхзвуковых течений. „Журн. вычислит, матем. и матем. физ.“, т. 12, N» 2, 1972.
2. Иванов М. Я., К р а й к о А. Н. Метод сквозного счета для двумерных и пространственных сверхзвуковых течений. „Журн. вы-числ. матем. и матем. физ.“, т. 12, № 3, 1972.
3. И в а н о в М. Я., К р а й к о А. Н. К расчету сверхзвукового обтекания конических тел. „Журн. вычислит, матем. и матем. физ.“, т. 13, № 6, 1973.
4. Бабенко К. И., Воскресенский Г. П., Любимов А. Н., Русанов В. В. Пространственное обтекание гладких тел идеальным газом. М., „Наука“, 1964.
5. Kutler P., Lomax Н. Shock-capturing finitedifference approach to supersonic flows. J. Spacecraft and Rockets, 1971, v. 8, No 12.
6. R a k i с h J. V., К u 11 e r P. Comparison of characteristics and shock-capturing methods with application to the space shuttle vehicle. A1AA Paper, No 72-191, 1972.
7. Kutler P., Lomax H., Warming R. F. Computation of space shuttle flow fields using noncentered finite-difference scht-mes. A1AA Paper, No 72-193, 1972.
Рукопись поступила 18/Х 1972
