Том ХЫН
УЧЕНЫЕ ЗАПИСКИ ЦАГИ 2012
№ 5
УДК 533.6.011.35; 629.7.024.36
ПОПЕРЕЧНОЕ ОБТЕКАНИЕ КРУГОВОГО ЦИЛИНДРА ТРАНСЗВУКОВЫМ (М = 0.8) ПОТОКОМ ПРИ БОЛЬШИХ ЧИСЛАХ РЕЙНОЛЬДСА
В. А. БАШКИН, И. В. ЕГОРОВ, И. В. ЕЖОВ, С. В. УТЮЖНИКОВ
На основе уравнений Навье — Стокса и Рейнольдса выполнено численное моделирование поперечного обтекания кругового цилиндра с теплоизолированной поверхностью трансзвуковым потоком совершенного газа при числах Маха М = 0.8 и Рейнольдса 103 < Яе < 107.
Показано влияние числа Рейнольдса на структуру поля течения, режим обтекания и аэродинамические характеристики цилиндра.
Ключевые слова: динамика вязкого газа; трансзвуковой поток; круговой цилиндр; численное моделирование.
ВВЕДЕНИЕ
В прикладной аэродинамике часто приходится иметь дело с поперечным обтеканием кругового цилиндра потоком совершенного газа, поэтому это тело простой конфигурации интенсивно исследовалось и исследуется теоретически и экспериментально в достаточно широком диапазоне изменения определяющих параметров задачи. В частности, информацию о характеристиках цилиндра при больших числах Рейнольдса в до- и трансзвуковом диапазоне скоростей можно найти в [1 — 6].
В [1, 2] приведены результаты экспериментального исследования поперечного обтекания кругового цилиндра до- и трансзвуковым потоком воздуха. В частности, для чисел Мда = 0.8 — 1.2 приведены результаты измерений распределения по поверхности цилиндра коэффициента давления и коэффициента сопротивления трения для двух чисел Рейнольдса Яе = 0.83 -105, 2.5 -105
БАШКИН Вячеслав Антонович
доктор физикоматематических наук, главный научный сотрудник ЦАГИ
ЕГОРОВ Иван Владимирович
доктор физикоматематических наук, член-корреспондент РАН, начальник отделения ЦАГИ
ЕЖОВ Иван Валерьевич
кандидат физикоматематических наук, научный сотрудник ЦАГИ
УТЮЖНИКОВ Сергей Владимирович
профессор, университет Манчестера
(Яе0 = 1.66-105, 5-105). При меньшем числе Яе реализуется ламинарное, а при большем — ламинарно-турбулентное обтекание цилиндра. Экспериментальные данные по обтеканию кругового цилиндра при числах Мо = 0.2 — 1.4 и Яед = 0.5 • 106 — 8.7 • 106 приведены в [3], а при числах Мо = 0.4, 0.55, 0.75, 0.85 и Яед я 105 — в [4].
Работа [5] посвящена численному исследованию течения сжимаемой жидкости около кругового цилиндра при числах Мо = 0.3, 0.8, 0.95, 0.98 и Яед = 103 и 5 • 106 на основе полных нестационарных двумерных уравнений Навье — Стокса, которые интегрировались модифицированным методом Бима — Уорминга — Стегера на неравномерных сетках 121 х 60 при Мо = 0.3, 0.8 и 181 х 69 при Моо = 0.95, 0.98. При расчете турбулентного течения использовалась алгебраическая модель турбулентности Болдуина — Ломакса.
В [6] численно исследовано течение сжимаемой жидкости около кругового цилиндра при числах Мо = 0.3 и Яед = 105, 3.7 • 105, 106, 3 • 106 и 7.83 • 106 на основе полных нестационарных двухмерных уравнений Навье — Стокса, которые интегрировались тем же методом на неравномерной сетке 481 х 120. Согласно экспериментальной классификации при первом значении числа Яе реализуется докритический режим обтекания цилиндра, при втором — критический режим, при третьем и четвертом — сверхкритический режим, при пятом — транскритический режим.
Приведенные в [1 — 6] результаты позволяют провести верификацию выполненного в данной работе численного моделирования.
Для теоретического исследования поперечного обтекания кругового цилиндра сверхзвуковым потоком совершенного газа в предположении о симметрии течения было выполнено численное моделирование на основе нестационарных двухмерных уравнений Навье — Стокса (ламинарное течение) [7, 8] и на основе нестационарных двухмерных уравнений Рейнольдса (ламинарно-турбулентное течение) [9, 10] с использованием дифференциальной д-ю модели турбулентности [11]. С помощью разработанного пакета программ были проведены параметрические расчеты сверхзвукового обтекания ряда тел простой конфигурации: круговые цилиндр и конус, сфера, сверхзвуковое течение в канале и др. Обширный расчетный материал позволил установить влияние определяющих параметров задачи на структуру поля течения и поведение аэродинамических характеристик рассмотренных тел в сверхзвуковом потоке.
Поскольку при трансзвуковых скоростях обтекание цилиндра может быть нестационарным и асимметричным, то при численном моделировании пришлось отказаться от предположения
о симметрии течения и рассматривать все поле течения. С учетом сказанного был модифицирован пакет программ, что позволило приступить к исследованию трансзвукового обтекания кругового цилиндра методами вычислительной аэродинамики.
Сначала была изучена эволюция структуры поля течения и аэродинамических характеристик при фиксированных значениях чисел М и Яе (Мш = 0.8, Яе = Ух Л/уш = 105, где Ух — скорость набегающего потока, Я — радиус цилиндра, — кинематический коэффициент вязкости),
когда обтекание цилиндра является нестационарным с образованием вихревой дорожки Кармана в ближнем следе. Также была проведена верификация численного моделирования путем сопоставления результатов расчетов с экспериментальными данными, которое показало в целом хорошее согласование их между собой [12].
Затем было исследовано влияние числа Маха (0.8 < Мш < 1.3) при фиксированном значении числа Рейнольдса (Яе = 105) на структуру поля течения и аэродинамические характеристики кругового цилиндра с теплоизолированной и изотермической (температурный фактор = 0.5) поверхностями и проведена дополнительная верификация численного моделирования [13]. В частности показано, что при числах Мш < М* с обтекаемой поверхности происходит периодический сход вихрей и течение около цилиндра носит нестационарный характер. При числах Мш >М * периодический сход вихрей не происходит и общая структура поля течения близка к симметричной; при этом в ближнем следе имеется узкая область течения в окрестности плоскости симметрии, в которой движение газа нестационарно. Согласно расчетам 0.9 < М* < 0.95 для теплоизолированного и 0.8 < М * < 0.9 для изотермического цилиндра.
Целью данной работы является исследование влияния числа Рейнольдса при фиксированном значении числа Маха (М«, = 0.8) на структуру поля течения и аэродинамические характери-
стики кругового цилиндра с теплоизолированной поверхностью. Она построена по следующему плану. В первом разделе кратко обсуждаются постановка задачи, метод численного моделирования и условия расчетов. Второй раздел посвящен анализу структуры поля течения около цилиндра и режимов его обтекания. В третьем разделе рассматриваются характеристики схода вихрей с обтекаемой поверхности цилиндра, а в четвертом разделе — интегральные аэродинамические характеристики цилиндра. Краткие итоги исследования подведены в заключении.
1. ПОСТАНОВКА ЗАДАЧИ И ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ.
Постановка задачи при численном моделировании на основе уравнений Навье — Стокса ламинарного поперечного обтекания кругового цилиндра трансзвуковым потоком совершенного газа описана в [12]. Поэтому ниже приводится краткое изложение постановки задачи при численном моделировании на основе нестационарных двухмерных уравнений Рейнольдса в предположении Буссинеска о рейнольдсовых напряжениях с использованием дифференциальной д-ю модели турбулентности [11]. Выбор этой модели обусловлен следующими причинами: 1) она протестирована ее разработчиками на ряде частных задач в широком диапазоне изменения параметров подобия; 2) большая устойчивость вычислительного процесса, в частности, при малых числах Яе.
Дифференциальные уравнения Рейнольдса. При теоретическом анализе различных режимов течения все большую роль играет численное моделирование на основе интегрирования осред-ненных по Рейнольдсу уравнений Навье — Стокса — уравнений Рейнольдса. Эти уравнения являются незамкнутыми, и для их замыкания используются различные модели турбулентности, как алгебраические, так и дифференциальные.
Осредненные по Рейнольдсу уравнения Навье — Стокса в произвольной криволинейной системе координат (^, п), где х = х(^, п), У = у(^, п) — декартовы координаты, можно записать в дивергентной форме
дО дЕ дС
дґ д£, дп
В .
(1)
Здесь О — вектор консервативных зависимых переменных; Е и С — векторы потоков в криволинейной системе координат; В — вектор источника. Векторы О, Е, С и В связаны с соответствующими векторами Ос, Ес, Сс и Вс в декартовой системе координат формулами
О = ^с, Е = ДЕс -ді + Сс ^),
С = ДЕс 5п + Сс 5п ),
В — 3Вс,
дх ду дх ду
в которых 3 — д(х,у)/д(^,п) — якобиан преобразования.
Декартовы компоненты векторов 0с, Ес, Сс и Вс для двухмерных осредненных по Рейнольдсу (методом Фавра) уравнений Навье — Стокса имеют вид:
рu Рv
р 2 і і 2 2 . рuv + т ху 0
Pu + P + “Р? + Тхх 0
рu 2 . . 2 2 .
рv Р™ + т ху РV + Р + -Р? + Т уу 0
Ос = Р(є + q 2) , Е с = тт 5 2 т рuH +—рuq + іх , С с = рл>Н + 5 рvq 2 + 1у , В с = 0
рq 3 3 Л1рю?
рю Ри? +1? РЩ +1? к2рю2
, т Ю рию + 1х р™ +ію
где р — плотность газа; и, V — декартовы компоненты вектора скорости V; р — давление; е = к - р/р + (и2 + v2)/2 — полная энергия на единицу массы; Н = к + (и2 + v2)/2 — полная энтальпия, к = срТ — статическая энтальпия; ср — удельная теплоемкость при постоянном давлении; т — тензор вязких напряжений с компонентами
. . і 2 Su
Txx = (^.т)| -dlvV - 2—
I — вектор теплового потока
, T xy =T yx т )
f
yx
Su Sv
Sy Sx
\
yy
= (.I. т )
-dlvV - 2—
,3 Sy.
\
I — — (А + А^га^Т) + тV,
ц и А — коэффициенты молекулярной вязкости и теплопроводности, Цг и А- — коэффициенты турбулентной вязкости и теплопроводности;
(
I
q
\
Pri)
grad(q), Iю=-
ь
Pr
grad(®).
2)
Вектор источника в уравнениях Рейнольдса для плоского (V = 0) и осесимметричного (V = 1) случаев имеет вид:
В = 3 ^ 0,0, ^ р + ■2divV - 2—^ ^ ,0, к1рю?, к2рю 2
В настоящей работе использована двухпараметрическая дифференциальная ?-ю модель турбулентности [11] с выражениями для турбулентной вязкости:
.т = С.f —, f =1 -expI -а
ю
Prwq
а = 0.02, С. = 0.09,
S 2 dlvV
h1 = С11(С.^^ "-------) - С12:
ю
S
S = -
3ю
Su Sv f Sv Sx Sy ІЄу
^2 = Сц(С. 2 - С23
ю
dlvV
) - С
Su Sv -------1------
Sy Sx
dlvV =
ю
Su Sv
------1------.
Sx Sy
где C11 = C12 = 0.5, C21 = 0.055 - 0.5f(q,rw,p,^), C22 = 0.833, C23 = 2.4, Pr = 2, Pr2 = 2, rw — расстояние от стенки, q и ю — скорость и частота турбулентности соответственно.
Для замыкания системы уравнений использованы: уравнение состояния совершенного газа p = pRgT /M , где Rg — универсальная газовая постоянная, М — молярный вес газа; степенная
зависимость молекулярного коэффициента вязкости от температуры ц/ц» = (T/T»)07 и условие постоянства чисел Прандтля Pr = цср/А = 0.7, PrT = ^,Tcp/AT = 0.9.
Граничные и начальные условия. При решении уравнений Рейнольдса (1) на газодинамические переменные накладываются те же граничные условия, что и для уравнений Навье — Стокса, поскольку сохраняется тип уравнений. Кроме того, должны быть заданы значения параметров турбулентности в набегающем потоке: q = qm, ю = юда.
На границе расчетной области, совпадающей с твердой поверхностью тела, ставятся граничные условия прилипания и непротекания (u = v = 0) и условие теплоизолированности ([ST / dn]w = 0) или изотермичности (Tw = const) обтекаемой поверхности. Кроме того, для уравнений, определяющих поведение параметров турбулентности, используются граничные условия на твердой поверхности: условие затухания турбулентных пульсаций (qw = 0) и условие частотной непроницаемости ([5ra/5n]w = 0).
На внешней границе расчетной области задаются условия излучения, записанные в инвариантах Римана и соответствующие расходящейся волне:
2a Se Se p Se Se
a1 =-------u-v— = const, a 2 = — = const, a3 = v—— u — = const,
y-1 Sx
2a
Sy
Sx Sy
аА =-
Y-1
де де
u-1 v— = const, а 5 = q = const, а 6 = ю = const,
Sx Sy
где а — скорость звука. При этом в каждой точке входной границы анализируются знаки собственных чисел
определяющих направление распространения возмущений относительно S, = const. При Xг- > 0 («входная граница») соответствующий инвариант на внешней границе вычисляется по значениям газодинамических переменных набегающего потока, а в случае Xг- < 0 («выходная граница») используется мягкая интерполяция вида Uk _ 2Uk _ 1 + Uk _ 2 = 0, где U — вектор инвариантов Римана.
На внутренней границе расчетной области, совпадающей с положительной осью абсцисс, ставятся условия периодичности.
В качестве начального приближения можно использовать условие однородного набегающего потока с последующим развитием поля течения в процессе решения нестационарной задачи. При проведении систематических расчетов по числам Маха и Рейнольдса в качестве начального приближения использовались ранее полученные варианты с наиболее близкими к необходимым значениями изменяющихся параметров.
Аппроксимация уравнений. Для численного решения система уравнений (1) приводится к безразмерному виду путем деления декартовых координат на радиус цилиндра R, компонентов
вектора скорости — на скорость V», давления — на удвоенный скоростной напор P»V» , времени — на характерное время пребывания жидкой частицы около тела t* = R/ V» ; остальные газодинамические переменные относятся к их значениям в набегающем потоке. Сформулированная выше начально-краевая задача решается численно на основе интегро-интерполяционного метода (метода конечного элемента). Его применение к уравнениям Рейнольдса (1) позволяет получить разностные аналоги законов сохранения.
Для монотонной разностной схемы вычисление потоков в полуцелых узлах осуществляется на основе решения задачи Римана о распаде произвольного разрыва. Математически эта задача сводится к решению нелинейной системы алгебраических уравнений. При аппроксимации конвективной составляющей векторов потоков E и G в полуцелых узлах используется монотонная схема типа Годунова [14, 15] и приближенный метод Роу [16] решения задачи Римана о распаде произвольного разрыва. Для повышения порядка аппроксимации (до второго) при интерполяции зависимых переменных на грань элементарной ячейки используется принцип минимальных производных (MUSCL) [17 — 19]. При аппроксимации диффузионной составляющей векторов потоков E и G на грани элементарной ячейки применена разностная схема типа центральных разностей второго порядка точности.
Условия расчетов. Для изучения влияния числа Re на структуру поля течения и аэродинамические характеристики кругового цилиндра с теплоизолированной поверхностью при числе Мш = 0.8 выполнено численное моделирование на основе уравнений Навье — Стокса при 103 < Re < 106 (ламинарное течение) и на основе уравнений Рейнольдса при 106 < Re < 107 (ламинарно-турбулентное течение). Расчеты проведены на ортогональной неравномерной сетке с числом узлов 201 х 401, которая строилась с помощью интегрального метода, основанного на численной реализации конформного преобразования. Сгущение узлов расчетной сетки к твердой поверхности для разрешения ламинарного и турбулентного пограничного слоя реализовано с применением квазиодномерного алгоритма, основанного на алгебраическом преобразовании. При сгущении
в окрестности твердой поверхности выделялись три зоны толщиной Re-1, 2Re-12 и 1.5Re-1/5, в каждой из которых содержалось соответственно 6, 15 и 20% от общего числа узлов в нормальном направлении. Расчетная область простиралась в горизонтальном направлении вверх и вниз по потоку на 50 калибров и в вертикальном направлении вверх и вниз от цилиндра на 100 калибров.
Численное интегрирование определяющих уравнений выполнено во временном интервале 0 < t = t*УХ /Я < 200 с постоянным шагом Дt = 0.01, при этом запоминание полей газодинамических переменных проводилось через интервал 5t = 0.1. Согласно численным экспериментам при получении стационарного решения для установления общей структуры поля течения требуется время t < 100, а для установления «тонкой» структуры — время t < 150. Таким образом, выбранный временной интервал вполне достаточен для выхода решения задачи на квазипериодический режим.
При анализе расчетного материала часто используется декартова система координат, начало которой помещено в центр цилиндра, а направление оси абсцисс совпадает с направлением скорости невозмущенного потока. На поверхности цилиндра имеются четыре характерные точки, расположенные в точках пересечения осей координат с окружностью единичного радиуса: -4(-1, 0), 5(0, 1), С(1, 0), Д0, -1).
2. СТРУКТУРА ПОЛЯ ТЕЧЕНИЯ И РЕЖИМЫ ОБТЕКАНИЯ
Информацию о режиме обтекания цилиндра и структуре поля течения обычно получают путем рассмотрения полей газодинамических переменных.
При всех рассмотренных числах Яе в качественном отношении около цилиндра реализуется однотипная нестационарная структура поля течения с образованием в ближнем следе вихревой дорожки Кармана. Изменение числа Яе влияет на количественные характеристики поля течения, т. е. на режим обтекания цилиндра.
Для того чтобы корректно оценить влияние числа Яе на поле завихренности или какой-либо другой газодинамической переменной, необходимо сопоставить поля, примерно соответствующие одной и той же стадии развития нестационарного течения. Такие базисные состояния удобно связать с экстремумами эволюционной зависимости первого нуля напряжения трения вниз по потоку от передней критической точки. В указанной точке происходит ветвление нулевой линии тока (особая точка типа седло), и с кинематической точки зрения ее можно интерпретировать как точку первичного отрыва. Мы ограничимся рассмотрением двух базисных состояний. Для них характерно наличие точки первичного отрыва на лобовой поверхности цилиндра, которая занимает предельно левое положение на верхней (0$1ь, первое базисное состояние) и нижней (0$2Ь, второе базисное состояние) сторонах цилиндра. Здесь 0 — центральный угол, отсчитываемый от передней критической точки по часовой стрелке.
Эволюция положения точек первичного отрыва на поверхности цилиндра для разных значений числа Рейнольдса показана на рис. 1. Согласно приведенным данным при всех числах Яе обтекание цилиндра носит нестационарный характер, но эта нестационарность проявляется по-разному в зависимости от числа Яе и режима течения.
При ламинарном режиме течения в зависимости от числа Яе реализуются различные режимы обтекания кругового цилиндра. При Яе — 103 эволюционные зависимости проявляют четко выраженный периодический характер поведения, поэтому этот режим обтекания часто называют регулярным (рис. 1, а). При последующем увеличении числа Яе наблюдается непрерывное понижение «качества» эволюционной зависимости. При Яе — 104 и 105 эволюционные зависимости имеют примерно постоянный период и переменную форму колебания, поэтому такие режимы обтекания называются нерегулярными (рис. 1, б, в). При Яе — 106 эволюционные зависимости осциллируют беспорядочным образом, поэтому этот режим обтекания часто называют хаотическим (рис. 1, г). При этом точка первичного отрыва в предельно левом положении располагается на лобовой поверхности цилиндра для всех расчетных вариантов, и с увеличением числа Рейнольдса глубина проникновения на лобовую поверхность сначала возрастает, а затем остается примерно постоянной. В качестве примера приведем данные для верхней стороны цилиндра 0$1ь ~ 80, 70, 70, 70° соответственно.
О турбулизации течения около цилиндра можно судить по картине поля параметра турбулентности д, а о положении начала ламинарно-турбулентного перехода в тонком пограничном слое на обтекаемой поверхности — по распределению параметра турбулентности ю, который на ней принимает максимальные значения в силу граничного условия частотной непроницаемости и критическое значение которого ю* оценивается по распределению местного коэффициента
ч ЛД Л Л, \ЛАЛЛ ДАЛЛЛ АЛЛА/1
•V V \] М \] V М М М \ $ X/ л; М л; V \1 М Л7
а /\ у\ л у\ /1 /1 /1 /1 л /\ Л /\ Л /\ У\ /\ /\ /\ А
VW1/V . УУ\У\ . . . . ЛЛЛгУ' . . . . АУУУ . . . .
О 50 100 150 г 200
Рис. 1. Эволюция положения точек первого нуля напряжения трения на поверхности цилиндра:
а — Яе — 103; б — Яе — 104; в — Яе — 105; г — Яе — 106 (ламинарное обтекание); д — Яе — 106; е — Яе — 107 (ламинарно-турбулентное обтекание)
3.35Е-01
а)
б)
Рис. 2. Поля параметра турбулентности д и вектора скорости около цилиндра (а) и распределение параметра турбулентности ю вдоль обтекаемого контура (б) при Яе = 106 в момент времени Г = 200
сопротивления трения, построенного в плоскости переменных С ^ л/Йё - 5 . Тогда режим течения в пограничном слое ламинарный при ю <и, и турбулентный при ю >и* .
Соответствующие данные для Яе = 106 приведены на рис. 2; при Яе = 107 в качественном отношении имеем схожие результаты. Различия между ними носят в основном количественный характер, так, например, максимальное значение параметра турбулентности ю « 110 в первом случае и ю « 1500 — во втором.
Согласно расчетным данным для обоих чисел Яe ламинарно-турбулентный переход наблюдается на лобовой поверхности цилиндра, и его начало соответствует значению ю « 40. Иными словами, область переходного течения начинается в сечении 0 « 20° в первом случае и в малой окрестности передней критической точки во втором случае.
Турбулизация течения около цилиндра приводит к подавлению осцилляций потока и появлению новых режимов обтекания. Так, например, для рассматриваемых условий эволюционные зависимости (рис. 1, д, е) соответствуют одному и тому же режиму обтекания, для которого характерно либо апериодическое поведение эволюционных зависимостей, либо квазипериодиче-ское поведение с переменной амплитудой и периодом, сравнимым по величине с полным временем движения. Из-за медленности изменения во времени состояния потока этот режим обычно называют квазистационарным.
Из приведенных результатов также следует, что при ламинарно-турбулентном обтекании цилиндра точка отрыва в крайнем левом положении располагается на его лобовой поверхности. Это означает, что при рассматриваемом числе Мш = 0.8 парадокс Прандтля — Эйфеля проявляется слабее, чем в несжимаемой жидкости.
Выше, при обсуждении режимов обтекания кругового цилиндра, использовалась их классификация согласно особенностям поведения эволюционной зависимости. Напомним, что наряду с нею применяется и другая, связанная с ламинарно-турбулентным переходом классификация режимов течения: докритический, критический, сверхкритический и транскритический. По этой
Рем=1031_АМ Ре„=10в1АМ
Ре„=ю6тикв 1^=107 ПЛ»
-50 -45 -40 -35 -30 -25 -20 -15 -10 -5 x0
б)
Рис. 3. Влияние числа Яе на распределение продольного компонента и вектора скорости (а) и коэффициента давления ср (б) вдоль оси абсцисс перед цилиндром в момент времени і = 200
0 5 10 15 20 25 30 35 40 45 х 50
а)
----- Кея = 103 ЬАМ
----- Ре„ = 106иш
..... Кеф = 10бТ1ЖВ
----- Иеш = 107ШВ
Рис. 4. Влияние числа Яе на распределение продольного компонента и вектора скорости (а) и коэффициента давления ср (б) вдоль оси абсцисс за цилиндром в момент времени і = 200
классификации в расчетных случаях имеют место следующие режимы течения: докритический при Яе = 103, 104 и 105, переходный к критическому при Яе = 106 и транскритический при Яе = 107.
Полезную дополнительную информацию можно получить в результате анализа распределений газодинамических переменных по оси абсцисс перед и за обтекаемым телом.
Согласно расчетам распределение газодинамических переменных вдоль оси абсцисс перед цилиндром практически не зависит от числа Яе в рассмотренном диапазоне его изменения. В качестве примера на рис. 3 показаны распределения коэффициента давления и продольного компонента вектора скорости. Можно видеть, что по мере удаления вверх по потоку от передней критической точки указанные газодинамические переменные асимптотически выходят на невязкие значения. Это говорит о корректности постановки задачи и достоверности получаемых результатов.
Распределения газодинамических переменных по оси абсцисс за цилиндром указывают на сильную нестационарность течения в ближнем следе, которая по мере продвижения вниз по потоку ослабевает, и на особенности выхода решения на дальний след. В качестве примера на рис. 4 приведены распределения некоторых газодинамических переменных по оси абсцисс за цилиндром. Можно видеть, что число Яе оказывает относительно небольшое влияние на искомые функции и газодинамические переменные, соответствующие различным числам Яе, образуют достаточно узкий коридор. Кроме того, пульсации компонентов скорости в ближнем следе настолько велики, что местные скорости становятся сверхзвуковыми. Также следует отметить, что на оси абсцисс дальнего следа относительная температура превышает единицу, т. е. дальний след является горячим.
3. СХОД ВИХРЕЙ
Как показано выше, в рассмотренном диапазоне чисел Яе обтекание кругового цилиндра носит нестационарный осциллирующий характер. Эта нестационарность наиболее сильно проявляется в ближнем следе; при этом в разных точках поля течения осцилляции происходят с различными амплитудами и частотами.
Рассмотрим частоту схода вихрей с обтекаемой поверхности цилиндра. Частотные характеристики этого явления определяются различными способами как экспериментально, так и теоретически. Для этой цели в работе использованы эволюционные зависимости местного коэффициента давления в характерных точках В и В, расположенных на обтекаемой поверхности в миде-левом сечении цилиндра сверху и снизу соответственно.
Эволюция коэффициента давления в характерных точках В и В имеет много общего с рассмотренной выше эволюцией положения точек первичного отрыва на поверхности цилиндра (см. рис. 1). В обоих случаях рассматриваемые величины осциллируют в противофазе, а сами зависимости имеют несколько различные форму и максимальную амплитуду; при этом в зависимости от числа Яе реализуется один и тот же режим обтекания цилиндра. В силу сказанного на рис. 5 показана эволюция коэффициента давления только в характерной точке В для граничных значений числа Яе.
Полученные эволюционные зависимости были подвергнуты Фурье-анализу с использованием последних п точек соответствующей кривой. Согласно [12] при ламинарном обтекании цилиндра независимость результатов осреднения от временного интервала, на котором оно проводится, наблюдается при условии п > 256. Для проверки этого условия расчеты были проведены для трех значений п = 256, 512 и 1024. Результаты Фурье-анализа обычно представляются в виде частотной характеристики Ат — Б, где Ат и Б = — безразмерные амплитуда
и частота. В качестве примера на рис. 6 показано влияние п на частотные характеристики при Яе = 103, когда реализуется регу-
Т аблица 1
Частотные характеристики коэффициента давления в миделевом сечении на поверхности цилиндра (п = 512)
ТочкаВ Точка В
Яе срт БИ Ат срт БИ Ат
Ламинарное обтекание
103 -1.2078 0.1953 0.2219 -1.2343 0.1953 0.2255
104 -1.235 0.1953 0.2375 -1.1432 0.1953 0.3036
105 -1.1638 0.1953 0.2543 -1.1962 0.1953 0.274
106 -1.2414 0.1953 0.1754 -1.0844 0.1953 0.2329
Ламинарно-турбулентное обтекание
106 -1.2746 0.0781 0.1976 -1.2887 0.0781
107 -1.1728 0.039 0.3121 -1.5527 0.039
-0.6
-0.9
-1.2
-1.5
-1.8
і
і і \ / / і, / „1_У„ \ / „и. \ 1 \ /
50
100
150
і 200
Рис. 5. Эволюция коэффициента давления ср в характерной точке В цилиндра:
гек;
(ламинарно-турбулентное обтекание)
а — Яе = 103; б — Яе = 106 (ламинарное обтекание); в — Яе = 106; г — Яе = 107
лярный режим обтекания цилиндра. Основные итоговые результаты анализа для п = 512 представлены в табл. 1.
Согласно рис. 6 частотная характеристика для п = 256 представляет собой узкополосный сплошной спектр, а частотные характеристики для п = 512 и 1024 — дискретный спектр, что свойственно для регулярного режима обтекания цилиндра, и совпадают между собой. При других числах Яе частотные характеристики однотипны, образуют узкополосный сплошной спектр и имеют аналогичный характер поведения по п. Таким образом, частотные характеристики при п = 256 определяются некорректно, а независимость результатов осреднения от временного интервала, на котором оно проводится, имеет место при условии п > 512.
Следует отметить, что в ламинарном потоке реализуются режимы обтекания со сходом вихрей с достаточно большой частотой, а используемая методика позволяет получать вполне
с
1
; ■ / 1. .
б)
0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 Р 1
А
0.21
в)
Р
0.2 0.3 0.4 0.5 0.0 0.7
0.9 Е 1
Рис. 6. Частотные характеристики коэффициента давления в характерной
3
точке В цилиндра при Яе = 10: а — п = 256; б — п = 512; в — п = 1024
надежные и достоверные данные. При ламинарно-турбулентном обтекании цилиндра имеют место низкочастотные (квазистационарные) режимы и возникает проблема достоверности полученных частотных характеристик. Несмотря на это, в табл. 1 ради любопытства и полноты информации приведены результаты анализа для обоих чисел Яе.
Сравнение осредненных значений коэффициента давления в симметричных точках В и В указывает на асимметрию осредненного течения. Она связана с неустойчивостью симметричного течения и реализацией около цилиндра нестационарного асимметричного течения. Если знак асимметрии определить как знак выражения А = Срт(В) - Срт(П), то согласно табл. 1 знак асимметрии может быть как положительным, так и отрицательным. Иными словами, выражение А является знакопеременной функцией, знак которой определяется случайными причинами — выбором начальных условий или возмущениями, возникающими в процессе счета.
При ламинарном обтекании частота пульсаций коэффициента давления не зависит от числа Яе, а амплитуда пульсаций изменяется по числу Яе немонотонно. Согласно приведенным данным турбулизация течения в ближнем следе приводит к уменьшению вдвое частоты пульсаций коэффициента давления при Яе = 106 и, следовательно, частоты схода вихрей. К сожалению, результаты расчетов не позволяют установить влияние числа Яе на частоту схода вихрей при ламинарно-турбулентном обтекании цилиндра, так как при Яе = 107 число Струхаля слишком мало и его значение нельзя надежно определить.
Далее отметим, что в [5] на основе уравнений Навье — Стокса смоделировано обтекание кругового цилиндра до- и трансзвуковым потоком совершенного газа и для рассматриваемого числа Маха приведены значения числа Струхаля схода вихрей и коэффициента сопротивления сх для двух значений числа Яе, после пересчета на нашу нормировку они примут вид: = 0.1,
сх = 1.9 для Яе = 5- 102 и = 0.202, сх = 1.7 для Яе = 2.5 • 106. Если при ламинарном обтекании согласно результатам наших расчетов число не зависит от числа Яе, то по данным [5] оно является переменной величиной и указанные его значения располагаются ниже и выше нашей
зависимости. При этом значение числа для второго числа Яе очень близко к нашему решению: их отношение 0.202/0.1953 = 1.034, что косвенно говорит о достоверности результатов численного моделирования. Ситуация с коэффициентом сопротивления обсуждается ниже в четвертом разделе статьи.
4. АЭРОДИНАМИЧЕСКИЕ ХАРАКТЕРИСТИКИ
В [12, 13] подробно изучена эволюция местных аэродинамических характеристик кругового цилиндра и проведена верификация использованного метода численного моделирования путем сравнения результатов расчетов с известными экспериментальными данными. Поэтому ниже основное внимание уделено рассмотрению поведения его интегральных аэродинамических характеристик--коэффициентов сопротивления Сх и подъемной силы Су.
Анализ аэродинамических характеристик цилиндра начнем с изучения поведения коэффициента давления в передней критической точке, где он принимает максимальное значение. В табл. 2 приведены значения Ср в зависимости от числа Яе, а также указаны его асимптотическое значение СрЕ при невязком обтекании (Яе ^ да) и величины отклонения вязкого решения от асимптотического 5 = (Ср - СрЕ)/СрЕ. Расчеты выполнены при больших числах Яе, когда справедливы предположения теории пограничного слоя, согласно которой коэффициент давления в передней критической точке не зависит от числа Яе и равен его асимптотическому значению. Из табл. 2 следует, что условие независимости коэффициента давления от числа Яе выполняется в пределах первых трех значащих цифр, но вместе с тем он примерно на 4% превышает асимптотическое значение. Эти результаты согласуются с данными [13], где изучено трансзвуковое обтекание цилиндра и, в частности, показано, что в передней критической точке как расчетные, так и экспериментальные значения коэффициента давления отличаются от соответствующего асимптотического значения; знак и модуль этого отличия зависят в основном от числа М набегающего потока.
Таблица 2
Коэффициент давления в передней критической точке цилиндра
Режим Ламинарный Ламинарно-турбулентный
Яе 103 104 105 106 106 107
Ср 1.216 1.2146 1.2148 1.2116 1.2184 1.2148
СрЕиЬЕЯ 1.1704 1.1704 1.1704 1.1704 1.1704 1.1704
5, % 3.8961 3.7765 3.7936 3.5202 4.1012 3.7936
Теперь перейдем к рассмотрению интегральных характеристик цилиндра. Коэффициент сопротивления определяется в основном проекциями нормальных напряжений на ось абсцисс, поэтому частота его осцилляций зависит от состояния потока в донной области и, следовательно, должна быть примерно вдвое больше частоты схода вихрей с обтекаемой поверхности. Коэффициент подъемной силы цилиндра определяется в основном проекциями нормальных напряжений на ось ординат, поэтому частота его осцилляций зависит от состояния потока в миделевом сечении тела и, следовательно, должна примерно совпадать с частотой схода вихрей с обтекаемой поверхности.
Этот качественный вывод подтверждается результатами расчетов и, в частности, эволюцией коэффициентов сопротивления и подъемной силы (рис. 7 и 8). Отсюда также видно, что коэффициент Сх осциллирует с малой амплитудой около большого осредненного значения, а коэффициент Су — с большой амплитудой около малого осредненного значения. В случае ламинарного течения увеличение числа Яе практически не влияет на частоту пульсаций и вызывает усложнение формы эволюционной зависимости. При фиксированном значении числа Яе турбулизация течения приводит к сильному уменьшению частоты и амплитуды колебаний; при последующем возрастании числа Яе происходит снижение рассматриваемых величин.
Фурье-анализ эволюционных зависимостей, выполненный для тех же трех значений п, что и выше, показал, что в большинстве случаев частотные характеристики представляют собой дискретный спектр (рис. 9), в котором частота, соответствующая наибольшей амплитуде, принимается
2.08
сх
2
1.92
1.84
1.76
200
б)
в)
г)
Рис. 7. Эволюция коэффициента сопротивления сх цилиндра:
а — Яе = 103; б — Яе = 106 (ламинарное обтекание); в — Яе = 106; г — Яе = 107 (ламинарно-турбулентное обтекание)
(ламинарно-турбулентное обтекание)
Рис. 9. Частотные характеристики коэффициентов сопротивления и подъемной силы цилиндра: а — Яе = 103; б — Яе = 106 (ламинарное обтекание); в — Яе = 106; г — Яе = 107 (ламинарно-турбулентное обтекание);
СХ? Су
в качестве основного тона колебаний рассматриваемой величины. Кроме того, поскольку рассматриваемые аэродинамические коэффициенты осциллируют с большей частотой, чем коэффициент давления в миделевом сечении, то для всех трех п получаются схожие результаты и практически для них справедлив принцип независимости от п. Основные результаты анализа приведены в табл. 3.
Т аблица 3
Частотные характеристики коэффициентов сопротивления и подъемной силы цилиндра (п = 512)
Режим Ламинарный Ламинарно-турбулентный
Яе 103 104 105 106 106 107
Сх 1.9005 1.8408 1.8309 1.7847 1.7898 1.7509
БИ 0.3906 0.3906 0.1171 0.3515 0.039 0.039
Ат 0.0449 0.047 0.053 0.0566 0.1406 0.0299
СУ -0.016 0.046 -0.0293 0.0955 -0.017 -0.1473
БИ 0.1953 0.1953 0.1953 0.1953 0.0781 0.039
Ат 0.2985 0.3065 0.3015 0.2397 0.2672 0.2349
Прокомментируем приведенные результаты расчетов. При ламинарном обтекании аэродинамические коэффициенты осциллируют с достаточно большой частотой, которая позволяет получать корректные результаты Фурье-анализа. В случае ламинарно-турбулентного обтекания цилиндра они осциллируют, как это можно видеть из эволюционных зависимостей, с очень малой частотой, что препятствует получению надежных частотных характеристик. Несмотря на это, в табл. 3 приведены результаты Фурье-анализа и для ламинарно-турбулентного обтекания цилиндра с целью полноты картины по влиянию числа Яе.
Согласно данным табл. 3, при ламинарном обтекании цилиндра его осредненный коэффициент сопротивления монотонно уменьшается по мере увеличения числа Яе. Турбулизация течения при фиксированном Яе приводит к уменьшению Сх, т. е. на рассматриваемом режиме обтекания цилиндра имеет место парадокс Прандтля — Эйфеля. Коэффициент подъемной силы является знакопеременной функцией по числу Яе, что объясняется следующим. Появление подъемной
Рис. 10. Влияние числа Re на коэффициент сопротивления сх цилиндра:
А — [5]; данная работа: о — уравнения Навье — Стокса, • — уравнения Рейнольдса
силы у симметричной конфигурации обусловлено неустойчивостью симметричного течения и асимметрией осредненного течения; знак асимметрии определяется случайными причинами — выбором начальных условий или возмущениями, возникающими в процессе счета.
Как отмечалось выше, в [5] приведены значения коэффициента сопротивления для двух чисел Рейнольдса. Если все эти результаты нанести на один график (рис. 10), то они прекрасно согласуются между собой и образуют единую зависимость, типичную для чисел М«, < 0.9. Это указывает на достоверность полученных результатов расчетов. Далее отметим, что результаты [5] для наибольшего числа Яе соответствуют критическому режиму обтекания цилиндра; в то же время численное моделирование на основе уравнений Рейнольдса при Яе = 106 дает для коэффициента сопротивления значение, которое больше критического. Отсюда следует, что переход от докритического режима обтекания к критическому является не скачкообразным, а эволюционным процессом, протекая хотя и в малом, но конечном интервале изменения числа Яе.
При ламинарном обтекании цилиндра частота осцилляций коэффициента сопротивления при Яе =103 и 104 вдвое превышает частоту схода вихрей, однако при последующем возрастании числа Яе она уменьшается. Это связано с тем, что при Яе = 103 и 104 частотная характеристика представляет собой дискретный спектр, а при последующих числах Яе — сплошной узкополосный спектр с двумя четко выраженными максимумами. В табл. 3 помещены данные, соответствующие первому максимуму; второму максимуму соответствует число 8Ь, примерно вдвое превосходящее его значение для схода вихрей. Частота осцилляций коэффициента подъемной силы, как и ожидалось, совпадает с частотой схода вихрей. При ламинарно-турбулентном обтекании частота осцилляций аэродинамических коэффициентов цилиндра одинакова и соответствует минимальному, отличному от нуля, значению числа 8Ь, которое выдает программа для заданных условий и которое не согласуется с эволюционными зависимостями соответствующих коэффициентов.
Заканчивая обсуждение расчетного материала, отметим следующее. В [20, 21] на частных примерах при различных числах Маха (0 < Мда < да) проведены исследования по эффективности численного моделирования турбулентных течений с использованием нескольких наиболее популярных моделей турбулентности. Сравнение расчетных и экспериментальных данных позволило сделать некоторые выводы общего характера. В частности, все рассмотренные модели турбулентности предназначены для расчета полностью развитых турбулентных течений при больших числах Яе, и все они выдают схожие результаты в том смысле, что различие между расчетными данными, соответствующими разным моделям, в общем много меньше разницы между расчетом и экспериментом.
ЗАКЛЮЧЕНИЕ
Численно на основе нестационарных двухмерных уравнений динамики вязкого совершенного газа смоделировано поперечное обтекание кругового цилиндра трансзвуковым потоком при
3 7
числе Маха Мш = 0.8 и числах Рейнольдса 10 < Re < 10 .
В исследованном диапазоне чисел Re обтекание цилиндра является нестационарным, а его изменение влияет на структуру поля течения и режим обтекания цилиндра. Для ламинарного обтекания цилиндра поле течения характеризуется наличием периодического схода вихрей с обтекаемой поверхности и формированием вихревой дорожки Кармана в ближнем следе; при этом частота схода вихрей не зависит от числа Рейнольдса и соответствует числу Струхала Sh = 0.1953. При наименьшем числе Re реализуется регулярный режим обтекания цилиндра согласно особенностям поведения рассматриваемой эволюционной зависимости. Последующее увеличение числа Re приводит к последовательной смене режимов обтекания — нерегулярным и хаотическим. Турбулизация течения обусловливает подавление осцилляций в потоке, в результате устанавливается квазистационарный режим обтекания.
Поведение коэффициента давления в передней критической точке согласуется с теорией пограничного слоя первого приближения: его значение не зависит от числа Рейнольдса и совпадает с его асимптотическим значением CpE. Результаты наших расчетов подтверждают эти закономерности, что свидетельствует также о достоверности информации, получаемой путем численного моделирования.
По соответствующим эволюционным зависимостям определены осредненные значения коэффициентов сопротивления и подъемной силы цилиндра. Результаты расчетов коэффициента сопротивления хорошо согласуются с расчетными данными [5], и в совокупности они дают полную картину поведения коэффициента сопротивления при переходе от докритического к транскритическому режиму обтекания цилиндра. Отсюда можно также видеть, что критический режим устанавливается не скачкообразно, а эволюционным путем по числу Рейнольдса на сравнительно небольшом интервале его изменения. Коэффициент подъемной силы цилиндра для рассмотренных условий обтекания является малой величиной, хотя его значение может изменяться на порядок при переходе от одного числа Рейнольдса к другому.
Работа выполнена при финансовой поддержке РФФИ (грант № 11-08-01099), ФЦП «Научные и научно-педагогические кадры инновационной России» (ГК № Р595), ФЦП «Исследование и разработки по приоритетным направлениям развития» (ГК № 16.513.11.3119), гранта Правительства РФ по постановлению № 220 «О мерах по привлечению ведущих ученых в российские образовательные учреждения высшего профессионального образования» по договору № 11.G34.31.0072, заключенному между Министерством образования и науки РФ, ведущим ученым Утюжнико-вым С. В. и Московским физико-техническим институтом (государственным университетом).
ЛИТЕРАТУРА
1. Мэрти В. С., Роуз В. К. Детальные измерения аэродинамических характеристик кругового цилиндра при поперечном обтекании // Ракетная техника и космонавтика. 1978.
Т. 16, № 6, с. 8 — 11.
2. Murthy V. S., Rose W. C. Form drag, skin friction and vortex shedding Frequencies for subsonic and transonic gross flows on circular cylinder // AIAA Paper. 1977. p. 77 — 687.
3. Macha J. M. Drag of circular cylinders at transonic Mach numbers // J. Aircraft. 1977.
V. 14, № 6, p. 605 — 607.
4. Rodriguez O. The circular cylinder in subsonic and transonic flow // AIAA. 1984.
V. 22, № 12, p. 1713 — 1718.
5. Ishii K., Kuwahara K. Computation of compressible flow around a circular cylinder // AIAA-84-1631. 1984, p. 11.
6. Ishii K., Kuwahara K., Ogawa S., Chyu W., Kawamura T. Computation of flow around a circular cylinder in a supercritical regime//AIAA-85-1660. 1985, p. 12.
7. Башкин В. А., Егоров И. В., Иванов Д. В. Применение метода Ньютона к расчету внутренних сверхзвуковых течений // ПМТФ. 1997. Т. 38, № 1, с. 30 — 42.
8. Бабаев И. Ю., Башкин В. А., Егоров И. В. Численное решение уравнений На-вье — Стокса с использованием итерационных методов вариационного типа // ЖВМиМФ.
1994. Т. 34, № 11, с. 1693 — 1703.
9. Башкин В. А., Егоров И. В., Егорова М. В., Иванов Д. В. Ламинарнотурбулентное обтекание кругового цилиндра сверхзвуковым потоком газа // Изв. РАН. МЖГ. 2000. № 5, с. 31 — 43.
10. Башкин В. А., Егоров И. В., Егорова М. В., Иванов Д. В. Развитие структуры поля течения около кругового цилиндра при наличии ламинарно-турбулентного перехода // ТВТ. 2000. Т. 38, № 5, с. 759 — 768.
11. Huang P. G., Coakley T. J. Turbulence modeling for high speed flows // AIAA Paper. 1993. № 92-0436.
12. Башкин В. А., Егоров И. В., Ежов И. В., Иванов Д. В. Круговой цилиндр в околозвуковом потоке вязкого совершенного газа // Ученые записки ЦАГИ. 2007. Т. XXXVIII, № 3 — 4, с. 1 — 13.
13. Башкин В. А., Ежов И. В. Круговой цилиндр в трансзвуковом потоке вязкого совершенного газа // Ученые записки ЦАГИ. 2011. Т. XLII, № 1, с. 12 — 30.
14. Годунов С. К. Конечно-разностный метод численного расчета разрывных решений уравнений газовой динамики // Мат. сб. 1959. Т. 47, с. 271 — 291.
15. Годунов С. К., Забродин А. В., Иванов М. Я., Крайко А. Н., Прокопов Г. П. Численное решение многомерных задач газовой динамики. — М.: Наука, 1976, с. 400.
16. R o e P. L. Approximate Riemann solvers, parameter vectors, and difference scheme // J. Comput. Phys. 1981. V. 43, p. 357 — 372.
17. Колган В. П. Применение принципа минимальных значений производной к построению конечно-разностных схем для расчета разрывных решений газовой динамики // Ученые записки ЦАГИ. 1972. Т. 3, № 6, с. 68 — 77.
18. Harten A. High resolution schemes for hyperbolic conservation laws // J. Computat. Phys. 1983. V. 49, p. 357 — 372.
19. Иванов М. Я., Крупа В. Г., Нигматуллин Р. З. Неявная схема С. К. Годунова повышенной точности для интегрирования уравнений Навье — Стокса // ЖВМиМФ. 1989. Т. 29, № 6, с. 888 — 901.
20. Rumsey C. L., Spalart P. R. Turbulence model behavior in low Reynolds number regions of aerodynamic flowfields//AIAA J. 2009. V. 47. № 4, p. 982 — 993
21.Rumsey C. L., Greenblatt D. Flow control predictions using unsteady Reynolds-averaged Navier — Stokes modeling: a parametric study // AIAA J. 2009. V. 47, N 9, р. 2259 — 2262.
Рукопись поступила 18/Х 2011 г.