CC BY
81
Том XLV
УЧЕНЫЕ ЗАПИСКИ ЦАГИ
2014
№ 6
УДК 533.6.011.6
ЛАМИНАРНО-ТУРБУЛЕНТНЫЙ ТЕПЛООБМЕН НА ПОВЕРХНОСТИ ПОЛУСФЕРЫ, ОБТЕКАЕМОЙ СВЕРХЗВУКОВЫМ ПОТОКОМ ВОЗДУХА
В. В. ГОРСКИЙ, М. А. ПУГАЧ
Сопоставлены результаты двух расчетных методов (метода эффективной длины и алгебраической модели турбулентности Себечи — Смита) с экспериментальными данными по распределению теплового потока по поверхности полусферы при числе Маха Мм = 5. Предложена модификация модели турбулентности Себечи — Смита.
Ключевые слова: теплообмен, пограничный слой, модель турбулентности, ламинарно-турбулентный переход.
ВВЕДЕНИЕ
Обгар тепловой защиты гиперзвукового летательного аппарата (ГЛА) во многих случаях оказывает заметное влияние на его аэродинамические характеристики [1, 2]. В этой связи значительный интерес вызывает описание турбулентного теплообмена на теплонапряженных элементах конструкции ГЛА, которое обычно базируется на использовании либо метода эффективной длины В. С. Авдуевского [3, 4], либо полуэмпирических (алгебраических или двухпараметриче-ских дифференциальных) моделей турбулентности [4].
Метод эффективной длины является одним из основных инструментов проектирования тепловой защиты в России [4]. Он основан на перенесении экспериментальных данных по конвективному теплообмену, полученных для пластины, на поверхность произвольной кривизны. Необходимо учитывать, что методы расчета турбулентного теплообмена по сути своей являются полуэмпирическими. Поэтому необходима апробация их в условиях, по возможности близких к натурным. В литературе основное внимание уделяется исследованию турбулентных течений на боковых поверхностях летательных аппаратов. В то же время полусферические конфигурации наиболее часто встречаются на практике при создании ГЛА.
Сказанное выше предопределяет повышенный интерес к анализу погрешностей, свойственных основным процедурам расчета турбулентного теплообмена на полусфере, и поиску путей
снижения этих погрешностей. Решению этих проблем и посвящена данная работа.
Доступные экспериментальные данные по распределению параметров ламинарно-турбу-лентного теплообмена по поверхности полусферы относятся лишь к малым числам Маха в набегающем потоке газа. Однако, на наш взгляд, это не препятствует использованию на практике результатов данной работы, так как в ней выполняются условия геометрического подобия и локального газодинамического моделирования исследуемого явления. К тому же, повсеместно используемые полуэмпирические модели турбулентности базируются на исполь-
ГОРСКИИ Валерий Владимирович
доктор технических наук,
профессор, главный научный сотрудник ЦАГИ
ПУГАЧ Михаил Александрович
инженер ЦАГИ, студент МФТИ
зовании понятия пути перемешивания Л. Прандтля, причем константы в этих моделях установлены на базе результатов исследований обтекания плоской пластины при малых числах М.
1. ФИЗИКО-МАТЕМАТИЧЕСКАЯ ПОСТАНОВКА ЗАДАЧИ
Объектом настоящих исследований являются экспериментальные данные Уидхопфа и Холла [5]. Расчет ламинарно-турбулентного перехода в пограничном слое проводился с использованием методики PANT [7], а в качестве математических моделей теплообмена рассматривались алгебраическая модель турбулентности Себечи — Смита [6] и метод эффективной длины В. С. Авдуевского [3].
Анализируемые экспериментальные данные [5] получены при числе Мда = 5 и температуре торможения порядка 700 К. Расчетные исследования выполнялись путем численного решения двумерных уравнений пограничного слоя, записанных для совершенного газа в переменных Лиза — Дородницына [8].
1.1. ОПРЕДЕЛЕНИЕ ПОЛОЖЕНИЯ ПЕРЕХОДНОЙ ОБЛАСТИ НА ПОВЕРХНОСТИ ТЕЛА
(МЕТОДИКА PANT)
Для определения положения переходной области на затуплении летательного аппарата на практике часто [4] используется модель из работы [7]. В рамках этой модели сохранение ламинарного режима в пограничном слое зависит от нескольких факторов:
отношения шероховатости «стенки» к толщине потери импульса в этом слое; уровня турбулентности в набегающем на тело газовом потоке; скорости вдува газа в пограничный слой через «стенку».
Для расчета критического значения числа Рейнольдса, которое соответствует границе ламинарного режима, используются формулы вида:
d = ГТ +
Ree,a1 = min ( 215/ d °7,300); RQeb t
c-i(-/rmi trnriQ
Ree ,»,1
(1)
0.1Gw +(1 + 0.25Gw )
Pw
_ peUebrou,trans . ^ _ peUe»
e,brou,trans .. ' e,» .. *
И e И e
Здесь Гт — уровень турбулентности в набегающем потоке; ре, Ue — плотность, коэффициент динамической вязкости и скорость газа на внешней границе пограничного слоя; pw — плотность газа при температуре «стенки»; Ree brou ^^ — число Рейнольдса, рассчитанное по параметрам течения на внешней границе пограничного слоя и высоте brou,trans шероховатости
«стенки», характерной для перехода от ламинарного режима течения к турбулентному течению; Ree » — число Рейнольдса, рассчитанное по параметрам течения на внешней границе пограничного слоя и толщине » потери импульса в ламинарном пограничном слое; Gw — массовая скорость вдува газа в пограничный слой, измеренная в долях от коэффициента теплообмена на непроницаемой «стенке».
Необходимым условием для возникновения турбулентного режима в пограничном слое, положение которого на поверхности тела определяется критерием (1), является выполнение неравенства
max (215.255(1 -ГТ/ d)
Ree,»,max > ^07 ,
где Ree» max — максимальное значение Ree» на поверхности тела.
Для расчета критического значения числа Яе(
'е,Ь„
, которое соответствует началу турбу-
лентного режима, используется зависимость, предложенная в работе [8]:
Л -1
Яе
е ,& ,2
= л/2Яе
е ,1
1
Здесь аей- — эффективный угол атаки, под которым понимается угол между вектором скорости летательного аппарата и плоскостью, касательной к поверхности тела.
Расчет коэффициента ц динамической вязкости газа производится с использованием коэффициента перемежаемости Г [10] по формуле вида:
Ц = Ць + ГЦт,
а величина этого коэффициента, так же как и в работе [9], рассчитывается с использованием кубического сплайна вида
Г = тах^ тт [ 1 е2 (3 - 2е)]}, е = (Я.^ - Я_ее,ад)/(Я^д - ).
Здесь Ць и — ламинарный и турбулентный коэффициенты динамической вязкости. Расчет коэффициента теплопроводности производится по аналогичной формуле. При этом турбулентное число Прандтля полагается равным единице.
1.2. АЛГЕБРАИЧЕСКАЯ МОДЕЛЬ ТУРБУЛЕНТНОСТИ СЕБЕЧИ — СМИТА
В рамках двухслойной полуэмпирической модели Себечи — Смита [6] во внутренней (пристеночной) и внешней областях течения используются различные формулы расчета составляющей коэффициента динамической вязкости, обусловленной турбулентными пульсациями газа.
В первой из указанных областей используется формула Прандтля цт = р12 |, в которой
масштаб турбулентности ¡т полагается пропорциональным удалению расчетного узла от «стенки» (т. е. координате у), а для плавного сопряжения течения газа в пристеночной области с течением в ламинарном подслое используется демпфирующая функция Ван — Дриста О(у). Таким образом, в этой области течения газа
Цт = Цт ,ш = Р [кУБ (У)]2 — |;
(2)
Б (у ) = 1 - ехр
У/У 26
У =-
Ц М
и =
Р
Здесь р — плотность газа; и — тангенциальная проекция вектора скорости газа; у — координата, отсчитываемая от «стенки» в направлении ее внешней нормали; к — константа Кармана, равная 0.41; у*, и* — характерные значения длины и скорости, используемые в качестве масштабов; х^ — напряжение трения на стенке. Индекс w относится к параметрам газа на «стенке».
Во второй (внешней) области для расчета Цт используется формула:
у ( \
Цт = Цт,си1 = 00168Р— П1 - — Ук Ф;
О I — J
Ук =
1 + в
У**
-1
(3)
а = 6; р = 5.5.
/ ** \
Предел интегрирования в этой формуле определяется из условия -= 0.995 . В качестве
ие
расчетного значения турбулентного коэффициента динамической вязкости принимается минимУм из Цт,ш и ЦТ,out •
2. ЭКСПЕРИМЕНТАЛЬНЫЕ ДАННЫЕ ПО КОНВЕКТИВНОМУ ТЕПЛООБМЕНУ
В работе [5] приведены результаты многочисленных экспериментальных исследований ла-минарно-турбулентного теплообмена, выполненных на трех различных режимах. Радиус Я сферического затупления модели во всех экспериментах был равен 0.0635 м. Информация об остальных условиях экспериментов приведена в табл. 1.
Здесь N — номер режима испытаний; Я — число Рейнольдса, рассчитанное по параметрам набегающего потока и радиусу сферического затупления модели; —удельный тепловой поток в критической точке сферического затупления модели, рассчитанный по формуле Фэя — Ридделла [11]; — отношение температуры «стенки» к температуре торможения.
Информации, представленной в табл. 1, явно недостаточно для проведения расчетных исследований. Недостающая информация может быть получена путем построения итерационного процесса определения температуры воздуха в набегающем потоке, соответствующей заданному значению удельного теплового потока в критической точке сферического затупления. Результаты решения указанной задачи приводятся в табл. 2.
Результаты экспериментальных исследований конвективного теплообмена на поверхности полусферы, полученные в работе [5], приводятся на рис. 1. Здесь — удельный тепловой поток, отнесенный к его значению в критической точке; 5 — длина криволинейной образующей сферы, измеренная от ее «критической точки» в долях от радиуса этой сферы Я. На рис. 1 приводятся также аналогичные результаты расчетных исследований, выполненных в различных постановках.
В экспериментальных данных наблюдается довольно большой разброс, особенно при
Таблица 1
Условия проведения экспериментов, указанные в работе [5]
N M« Re«, r qw,0, Вт/м2 T L w
I 5 107 0.41 ■ 106 0.25
II 5 4 ■ 106 0.25 ■ 106 0.26
III 5 2.5 ■ 106 0.19 ■ 106 0.23
Т аблица 2 Установленные условия проведения экспериментов
N T«, K p«, кг/м3 V«, м/с T0, K p0, МПа Tw, K
I 82.3 0.9895 909.3 493.8 0.762 123.45
II 81.4 0.3872 904.3 488.4 0.296 126.984
III 79.9 0.2310 896.0 479.4 0.173 110.26
Рис. 1. Распределение удельного теплового потока по поверхности полусферы:
режим испытаний а — I, б — II, в — III; 1 — экспериментальные данные; 2 — расчет по модифицированной модели Себечи — Смита; 3 — расчет по стандартной модели Себечи — Смита; 4 — расчет по методу эффективной длины
больших числах Рейнольдса, что характерно для большинства тепловых измерений рассматриваемого класса. При практическом решении тепловых задач результаты измерений обычно аппроксимируют по верхней границе коридора разброса экспериментальных данных, что обеспечивает некоторый запас по величине прогрева конструкции. Однако при решении абляционных задач необходимо ориентироваться на средние значения экспериментальных данных по теплообмену, что и сделано при проведении данных исследований.
3. РЕЗУЛЬТАТЫ РАСЧЕТНЫХ ИССЛЕДОВАНИЙ
Во всех расчетах газ считался совершенным, а условия на внешней границе пограничного слоя формировались на базе решений уравнений Эйлера с последующим применением процедуры сплайновой аппроксимации [12].
Как видно из представленных на рис. 1 экспериментальных данных, на всех режимах испытаний четко прослеживается положение переходной области. Поэтому первой задачей, решенной в данной работе, было определение величины шероховатости стенки по положению переходной области для каждого режима испытаний. Результаты решения этой задачи представлены в табл. 3.
Как видно из приведенных данных, поверхность всех испытанных в работе [5] моделей характеризуется шероховатостью «стенки» порядка 5 мкм. Из этого следует, что усиления теплообмена за счет шероховатости «стенки» в этих экспериментах не происходило, так как в рассматриваемых условиях оно должно происходить при существенно большей шероховатости «стенки» [4].
Экспериментальные данные по распределению конвективного теплообмена по поверхности полусферы сопоставлены на рис. 1 с результатами расчетов в рамках метода эффективной длины В. С. Авдуевского и стандартной модели турбулентности Себечи — Смита. Модель Себечи — Смита качественно правильно описывает экспериментальные данные, но дает завышенную тепловую нагрузку по сравнению с экспериментом. Эта тенденция становится существенной по мере увеличения числа (вплоть до увеличения теплового потока приблизительно в полтора раза при Ке,»^ ~ 10 ). Метод эффективной длины дает аналогичные результаты. Необходимо учитывать, что именно при больших числах Re турбулентный теплообмен оказывает наибольшее влияние на обгар тепловой защиты. Поэтому представляет интерес модификация расчетной модели с целью существенного улучшения описания экспериментальных данных работы [5].
Решение этой задачи, рассматриваемое ниже, базируется на оптимизации ряда констант стандартной модели Себечи — Смита с целью обеспечения минимума среднеквадратичного рассогласования между расчетными и экспериментальными данными.
Оптимизировались следующие константы: константа Кармана к, входящая в масштаб турбулентности во внутренней части пограничного слоя, рассчитываемый по формуле (2), и константы а и р, входящие в функцию Клебанова ук , используемую во внешней части пограничного слоя и рассчитываемую по формуле (3).
Оптимизационная задача решается с помощью одного из вариантов метода Хука — Дживса [13], который служит для поиска безусловного локального экстремума функции и опирается непосредственно на значения функции. Результаты решения показали, что трехпараметрическая оптимизация не имеет преимущества перед оптимизацией по одному параметру — константе Кармана, оптимальное значение которой равно 0.252. Результаты оптимизации приведены на рис. 1.
Необходимо отметить, что для расчетов тепловых потоков на боковой поверхности тел (цилиндра, конуса) стандартная модель с константой к = 0.41дает хорошие результаты. Она также хорошо описывает экспериментальные данные на полусфере для точек с центральным углом
9 > 70° (см. рис. 1). Поэтому для задач,
Таблица 3
Параметры, характеризующие положение переходной области на поверхности модели
N Ьти,и, мкм Кее,э,1 Кее,8,2 $2
I 3 103 0.25 203 0.48
II 4.5 118 0.45 233 0.97
III 9 75 0.36 149 0.77
в которых затупление плавно сопряжено с боковой поверхностью, можно предположить, что к = к(я) и что эта зависимость подбирается таким образом, чтобы обеспечить изменение константы Кармана к от к0 (при 9 = 0) до к = 0.41(при 9 = 90°), где к0 определяется из решения оптимизационной задачи. Из постановки ясно, что в качестве аргумента для такой
зависимости необходимо использовать функцию, изменяющуюся от нуля до 1. В качестве возможного варианта была предложена следующая зависимость (кубический сплайн от sin(s)):
k(s) = к0 + (0.41 - к0 ) sin2 (s)(3 " 2 sin (s)).
Оптимизационная задача, как было сказано выше, решалась для константы К0. Результаты расчетов (рис. 2) показали, что использование такой зависимости позволяет существенно повысить точность описания экспериментальных данных во всем рассмотренном интервале изменения Re«^. При этом оптимальное значение к0 = 0.04.
ВЫВОДЫ
1. Показано, что использование методики PANT позволяет удовлетворительно описать расположение переходной области на поверхности полусферы с шероховатостью порядка 5 мкм во всех рассмотренных экспериментах.
2. Отмечено, что при столь малой шероховатости «стенки» можно пренебречь усилением турбулентного теплообмена в анализируемых экспериментах за счет шероховатости.
3. Показано, что стандартная модель турбулентности Себечи — Смита и метод эффективной длины Авдуевского качественно правильно описывают рассмотренные экспериментальные данные, но характеризуются завышением удельного теплового потока примерно в полтора раза при больших числах Рейнольдса.
4. Предложена модифицированная модель турбулентности Себечи — Смита, использование которой позволяет удовлетворительно описать весь круг рассмотренных экспериментальных данных.
Работа выполнена в МФТИ при финансовой поддержке РНФ (грант №14-19-00821).
Рис. 2. Распределение удельного теплового потока по поверхности полусферы:
режим испытаний а — I, б — II, в — III; 1 — экспериментальные данные, 2 — расчет по модифицированной модели Себечи — Смита с к = k(s)
ЛИТЕРАТУРА
1. Егоров И. В., Горский В. В., Лунев В. В. и др. Гиперзвуковая аэродинамика и тепломассообмен спускаемых космических аппаратов и планетных зондов. — М.: Физмат-лит, 2011, 548 с.
2. Бабылев А. В., Ваганов А. В., Дмитриев В. Г. и др. Разработка аэродинамической компоновки и исследования аэротермодинамических характеристик малоразмерного крылатого возвращаемого аппарата // Ученые записки ЦАГИ. 2009. Т. XL, № 3, с. 3 — 15.
3. Авдуевский В. С., Галицейский Б. М., Глебов Г. А. и д р.Основы теплопередачи в авиационной и ракетно-космической технике / Под общей ред. В. К. Кошкина. — М.: Машиностроение, 1975, 623 с.
4. Руководство для конструкторов. Конвективный теплообмен изделий РКТ. — г. Королев, МО: ФГУП «ЦНИИмаш», 2010, 397 с.
5. W i d h o p f G. F., H o 11 R. Transitional and turbulent heat-transfer measurements on a yawed blunt conical nosetip // AIAA J. 1972. V. 10, N 10, p. 125 — 132.
6. Cebeci T., Smith A. M. O. Analysis of turbulent boundary layers. — New York: Academic Press, 1974.
7. Лойцянский Л. Г. Механика жидкости и газа. — М.: Дрофа, 2003, 840 с.
8. Anderson A. D. Surface roughness effect. Boundary layer transition data correlation and analysis // PANT. 1974. Part III, SAMSO TR-74-86.
9. Горский В. В., Носатенко П. Я. Математическое моделирование процессов тепло- и массообмена при аэротермохимическом разрушении композиционных теплозащитных материалов на кремнеземной основе. — М.: Научный мир, 2008, 255 с.
10. Сафиуллин Р. А. Теплообмен в области перехода ламинарного пограничного слоя в турбулентный // Изв. АН СССР. МЖГ. 1971, № 6, с. 102 — 108.
11. Фэй Дж., Ридделл Ф. К. Теоретический анализ теплообмена в лобовой точке, омываемой диссоциированным воздухом / В кн. : Проблемы движения головных частей ракет дальнего действия. — М.: Изд. иностр. лит., 1959.
12. Горский В. В. Метод сплайновой аппроксимации // ЖВМиВФ. 2007. Т. 47, № 6, с. 939 — 943.
13. Аоки М. Введение в методы оптимизации. Основы и приложения нелинейного программирования. — М.: Наука, 1977, 343 с.
Рукопись поступила 17/III2014 г.
