Том XXXVI
УЧЕНЫЕ ЗАПИСКИ ЦАГИ 2 00 5
№ 3 — 4
УДК 629.735. 45.015.3.035.62
СВОБОДНЫЙ ВИХРЕВОЙ СЛЕД В НЕЛИНЕЙНОЙ ДИСКОВОЙ ТЕОРИИ НЕСУЩЕГО ВИНТА В КОСОМ ПОТОКЕ
В. М. ЩЕГЛОВА
Для винта с бесконечным числом лопастей с использованием понятия о плоском движении части поверхности свободных вихрей проведены расчеты для приближенной оценки качественной картины состояния потока в следе несущего винта. Учитывалась толщина сходяшдх с диска свободных вихрей для двух случаев. В первом из них толщина вихрей, будучи заданной, не изменялась по времени, во втором — была учтена простейшая диффузия. Приведены результаты расчетов по данному методу формы свернувшихся жгутов для режимов полета ц = 0.0778 и ц = 0.326 с учетом как переменного по времени, так и постоянного параметра толщины вихрей.
Экспериментальные, а также расчетные исследования поведения вихревой поверхности представляют определенный интерес для выяснения качественной физической картины обтекания не только лопастей винта, но и остальных частей вертолета. В следе за несущим винтом, или вблизи него постоянно находятся конструктивные части вертолета, попадающие
в неблагоприятные условия обтекания. Поэтому важное значение имеет знание характеристик поведения свернувшейся вихревой поверхности в пространстве. В приложении к несущему винту представляется интересным определение местоположения вихревых шнуров и темпов наступления конечной стадии их развития в зависимости от нагрузки на диск несущего винта.
По лопастной вихревой теории уравнения, определяющие движение имеющихся в начальный момент в жидкости вихрей и процесс создания новых, приведены ниже [1], [2]. Для индуктивных скоростей от всех вихрей в точке с(г, т, t) на пелене Р (рис. 1) имеем:
таким
образом
(1)
где I = с (г, т, ґ)-с (гт', ґ), Ї1 = с (г, т, ґ)-В (гиґ),
параметры со штрихами относятся к точкам с расположением вихрей, от которых считаются
индуктивные скорости, без штрихов — к точкам, где эти скорости считаются.
Для определения вектора с имеем:
Рис. 1. Вихревая пелена за диском несущего винта
Стт= = в(г,и, t),
(3)
в выражении (2) V — скорость невозмущенного потока.
В этих уравнениях система параметров г и т характеризует положение точек поверхности Р в момент времени t, при этом считается, что параметр т есть время схода каждой жидкой точки с задней кромки. На лопасти положение точки задается величинами г и и. Циркуляция пелены определится условием
Г(г, т, t) = Г(г, и, t). (4)
Величина Г(г, и, t) задается в точках лопасти с векторами I (г, и, t), а движение поверхности Ь характеризуется радиусом-вектором ее точек В (г, и, t). Радиус-вектор точек пелены Р обозначен через вектор с(г, т, t).
Выражения (1) — (4) сохраняют силу при любых деформациях пелены вихрей Р. Например, несущий винт с бесконечным числом лопастей может быть рассмотрен как частный случай рассмотренной выше несущей поверхности. Для этого нужно взять поверхность Ь в виде диска, покрытого радиально расходящимися присоединенными вихрями, и вращающегося относительно перпендикуляра, проходящего через его центр (см. рис. 1). Центр же этого диска движется со скоростью полета вертолета.
Между конечным множеством точек с фиксированными г и т, находящимися на форме Р, примем какой-либо закон для аппроксимации, например, метод Эйлера. Тогда совокупность уравнений (1) и (2), написанных для каждой такой точки, образует систему связанных обыкновенных дифференциальных уравнений с начальными условиями типа Коши, решение которой может быть получено численными методами.
Вследствие трудности построения методов численного решения, а также громоздкости такого решения при определении формы пелены вихрей за несущим винтом вертолета в целом приходится вводить упрощения, направленные на ускорение вычислений, которые обычно сводятся к построению той или иной модели непрерывного распределения вихрей.
В настоящей работе рассматривается частный случай исследования системы уравнений (1) — (4), когда движение части поверхности свободных вихрей можно трактовать как плоское. Хотя такой подход нельзя использовать для получения данных о форме следа в целом, однако его можно применить для изучения деформаций частей пелены в местах сильного сгущения вихрей. Поскольку на интересующее здесь сворачивание пелены вихрей влияют скорости, вызываемые каждой вихревой нитью только в малой ее окрестности, указанная замена отрезков вихрей бесконечными прямолинейными вихрями представляется допустимой. Замена пространственного движения плоским в случае винта вертолета означает, что криволинейные вихри лопасти заменяются прямолинейными. Эта система вихрей переходит как бы в систему вихрей крыльев. Таким образом задача определения движения вихрей сводится к определению движения цепочки плоских вихрей в некоторой плоскости 2оу под действием скоростей от них самих. Подобный прием уже использовался ранее применительно к исследованию сворачивания пелены вихрей за крылом и лопастью [1].
Большинство прикладных методов расчета обтекания лопастей винтов базируется на основе понятия о бесконечной тонкости вихревой трубки, что объясняется относительной простотой расчетов вызываемых ими скоростей. Это приводит в ряде случаев к нежелательным последствиям. Появляющиеся трудности в виде возникновения разного рода особенностей в проведении расчетов полей индуктивных скоростей преодолеваются введением в расчеты вихрей с конечными поперечными размерами, что и проделано в данной работе. Преимущества простоты вычисления скоростей от тонких нитей сохраняются и исчезают недостатки типа особенностей.
Постановка задачи. Рассматривается ^-лопастный несущий винт радиуса Я, который вращается с постоянной угловой скоростью вращения ю. Винт работает в поступательном потоке, скорость которого V и его плоскость вращения составляют со скоростью полета угол а.
Рис. 2. Схема, используемая при выводе выражения для определения интенсивности вихорьков, сходящих с контура диска несущего винта
0X2 — плоскость вращения и у направления оси ох (рис. 2). На винт набегает поток со скоростью
Для упрощения анализа движения системы вихрей за винтом принимаются некоторые допущения.
Пространство заполнено несжимаемой жидкостью. При движении несущего винта вперед в качестве рабочей модели несущего винта рассматриваем вихревую систему за диском винта в виде наклонного вихревого цилиндра, уходящего в бесконечность, основанием которому служит диск несущего винта. Образующий винтовую поверхность концевой вихрь интенсивности О, сходящий с конца вращающейся лопасти, лежит на по-верхности этого эллиптического цилиндра. Ось цилиндра совпадает с направлением скорости набегающего потока. Винт с конечным числом лопастей заменен винтом, эквивалентным ему, с бесконечно большим числом бесконечно узких лопастей. Поперечные свободные вихри во внимание не принимаются.
Определим погонную интенсивность
цилиндрического слоя вихрей, сходящих с этого бесконечно большого числа лопастей и направленных по образующей цилиндра.
Пусть оху2 связанные со втулкой оси, причем азимут лопасти, отсчитываемый от отрицательного
V = -iVx - jV'y, Vx = V cos a, V'y = V sin a + V
(5)
Уравнения вихревой нити, сходящей с сечения радиуса r лопасти на азимуте у, будут
x = —r cos 0cos a — ra0 sin a — Vx (у — 0), y = -V'y (у — 0),
(6)
z = r sin a,
где -да < 0 < у, причем у — азимут, в котором лопасть находится в данный момент; а — угол атаки несущего винта, линейные размеры отнесены к Я, скорости — к юЛ.
Элемент нити А/ (см. рис. 2), соответствующий приросту 0 на Л0 имеет длину
л/ = i — + j — + к — или л/ = {i (r sin Є cos a + Vx )+ jV'. + kr cos ЄІЛЄ.
ае ае ае н x> y >
(7)
Проекция этого элемента, обозначенного через А/1, на образующую косого вихревого
- V+V
цилиндра, направленную вдоль единичного вектора / =
л// =лТ -1 =
V
(v =^Vx2 + v;2) будет
_V - V2 ^
V + r -ё- sin е cos a + 2sin a-V + ^т-
v V v2 У
ле.
(8)
а)
у
1,5
Рис. 3. Распределение циркуляции вихрей по контуру диска (а) и начальные условия для линии их расположения (б) на режимах ц = 0.0778 и 0.326
При переходе к бесконечно большому числу лопастей элемент завихренности -GЛ/l
распределяется на участке образующей косого цилиндра длины
( 2пУЛ
V к у
что соответствует
направленности вдоль этой образующей прямолинейной нити циркуляции [3]
Ч = -
А/гОк Ла
-А0 или
кв_
2п
2пУ
^ г sin0cos3 а 2sinаУі Уі2 ^
1 +-------------+ 1 ' 1
ц
У
А0.
(9)
Здесь О = -С— — отнесенная к юЯ2 интенсивность вихря, сбегающего с конца лопасти и Х-к
т7 2— „
повторяющегося через расстояние V —; к — число лопастей несущего винта; а — угол атаки
к
несущего винта; — коэффициент силы тяги; х — коэффициент, учитывающий влияние
концевого эффекта; Уі = -
средняя по диску индуктивная скорость.
0.313-12-ц
Для исследования движения пелены в настоящей работе использована плоская модель явления, где пространственные участки пелены изображаются в виде плоского движения тангенциального разрыва скоростей. В работе исследуется деформация линии Ь, представляющей собой
с
Ґ
в начальный момент времени эллипс, получающийся при сечении цилиндра плоскостью, перпендикулярной направлению скорости набегающего потока. Проследив за изменениями во времени соответствующей линии Ь на вихревой поверхности (рис. 3), можно составить заключение
о форме этой поверхности.
Вопрос о движении пелены вихрей и быстроте ее сворачивания исследовался на основе [1], где эта задача сводится к определению движения цепочки плоских вихрей в некоторой плоскости 2ву в соответствии с полем скоростей, вызываемом ими самими. Вихревые нити считаются параллельными друг другу и продолженными в обе стороны от плоскости 2ву до бесконечности.
Уравнения движения цепочки свободных плоских вихрей заданной циркуляции запишутся в виде
точка, где расположены элементы вихрей, от которых определяются скорости.
При численном решении система интегродифференциальных уравнений (10) заменялась следующей аппроксимирующей ее системой обыкновенных дифференциальных уравнений относительно искомых функций:
В уравнениях (11) время ґ заменено на новую безразмерную переменную интегрирования у = юґ, имеющую смысл азимутального угла поворота лопасти за соответствующее время ґ. Сделан переход к безразмерным величинам координат г и у, отнеся их к радиусу винта г = г/ Я,
Нелинейную систему обыкновенных дифференциальных уравнений (11) можно решить только численно. Для ее решения применяется модифицированный метод Эйлера с двойным просчетом на каждом шаге по времени у. В первом приближении индуктивные скорости при переходе от времени у к следующему шагу расчета со временем у + Ду вычисляются по перемещениям
и — точка, где определяется скорость; и'
і max
(11)
где Хі = х (и і, і), х} = х (иу, ґ), Уі = у (и, ґ), у} = у (и}, ґ)
у = у/Я. Введено обозначение
гу+ Лу гу + Лу' гу, уІ+Ау = уу+ЛУ' уу
(12а)
Во втором же приближении г и у будут равны
гу+ Лу гу + Лу' гср, уї+А£= уу+ЛУ' уср,
(12б)
где
гу + гу+Ду уу + уу +Ду
гср - 2 , Уср- 2 .
Полученные во втором приближении значения перемещений и скоростей принимаются за окончательные на данном шаге и являются исходными при расчете г и у в последующий
момент времени + 2Д .
Для знаменателей, входящих в уравнения для определения скоростей (11), введем обозначение
(у -у )2 +(ъ-^)2 -12. (13)
Появление особенностей и разрывов при подсчете индуктивных скоростей может быть устранено тем, что при их вычислении в выражение (13) к 12 добавляется еще величина в2, где в — малая положительная величина. Физический смысл добавки в том, что тем самым элементарные вихри приобретают толщину порядка в. Если предположить, что сходящие вихри диффундирующие, то, сделав в функцией времени жизни элемента свободного вихря можно приближенно учесть влияние диффузии вихрей за счет влияния вязкости. В данном случае для практических вычислений сделана попытка аппроксимировать в более простыми зависимостями.
При отсутствии диффузии параметр толщины слоя свободных вихрей В0 должен быть порядка толщины пограничного слоя. Выражая последнюю через коэффициент профильного сопротивления схр [4], получим
V
в0 -1.7216,—, (14)
0 XV
* ГО
где Ь — хорда профиля; Vcю — скорость обтекания профиля; V — коэффициент кинематической вязкости, который следовало бы вводить с поправкой на учет турбулентных пульсаций внутри вихря.
Для одновременного учета начальной толщины слоя и его диффузии все время увеличивается время жизни каждого вихря т, тогда
в2 =в0 + 4vт. (15)
Турбулентный характер движения воздуха в следе можно приближенно учесть тем, что коэффициент вязкости V берут гораздо большим, чем коэффициент вязкости для ламинарного течения воздуха, который равен V = 0.0000133 м /с. Эту величину по [5] следует брать, равной 105—106, но предыдущая расчетная практика по определению индуктивных скоростей и деформаций пелены несущего винта показала, что вряд ли следует увеличивать V больше, чем в 1000 раз, т. е. V = 0.0133 м2/с.
Результаты расчетов. Рассматривались режимы горизонтального полета с ц = 0.0778 и ц = 0.326. Им соответствовали следующие характеристики несущего винта: для меньшей
скорости — а = -3.7°, а0 -4.24°, с{ -0.01134 и для большей — а = -5.5°, а0 -4.24°,
с( - 0.01134. Винт пятилопастный. В начальный момент времени задана кривая Ь в виде эллипса с расположенными по ней вихрями интенсивности q (рис. 3, а). В практических расчетах непрерывное распределение свободных вихрей заменяется конечным числом точек с сошедшими с них вихрями. В данном работе эллипс поделен на 100 участков.
Задача состоит в отыскании кривой Ь в последующие моменты времени. Рассматривается движение каждого вихря в отдельности, а задача об определении формы вихревой поверхности сводится к задаче определения движения дискретных вихрей.
На рис. 3, б приведены начальные условия, от которых начинается движение линии Ь и которые выражаются соотношениями
у (і) = -cos у sin а + а0 cos а, г (і) = sin у,
(16)
где а — угол атаки несущего винта, а0 — угол конусности лопастей, у — азимутальным угол расположения точки с индексом /.
Одна группа расчетов производилась для случая, когда величина в2, характеризующая толщину вихрей, сходящих с линии Ь, была одинаковой для всех вихрей и постоянной по всей длине вихря на протяжении всего времени счета. Эта величина выбиралась не совсем произвольно,
учитывались рекомендации работы [2]. В табл. 1 приведены значения величин в2, употребленные в расчетах, значения индуктивных скоростей V по оси г и перемещений по этой же оси для
точек линии Ь при у = 0 и у = 25°. Расчеты показали, что дальнейшее уменьшение в после значения 0.001 практически не влияет на изменение скоростей и перемещений.
-1.5
6)
Рис. 4. Поле индуктивных скоростей по оси у и сдеформированная пелена в конце счета для не меняющихся по времени различных значений пара-
Таблица 1
^ (0) Кг (25°) г (0) г(25°)
0.1 -0.1335 -0.1106 -0.77184 -0.70403
0.001 -0.1413 -0.1046 -0.91833 -0.84417
0.00001 -0.14135 -0.1037 -0.92456 -0.84961
1Е-07 -0.1414 -0.1036 -0.92462 -0.84962
1Е-09 -0.1414 -0.1034 -0.92462 -0.84965
метра толщины вихрей є2:
а — индуктивные скорости в точках линии Ь; б — деформация пелены в конце пятого
оборота
Как показывают графики на рис. 4, а для индуктивных скоростей по оси у в зависимости от у (линия Ь на графике развернута в прямую), величины амплитуд скоростей в местах
свертывания жгутов примерно одного и того же порядка для значений величин є < 0.001.
Величина Уу при є2 = 0.1 является как бы средней для значений Уу при є2 < 0.001. Графики
построены для времени 0 = 5-2п. На рис. 4, б нанесены свернувшиеся следы для є2, равные 0.0000001, 0.001, 0.1. Для всех случаев ядро свернувшихся жгутов визуально находится на одном
месте, почти одинакового размера и имеет один центр, параметр толщины вихря при є < 0.001 мало влияет и на деформации.
В данной работе сделана попытка приближенно оценить влияние изменений толщины вихрей, вследствие вязкости окружающего воздуха, на результаты расчетов для разных коэффициентов кинематической вязкости V. В момент схода вихрей с линии L, в момент времени 0 = 0, тол-
Рис. 5. Изменение толщины сходящих с контура диска вихрей в течение времени счета 0
Т аблица 2
V 0 0.013 0.13 0.565
0 0.0004 0.0039 0.0169
6.23 0.0006 0.0057 0.024
12.53 0.0007 0.0073 0.0311
18.83 0.0009 0.0089 0.0382
25.13 0.0011 0.0106 0.0453
31.34 0.0012 0.0123 0.0524
Ц = 0.0778
-0.013. —0.13, -*-0.23. —<—0.34, —I—0.565
Рис. 6. Индуктивные скорости по оси у в следе в конце пятого и первого оборотов для различных коэффициентов вязкости: а — конец пятого оборота; б — конец первого оборота
щина всех вихрей одинакова и подсчитана по выраже-
нию (14). По мере движения эллиптической поверхности, окаймленной кривой Ь, при 0 >0 толщина сходящих вихрей начинает постепенно увеличиваться согласно формуле (15). В табл. 2 и на рис. 5 для разных коэффициентов V даются расчетные значения толщин вихрей на всем протяжении времени счета 0 < 0 < 5 • 2п (0 в радианах).
В случае, когда величина в2 подсчитывается на каждом шаге интегрирования и вследствие этого увеличивается со временем, в конце расчета наблюдаются менее резкие изменения в величинах индуктивных скоростей для разных коэффициентов вязкости (рис. 6). Уменьшение
параметра в2 до значений близких к реальным, например, при диапазоне изменения в2 от 0.00039 до 0.00122 при V = 0.0133 м /с, приводит, судя по поведению индуктивных скоростей, и к более резким деформациям и, как следствие, более сумбурному вхождению отдельных вихрей в жгуты (рис. 6, а, 7).
Влияние диффузии на индуктивные скорости лучше всего проследить при большом значении коэффициента вязкости V = 0.565 м2/с. Толщина вихрей за время счета изменяется в
пределах в = 0.0169 0.0524, что стыкуется со значениями толщины вихрей для малых скоростей ра-
боты [2]. Значения индуктивных скоростей здесь примерно в два раза меньше (рис. 6, б), чем при
V = 0.0133 м2/с. На рис. 8 для этой вязкости представлена свернувшаяся пелена вначале и конце счета, соответствующая скоростям рис. 6, б. Из графиков видно, как точки спокойно одна за другой втягиваются в жгуты.
Рис. 7. Деформированная пелена на режиме ц = 0.0778 в начале и конце счета с коэффициентом вязкости V = 0.013 м2/с:
а — деформации в конце первого оборота; б — деформации в конце пятого оборота
Рис. 8. То же, что на рис. 7, но V = 0.565 м2/с
По результатам расчетов по свертыванию эллиптической поверхности были построены графики рис. 9, а, 10, а, 11, а. Как и для крыла, свободные вихри на некотором расстоянии от несущего винта сворачиваются в два жгута. Процесс формирования жгутов заканчивается за время, в течение которого винт повернется на угол 0 « 8 рад для ц = 0.0778 и на угол 0 « 12.6 рад для ц = 0.326. За это время суммарная интенсивность вихрей, входящих в жгуты, достигает от максимальной примерно 60% для правого и 28% для левого жгутов, при этом винт проходит в поступательном движении ~0.7Л для ц = 0.0778 и ~4.06Л для ц = 0.326. Таким образом, точка свертывания будет ближе к диску по мере уменьшения характеристики режима ц и по мере роста интенсивности вихрей. Устойчивых свернувшихся шнуров, сходящих с центральных частей диска, не обнаружено. Вихри, не вошедшие в вихреобразования, остаются в виде пелены.
ц = 0.0778; V = 0.13 У
-------в*-,----------------
б)
Рис. 9. Сравнение деформированных следов в плоскости гоу для режима ц = 0.0778, полученных по дисковой и нелинейной вихревым теориям, в диапазоне изменения параметра в2 от 0.0039 в начале счета до 0.0123 в конце его (V = 0,13 м2/с):
а — деформация по дисковой теории; б — деформация по нелинейной теории
После отделения от летательного аппарата пелена вихрей, индуцируя индуктивные скорости, постепенно сама смещается вниз, сохраняя при этом примерно постоянное расстояние между жгутами после того, как они сформировались. Для большей скорости поступательного полета это смещение меньше.
На рис. 9, б, 10, б и 11, б представлены для сравнения деформированные следы, полученные по лопастной теории. Сравнение графиков на этих рисунках с графиками рис. 9, а, 10, а, 11, а показывает неплохое качественное совпадение в целом. Поведение пелены за несущим винтом, определенной по дисковой вихревой теории, сохраняет тенденции поведения деформированной пелены, полученной по лопастной вихревой теории. Границы распространения обоих следов находятся примерно в одних пределах. Поджатие струй в начале счета, размеры вихреобразований и их ядер, полученные по этим двум теориям, близки друг к другу. Имеющиеся различия, например в длине следа (см. рис. 11) и смещении пелены вниз (см. рис. 9, 10), могут быть объяснены наличием составляющей по оси х в следе по лопастной теории и ее отсутствием в следе по дисковой теории.
То, что за боковыми кромками диска несущего винта при косом обтекании имеется спиральное движение части воздуха, аналогичное картине обтекания крыла малого удлинения, подтверждается и визуальными наблюдениями. В качестве примера можно привести работы [6] и [7] и рис. 12 из работы [7].
Рис. 10. То же, что на рис. 9, но ц = 0.326
Рис. 11. Сравнение деформированных следов в плоскости хох для режима ц = 0.0778, полученных по дисковой и нелинейной вихревым теориям, в диапазоне изменения параметра в2 от 0.0039 в начале счета до 0.0123 в конце его (V = 0,13 м2/с):
а — деформация по дисковой теории; б — деформация по нелинейной
теории
Рис. 12. Сечения вихревой пелены за прямоугольным крылом с углом атаки 9°. След показан на расстояниях от задней кромки, равных соответственно 1.0, 1.6, 2.9, 5.5,
11.2 длинам хорды
Данная программа по предложенной методике удобна для применения в случае, когда надо быстро оценить физическую картину обтекания несущего винта, так как расчет по ней может быть проведен достаточно быстро.
ЛИТЕРАТУРА
1. Баскин В. Э., Щеглова В. М. О деформациях пелены вихрей винта в косом потоке // Труды ЦАГИ. — 1968. Вып. 1104.
2. Баскин В. Э., Королева К. К. К нелинейной вихревой теории винта вертолета в горизонтальном полете // Труды ЦАГИ. — 1966. Вып. 1013.
3. Ван Ши-цунь. Обобщенная вихревая теория несущего винта вертолета // В сб.: Вопросы аэродинамики несущих винтов вертолетов // Труды МАИ. — 1961.
4. Прандтль Л., Титьенс О. Гидроаэромеханика // ОНТИ, НКТП СССР. —
1935.
5. К о чин Н. Е., Кибель И. А., Розе Н. В. Теоретическая гидромеханика. Ч. 2. — М.: Физматгиз. — 1963.
6. Миль М. Л., Сперанский М. К. Исследование поля скоростей вокруг ротора геликоптера при осевом и косом обтекании // Труды ЦАГИ. — 1949. Вып. 966.
7. Ван - Дайк М. Альбом течений жидкости и газа. — М.: Мир. — 1986.
Рукопись поступила 3/ХІ2003 г.