CC BY
350
Вычислительные технологии
Том 18, № 4, 2013
Метод расчёта осевых и радиальных нагрузок
на рабочее колесо гидротурбины
>к
в нестационарном потоке*
А. Ю. АвдюшЕнко, С. Г. Чёрный, Д. В. Чирков Институт вычислительных технологий СО РАН, Новосибирск, Россия e-mail: ovalur@gmail.com, cher@ict.nsc.ru, Dchirkov@ngs.ru
Предлагается методика определения осевых и радиальных нагрузок на рабочее колесо гидротурбины, вызванных нестационарным течением рабочей жидкости в его межлопастных каналах, а также протечками в лабиринтных уплотнениях, полостях и разгрузочных отверстиях. Найдены расходы жидкости через указанные области течений. Данная методика превосходит по точности инженерно-эмпирические подходы, так как основана на моделях трёхмерных нестационарных турбулентных течений и позволяет учитывать окружную неравномерность потока, прецессию вихревого жгута за рабочим колесом и другие особенности трёхмерных нестационарных течений. Проведено усовершенствование инженерно-эмпирического метода расчёта радиальных нагрузок в лабиринтных уплотнениях, вызванных смещением оси вращения рабочего колеса. При этом учтены влияние вращения ротора на коэффициент сопротивления узкой части лабиринта, сопротивления ячеек расширения, зависимости коэффициента сопротивления узкой части лабиринта и ячеек расширения от переменного зазора между статором и ротором.
Ключевые слова: гидротурбина, нестационарные течения, лабиринтные уплотнения, радиальные и осевые нагрузки.
Введение
Между вращающимися и неподвижными частями гидротурбины (ГТ) имеются зазоры протечек жидкости (рис. 1). Для снижения объёмных потерь в зазорах устанавливаются лабиринтные уплотнения (ЛУ). Кольца лабиринтных уплотнений крепятся на верхнем и нижнем ободах или на камере рабочего колеса (РК). Ступица РК имеет разгрузочные отверстия, расположенные за выходными кромками лопастей, через которые полость над РК сообщается с его основным проточным трактом (ПТ). На рис. 2 приведена схема течения жидкости в проточном тракте ГТ, лабиринтных уплотнениях, полости между ступицей и крышкой турбины, разгрузочных отверстиях, полости между ободом и камерой РК.
Кроме уменьшения потерь мощности в турбине за счёт ограничения величин протечек воды между вращающимися и неподвижными частями, уплотнения обеспечивают также снижение осевых и радиальных нагрузок, действующих на РК [1].
* Работа выполнена при финансовой поддержке РФФИ (грант № 11-01-00475-а) и Интеграционного проекта СО РАН № 130.
Рис. 1. Рабочее колесо с зазорами, полостями Рис. 2. Основной проточный тракт и обла-и разгрузочными отверстиями сти протечек за ступицей и ободом рабочего
колеса гидротурбины
Полная гидравлическая нагрузка Е, действующая на РК, складывается из трёх составляющих:
Е — Ех + Е2 + Е3,
(1)
где Ех — нагрузка от основного потока на поверхность проточной части РК, Е2 — нагрузка от течения через верхнюю область протечки, Е3 — нагрузка от течения через нижнюю область протечки (см. рис. 2). В полной гидравлической нагрузке Е — (^Х, ¥у, ^) выделяются осевые ^ и радиальные — (¥х, ¥у) усилия. Индексами х, у, г обозначены компоненты сил, действующие по соответствующим осям декартовой системы координат. Ось О г совпадает с осью РК.
Основными источниками осевых и радиальных сил в нагрузке Ех являются неравномерность потока на выходе из спиральной камеры, ротор-статор взаимодействие между лопастями РК и лопатками направляющего аппарата (НА), прецессирующий вихревой жгут в конусе отсасывающей трубы (ОТ), гидроудар.
Существенный вклад в осевые и радиальные усилия вносят также составляющие Е2 и Е3, вызванные течениями в ЛУ, полостях и разгрузочных отверстиях.
Таким образом, для наиболее точного определения осевых и радиальных нагрузок (ОРН) следует проводить совместные расчёты полей во всех указанных областях течений основного ПТ и протечек за ступицей и ободом РК и определять по ним нагрузку Е (1). Однако такой расчёт потребует огромных затрат вычислительных ресурсов из-за необходимости его совместного проведения в областях, имеющих значительно различающиеся характерные размеры. Так, диаметр РК натурной гидротурбины составляет Бх ~ 6 м, а ширина зазоров d при этом варьируется в пределах 2-3 мм.
В настоящее время для предсказания ОРН используются приближённые методики, которые можно разделить на инженерно-эмпирические [2-4] и основанные на методах вычислительной гидродинамики [5-7].
1. Обзор существующих методик
1.1. Инженерно-эмпирическая методика определения осевых нагрузок
В этой методике расходы жидкости через верхнее и нижнее уплотнения находятся по формулам
Qs = /i(£b , Ahs ),
qh = /2 (în, AhH),
где îs ,îh — коэффициенты гидравлического сопротивления верхнего и нижнего уплотнений, приведённые к входному сечению уплотнения, A h s, AhH — напоры, срабатываемые в уплотнениях. Конкретный вид зависимостей /1, /2 приведён в [3].
Коэффициент гидравлического сопротивления для уплотнения ступицы складывается из коэффициентов сопротивления лабиринта и разгрузочных устройств
S2 S2 S2
îs = Свх SI + £лу + îp c 2 + Si,
SBX (NpSp) Sщ
где Si — площадь сечения узкой части зазора ЛУ; îBX — коэффициент гидравлического сопротивления щели входа в уплотнительный тракт ступицы в случае установки уплотнения не на периферии; £лу — коэффициент гидравлического сопротивления ЛУ; îp — коэффициент гидравлического сопротивления разгрузочных отверстий; Np — число разгрузочных отверстий; Sp — площадь сечения разгрузочного отверстия, — коэффициент гидравлического сопротивления выходной щели.
Коэффициент гидравлического сопротивления ЛУ вычисляется по формуле
Слу = îлу(А,îo),
где А = X(u,q), ш — частота вращения РК; î0 — коэффициент гидравлического сопротивления ячейки расширения, зависящий от геометрии ячеек и числа Рейнольдса.
Коэффициент гидравлического сопротивления вращающихся разгрузочных отверстий определяется соотношением
îp = îp(Wo,Vp,S'p ),
где W0 — относительная скорость жидкости перед разгрузочным отверстием, Vp — сред-нерасходная скорость воды в разгрузочном отверстии, Slp — относительная ширина полости над отверстием.
Коэффициент гидравлического сопротивления уплотнения обода определяется сопротивлением ЛУ, входа в ЛУ, если таковой имеется, и выхода из полости. Напор в уплотнении
A h hin hout hui,
где hin — напор на входе в уплотнение, hout — напор на выходе из уплотнения, hM — напор, создаваемый вращающимся объёмом жидкости, заключённым между ступицей и ободом и неподвижными крышкой турбины и камерой РК, hu — функция частоты вращения РК ш и геометрических размеров проточного тракта.
Так как коэффициенты гидравлического сопротивления лабиринтных уплотнений Слу и разгрузки îp определяются расходом жидкости через уплотнения, который в свою
очередь зависит от коэффициентов £в и £н, то расчёт протечек ведётся итерационным методом последовательных приближений.
Полная гидравлическая нагрузка на РК рассчитывается по формуле (1). Нагрузка воздействия потока на внутреннюю полость РК вычисляется по инженерно-эмпирическим формулам. Нагрузка , действующая на внешнюю часть ступицы РК, в общем случае складывается из трёх составляющих
Г2,г = Г2ст + Г2л + ^
где ^2ст — сила давления жидкости на часть поверхности ступицы, расположенную между входом в полость "ступица — крышка" и ЛУ (если ЛУ находится на периферии ступицы, то Гст = 0); Г2л — сила давления жидкости на ЛУ; — сила давления жидкости на часть поверхности ступицы, расположенную между ЛУ и уплотнением вала. Для однощелевых уплотнений с односторонней и двусторонней ячейками расширения Г2л = 0, а для уплотнений елочного и гребенчатого типов эта сила определяется по формуле [3]
Г2л = I (£в ,яв).
Силы Г2ст и Г2Ст находятся интегрированием пьезометрического напора по поверхности ступицы. Приближённо можно положить
Гст = 4 прв Мх - °2)
- £•*< й2 - С (в2х -
-,2
1 ( 2 2 шш (2 2 \
Г2ст = 4(^вых - ^ва^ "-вых - (^вых - ^вал)
Ш
16д
где Двх — диаметр входа в полость; Дл — диаметр ЛУ; Двал — диаметр уплотнения вала; Двых — диаметр выхода из уплотнения; Ш = ш(ш) — угловая скорость вращения жидкости; квх, квых — пьезометрический напор на входе и выходе.
Аналогично по инженерно-эмпирическим формулам определяется нагрузка на внешнюю часть обода РК Г3,г.
Отметим, что в силу предположения об осесимметричности потока в уплотнениях и повторяемости течения в каналах РК определение радиальных нагрузок по данной методике невозможно.
1.2. Инженерно-эмпирическая методика определения радиальных нагрузок
Для оценки радиальных сил, вызванных несоосностью ротора и статора для однощеле-вых уплотнений без ячеек расширения и без учёта вращения смещённого ротора, в [2] предложена инженерно-эмпирическая методика. В ней, как ив [8], предполагается, что падение давления на входе в уплотнение определяется выражением
VI (в)
^Рт = РА - Ргп = (1 + £гп)
2д
а изменение давления по всей длине уплотнения равно сумме потерь на входе и потерь от гидравлического трения по всей длине
рА - роиг = рА - ргп + ргп - роЫ = Аргп + ргп - роЫ =
. ^ \ VI (в)
1 + £гп + \ гп У 1
21* (в)) 2д
где £гП — 0.5 — коэффициент гидравлического сопротивления входа; Ут(в) — скорость на входе в ЛУ, зависящая от окружной координаты в; X — коэффициент гидравлического трения; I — длина уплотнения; d* (в) — переменный зазор. Авторы [2] не учитывали в (3) зависимости X от частоты вращения ротора ш и наличия ячеек расширения с их гидравлическим сопротивлением. Эти упрощения были сделаны с целью существенного сокращения времени расчёта нагрузки. При известном перепаде давления скорость У^п выражается из уравнения (3) в виде
Угп (0) =
\ 1 + А,. +
2d* (0)
и подставляется в уравнение (2). Следовательно, перепад давлений на входе
2g(pa - Pout)
Ml
1 + +
^Pin(O) =-1 + — (pa - Pout)
1 + ^in + 2d4T)
меняется по окружности. Принимается, что давление p по длине щели z изменяется по линейному закону
z
p(z, в) = PA - Apin - (PA - Pout - APin) j.
Наличие в уплотнении изменяющегося по окружности давления вызывает силу, действующую на ротор в радиальном направлении. Эта сила вследствие симметрии величины зазоров будет направлена по линии смещения центра ротора от оси уплотнитель-ного кольца и равна
l 2п
F = -j J P (z, в) Rr cos edzde, 0 0
где Rr — радиус ротора. При определении F не учитывается возможность частичного выравнивания неравномерности распределения давления p(z, в) по окружности из-за перетекания жидкости.
В настоящей работе данная методика усовершенствована путём учёта вращения ротора и ячеек расширения лабиринтных уплотнений (см. ниже раздел 4).
Существенный недостаток инженерно-эмпирических методик — грубое представление в них потока жидкости в основном проточном тракте, не учитывающее нестационарность потока и форму лопасти, следствием чего становится неточное задание давлений на входе и выходе зазоров, приводящее к неправильному определению осевой силы, действующей на РК.
В связи с этим в последнее время стали появляться работы, в которых моделирование течений в зазорах и полостях ГТ осуществляется на основе более совершенных методов вычислительной гидродинамики.
1.3. Подходы, основанные на методах вычислительной гидродинамики
В работе [5] в полной постановке рассчитывается трёхмерное течение в области НА — РК — ЛУ — полость над ступицей — разгрузочные отверстия — полость у обода. Во входном сечении направляющего аппарата задается расход, в выходном сечении РК — давление. Области РК, полости за ступицей и разгрузочных отверстий рассчитываются во
вращающейся системе координат, а область НА — в неподвижной. На стыке вращающихся и неподвижных областей осуществляется обмен всеми параметрами течения.
В рассматриваемых областях решаются нестационарные уравнения Рейнольдса с замыканием 88Т-моделью турбулентности. 88Т-модель турбулентности является синтезом моделей к — £ и к — ш. Она основана на том, что модели типа к — £ лучше описывают свойства свободных сдвиговых течений, а модель к — ш имеет преимущество при моделировании пристеночных течений. Плавный переход от к — ш-модели в пристеночной области к к — £-модели вдали от твёрдых стенок обеспечивается введением весовой эмпирической функции. Средняя величина нормализованного расстояния до стенки в расчётах [5] равна у+ — 70.
Приводятся распределения давления в ЛУ и по верхней части ступицы для трёх режимов работы ГТ: неполной загрузки, оптимального КПД и максимальной мощности, и распределения давления по верхней части ступицы при наличии разгрузочного устройства для тех же режимов. Полученные зависимости качественно адекватно отражают моделируемые явления.
В [6] исследуется применимость различных вычислительных моделей и сеток для расчёта течения в лабиринтных уплотнениях радиально-осевых турбин. Осевое число Рейнольдса для рассмотренных в [6] конфигураций ЛУ и режимов течения варьируется
в диапазоне Иеах1а1 — ах1а1 — 2 • 103 ^6.4• 103, число Тейлора Та — ^ 1 — 300^962,
V V
в связи с чем на основании рис. 3, взятого из работы [6], делается вывод, что течение в зазоре является турбулентным. Утверждается также [6], что при ширине зазора лабиринтного уплотнения d много меньшей радиуса рабочего колеса Ях вихри Тейлора — источник нестационарности — в нём не возникают. Поэтому течение в зазоре ЛУ моделируется стационарными уравнениями Рейнольдса в двумерном осесимметричном приближении. На входе в ЛУ задается либо расход, либо полное давление, на выходе — статическое давление, которые берутся из эксперимента. Для замыкания уравнений ис-
4000
3500
^ 3000 с*
ч 2500 ч
2000
<и рц
§ 1500 к
V
§ юоо-«
<и
о 500
А [>
□
х □
> <1
А д
X >
V
V
НУ V
□ < □
У
□А
V V
Зазор ЛУ
Ж
у Турбулентное течение
ж
ж
* Ж X вихревое течение Ь-оШг
?еТеГРН°7ЛаМИНа^е
_^—вихревое течешь
Турбулентное
Ж — 0.006" - Длина ЛУ 2", 2.5" + —0.010" —1 -х—0.015" А —0.020" <1 —0.030" _1
Длина ЛУ 2" Длина ЛУ 2.5"
50 100 150 200 250 300 Число Тейлора Та
350 400
Рис. 3. Экспериментальные данные о зависимости типа течения в прямом ЛУ от величины осевого числа Рейнольдса и числа Тейлора [6]
пользуется Ьс^-Ке-ЯЯТ-модель турбулентности без пристеночных функций с использованием подробной сетки, разрешающей пограничный слой (у+ ~ 1). При низких числах Рейнольдса толщина логарифмического пограничного слоя сопоставима с шириной зазора 1. Поэтому метод пристеночных функций, требующий расположения первого слоя сетки на расстоянии у+ ~ 30 ^ 300 от стенки, должен применяться с осторожностью. В [6] исследуется влияние сгущения сетки и режима течения, подбираются сетка, модель турбулентности и краевые условия для адекватного количественного описания течения в ЛУ. Проводится сравнение результатов расчётов с экспериментом и результатами, полученными по эмпирической модели.
Использование подходов работ [5, 6] позволяет построить вычислительные модели для расчёта объёмных потерь через уплотнения, осевых и радиальных нагрузок на РК при условии совпадения оси вращения колеса и оси симметрии неподвижной камеры РК. Однако эти подходы не годятся для моделирования радиальных сил, связанных с эксцентричностью рабочего колеса и камеры РК.
1.4. Моделирование течения в зазорах между вращающимися
и неподвижными частями гидротурбины с учётом биения РК
При наличии эксцентриситета оси РК (рис. 4) область течения в зазорах на ступице и ободе становится областью с подвижными границами, что при расчёте этих течений требует использования подвижных сеток, перестраивающихся на каждом шаге по времени. Это обстоятельство усложняет расчёт осевых и радиальных нагрузок. Однако в ряде случаев постановка задачи может быть упрощена.
В работе [7] численно моделируются течения в зазорах ЛУ паровой турбины. При этом рассчитываются радиальные силы, действующие на ротор со стороны потока при малых биениях ротора в радиальном направлении (в ^ 1). Предполагается, что биение оси ротора периодическое по времени и вместе с тем она перемещается в окружном направлении (прецессирует) с угловой скоростью ^ (см. рис. 4). В этом случае течение приближённо представимо в виде суперпозиции стационарного осесимметричного течения при фиксированном положении ротора в центре и малого периодически нестационарного несимметричного возмущения, порождаемого движением центра ротора относительно центра полости
в
Ф(х,Т,в,г) = Ф0(Х,Г) + - (Ф1с(х,Г)еОБ(Ш-в) - ФЬ(х,г)8т(Ш-в)) .
1
Рис. 4. Кольцевая щель лабиринтного уплотнения е осью вращения РК, смещённой относительно центра полости
Для нахождения первой, осесимметричной, составляющей течения Ф0(х,г) используются уравнения Рейнольдса, замкнутые стандартной к — е-моделью турбулентности. Для определения второй, несимметричной, составляющей применяются уравнения Рей-нольдса для возмущённых течений. На основе полученной несимметричной составляющей поля давления вычисляется радиальная сила, действующая на ротор. Далее найденная сила используется для вычисления коэффициентов динамической устойчивости ротора, которые показывают, насколько лабиринтное уплотнение демпфирует колебания вала в радиальном направлении. Расщепление на осесимметричную невозмущенную и несимметричную возмущенную части позволило в [7] находить эти составляющие решением двумерных осесимметричных задач.
Подход [7] может быть применён для вычисления радиальных гидродинамических сил, действующих на РК при наличии эксцентриситета и неподвижной поверхности уплотнений. В более общем случае можно использовать трёхмерное моделирование течения во всей кольцевой щели переменной ширины в неинерциальной системе отсчёта, вращающейся со скоростью прецессии вала П относительно центра камеры РК. При этом предполагается, что разгрузочные отверстия заменяются на кольцевую щель. В данном случае не требуется перестроения сетки.
2. Основные уравнения и численный метод их решения
Для моделирования течения жидкости в проточных частях гидротурбины используются нестационарные трёхмерные уравнения Навье — Стокса, осреднённые по Рейнольдсу, имеющие в декартовой системе координат (х^х2,х3) следующий вид:
диц
д'Шг д'Шг'Ш^
дЬ
дхп
+
др
д
дхг дхп
дх.
Vей
= 0, д'Шг
+
д™3
дхп дх
+ ¡г, г = 1, 2, 3,
4)
(5)
где (т1,т2,т3) — компоненты вектора скорости (м/с), р = — +— к, рс — гидростати-
р 3
ческое давление (Па), к — кинетическая энергия турбулентных пульсаций (м2/с2). При расчёте течения во вращающейся с угловой скоростью ш вокруг вертикальной оси Ох3 системе координат вектор массовых сил / = (/1, /2, /3) помимо силы тяжести включает также центробежную и кориолисову силы: / = (х1ш2 + 2и>2ш, х2ш2 — 2,ш1ш,д). Координату х3 будем обозначать буквой г, а давление р — делить на д и измерять в м вод. ст. Величина иед есть сумма молекулярной V и турбулентной VI вязкости
Veff = V + VI.
(6)
к2
Для определения турбулентной вязкости ^ = Сп— система (4)-(5) замыкается стане
дартной к — е-моделью турбулентности.
Численный алгоритм решения уравнений (4)-(5) основан на методе искусственной сжимаемости, который состоит во введении производной по псевдовремени т от давления в уравнение неразрывности и производных по т от соответствующих компонент скорости в уравнения сохранения количества движения. Полученные уравнения дискре-тизуются неявным методом конечных объёмов с использованием МиБСЬ-схемы третьего порядка аппроксимации для конвективных членов и центрально-разностной схемы
второго порядка для вязких членов. Производные по физическому времени Ь аппроксимируются направленными разностями второго порядка, а введённые производные по псевдовремени — направленными разностями первого порядка. При этом на каждом шаге по Ь проводится установление решения по т. Полученная в результате дискретизации система уравнений линеаризуется методом Ньютона. Система линейных алгебраических уравнений на приращение неизвестных при переходе со слоя в на слой в + 1 по псевдовремени приближённо Ш-факторизуется и решается бегущим счётом. Детали данного численного алгоритма изложены в [9].
3. Методика определения полной гидравлической нагрузки на рабочее колесо
3.1. Определение нагрузки Ех от основного потока
Для определения гидравлической нагрузки Ех, действующей на РК, рассчитывается нестационарное течение в основной проточной части гидротурбины, включающей спиральную камеру, все каналы направляющего аппарата, все межлопастные каналы РК и отсасывающую трубу. В качестве граничного условия во входном сечении спиральной камеры задаётся равномерное распределение скорости, перпендикулярной сечению
IV* | — О-,
где Q — фиксированный расход жидкости, Бгп — площадь входного сечения.
Для определения абсолютных значений давления в точках А, В, С, О (см. рис. 2) рабочего колеса, требуемых для постановки граничных условий при расчёте течения жидкости через области протечки, нужно задать абсолютное давление в выходном сечении ОТ
Р(г) — Рагш + ^ — ¿н.б), (7)
где ра1т — атмосферное давление, гн.б — уровень нижнего бьефа. Полная нагрузка Ех, действующая на рабочее колесо, не зависит от константы, входящей в распределение давления (атмосферное давление). Но мы не можем посчитать вклад в Ех атмосферного давления, действующего сверху на вал гидротурбины. Поэтому для корректного расчёта силы Ех необходимо проводить вычисления при р^т — 0. Этот же подход использован в работе [3].
3.2. Определение осевой составляющей нагрузки Е2
Область протечки над ступицей включает (рис. 5) верхнее лабиринтное уплотнение, полость между крышкой турбины и внешней поверхностью ступицы, разгрузочные отверстия (РО), кольцевую полость и щель, через которую жидкость попадает за лопасти РК. Разгрузочные отверстия служат для снижения давления на крышку турбины.
Из-за разгрузочных отверстий течение в полости не является осесимметричным. Поэтому все течение в верхней области протечки разбивается на два. Течение в самом ЛУ считается в осесимметричной постановке. Область ЛУ содержит один слой сетки по окружному направлению. Трёхмерное течение в остальной части предполагается циклически повторяющимся по числу разгрузочных отверстий и рассчитывается в секторе полости между крышкой и ступицей с одним РО и в секторе кольцевой полости.
Рис. 5. Верхняя область протечки: а — расчётная область периодического трёхмерного течения в полости и осесимметричного течения в ЛУ; б — сетка в одной ячейке ЛУ для расчёта в ней осесимметричного течения; в — г г-проекция верхней области протечки
При передаче параметров потока между этой областью и областью ЛУ проводится их усреднение по окружному направлению.
Расход воды через области протечки априори не известен и должен быть найден в ходе решения задачи.
Поскольку в предлагаемом методе определения нагрузки на РК мы отказались от одновременного совместного расчёта течений в основном проточном тракте и в областях протечек, то возникает проблема обмена между ними параметрами течений в точке А, расположенной у входного сечения Бгп в ЛУ, и в точке В у выходного сечения из кольцевой щели (см. рис. 5, в). В этом обмене исключено влияние течения в зазорах на течение в основном ПТ. В построенной методике течение в зазорах определяется течением в основном ПТ. Для сопряжения этих течений необходимо учесть входные потери энергии в точке А и выходные в точке В. Согласно [8], потери полной энергии на входе в лабиринтное уплотнение могут быть найдены по формуле
(ДЕ)ш = ЕА - Егп = £гп I 2; I ,
(8)
где £гп — коэффициент гидравлического сопротивления входа в ЛУ полагается равным £гп = 0.5. Поскольку полная энергия
Е = р - г +
И
2 9
(9)
то из (8) следует
Ра - гл +
2
2 9
ргп ггп +
|Угп I2 2 9
= &
|Угп I2 2 9
(10)
а
в
Считается [3], что течение в лабиринтном уплотнении не зависит от скорости жидкости в основной части ГТ, а зависит только от давления. Поэтому в (10) полагаем Уа = 0,
а также считаем, что zA = zin. Тогда
PA - Pin--2>f~ = 2^' ( )
откуда получаем
'^'Чт® (12)
Таким образом, во входном сечении Sin задаётся условие на модуль скорости, из которого находится расходная z-я составляющая, при этом тангенциальные составляющие скорости полагаются нулевыми. Давление pin в соотношении (10) находится в процессе решения задачи.
Течение в выходной щели K (рис. 6) не рассчитывается. Поэтому для задания условия в сечении Sout необходимо связать параметры потока в этом сечении с параметрами потока в точке B. Потери энергии на выходе из кольцевой полости (сечение Sout) складываются из потерь на входе в кольцевую щель
IV J2
771 771 £ ' out' iTQ\
Eout — EK — Qin-^--(13)
2g
и потерь на выходе из неё в область основного ПТ (см. рис. 6)
|voui |
2
Ek — EB — , (14)
2g
где, согласно [8], коэффициент гидравлического сопротивления входа £in — 0.5, коэффициент гидравлического сопротивления на выходе из кольцевой щели £out — 1. Тогда, складывая соотношения (13), (14) и подставляя выражения для полной энергии, получим
, Vout'2 ( JVb|Л u ^с , Vout'2 п^
Pout — Zout --( PB — ZB + ^^ j — (£in + ^ut) 2g • (15)
Рис. 6. К заданию давления в выходном сечении
Положив в (15) VB = 0 и zB = zout, находим
Pout = PB +
IVout I
2 g
(16)
Таким образом, скорость "1/Ои находится из основных уравнений в процессе их решения, а давление вычисляется по формуле (16), где рв предварительно определено из расчёта потока в основной части ГТ.
3.3. Определение радиальных составляющих нагрузок ¥2 и Е3
Силы Е2,# и Е3,# возникают главным образом вследствие окружной неравномерности потока в лабиринтных уплотнениях при смещении центра вращения РК относительно геометрического центра камеры рабочего колеса.
В общем случае полагаем, что имеются смещение оси РК на эксцентриситет е и её прецессия вокруг центра статорной части с угловой скоростью ^ (см. рис. 4). Течение в кольцевом зазоре переменной ширины приводит к появлению радиальных сил, действующих со стороны протекающей в зазоре жидкости. Используются две методики расчёта радиальных составляющих нагрузок Е2 и Е3 в модели с эксцентриситетом. Первая основана на расчёте трёхмерного течения в несимметричной щели, вторая — усовершенствованный инженерно-эмпирический подход [2].
Расчёт трёхмерного течения в ЛУ с эксцентриситетом. Если для описания течения в кольцевой щели лабиринтного уплотнения с осью вращения РК, смещённой относительно центра полости, перейти во вращающуюся вокруг центра полости систему координат (х', у') (см. рис. 4), то относительно этой системы течение будет стационарным. В самом деле, форма области зазора в координатах (х', у') не меняется, сетку можно считать неподвижной и необходимо учесть только вращение ротора с угловой скоростью ^ — и и статора — со скоростью Компоненты скорости движения точки (хрк,Урк) на роторе будут следующими:
u
= yPK (Q — и), v = — (x'PK — e) (Q — и), w = 0
а точки (хСт, уСт) на статоре — следующими:
u = Уст Q, v = — х'ст Q, w = 0.
(17)
(18)
Рис. 7. К определению давления рои
2
При этом в уравнениях количества движения (4)-(5) в источниковых членах величину и следует заменить на О. Условия для входного и выходного сечений выводятся с учётом входных и выходных потерь в полном соответствии с данными раздела 3.2. Связь скорости и давления на входе задаёт соотношение (12), из которого находится расходная г-компонента скорости , остальные компоненты скорости полагаются нулевыми.
Выражение для давления в выходном сечении (рис. 7) выводится из потерь энергии
С учётом того что VL = 0, zL = zout и £out = 1, находим pout = pL. Давление pL берется из расчёта течения в верхней области протечки (см. раздел 3.2).
4. Улучшенная инженерно-эмпирическая методика определения радиальных составляющих нагрузок F2 и F3
Методика [2] позволяет вычислять действующие на ось РК радиальные нагрузки, вызванные наличием эксцентриситета (см. рис. 4). Она применяется только для щелевых уплотнений (без ячеек расширения) для неподвижного РК и нулевой скорости прецессии вала Q = 0. Ниже этот подход распространяется на случай лабиринтных уплотнений с ячейками расширения и вращающимся рабочим колесом.
Как ив [2], предполагается, что падение давления на входе в ЛУ (рис. 8) определяется выражением
ТТЛ ТТЛ _ £ 1 Voutl
E out — EL = £out-~-
(19)
(20)
•P.
Рис. 8. Давления pa, Pin, Pout и коэффициенты сопротивления верхнего ЛУ
Для изменения давления по всей длине ЛУ вместо (3) используем формулу из [8], в которой А зависит от частоты вращения ротора и, учитываются ячейки расширения и коэффициенты их гидравлического сопротивления. Тогда
Рл — Рсаг = Рл — Рт + Рт — Роаг =
А (в, и) 1уз , ,„Л V? (0)
= Apin + pm — Pout = + + 2d* (в) + nCc (в) J ^n^, (21)
где £in = 0.5 — коэффициент гидравлического сопротивления входа; V¿n(0) — скорость на входе в ЛУ, зависящая от окружной координаты в; X (в,ш) — коэффициент гидравлического трения; 1уз — суммарная длина узкой части (щели без ячеек) уплотнения; d* (в) — переменный зазор; n — число ячеек расширения; £0 (в) — коэффициент гидравлического сопротивления одной ячейки.
Найдём зависимость величины зазора от окружного угла d* (в). Пусть Rs — радиус внутренней стенки статора, Rr — радиус ротора. Тогда (см. рис. 4)
Rs = e cos (в — y) + (Rr + d*) cos y.
При эксцентриситете e << Rr угол y << 1, поэтому можно считать
Rs = e cos в + R r + d* .
Тогда
d* = Rs — Rr — e cos в = d — e cos в. Окончательно переменный зазор ЛУ
d* (в) = d (1 — е cos в),
где d = Rs — Rr — радиальный зазор, е = e/d — относительный эксцентриситет, отсчёт угла в идет от самой узкой части зазора. Из соотношения (21) находим
ЛЛ (а\ _(pA — pout,l)__, ,
Vin (в) =-Vm—w-. (22)
X [в,Ш) lV3 „ ч
1 + ^ + -bFW + ^(в)
В уравнении (22) коэффициент гидравлического трения кольцевой щели X (в,ш) определяется по формуле [8]
X (в,") = к°5тву
где
2-i 0.375
i + . Reu (в,ш)х + 1 Re (в)
(23)
Re = К, (O)2^, Re„ = ^, (24)
v — коэффициент кинематической вязкости.
Коэффициент гидравлического сопротивления ячеек расширения лабиринта £0 (O) определяется по формуле [8]
Со = a (O) + 0.5b (O), (25)
где а — коэффициент потерь энергии постоянной массы в начале и в конце ячейки, Ь — коэффицент потерь энергии от внезапного сжатия ядра постоянной массы при втекании в щель. Значения этих коэффициентов определяются в зависимости от соотношения величин
¿Тр = 6* (9) + 0.24£ (26)
и глубины камеры ячейки Ни (см. рис. 8). Если Ни > 6гр, то а и Ь находятся из эмпирических зависимостей (рис. 9), если Ни < 6гр, то по формулам
а = 1 - ^, Ь = (\ - ^ Нк V Нк
(27)
Далее, так как коэффициент гидравлического трения кольцевой щели Л (9, и) зависит от расходной составляющей скорости в щели —п (9) (23)-(24), которая в свою очередь зависит от коэффициента гидравлического трения Л (9, и), то расчёт этих двух величин ведётся итерационным методом последовательных приближений.
Таким образом, в соответствии с (21) распределение давления р(г,9) по длине ЛУ
г Е [0,1] состоит из потерь давления на трение
Л (9, и) /УЛ V? (9)
26* (9)
ке узкой части уплотнения и скачкообразных потерь (9)
2д
-Ш
2д
на каждом участ-
на входе и выходе
каждой ячейки расширения. На рис. 10 приведена зависимость р(г, 0) для значения относительного эксцентриситета £ = 0.5. Конфигурация ячеек расширения показана над осью абсцисс. При построении зависимости скачкообразные потери на каждой ячейке расширения разделены на входные и выходные составляющие
а(9)^^ и 0.5Ь(9)-п(9)
2д
2д
и отнесены к началу и концу ячейки соответственно. Поэтому на внутренней части ячейки давление не меняется.
0 10 30 50 70 90 Ш*(6) Рис. 9. Определение коэффициентов а и Ь (25) при Ки > ^гр
2
Рис. 10. Распределение давления вдоль длины ЛУ
Радиальную составляющую нагрузки F2 находим численным интегрированием найденного распределения давления p(z, в) по поверхности ротора:
l 2п
\F2,R| = - J Jp (z, в) Rr cos edzde. (28)
0 0
Эта радиальная сила при ш = 0 в силу симметрии зазоров будет направлена по линии смещения центра ротора к оси статора.
Радиальная составляющая нагрузки F3 определяется аналогично, отличаются только геометрические размеры и число ячеек расширения нижнего лабиринтного уплотнения.
В случае прецессии вала РК с ненулевой скоростью Q этот подход определения радиальной силы не применим.
5. Результаты расчётов
Разработанный метод использован для расчёта основных радиальных нагрузок на РК в двух режимах работы гидротурбины, приведённых в табл. 1. Для режима частичной нагрузки характерно наличие прецессирующего вихревого жгута в конусе отсасывающей трубы, влияние которого на гидродинамические силы в проточном тракте представляет практический интерес. Расчёты проведены в натурных параметрах, шаг по физическому времени Д£ = 0.003126 с соответствует повороту РК на угол в = 3.75°.
Таблица 1. Режимы работы гидротурбины, в которых рассчитывались основные радиальные нагрузки
Режим Открытие НА а0, мм Частота вращения РК п, об/мин Напор Н, м Расход Я, м3/с
Максимальная мощность 219.13 200 73.5 87.980
Частичная нагрузка 123.26 200 73.5 50.615
Расчётная сетка содержит суммарно по рассмотренным спиральной камере, НА, РК и ОТ 2 615 925 внутренних ячеек.
5.1. Режим максимальной мощности 5.1.1. Сила
На рис. 11 представлены рассчитанные в полной постановке зависимости компонент силы Ех от времени. Период пульсаций Т = 0.01875 с, их частота /ь = 1/Т = 53.28 Гц. При этом частота вращения РК /п = 3.33 Гц. Отношение ///п = 16 совпадает с количеством лопастей РК, т. е. наблюдаемые пульсации обусловлены лопастной системой РК и имеют в своей основе продольную (вдоль потока) природу. Та же самая частота наблюдается для пульсаций радиальных компонент и нагрузки Ех. Среднее
значение модуля радиальной силы |Е1?д| = у + = 3077 кГс. Амплитуда колебаний в данном режиме пренебрежимо мала и составляет 1 %.
5.1.2. Осевая составляющая
Расчёт осевой составляющей проводился в соответствии с постановкой, описанной в разделе 3.2. В частности, давления рА перед входной границей расчётной области и рв за выходной границей (см. рис. 5) взяты из результатов расчёта потока в полной поста-
Рис. 11. Пульсации компонент нагрузки Ех в режиме максимальной мощности: 1 — , 2
, 3 —
новке. С целью сравнения с результатами, полученными по инженерно-эмпирической методике [3], проведён расчёт с давлениями рА и рв, взятыми из этой работы (табл. 2).
Расход и осевая сила , найденные в расчётах, сравниваются с данными работы [3] в табл. 3. Видно, что в нашем расчёте при одинаковых с [3] давлениях рА и рв расход получается на 12% ниже, а сила хорошо совпадает с результатом [3]. На рис. 12, а показано распределение давления во всей расчётной области, б — представлена картина течения в полости над ступицей и в разгрузочных отверстиях. На рис. 13 приведены распределения давления (а) и отношения окружной компоненты скорости си к шЯг (б)
Таблица 2. Давления ра и рв для расчёта течения в верхней области протечки
Давление Инженерно-эмпирическая Расчёт
методика [3] в полной постановке
Ра 54.570 52.465
Рв -4.100 -0.535
Ра - Рв 58.670 53.000
Таблица 3. Характеристики течения в верхней области протечки
Расчёт с ра и рв,
Параметр Данные [3] Расчёт с ра и рв из [3] полученными в полной постановке
Расход через уплотне- 0.295 0.258 0.239
ние Q, м3/с
Осевая сила , кГс 34546 33113 47694
Рис. 12. Распределение давления в верхней области протечки (а) и линии тока вблизи разгрузочного отверстия в относительном движении (б)
a
б
14
1.2
2
1
0.2
0.......
0.6 0.7 0.8 0.9
1.1 1.2 1.3П 1.4
R,m
Рис. 13. Давление (а) и распределение отношения cu/uRr (б) в полости между ступицей и крышкой турбины в режиме максимальной мощности: 1 — расчёт при pa и рв, взятых из [3]; 2 — расчёт при pa и рв из полной постановки; □ и серый цвет — данные [3]
вдоль средней линии полости. Результаты, полученные в расчёте, сравниваются с данными, полученными по [3]. Видно, что предположение о квадратичной зависимости давления от расстояния до оси РК, принятое в инженерно-эмпирической методике [3], не выполняется.
5.1.3. Радиальные составляющие F2,R и F3,R
Рассмотрена несоосность статора и ротора без прецессии (Q = 0) и с прецессией (с заданной частотой Q) вала ротора. В последнем случае исследован наиболее вероятный вариант — со скоростью вращения рабочего колеса Q = и, соответствующий изгибу вала ротора.
При расчёте радиальной составляющей F2,R проведено сравнение подхода (см. раздел 3.3), основанного на трёхмерном расчёте всего ЛУ, и подхода (см. раздел 4), реализующего улучшенную инженерно-эмпирическую методику. Результаты расчётов F2,R в случае непрецессирующего ротора (Q = 0) представлены на рис. 14. Видно, что модифицированная методика и трёхмерная модель течения жидкости дают близкие между собой результаты, отличающиеся от данных расчёта по методике [2].
Давления перед входом в верхнее ЛУ pA = 54.57 м, за выходом из ЛУ pout = 12.16 м взяты из расчёта в разделе 5.1.2. Дополнительные входные данные, необходимые для определения радиальной нагрузки, следующие: радиус ротора Rr = 1.38225 м, постоянный зазор щели ЛУ d = 0.002 м, длина ЛУ l = 0.11 м, число ячеек расширения n = 4.
Для трёхмерного расчёта требуется 2-3 ч, для улучшенной методики — 2-3 мин, что позволяет рекомендовать последнюю для практической оценки модуля радиальной силы.
Результаты трёхмерных расчётов показали, что вращение вала РК и скорость его прецессии оказывают существенное влияние на направление радиальной силы. Этот
факт находится в соответствии с результатами работы [10]. Если РК не вращается (и = 0), то = (ЕХ, 0). Однако при вращении РК с угловой скоростью и, согласно [10], в уплотнении помимо компоненты ¥Х радиальной силы, пропорциональной перепаду давления, возникает дополнительная гидродинамическая сила, перпендикулярная
Рис. 14. Зависимости модуля радиальной силы |Р2,#| от относительного эксцентриситета £ в случае непрецессирующего ротора с невращающимся (а) и вращающимся (б) РК: о — инженерно-эмпирическая методика [2], Д — улучшенная инженерно-эмпирическая методика, □ — трёхмерный расчёт
Рис. 15. Зависимости —^2,Х и от £ в расчёте трёхмерного течения в ЛУ: Д — неподвижное РК (и = 0, ^ = 0), □ — РК, вращающееся со скоростью и без прецессии вала (^ = 0), о — РК, вращающееся со скоростью и с прецессией вала (^ = и)
Рис. 16. Пульсации осевой и радиальных Рис. 17. Пульсации давления в точках А и В , составляющих нагрузки Е1
оси Ох и пропорциональная угловой скорости вращения РК (компонента ¥2у, рис. 15). Кроме того, вращение РК приводит к появлению радиальной силы, направленной в сторону увеличения эксцентриситета, уменьшающей компоненту ¥2х радиальной силы. Наличие прецессии ^ существенно меняет направление радиальной силы Е2^. За счёт возникающего при ^ = и вязкого сопротивления поступательному перемещению центра вала компонента ¥2у меняет знак (см. рис. 15).
5.2. Режим частичной нагрузки
Для данного режима характерно наличие прецессирующего вихревого жгута в конусе отсасывающей трубы, оказывающего влияние на ГТ вверх по потоку. Частота прецессии вихря / = 0.74 Гц явно выделяется в пульсациях осевой и радиальной нагрузок, действующих на РК (рис. 16). В режиме неполной нагрузки присутствует сильная динамическая составляющая сил , Е^^, вызванная прецессией вихревого жгута и имеющая для величину порядка 10 %, для Е^^ — 100 % от среднего значения. Частота динамических составляющих равна 0.222/. На рис. 17 показаны пульсации давления в точках А и В в режиме частичной нагрузки. При расчёте нагрузок Е2 и Е3 давления рА и рв усредняются по времени.
Заключение
Разработан метод определения нестационарных осевых и радиальных нагрузок на рабочее колесо, основанный на расчётах трёхмерного течения несжимаемой жидкости в основной проточной части и в областях протечек гидротурбины, позволяющий рассчитывать осевые и радиальные нагрузки, вызванные следующими факторами: неравномерностью потока в спиральной камере и статоре, ротор-статор взаимодействием, неравномерностью потока за РК вследствие нестационарности потока в конусе ОТ (вихревой жгут и т.д.), несоосностью статора и ротора, изгибом вала ротора. Кроме того,
метод позволяет рассчитывать расход жидкости через зазоры между вращающимися и неподвижными частями гидротурбины. Потери давления в лабиринтных уплотнениях, щелях и разгрузочных отверстиях вычисляются непосредственно с использованием турбулентной модели, тем самым не требуется эмпирическая информация о гидравлических сопротивлениях этих элементов.
Показаны хорошее соответствие рассчитанных ОРН и расходов через области протечки с инженерно-эмпирическими данными и преимущества предложенного метода.
Усовершенствована инженерно-эмпирическая методика [2] для расчёта радиальных нагрузок, действующих на лабиринтное уплотнение вследствие смещения оси вращения ротора. В предложенной методике учтены влияние вращения ротора на коэффициент сопротивления узкой части лабиринта, сопротивление ячеек расширения, зависимость коэффициента сопротивления узкой части и ячеек расширения от переменного зазора между статором и ротором. Усовершенствованная методика применима для расчёта радиальных нагрузок в случае несоосности статора и ротора и отсутствия прецесии последнего. Результаты расчётов радиальной силы по предложенной методике хорошо совпадают с данными трёхмерных расчётов течения во всем ЛУ для диапазона рассмотренных значений относительного эксцентриситета е от 0 до 1.
Трёхмерный расчёт всего лабиринтного уплотнения позволяет также определить радиальную силу при ненулевой скорости прецессии ротора П, в частности, при П = ш, соответствующей изгибу вала. Показано, что при е < 0.5 прецессия ротора не оказывает заметного влияния на модуль радиальной силы. Влияние прецессии растёт по мере увеличения эксцентриситета. Так, при е = 0.9 в случае П = ш радиальная сила на 30 % больше, чем в случае простой несоосности ротора и статора.
Показано, что при относительном эксцентриситете е > 0.6 радиальные силы, действующие на лабиринтные уплотнения, дают более 50 % итоговой радиальной нагрузки, действующей на всё рабочее колесо гидротурбины.
Список литературы
[1] Грлновский С.А., Малышев В.М., Орго В.М., Смоляров Л.Г. Конструкции и расчёт гидротурбин. Л.: Машиностроение, 1974. 408 с.
[2] Бедчер Ф.С., Ломакин А.А. Определение критического числа оборотов ротора насоса с учётом сил, возникающих в уплотнениях // Паро- и газотурбостроение. Тр. ЛМЗ. Вып. 5. М., Л.: Машгиз, 1957. C. 249-269.
[3] Макаров В.В., Пылёв И.М., Пьянов В.И. Объёмные и дисковые потери в радиально-осевых гидротурбинах // Энергомашиностроение. 1982. № 1. С. 11-15.
[4] Кузьминский С.С., Фёдорова И.И., Пылёв И.М. Методические вопросы определения давлений и осевой силы в радиально-осевых гидротурбинах // Там же. 1979. № 3. С. 6-9.
[5] Staübli T., Sallaberger M., Senn F. Parameters to adjust axial thrust // Proc. of Intern. Conf. HYDRO2006. Porto Carras, Greece, 2006.
[6] Le Roy V., Güibaült F., Vu T. Validation of CFD model for hydraulic seals // Intern. J. of Fluid Machinery and Systems. 2009. Vol. 2, No. 4. Р. 400-408.
[7] Xi J., Rhode D.L. Rotordynamics of turbine labyrinth seals with rotor axial shifting // Intern. J. of Rotating Machinery. 2006. Vol. 2006. P. 1-11.
[8] Идельчик И.Е. Справочник по гидравлическим сопротивлениям. М.: Машиностроение, 1992. 672 с.
[9] Чёрный С.Г., Чирков Д.В., Ллпин В.Н. и др. Численное моделирование течений в турбомашинах. Новосибирск: Наука, 2006. 202 с.
[10] Млрцинковский В.А. Гидродинамика и прочность центробежных насосов. М.: Машиностроение, 1970. 296 с.
Поступила в 'редакцию 7 мая 2013 г.
