Моделирование термонапряжённого состояния ротора компрессора выского давления с учётом вторичных течений в полостях Текст научной статьи по специальности «Физика»

УДК 621.431.75

МОДЕЛИРОВАНИЕ ТЕРМОНАПРЯЖЁННОГО СОСТОЯНИЯ РОТОРА КОМПРЕССОРА ВЫСКОГО ДАВЛЕНИЯ С УЧЁТОМ ВТОРИЧНЫХ ТЕЧЕНИЙ В ПОЛОСТЯХ

© 2011 Ю. М. Темис, А. В. Селиванов, Г. Г. Юрченко ФГУП «ЦИАМ им. П. И. Баранова», г. Москва

На основе объединения одномерных и двумерных моделей разработана методика расчёта течения в полостях КВД, что позволило получить более точное распределение коэффициента теплоотдачи и других параметров газа во вторичном тракте. Полученные распределения использованы в качестве граничных условий при определении термонапряжённого состояния ротора КВД.

Термонапряжённое состояние, вторичный тракт, гидравлическая сеть, коэффициент теплоотдачи.

Введение

Определение теплового и напряжённо-деформированного состояния роторов и корпусов газотурбинных двигателей (ГТД) относится к сложным проблемам многодисциплинарного математического моделирования. Проведение таких расчётов необходимо на стадии проектирования компрессора выского давления (КВД), в том числе для определения прочности и долговечности конструкции и расчёта кинетики радиальных за-

зоров (рис. 1).

Моделирование термонапряжённого состояния предполагает предварительный расчёт течения газа во вторичном тракте для получения распределений давления, температуры и коэффициента теплоотдачи, которые используются в качестве граничных условий. Известны методики такого расчёта, основанные на построении гидравлической сети и одномерных моделей течения [1]. Однако в междисковых полостях КВД реализу-

Рис. 1. Алгоритм расчёта кинетики радиальных зазоров

ется циркуляционное течение, для которого применение одномерных моделей не даёт удовлетворительных результатов. В настоящей работе в единую расчётную схему объединены одномерные и двумерные модели, что позволило получить достаточно точные результаты без чрезмерных вычислительных затрат.

Постановка задачи

Распределение температуры в деталях ротора и статора зависит от их формы и теплофизических свойств материалов и определяется, в основном, конвективным теплообменом с потоками газа в тракте, полостях и каналах охлаждения и контактным теплообменом сопряжённых деталей:

Т = Т(а ср ,Т

ср

а к

Л, РЛБ д)

где а ср - коэффициенты теплоотдачи газа на

поверхностях ротора и статора; Тср - температуры газа, омывающего поверхности ротора и статора; а* - коэффициенты контактного теплообмена; Бд - вектор геометрических параметров, определяющих форму деталей; Cp (Т), А(Т) и р(Т) - теплоёмкость, теплопроводность и плотность деталей ротора и статора соответственно.

При решении задачи теплопроводности для ротора и статора КВД в первом приближении можно пренебречь лучистым теплообменом, а в качестве граничных условий использовать соотношения:

Яср = а ср(Т„ - Тср);

Чк =а к (Т*2 - Т,^

где ^р , Я к, Tw, тк1, Тк2 - удельные тепловые потоки и температуры металла на поверхностях конвективного и контактного теплообмена, соответственно.

В свою очередь, параметры конвективного теплообмена ( а ср, Тср ) являются функциями расходов G в системе вторичных потоков и параметров течения газа (скоростей, плотности, теплофизических свойств и др.).

Они зависят от формы каналов системы вторичных потоков F о, текущих зазоров в каналах и уплотнениях А, коэффициентов гидравлического сопротивления X , угловой ско-рости ротора о и параметров газа на входе (Ротб и Тотб ) и выходе (Рвых ) из системы:

а = а( а , х ,р отб ,т<

,й ,1^ т ) .

Таким образом, для расчёта термонапряжённого состояния ротора КВД необходимо предварительно определить параметры течения во вторичном тракте. Для этого можно использовать подход, основанный на построении эквивалентной гидравлической сети [1]. Идея этого метода заключается в разбиении вторичного тракта на типовые элементы: кольцевые и радиальные каналы, радиальные полости, лабиринтные уплотнения, элементарные сопротивления, - для каждого из которых строятся одномерные модели течения.

Одномерные модели течения газа в типовых элементах

Кольцевой канал

Для описания течения в кольцевом канале используются уравнения сохранения массы, импульса и энергии в цилиндрической системе координат:

dVz ЛР

Р^—Г = -Т -(і~г dz dz

= -т

Е <р ;

d pVzd

dz

= 0

(1)

(2)

(3)

(4)

где тЕг и тЕ(р - суммарные напряжения трения в осевом и окружном направлениях соответственно; Ь - работа, совершаемая вращающейся поверхностью; Q - поток тепла

за счёт теплообмена; е - внутренняя энергия; где а - поправочный коэффициент на на-d - зазор в канале; Я - внутренний радиус канала;

L =

Яо

(У, -Яо)^у; + (У, - Яо)

2р 7 а;

Q = а^Т +

(У, - оЯ)2 + У,2

2с,

- Т) +

22

+ а2к2 (Т +

(Уф-тЯУ + У,

2с,

- Т2);

чальный участок; /1 и /2 - безразмерные коэффициенты трения на поверхности ротора и статора, соответственно [3]:

/1,2 =

0.3164

"п 025 Яе1,2

Ь:

Яе1 =

22

(V, -Яо)2 + У,

ТЕг к1Т1г + к2Т2г ;

где т1 г ,т1ф ,т2г и т2ф - проекции напряжения

трения в осевом и окружном направлениях на поверхности ротора и статора, соответственно [2]; к1 и к2 - отношение площади

действительной поверхности ротора к кольцевой площади на участке интегрирования. Напряжения трения вычисляются по формулам:

/ У^У,2 + (К, - Яо)2 ^

—ар--------------------;

т1ф =

/1 (Уф - Яо)^У,2 + (Уф - Яо)2

ар—----------------------у— -;

Яе2 =

т

где Ь = 1.5 - поправочный коэффициент для течения в кольцевом канале [3].

В формуле для вычисления потока тепла Q присутствует коэффициент теплообмена а = Ши /1, где число Нуссельта определяется по эмпирической зависимости через число Рейнольдса Яв и Прандтля Рг : [1]

0.43 г>„0.8

Ш = 0.021є,РгЯе

(5)

Радиальный канал и радиальная полость

Уравнения, определяющие течение для радиального канала и радиальной полости, схожи по структуре, но отличаются значениями массовых сил ¥г и ^, а также значениями работы Ь и потока тепла Q:

_ Уф=Р_ 1 &

г dr г г р dr ’

dУ У

Уг—ф + Угф = ^

dr г

(6)

(7)

/2ар ,

d (ерУ^ ) dr

+ L + Q = 0;

(8)

т„ = /2 арЩ^-

4 2

' 2ф

d (рУгБ )

dr

=0.

(9)

2

1г

Для радиального канала

^ =-р^ + Го2; Fф = 2оУг

V2 + У2

L = 28о2УгГ ; Q = 2пЯа(Т+ ^—^-Т,),

2с

где а определяется через число Ыи (5), £ - площадь сечения канала.

Для радиальной полости

F 1 F 1

г =-----------тЕ • г ф =-----------т, •

г рd Ег ; ф рd Еф;

нена модель двух контрольных объёмов (рис. 2).

Рис. 2. Модель лабиринтного уплотнения

L = о1Г Тф1 + о 2Г Тф2';

В осесимметричном случае определяющее уравнение этой модели примет вид

Q = «1(Т +

(Уф - о/)2 + Уг2

2

22

- Т1) +

, (Т, (Уф о2гу+у; Т ч

+ а2(Т + —-----2-------T2),

дрА- + РА К

ді

+

+ т0т-+1Нк+1 VЯТ (Рк Р-+1)

- ^0 т-Нку1ЯТ(Рк 1 - Рк ) = 0,

где тЕг - суммарное касательное напряжение в радиальном направлении. Коэффициенты трения определяются по эмпирическим зависимостям, аналогичным зависимостям для кольцевого канала. Вычисление коэффициента теплоотдачи зависит от того, вращается стенка или она неподвижна [4]. Для вращающейся стенки:

Ыи = 0.0196 11 - Ь |02 Яв0 8К , где Ь - закрутка потока. Для неподвижной стенки:

Ыи = 0.037Яв0 8Рг033.

Лабиринтное уплотнение

Исследования течений в лабиринтных уплотнениях показывают, что структура потока в камере и в кольцевом зазоре над ней отличаются друг от друга. Течение в камере характеризуется наличием вихревой структуры, а течение в зазоре имеет ярко выраженную осевую направленность. Однако давление в этих областях практически одинаково. По этим причинам для моделирования течения через лабиринтное уплотнение приме-

где коэффициент расхода тк определяется в соответствии с формулой Чаплыгина [2]:

п

п + 2-5Д. + 2£2 ,

в которой

= -1 +

Р

N(7-1)/7

к-1

V Рк 0

а коэффициент т0 может быть вычислен по формуле:

тс =(1 -у)

•0.5

где У = 8 52Нк /[ьк - Р + 7 23Нк ].

Элементарное сопротивление

Типовое элементарное сопротивление вводится, например, при внезапном расширении канала. Коэффициент сопротивления

в случае турбулентного течения и равномерного распределения скоростей по сечению канала перед его расширением зависит только от отношения площадей узкого и широкого сечений т и вычисляется по формуле Борда-Карно [3]:

2 АР )2

S = PV2 = (1 - ПІ) ■

где DP - потеря давления на участке.

Для всех типовых элементов созданы отдельные расчётные модули, реализованные в программной среде Maple. Эти модули могут быть легко интегрированы в единую расчётную схему, что позволяет автоматизировать процесс построения гидравлической сети. Для решения нелинейных дифференциальных уравнений применён метод Рунге-Кутты четвёртого порядка.

Верификация типовых элементов

Верификация одномерных моделей типовых элементов выполнена с помощью программного комплекса Fluent. Для каждого типового элемента построены модели, отражающие топологические особенности течения. В одномерных и поверочных двумерных моделях были заданы одинаковые граничные условия (давление, температура, закрутка). Для двумерных моделей применялась к - e модель турбулентности и расширенная пристеночная функция для разрешения пограничного слоя.

Погрешность вычисления давления относительно перепада давлений между входом и выходом составила менее 4.6 %. Абсолютная погрешность для выходной температуры - менее 5 К, для скоростей - порядка 4 м/с. Полученные результаты свидетельствуют о том, что одномерные модели достаточно хорошо описывают течения во вторичном тракте с точки зрения инженерной практики.

Двумерные модели полостей

Радиальная полость, в которой преобладает сквозное течение, может быть описана при помощи одномерной модели, однако для междисковых полостей КВД с характерным циркуляционным течением необходимо

построение двумерных математических моделей (рис. 3). Такой подход был реализован в работе [7], однако полости рассматривались в виде изолированных областей, не связанных с остальными элементами вторичного тракта. В настоящей работе 2Б модели полостей и Ш модели гидравлической сети были интегрированы в единую расчётную схему.

Рис. 3. Модель междисковой полости КВД

Для моделирования 2Б течения в полостях были построены расчётные сетки со сгущением возле стенок, обеспечивающим

значения характерного параметра У + менее

3 (пограничный слой содержит достаточное количество ячеек). Сетки с такими параметрами позволяют достаточно точно рассчитывать пограничный слой и находить распределение коэффициента теплоотдачи по поверхности полости.

Согласование одномерной гидравлической сети и двумерных моделей полостей (по давлению, температуре и скорости) проводится следующим образом. Для полости строится расчётная сетка, содержащая область циркуляционного течения I и область II, прилегающую к валу (рис. 3). Область II одновременно является ветвью гидравлической сети (кольцевой канал с нулевым трением на одной стенке), и поэтому 2Б модель полости и Ш модель канала будут иметь одинаковые граничные условия. После начального расчёта гидравлической сети определяются параметры газа в её узлах, которые пе-

редаются в виде граничных условий в двумерную модель полости. По результатам двумерного моделирования, если в этом есть необходимость, проводится коррекция одномерной модели канала II. Затем производится повторный расчёт гидравлической сети.

Необходимо отметить, что в междис-ковых полостях экспериментально наблюдается трёхмерная структура течения, вызванная “ эффектом плавучести газа” (области газа с различными температурами начинают перемещаться в поле массовых сил) [5]. Однако для инженерных расчётов использование трёхмерных моделей для учёта этого эффекта является нецелесообразным с точки зрения вычислительных затрат. Возможным решением является применение двумерных моделей с модифицированными коэффициентами для молекулярной вязкости

ц' = СЯаПц и теплопроводности газа 1 = СЯа^Х [6], где С, п - эмпирические константы; Яа{ - локальное число Рэлея, которое вычисляется по формуле:

Щ - у = °.

(11)

Я Р р2о2ЬА Яа, = Рг-------- —тах

т

(г \2 л Л

ог | г ар 0

р Ф ’

Построение гидравлической сети КВД

При построении гидравлической сети вторичный тракт разбивается на типовые элементы, которые образуют ветви сети, связанные между собой узлами.

Для каждой ветви выполняется уравнение, которое связывает значения давлений на её концах [1, 8]:

р(11) р(1,)

О. .

= хе 1 -

- + АР

2р-Б2 цс

(10)

где р(1) = Р / Р*; X - коэффициент гидравлического сопротивления; £ - площадь проходного сечения; р - плотность воздуха;

М - число Маха; G - расход воздуха на ветви; 0 - функция учёта сжимаемости [4].

Также в каждом внутреннем узле выполняется второй закон Кирхгофа для расходов:

В результате получаем нелинейную систему уравнений относительно неизвестных давлений в узлах сети, решение которой выполняем методом последовательных приближений. На первой итерации при п(А)=1 и

АР*с =0 получаем значения давления в узлах

сети, а также распределение расходов по ветвям. Далее для каждого типового элемента решаем систему уравнений, описывающую течение в данном элементе, после чего получаем значения давления, температуры и закрутки на концах ветвей. В точках слияния потоков проводится осреднение (по расходу) закрутки и полной температуры потоков. Затем при тех же значениях давлений в узлах проводится пересчёт типовых элементов. Процедура выполняется до тех пор, пока значения закрутки и полной температуры не установятся (разность значений на двух соседних итерациях будет меньше заданной точности). Затем определяется вектор значений (п(А)} и {АРЦс}, решается система уравнений (10)-(11) с новыми коэффициентами и т.д. до установления значений {п(А)} и

{АР,*}.

Расчёт термонапряжённого состояния КВД

При помощи разработанной методики осуществлён расчёт коэффициента теплоотдачи и параметров газа во вторичном тракте КВД, включая междисковые полости (рис. 4). Параметры газа в основном тракте считались известными. Полученные данные приложены в качестве граничных условий (третьего рода) к конечно-элементной модели ротора КВД, которая состоит из осесимметричных элементов на дисках и плосконапряженных элементов на лопатках.

После расчёта теплового состояния ротора КВД проведён расчёт напряжённо-деформируемого состояния (НДС) с учётом тепловых и центробежных нагрузок по полётному циклу. На рис. 5 показано перемещение ротора на взлётном режиме.

На основе решения задачи в квазиста-ционарной постановке можно оценить зна-

чения радиальных зазоров на различных режимах. Интерполируя значения коэффициента теплоотдачи, температуры и давления между стационарными точками, можно рассчитать кинетику радиальных зазоров по полётному циклу, что позволит предпринять необходимые конструктивные решения по минимизации зазоров ещё на стадии проектирования.

Заключение

Разработана методика расчёта термонапряжённого состояния ротора КВД с учётом течения газа во вторичном тракте. Расчёт вторичных потоков основан на построении эквивалентной одномерной гидравлической сети и двумерных моделей междис-ковых полостей, что позволило получить

более точное распределение температуры и коэффициента теплоотдачи. Одномерные и двумерные модели согласованы по входным и выходным параметрам. При помощи разработанной методики рассчитано тепловое и напряжённо-деформированное состояние ротора КВД с уточнённым распределением температуры и коэффициента теплоотдачи в междисковых полостях.

Библиографический список

1. Слитенко, А.Ф. Конструкция и расчёт систем охлаждения ГТД [Текст] / А.Ф. Слитенко, С.З. Копелев. - Харьков: Издательство “Основа” при Харьковском государственном университете, 1994. - 239 с.

2. Childs D.W. Turbomachinery Rotor-dynamics: phenomena, modeling, and analysis

[Текст] / D.W. Childs // John Wiley & Sons Inc., 1993, - 476 p.

3. Идельчик, И.Е. Справочник по гидравлическим сопротивлениям [Текст] / И.Е. Идельчик. - М.: Гос. Энерг. Издат., 1960.

- 464 с.

4. Дорфман, Л. А. Гидродинамическое сопротивление и теплоотдача вращающихся тел [Текст] / Л. А. Дорфман. - М.: Физмат-гиз, 1960. - 260 с.

5. Bohn D., Deutsch G., Simon B., Burkhardt C. Flow visualization in a rotating cavity with axial throughflow [ Текст] / D. Bohn, G. Deutsch, B. Simon, C. Burkhardt // ASME Turbo Expo 2005, paper № 2000-GT-280.

6. Kilfoil A.S.R., Chew J.W. Modeling of Buoyancy-Affected Flow in Co-rotating Disk

Cavities [Текст] / A.S.R. Kilfoil, J.W. Chew // Proceedings of ASME Turbo Expo 2009: Power for Land, Sea and Air GT2009, June 8-12, 2009, Orlando, Florida, USA / GT2009-59214. - 10 p.

7. Каданер, Я.С. Течение воздуха и конвективный теплообмен в междисковых полостях ротора компрессора [Текст] / Я.С. Каданер, Ю.М. Темис, Л.Ф. Гришина // Вопр. авиац. науки и техники. Авиац. двига-телестроение: ЦИАМ, 2001. - Вып. 1 (1320).

- С. 24-36.

8. Харьковский, С.В. Расчётное определение параметров среды в разветвленной системе воздухоподвода турбины и граничных условий теплообмена на поверхности диска [Текст] / С.В. Харьковский, В.К. Кос-теж // Труды ЦИАМ, 1269. - 1990. - С. 116128.

NUMERICAL SIMULATION OF THE THERMAL-STRESS STATE

OF AN HIGHPRESSURE COMPRESSOR ROTOR WITH CONSIDERATION

FOR SECONDARY FLOWS

© 2011 Yu. M. Temis, A. V. Selivanov, G. G. Yurchenko Central Institute of Aviation Motors

A method of calculating secondary flow on the basis of combining one-dimensional and two-dimensional flow models is developed. This method permits us to obtain more accurate data about the distribution of the heat-transfer coefficient and other secondary flow parameters. The data are used as boundary conditions for the calculation of the thermal-stress state of an HPC rotor.

Thermal-stress state, secondary flow, hydraulic net, heat-transfer coefficient.

Информация об авторах Темис Юрий Моисеевич, доктор технических наук, профессор, начальник отдела математического моделирования, ФГУП «ЦИАМ им. П.И. Баранова». E-mail: [email protected]. Область научных интересов: многодисциплинарное математическое моделирование; оптимизация, динамика и прочность конструкций; проблемы ресурса; численные методы; нелинейные проблемы пластичности, ползучести.

Селиванов Алексей Валерьевич, начальник сектора отдела математического моделирования, ФГУП «ЦИАМ им. П.И. Баранова». E-mail: [email protected]. Область научных интересов: многодисциплинарное математическое моделирование; аэроупругие колебания; газовая динамика; теплообмен.

Юрченко Геннадий Геннадиевич, аспирант Московского физико-технического института (государственного университета), инженер ФГУП «ЦИАМ им. П.И. Баранова». Е-mail: [email protected]. Область научных интересов: многодисциплинарное математическое моделирование; газовая динамика; теплообмен.

Temis Yuriy Moiseyevitch, doctor of technical sciences, professor, head of the department of mathematical simulation, Central Institute of Aviation Motors. E-mail: [email protected]. Area of research: multidisciplinary mathematical simulation, optimization, dynamics and strength of constructions, durability, numerical methods, nonlinear problems of plasticity and creep.

Selivanov Alexey Valeryevitch, sector leader of the department of mathematical simulation, Central Institute of Aviation Motors. E-mail: [email protected]. Area of research: multidisciplinary mathematical simulation, aeroelastic vibrations, gas dynamics, heat exchange.

Yurchenko Gennadiy Gennadyevitch, post-graduate student, Moscow Institute of Physics and Technology, engineer of the department of mathematical simulation, Central Institute of Aviation Motors. E-mail: [email protected]. Area of research: multidisciplinary mathematical simulation, gas dynamics, heat exchange.