________УЧЕНЫЕ ЗАПИСКИ Ц А Г И
Т о м XX 198 9
№ 6
УДК 629.7.015.3.036 :533.697.2
РАСЧЕТ ПРОСТРАНСТВЕННОГО ОБТЕКАНИЯ ПЛОСКОГО СВЕРХЗВУКОВОГО ВОЗДУХОЗАБОРНИКА ПРИ НАЛИЧИИ УГЛОВ АТАКИ И СКОЛЬЖЕНИЯ
С. М. Босняков, С. В. Михайлов, Н. С. Яцкевич
В рамках предложенного в работах [1—3] алгоритма проведен расчет пространственного обтекания плоского сверхзвукового воздухозаборника с толстыми боковыми щеками при наличии углов атаки и скольжения. Учтено вытесняющее действие пограничного слоя. Головной скачок уплотнения выделен явным образом. Полученные результаты сопоставлены с результатами расчета методом [4—6] и экспериментальными данными.
Расчет пространственного обтекания плоского сверхзвукового воздухозаборника достаточно важен с точки зрения практических приложений. Так, при проектировании летательного аппарата процесс разбивается на два этапа: выбор места установки воздухозаборника на планере летательного аппарата и выбор воздухозаборника. В настоящее время выбор воздухозаборника осуществляется на основании экспериментальных и расчетных исследований. Расчет воздухозаборника выполняется методом [4—6], который в дальнейшем будет называться методом I. В основу метода положена конечно-разностная схема [7] повышенного порядка аппроксимации, реализованная на прямоугольной декартовой сетке, не связанной с поверхностью тела и особенностями течения [8]. У поверхности тела образуются дробные ячейки нестандартной формы, в которых разностная схема записывается специальным образом. Скачки уплотнения в решении не выделяются. Другой подход к решению задачи пространственного обтекания воздухозаборника предложен в работе [9]. Он предполагает использование конечноразностной схемы С. К- Годунова [10], построение адаптированной расчетной сетки и выделение головного скачка уплотнения. Однако практического применения указанный подход не получил. В данной работе развиваются идеи второго подхода на базе комплекса программ, описанного в работах [1—3]. В работе [1] уже приведен пример расчета плоского воздухозаборника с выделением трех скачков уплотнения от клиньев торможения, который демонстрирует возможности метода. Но программа, реализующая многоскачковый вариант, не предусматривает возможности пересечения скачков. Поэтому время счета значительно возрастает при сближении скачков друг с другом, что характерно для воздухозаборника.
В настоящей работе рассмотрены методические вопросы расчета пространственного обтекания плоского воздухозаборника в рамках схемы первого порядка аппроксимации с выделением головного скачка уплотнения и «размазыванием» внутренних разрывов. Это исследование проведено путем сопоставления двух подходов — первого и второго.
1. Решается задача пространственного обтекания плоского сверхзвукового воздухозаборника с целью определения его характеристик. Воздухозаборник расположен в декартовой системе координат таким образом, что его передняя кромка совпадает с плоскостью Х=0, плоскость симметрии с плоскостью 2=0, а строительная горизонталь— с плоскостью У=const. Клинья торможения и боковые щеки развернуты вверх. Положительным углом атаки считается угол, при котором увеличивается степень сжатия потока на первой ступени клина торможения (рис. 1). Все линейные размеры за-
дачи отнесены к полуширине клина торможения воздухозаборника в сечении Х=0. Расчет проводится при условии, что во всей возмущенной области выполняется условие М*> 1, где Мх — проекция вектора числа М на ось ОХ. В результате расчета определяются значения коэффициентов: расхода /, сопротивления по жидкому контуру сх ж, восстановления полного давления V, а также средние значения числа М, относительного статического давления р/р<*>, углов скоса аск и рск в плоскости Х=Хвх входа в канал воздухозаборника. Углы скоса определены как
I V \ /№\ -
аск = агс1д1 —I и Эск=ап^ 1“1 где и, V, — проекции вектора скорости V
на оси декартовой системы координат.
2. Метод 1 реализует конечно-разностную схему повышенного порядка аппроксимации [7]. Расчетная сетка строится в декартовой системе координат таким образом, что все пространство оказывается заполненным ячейками двух типов — стандартными и нестандартными [8]. Стандартные ячейки — это параллелепипеды, не пересекающиеся с поверхностью тела. Нестандартные ячейки являются всевозможными усеченными пирамидами, образующимися вследствие пересечения прямоугольной сетки с поверхностью тела. Объемы нестандартных ячеек дополняются до объема стандартных ячеек и образующиеся ячейки называются эквивалентными.
Размеры расчетной области определяются ОЗУ (оперативное запоминающее устройство) используемой ЭВМ, так как расчетная сетка создается заранее и не связана с физическими возмущениями. В расчете участвуют только те ячейки, которые попадают в возмущенную область, так как ее граница приблизительно отслеживается. В тех случаях, когда возмущенная область превышает размеры расчетной, на границах реализуется алгоритм экстраполяции решения, который в ряде случаев позволяет продолжить расчет. Ячейки, попадающие внутрь тела, в расчете не участвуют, но место в ОЗУ занимают. Использование декартовой системы координат практически обеспечивает консервативность метода, которая нарушается только в эквивалентных ячейках. Граничные условия в программе [8] выполнены в рамках схемы С. К. Годунова, т. е. с первым порядком аппроксимации. В работе [4] этот недостаток исправлен путем использования градиентов параметров, взятых из одностороннего шаблона со стороны потока, см. [7].
В качестве начальных данных при расчете воздухозаборника в сечении Х=0 задаются параметры невозмущенного потока. Специального установления течения в областях коничного течения не производится. Боковые щеки воздухозаборника предполагаются тонкими, а углы скольжения равными нулю. Это является следствием ограниченности ресурса ОЗУ ЭВМ. Расчетная сетка и форма воздухозаборника в сечении Х=сопзг представлены на рис. 1 слева. Границы сетки совпадают с границами расчетной области. Угловые точки контура воздухозаборника с узлами сетки в общем случае не связаны. Число ячеек, попадающих в плоскость входа воздухозаборника, зависит от нескольких факторов; один из существенных — отношение площади входа воздухозаборника /ч к площади его миделя Г0. Площадь входа на рис. 1 заштрихована. В ряде случаев «с/-’о, что приводит к существенным трудностям в использовании метода 1. Для надежной аппроксимации решения по входу воздухозаборника требуется расположить достаточное количество ячеек. Их число оценивается величиной порядка 80- 120. Это накладывает ограничение на размер ячеек, определяемый шагом сетки ку и /г* (см. рис. 1).
Значительное количество «маленьких» ячеек попадает «внутрь» клина торможения воздухозаборника и из счета выбывает. На практике встречаются ситуации, когда число «бесполезных» ячеек в 8—11 раз превышает число «полезных». Уменьшение шага расчетной сетки приводит также к сокращению размеров расчетной области, определяемых как Hy=hyXm и Hz=hzXm (см. рис. 1), где тип число ячеек в вертикальном и горизонтальном направлениях соответственно. При расчете достаточно длинных воздухозаборников физические возмущения достигают границ и в ряде случаев дают отраженные волны, которые попадают на вход воздухозаборника. Попытка расчета воздухозборников с толстыми щеками и при наличии углов скольжения [5 еще больше осложняет указанные проблемы, так как относительней1 объем нерационально используемой памяти продолжает возрастать. Это неприемлемо в условиях эксплуатации ЭВМ средней мощности.
3. Метод 2 [1—3] реализует конечно-разностную схему первого порядка аппроксимации С. К. Годунова [10]. Расчетная область разбивается на подобласти (см. рис. 1 справа). Подобластью называется криволинейный четырёхугольник, границами которого являются либо физические возмущения, либо поверхности непротекания, либо эйлеровы поверхности. Эйлеровой считается поверхность, на которой сеточная функция теряет свойство регулярности, а решение остается непрерывным с точностью метода [10]. В каждой подобласти строится регулярная сетка путем прямого и обратного преобразования, переводящего подобласть в элементарный прямоугольник. Разностная схема расписывается на криволинейной сетке в декартовой системе координат консервативным образом. При образовании в процессе счета мелких ячеек они автоматически объединяются в крупные по заранее заданному критерию. Процесс объединения строго консервативен. Головной скачок уплотнения выделяется явным образом с использованием принципа [9].
На рис. 1 (справа) представлено сечение воздухозаборника плоскостью Х=const. Расчетная область совпадает с возмущенной. Границей возмущенной области может быть как скачок уплотнения, так и поверхность слабого разрыва. Расчет проведен с использованием четырех подобластей. Эйлеровы поверхности (штриховые линии на рис. 1) привязаны к угловым точкам контура, что облегчает их выделение в процессе расчета. Одна из подобластей совпадает с плоскостью входа воздухозаборника. Расчет начинается в сечении X—Dx, где Dx = 0,1. В начальном сечении приближенно рассчитываются положение головного скачка уплотнения и поле течения за ним. В приближенном расчете используются двумерные соотношения и гипотеза плоских сечений. Коническое течение на начальном участке воздухозаборника специальным образом не устанавливается, но контроль за установлением ведется. Как правило, при Х=2 точность установления характеризуется величиной 4-f-6% для полей р и М. При расчете обтекания воздухозаборника с углом скольжения р число подобластей удваивается, причем в плоскости симметрии располагается эйлерова поверхность.
4. Оба рассматриваемых метода предполагают учет вытесняющего действия пограничного слоя, который проводится с использованием соотношений работ [11] и [12]. Расчетный шаг разбивается на два этапа: вначале определение интегральных толщин пограничного слоя и поправка контура клина торможения и боковых щек на величину 6*, а затем перестроение сетки и расчет параметров в ее ячейках. На рис. 1 изображены как основной контур воздухозаборника, так и исправленный. При расчете поправок используется гипотеза плоских сечений, причем плоскости выбираются с учетом локальных направлений линий тока у поверхности тела. Глобальные итерации с целью учета вязко-невязкого взаимодействия не производятся вследствие малого значения получаемой поправки по сравнению с точностью используемого метода. Возможные отрывы пограничного слоя не учитываются. В углах применяется принцип суперпозиции решений, полученных для каждой поверхности отдельно. Граница перехода ламинарного пограничного слоя в турбулентный определяется по эмпирической зависимости [13] критических значений числа Re, полученной для плоской пластины. Несмотря на приближенный характер описанного подхода, в работе [6] показана его хорошая работоспособность, что делает его приемлемым для прикладных исследований.
5. Тестирование программы расчета воздухозаборника методом 2 проведено путем сопоставления результатов расчета с соответствующими результатами, полученными методом 1. Для примера выбран воздухозаборник, описанный в работах [4—6]. Его основные параметры: расчетное число Мр=4; две ступени имеют углы наклона Si= 12°, 62 = 22,8°; отношение высоты к ширине h/b= 1,5; боковая щека тонкая и срезана по второму скачку уплотнения на расчетном режиме. На рис. 2 приведено поле статических давлений, построенное в сечении Х= 1,255. Рисунок подготовлен путем монтажа. На его левой части сплошными линиями нанесены изолинии из работы [5], а точками — соответствующие изолинии, построенные в данной работе (метод 2). Точками представлено также положение головного скачка уплотнения. В выбранном сечении Х= 1,255 первый и второй скачки уплотнения вышли за пределы кромок боковых щек, а третий скачок уплотнения (см. [5]), образующийся при обтекании боковой щеки скошенным потоком, еще не дошел до плоскости симметрии. Расчет ме-
N = 2,50 , «.=¿7, р=0, Х-1, 255
тодом 1 проведен на сетке, содержащей 3600 ячеек (от=60, я = 60), методом 2 — 4096 ячеек (от=32, п=32 в каждой подобласти); Изолинии достаточно хорошо совпадают друг с другом. Наибольшее отличие наблюдается в окрестности скачка уплотнения, который «размывается» в методе 1. На правой стороне рис. 2 сопоставлены решения, полученные методом 2 на разных сетках: сплошная линия — 4096 ячеек (т=32, п=32 в каждой подобласти), крестики— 1024 ячейки (от =16, п= 16 в каждой подобласти). Построение изолиний организовано следующим образом: во всем расчетном поле находятся минимальное и максимальное значения параметра, их разность делится на целое число, что определяет шаг б для построения. Изолинии строятся графопостроителем, начиная с минимального значения. Совпадение минимальных значений параметров и величин шага, рассчитанных на разных сетках, с точностью порядка 0,001 указывает на достаточно высокое качество расчета уже на крупной сетке. Заметно увеличивается лишь «размывание» второго скачка уплотнения от клина торможения (на рис. 2 нанесены только внешняя и внутренняя границы области «размывания» скачка на сетке в 1024 ячейки). Область сильного градиента параметров, соответствующая внутреннему скачку уплотнения (третьему), остается, практически неизменной. Возможно, это следствие ошибочности утверждения работы [5] о наличии третьего скачка уплотнения. В условиях неравномерного потока, набегающего на боковую щеку, может сформироваться интенсивная волна сжатия, размеры которой и определяет область сгущения изолиний. На это указывает факт соответствия с точностью порядка 10% размеров области «размывания» волны (скачка), рассчитанной обоими методами. (Скачок уплотнения, полученный методом второго порядка, было «размыт» слабее, чем полученный методом первого порядка.)
На рис. 3 сопоставлены значения коэффициентов расхода / и сопротивления по жидкому контуру сх ж, рассчитанные методами 1 и 2 как с учетом вытесняющего действия пограничного слоя (темные кружки и штриховые линии), так и без его учета (звездочки и сплошные линии). Приведены также экспериментальные значения [4], полученные в сверхзвуковой аэродинамической трубе (светлые кружки). Видно, что результаты расчета практически точно соответствуют друг другу. Отличие наблюдается только при М=4. Вследствие «размытия» скачков уплотнения методом 1 на рас-
четном режиме, как отмечено в работе [4], точность счета снижается на 3%. Ошибка устраняется линейной экстраполяцией на нулевой шаг расчетной сетки (экстраполяция по Ричардсону). Точка (звездочка на рис. 3 при М=4), взятая из работы [4], получена при помощи указанной экстраполяции. Сплошные кривые (метод 2) рассчитаны на сетке в 1024 ячейки без каких-либо поправок. Дополнительное исследование показывает, что для определения характеристик воздухозаборника не требуется задания подробной сетки. На рис. 3 (справа) приведена зависимость ошибки расчета (бг=-------двуМ , Удвум = 0,838 при М=4)
^двум
от числа ячеек расчетной сетки, из которой видно, что точность расчета заметно падает только при использовании очень разреженных сеток. Известно, что из трех характеристик воздухозаборника (}, с* ж, \) коэффициент восстановления полного давления V имеет наименьшую точность. Это замечание позволяет экономить расчетное время при проведении массовых расчетов. Так, в большинстве случаев время расчета одного варианта при р = 0 не превышает 10 минут. Это не означает отсутствия проблемы «размытых» скачков уплотнения при расчете обтекания плоских воздухозаборников. Решение этой проблемы связано с повышением порядка аппроксимации разностной схемы и разработкой надежного алгоритма пересечения выделенных скачков уплотнения.
Значения / и сх ж (штриховые линии на рис. 3 соответствуют методу 2, темные кружки—работе [6]), рассчитанные с учетом вытесняющего действия пограничного слоя, находятся в хорошем соответствии с экспериментальными данными (светлые кружки на рис. 3). Этот вывод подтвержден опытом многочисленных расчетов, проведенных авторами в различное время. Заметное расхождение результатов расчета и эксперимента возникает только на режимах течения в воздухозаборнике с сильными отрывами пограничного слоя. В настоящее время отсутствуют надежные методы расчета отрывных течений, возникающих на внутренних поверхностях воздухозаборников, что делает трудноразрешимой задачу расчета воздухозаборника с сильным сжатием потока или с прикрытым дросселем, когда ударные волны вызывают интенсивные отрывы потока. В частности, открытым остается вопрос определения границы помпажа, имеющий важное прикладное значение. Практически, рассматриваемые методы 1 и
2 надежно дают только одну (угловую) точку на дроссельной характеристике воздухозаборника, соответствующую максимальным значениям [ и V. Расчет этой точки важен, так как позволяет провести приблизительное согласование воздухозаборника с двигателем. Окончательное согласование может быть проведено только на основании полномасштабного эксперимента.
Реальный воздухозаборник имеет толстые боковые щеки, которые в ряде случаев оказывают влияние на его характеристики. Например, на расчетном режиме обтекание толстой щеки может привести к искажению головного скачка уплотнения, его отходу от кромок щек и, как следствие, перепуску части расхода воздуха. Чаще всего воздухозаборник работает на режимах с углами атаки и скольжения. Рассмотрим одну из типичных задач — определение характеристик воздухозаборника, установленного на фюзеляже самолета. Возможное решение состоит в определении средних параметров потока перед воздухозаборником и последующий расчет изолированного воздухозаборника. В большинстве случаев воздухозаборник у фюзеляжа обтекается потоком с углами атаки и скольжения.
На рис. 4 приведен пример расчета изолированного воздухозаборника с толстыми боковыми щеками. Геометрия воздухозаборника полностью совпадает с геометрией, описанной выше. Боковые щеки имеют толщину, равную бщ = 0,2 и срезаны по первому скачку уплотнения на расчетном режиме [4]. Кроме боковой щеки выполнена путем снятия острой фаски постоянной высоты Лф = 0,6, что полностью соответствует экспериментальной модели, испытанной в работах [4—6]. На рис. 4 представлено поле линий равных значений статического давления, построенных в сечении плоскости входа в воздухозаборник.
Число М набегающего потока равно М = 3, угол атаки а=5°, угол скольжения ¡5=10°. В рассматриваемом сечении второй («размытый») скачок уплотнения пересекается с первым (выделенным). В результате взаимодействия образуется сильный ска-
а = 0
М-3,0 , а=.Г; /}=10х = 6,6&
/
.V
- 0,80
0,23
о,го
0,19
9,21
015
О
Г 80 12° $
Рис. 4
Рис. 5
чок уплотнения, который выделен как линия на рис. 4, и слабый, который решением не улавливается.
Боковые щеки воздухозаборника обтекаются существенно асимметричным образом. С внешней стороны наветренной боковой щеки (справа) создается довольно интенсивное сжатие. Это препятствует боковому вытеканию потока из области высокого давления за головным скачком уплотнения в указанном направлении. Основное вытекание происходит в сторону подветренной поверхности. С внутренней стороны боковые щеки действуют на поток противоположным образом. Левая щека становится наветренной, а правая — подветренной. Волна сжатия от внутренне-наветренной боковой щеки формируется вследствие наличия угла скоса потока (также, как и в примере на рис. 2) и взаимодействует с волной разрежения от противоположной щеки. Область максимального давления р/р„О = 6,0-н6,1 смещается в сторону от середины воздухозаборника к внутренне-подветренной боковой щеке.
Угол скольжения оказывает существенное влияние на характеристики воздухозаборника. На рис. 5 приведены зависимости значений коэффициентов [, с* ж, V от угла (5 при а=5°, М=3. Расчет проведен как с учетом пограничного слоя (штриховые линии), так и без учета (сплошные линии). С ростом угла Р коэффициент расхода сначала возрастает приблизительно на 1,5% при Р«8-;-Ш0, а потом уменьшается. Рост /' обусловлен уменьшением бокового вытекания потока со стороны наветренной боковой щеки и дополнительным сжатием потока в образующейся внутренней волне. Последующее уменьшение { (при Р> 10°) вызвано уменьшением видимого со стороны потока нроходного сечения входа воздухозаборника. Действительно, у воздухозаборника, поставленного перпендикулярно потоку (Р=90°), коэффициент расхода ^О. Коэффициент сопротивления по жидкой линии тока с увеличением р монотонно возрастает. Это обусловлено смещением забираемой трубки (по линиям тока) в область высокого давления у наветренной боковой щеки и, как следствие, ростом интеграла избыточного давления по ее поверхности. Коэффициент восстановления полного давления V с ростом р уменьшается. При умеренных значениях Р<8° уменьшение обусловлено физическими факторами (волной сжатия), но с ростом значительную роль начинает играть энтропийная ошибка, возникающая в процессе счета. Дополнительные оценки показывают, что при Р=0 точность расчета V вполне приемлема и ошибка равна приблизительно 1—3% (см. рис. 3), при Р«10° ошибка может возрастать до 10—15%. Снижение точности расчета обусловлено наличием нтенсивной волны разрежения, образующейся при обтекании внутренне-подветренной поверхности боковой щеки. Возникающая при расчете волны энтропийная ошибка значительно занижает значение V. Уменьшение шага расчетной сетки не позволяет кардинально исправить ситуацию. Повышение точности расчета V на 2—3% достигается путем восьми — десятикратного увеличения времени счета.
Так как на практике углы скольжения потока в местах установки воздухозаборников редко превышают 4°—5°, для расчета обтекания плоских воздухозаборных устройств целесообразно применегие метода 2.
ЛИТЕРАТУРА
1. Босняков С. М., Михайлов С. В. Программа расчета трехмерных сверхзвуковых течений идеального газа. Труды IX конференции молодых ученых ФАЛТ, Долгопрудный. МФТИ, 1984 (Рукопись депонирована в ВИНИТИ 13.09.84 № 6242—84 дсп).
2. Босняков С. М., Карпов Е. В., Михайлов С. В. Расчет сверхзвукового обтекания крыла и фюзеляжа с выделением скачка уплотнения от крыла. — Ученые записки ЦАГИ, 1988, т. 19, № 5.
3. Босняков С. М., Михайлов С. В., Коваленко В. В., Рем ее в Н. X. Численное решение задачи обтекания трапециевидного клина сверхзвуковым потоком идеального газа. — Ученые записки ЦАГИ, 1989, т. 20, № 1.
4. Босняков С. М., Р е м е е в Н. X. Исследование пространственного обтекания плоского воздухозаборника с боковыми щеками сверхзвуковым потоком газа. — Ученые записки ЦАГИ, 1980, т. 11, № 5.
5. Босняков С. М., Быкова С. А., Ремеев Н. X. Исследование пространственного обтекания и аэродинамических характеристик плоских воздухозаборников с различной формой входа и размерами боковых щек. — Ученые записки ЦАГИ, 1983, т. 14, № 3.
6. К а р п о в Е. В., Ремеев Н. X. Расчет пространственного обтекания плоского сверхзвукового воздухозаборника с учетом вытесняющего действия пограничного слоя. — Труды ЦАГИ, вып. 2274, 1985.
7. К о л г а н В: П. Применение принципа минимальных значений производной к построению конечно-разностных схем для расчета разрывных решений газовой динамики.— Ученые записки ЦАГИ, 1972, т. 3, № 6.
8. Косых А. П., М и н а й л о с А. Н. Расчет сверхзвукового невязкого течения у пирамидального тела, моделирующего дельтовидный летательный аппарат. — Изв. АН СССР, МЖ1, 1975, № 3.
9. Макаров В. Е. К выделению поверхностей разрывов при численном расчете сверхзвуковых конических течений. — Ж. вычисл. матем. и матем. физ., 1982, т. 22, № 5.
10. Иванов М. Я., Крайко А. Н., Михайлов Н. В. Метод сквозного счета для двумерных и пространственных сверхзвуковых течений,—Ж. вычисл. матем. и матем. физ., 1972, т. 12, № 2.
11. Репик Е. У. Приближенный расчет турбулентного пограничного слоя в сжимаемой жидкости при наличии градиента давления.—Технические отчеты ЦАГИ, вып. 167, 1960.
12. Старухин В. П., Тарышкин А. Г. Экспериментальное исследование турбулентного пограничного слоя в сверхзвуковом потоке на плоских тормозящих поверхностях с изломами образующих.—Ученые записки ЦАГИ, 1975, т. 6, № 4.
13. Босняков С. М., Ремеев Н. X. Расчетное и экспериментальное исследование интегральных параметров пограничного слоя на ступенчатых клиньях плоских сверхзвуковых воздухозаборников.—Труды ЦАГИ, 1981, вып. 2098.
Рукопись поступила 1/ХП 1988