CC BY
40
Том XVI
УЧЕНЫЕ ЗАПИСКИ Ц А Г И
198 5
№ 3
УДК 629.7.015.3
РЕШЕНИЕ ПАРАБОЛИЗИРОВАННЫХ УРАВНЕНИЙ НАВЬЕ—СТОКСА ПРИМЕНИТЕЛЬНО К ЗАДАЧЕ О ВЗАИМОДЕЙСТВИИ ПОГРАНИЧНОГО СЛОЯ СО СВЕРХЗВУКОВЫМ потоком
О. М. Колесников
Изложен маршевый метод численного интегрирования параболизи-рованных уравнений Навье—Стокеа для исследования при больших числах Ие безотрывного взаимодействия двумерного пограничного слоя со сверхзвуковым потоком. В качестве примеров рассмотрены задачи о гиперзвуковом обтекании заостренной пластины в режимах сильного и слабого взаимодействия, а также об интенсивном отсосе турбулентного пограничного слоя с поверхности пластины.
Расчеты по предложенному методу сравниваются с известными результатами теории сильного взаимодействия и результатами решения полной системы уравнений Навье—Стокса.
В последнее время для теоретического исследования сверхзвуковых вязких течений при больших числах Рейнольдса широкое распространение получили численные методы, основанные на решении пара-болизированных уравнений Навье—Стокса [1—8]. Эта система уравнений справедлива как во внешнем потоке, так и в пограничном слое. Поэтому ее удобно привлекать для решения тех задач, где неприменимо разбиение расчетной области на область невязкого потока и пограничного слоя. Упрощение основано на исключении из полных уравнений Навье—Стокса вязких членов, имеющих порядок 0(1/Яе) и высший, а также тех или иных вязких членов порядка О (1/1/Не) [5, 7]. При сверхзвуковых скоростях система упрощенных уравнений относится к смешанному гиперболическо-параболическому типу, т. е. для ее решения можно разработать очень эффективные с точки зрения реализации на вычислительной машине маршевые методы. Однако в дозвуковой части пограничного слоя система уравнений эллиптична и постановка задачи Коши принадлежит к некорректным задачам. Успех применения маршевых методов возможен только в том случае, если удается разработать подходящую процедуру регуляризации.
В настоящей работе излагается маршевый метод решения параболи-зированных уравнений Навье—Стокса применительно к взаимодействию двумерного пограничного слоя со сверхзвуковым потоком, развитый на основе метода, предложенного ранее в работе [6]. Рассмотрены две за-
дачи: гиперзвуковое обтекание заостренной теплоизолированной пластины и интенсивный отсос турбулентного пограничного слоя с поверхности пластины.
В первой задаче анализируется величина возможных ошибок, вызванных допущениями, лежащими в основе используемой процедуры регуляризации, а во второй — способность метода «проходить» через области с большими продольными градиентами давления (начало и конец зоны отсоса). Применение любых маршевых методов в этих областях всегда сопряжено с определенными затруднениями, что является следствием некоторого несоответствия математической модели физике явления: математическая модель не допускает возможности распространения возмущений по дозвуковому пристеночному слою вверх по потоку.
1. В двумерной постановке используемая система безразмерных параболизированных уравнений Навье—Стокса записывается в следующем виде:
дри . дру_р дри2 . дрьи _ др . 1 д / ди
дриу , дри'2_ др , 4 1 д [ ду\
!А1 .77 >
дх ду ду 3 ду \ ду ]
^ + ^ = + _1).ма.
дх ду \ дх дх 1 Ке \ду /
Не ду \ ду
где + ргИе; ¡а2 = [а/Ргл + реИе/Ргт; Ргл, Ргт — числа Прандтля
(ламинарное и турбулентное); х — отношение удельных теплоемко-стей в набегающем потоке.
Компоненты скорости и, плотность р, температура Г, коэффициент ламинарной вязкости ^ отнесены к соответствующим значениям для невозмущенного потока, а коэффициент турбулентной вязкости е, статическое давление р — соответственно к и рс» н^ .
В уравнении поперечной компоненты количества движения внепо-рядковый вязкий член сохранен для упрощения вычислительного алгоритма и повышения устойчивости конечно-разностной схемы.
На стенке помимо условия 4^- = 0 ставятся либо условия прилипания и — V = 0, либо и = О, v = vw. В конечно-разностной постановке используется еще одно условие:
^ = 0, ду
являющееся упрощенным вариантом уравнения сохранения поперечной компоненты количества движения, записанного с учетом условий прилипания:
др__3_ д / ^ ду
ду 4 ду \Ие ду,
Задача о гиперзвуковом обтекании заостренной пластины решалась с привлечением обоих уравнений. Результаты практически не отличались друг от друга, что является лишним подтверждением внепоряд-
4 д [ ц, ди)
ковости вязкого члена--— —
3 ду \Ие ду
При решении задачи об отсосе используется первое условие:
^ = 0.
ду
На внешней границе ставятся граничные условия двух типов: без выделения головного скачка уплотнения и с выделением. В первом случае предполагается, что граница располагается в невязкой области, где выполняются условия «простой» волны, т. е. вдоль характеристики, выпущенной из последнего узла назад, плотность, температура и компоненты скорости сохраняются постоянными. Внешняя граница располагается вблизи границы пограничного слоя и, в отличие от условий невозмущенного потока, нет никакой необходимости следить за тем, чтобы возникающие скачки уплотнения всегда находились внутри области интегрирования. Однако, если головной скачок обладает достаточной интенсивностью, его полезно выделить с помощью постановки на нем граничных условий Рэнкина—Гюгонио. После введения новой независимой переменной г\=у/5, где 5 — расстояние от поверхности пластины до головного скачка уплотнения, они записываются в следующем виде:
(1— + ®* = 0,
дя\2 I д& \2 . Г, . /¿5\2
1
хМ^
=0,
0, ----!-г = 0,
1 дх 2 (*-!>*£ 2 (х-1)Л£
где индекс е обозначает условия непосредственно за головным скачком уплотнения.
В работе [7] положение ударной волны 5 определяется в результате итерационного процесса, который прекращается, когда расход за головным скачком совпадает с расходом, вычисленным по параметрам набегающего потока. Такой способ прост, но не экономичен. Здесь использован более сложный, зато не требующий итераций способ, предложенный в работе [2]. Исходная система уравнений дополняется формальным уравнением:
-^ = 0,
а соотношения Рэнкина—Гюгонио — уравнением расхода, записанным для узлов, примыкающих к внутренней стороне скачка.
Численная процедура решения исходной системы уравнений, как и в работе [6], основывается на ее квазилинеаризации с помощью метода Ньютона и последующем неявном численном интегрировании методом векторной прогонки. Для аппроксимации конвективных членов используются центральные разности. В итоге схема имеет в поперечном направлении второй порядок точности, а в продольном — первый. Схема устойчива, если выполняется известное ограничение на сеточное число Рейнольдса, которое применительно к приведенной выше системе уравнений можно записать в виде:
рг Ду
р-
<2.
Это ограничение в работе [6] оказывается необременительным: V относительно невелико, а турбулентная вязкость значительна. Тем не
менее оно приводит к необходимости сохранить внепорядковый вязкий член в уравнении сохранения количества движения в поперечном направлении. В приведенных ниже примерах удовлетворить этому ограничению сложно: потребовалось бы чрезмерное измельчение сетки во внешней части пограничного слоя, характеризующейся малой вязкостью. Поэтому в этих областях в настоящей работе используется схема с направленными разностями, которая уже свободна от ограничений на сеточное число Рейнольдса.
Отличия заключаются в разностном представлении конвективных членов, например:
ц(л) _ ц оо
дУ/_1
_ДУ/_1 (П-1) -1- 2 иУУ
диУп)
ду/1 '
д \(л) „(Л) _„(") д ди _ "¿+1 '_ дУг (л-1
при V > О,
при ^ <0,
где / — номер узла в ¿/-направлении, п — номер итерации (в настоящих расчетах проводилась только одна итерация). На равномерной сетке эта схема в поперечном направлении также имеет второй порядок точности.
Как уже отмечалось, задача Коши для системы параболизирован-ных уравнений Навье—Стокса в дозвуковой области является некорректной. Это обусловливается наличием в уравнении продольной компоненты количества движения члена др/дх, который описывает возмущения, распространяющиеся в дозвуковом слое вверх по потоку. Поэтому, если в этом слое исключить из исходной системы уравнений др/дх, эллиптичность устраняется, и исходную систему уравнений можно решать маршевыми методами (см., например, [2]). Конечно, такого рода упрощение исходных уравнений при определенных обстоятельствах может привести к ошибкам, поэтому в ряде работ [1—4, 8] были предприняты попытки подавить неустойчивость, вызванную некорректностью постановки задачи Коши, иными способами. Так, в работе [2] удалось решить без каких-либо допущений задачу о гиперзвуковом обтекании конуса под углом атаки. В этой задаче продольный градиент давления достаточно мал и успех был достигнут только за счет выбора размера шага в продольном направлении. Оказалось, что шаг не может быть произвольно уменьшен: с уменьшением шага уменьшается схемное демпфирование и, начиная с некоторого момента, его уже не хватает для подавления неустойчивости. Очевидно, такой подход не пригоден для решения задач, где требуется достаточно малый шаг, например, в задачах с неравновесными процессами. В работах [3, 4] неустойчивость подавлялась путем решения уравнений, записанных в нестационарной форме. В каждом новом сечении стационарное решение находилось методом установления. Недостатком такого подхода является низкая экономичность: для установления требуется несколько десятков шагов по времени.
Другой подход, получивший также широкое распространение [8], заключается в вычислении члена др/дх за пределами дозвукового слоя.
В связи с вышеизложенным необходимо отметить, что в областях с большими продольными градиентами давления подавление неустойчивости только за счет одного схемного демпфирования далеко не всегда возможно. Для плоских течений по сравнению с трехмерными соответствующие трудности возрастают. Поэтому в настоящей работе используется процедура регуляризации, основанная на исключении чле-
на др/дх в дозвуковом пристеночном слое. Возможные ошибки, вносимые этим упрощением исходной системы уравнений, можно выявить с помощью глобальных итераций. В этом случае вся область интегрирования просчитывается несколько раз, причем член др/дх в дозвуковом слое аппроксимируется явно на основании данных, полученных в предыдущей итерации.
2. Рассмотрим задачу о гиперзвуковом обтекании заостренной пластины в режимах как сильного, так и слабого взаимодействия. Непосредственно вблизи кромки пластины находится область молекулярного течения, не подчиняющаяся уравнениям механики сплошной среды. Затем идет область слившегося с пограничным слоем головного скачка уплотнения, а на стенке выполняются условия скольжения. Далее происходят выделение скачка и его отход от пограничного слоя. В этой области на скачке уплотнения выполняются соотношения Рэн-кина—Гюгонио, а на стенке — условия прилипания.
Известно, что эффектами скольжения можно пренебрегать на таких расстояниях от передней кромки, когда параметр скольжения V (У=М00Ус1^х> с■=■ {рТсо^юТ)) становится меньше 0,1. В представленном численном примере рассматривается обтекание теплоизолированной пластины потоком воздуха (з< = 1,4) при числе Маха, равном 16,35 (том же, что и в работе [9]). Число Рейнольдса, вычисленное по расстоянию х0 от передней кромки до начального сечения, равно 3,3 • 104. Константа с, входящая в зависимость коэффициента ламинарной вязкости от температуры, выбирается равной единице. При этих
венно равны 24 и 0,09, т. е. начальное сечение располагается в области сильного взаимодействия. Необходимо отметить, что начальное сечение можно было бы взять и выше по потоку в области, где пограничный слой и головной скачок уплотнения еще не отделены друг от друга. Однако там нужно было бы использовать уже другие граничные условия.
Задание начальных профилей и, р, Т, а также поперечной координаты головного скачка уплотнения осуществляется из общих соображений и с помощью данных, приведенных в работах [9—10]. Как видно из рис. 1, где показано распределение вдоль стенки статического давления и коэффициента трения, давление на стенке в начальном сечении превышает давление в набегающем потоке в 11 раз. Ниже по потоку давление падает и на расстоянии х/хо= 100 равно уже только 2. Согласование настоящих расчетов с результатами теории сильного взаимодействия [10, 11] очень хорошее, за исключением области, примыкающей к начальному сечению (х<2х0). В этой области степень согласования зависит от точности задания начальных данных. Так, например, цифрами 1 и 2 отмечены результаты, полученные при различных начальных профилях давления (рис. 2,а).
О нарастании толщины ударного слоя и его вязкой части — пограничного слоя — можно судить по рис. 1,6. Как известно, в режиме сильного взаимодействия статическое давление переменно поперек ударного слоя. Это хорошо иллюстрируется расчетами, показанными на рис. 2, а, где поперечная координата у отнесена к расстоянию в от стенки до головного скачка уплотнения в каждом сечении. Непосредственно за головным скачком уплотнения в режиме сильного взаимодействия давление в полтора — два раза превышает давление на стенке. Весь перепад давления приходится на невязкую часть ударного слоя, в то время как поперек пограничного слоя давление постоянно.
условиях параметр подобия
и параметр V соответст-
3—«Ученые записки ЦАГИ» № 3
33
На рис. 2, а приведены два варианта начальных профилей давления (кривые- 1 и 2), используемые в расчетах. Сравнивая графики в различных, сечениях, можно отметить, что первый вариант менее удачен, это и явилось причиной некоторого повышения давления на стенке
вблизи начального сечения (см. рис. 1,а). В отличие от распределения давления величина продольной компоненты скорости в невязкой части Ударного слоя постоянна, а ее профиль в пограничном слое линеен (рис. 2,6).
Известно, что чувствительность пограничного слоя к продольному градиенту давления падает по мере уменьшения температуры стенки. Соответствующее экспериментальное исследование сильного взаимодействия представлено в работе [12]. По-видимому, для холодного пограничного слоя использование описанной выше процедуры регуляризации не приведет к заметным ошибкам. Однако в рассматриваемом примере стенка теплоизолирована, поэтому были проведены расчеты с помощью глобальных итераций, что позволяет учесть продольный градиент давления в дозвуковой части пограничного слоя. Как видно из рис. 1 , а, изменения в распределении статического давления незначительны (штрихпунктирная кривая), т. е. и в этом случае вклад дозвуковой части пограничного слоя в общее изменение толщины вытеснения также невелик.
3. Схематично картина течения при отсосе пограничного слоя с поверхности пластины показана на рис. 3, а: начало отсоса сопровождается падением давления, его окончание — повышением; соответственно во внешнем потоке формируется веер волн разрежения и сжатия. Если величина отсоса достаточно велика, то поперечный градиент давления становится отличным от нуля и, как следствие, методы, основанные на
^ Отсос
а)
Рис. 3
теории пограничного слоя, малопригодными. В представленном здесь примере расчеты проводятся при тех же параметрах набегающего потока, что и в работе [13], где решена полная система уравнений Навье—
Стокса: Мта= 1,5; Г»—100.К, Р=о= 6-106 м^', р°°"°°8° =2,4-Ю4,
$о — толщина пограничного слоя непосредственно перед началом области отсоса. Область отсоса располагается на таких расстояниях от передней кромки пластины, где пограничный слой турбулентный.
Рис. 4
Для описания коэффициента турбулентной вязкости используется известная двухслойная модель Себеси—Смита [14]. Расчеты проведены при двух значениях параметра отсоса (рг* )и,/(ри)ос>: 0,006 и 0,05. Соответствующие графики по распределению давления вдоль стенки показаны на рис. 3, б. Интенсивный отсос сопровождается быстрым уменьшением толщины пограничного слоя, что в свою очередь согласно используемой модели турбулентности приводит и к быстрому уменьшению коэффициента турбулентной вязкости. Равновесная модель Себеси—Смита не пригодна для описания быстрого процесса ламинариза-ции потока. Поэтому при (ри)ю/(ри)оо = 0,05 длина области отсоса была уменьшена (как это видно из рис. 3, б) для того, чтобы коэффициент турбулентной вязкости не стал чрезмерно малым. По-видимому, в этом случае предпочтительней привлекать двухпараметрические модели турбулентности, однако вопрос о пригодности той или иной модели турбулентности уже выходит за рамки данной работы.
Как видно из рис. 3, б, в окрестности начала и конца зоны отсоса градиенты давления очень велики. В областях с большими продольными градиентами давления достоверность результатов, полученных с помощью любых маршевых методов, сомнительна, они не в состоянии учесть процесс передачи давления вверх по потоку. Однако вклад этих областей в общий интеграл давления мал и от маршевых методов достаточно потребовать, чтобы в этих областях не возникало трудностей с получением решения, а его точность ниже по потоку восстанавливалась. Первые попытки получить решение при (ри)№/(ры)оо = 0,05 по изложенной выше методике оказались неудачными. В окрестности звуковой линии возникали усиливающиеся колебания давления, приводящие к быстрому разрушению решения. Их появление было вызвано тем, что по одну сторону от звуковой линии продольный градиент давления учитывался, а по другую нет. Успех был достигнут с помощью сглаживания профиля давления в примыкающих к звуковой линии узлах расчетной сетки (как показал опыт, эта процедура сглаживания оказалась полезной и в других примерах и не только в областях с сильным продольным градиентом давления).
На рис. 4,а показано распределение давления в нескольких сечениях, располагающихся ниже по потоку за областью отсоса пограничного слоя. Как видно из этого рисунка, прекращение отсоса сопровождается образованием замыкающего скачка уплотнения. При сильном отсосе перепад давления на нем составляет 25% (рис. 4,а), при слабом— 6% (рис. 4,6). Согласование с аналогичными данными, полученными в [13] на основе решения полных уравнений Навье—Стокса, хорошее. На рис. 4, б данные [13] приведены в сечении, располагающемся в точке окончания отсоса, а полученные здесь — несколько ниже по потоку, где скачок уже отошел от стенки.
Таким образом, представленные здесь примеры свидетельствуют, что предложенный маршевый метод вполне пригоден для решения задач о безотрывном взаимодействии двумерного пограничного слоя со сверхзвуковым потоком. Этот метод достаточно точен в областях с умеренным продольным градиентом давления и способен «проходить» через области с большим градиентом без потери точности ниже по потоку.
ЛИТЕРАТУРА
1. К овен я В. М., Черный С. Г. Решение упрощенных уравнений вязкого газа маршевым методом. — В сб. Численные методы механики сплошной среды, 1979, т. 10, № 1.
2. JI у б а р д С. К., X е л л и в е л В. С. Расчет обтекания конуса под большим углом атаки. — РТК, 1974, № 8, т. 12.
3. Головачев Ю. П., Ф у р с е н к о А. А. Маршевый метод расчета течений вязкого газа. — Ж. вычислительной математики и математической физики, 1981, т. 21, № 6.
4. Н u n g С. М., ChausseeD. S. Computation of supersonic turbulent flows over an inclined ogive-cylinder—flare. — AIAA J., 1981, vol. 19, N 9.
5. Численное исследование современных задач газовой динамики. — Сб. статей./Под редакцией О. М. Белоцерковского. — М.: Наука, 1974.
6. Колесников О. М. Расчет повышения давления при горении плоской сверхзвуковой струи водорода в сверхзвуковом спутном потоке. — Ученые записки ЦАГИ, 1982, т. 13, № 6.
7. Кокошинская Н. С., Павлов Б. М., Пасконов В. М. Численное исследование сверхзвукового обтекания тел вязким газом. — М.: МГУ, 1980.
8. Schiff L. В., Steger J. L. Numerical simulation of steady supersonic viscous flow. — AIAA Paper, 1979, N 0130.
9. G a r v i n e R. W. Upstrem influence in viscous interaction problems. — The Physics of Fluids, 1968, vol. 11, N 7.
10. Тендер сон А. Гиперэвуковые вязкие течения.— В сб. Современные проблемы газовой динамики./Под ред. Лоха У.Х.Т. — М.: Мир,
1971.
11. Дорренс У. X. Гиперзвуковые течения вязкого газа. — М.: Мир, 1966.
12. Королев А. С. Экспериментальное исследование теплопередачи и распределения давления на плоской заостренной пластине в режиме сильного взаимодействия течения с пограничным слоем. — Труды ЦАГИ,
1972, вып. 1374.
13. Rose W. С. Comparison of solutions to boundary-layer and Navier—Stokes equations for the case of mass removal from a turbulent boundary layer.— AIAA Paper, Г976, N 150.
14. Cebeci Т., Smith A. M. O. Analysis of turbulent boundary layer. — Academic Press, New Jork — London, 1974.
Рукопись поступила 27/IX 1983 г.
