УДК 517.9
В. И. В а н ь к о, И. К. М а р ч е в с к и й, Г. А. Щеглов
ЧИСЛЕННО-АНАЛИТИЧЕСКИЙ МЕТОД ИССЛЕДОВАНИЯ УСТОЙЧИВОСТИ ПОЛОЖЕНИЙ РАВНОВЕСИЯ ПРОФИЛЯ В ПОТОКЕ
Предложен эффективный метод исследования устойчивости положений равновесия профиля, находящегося в потоке среды. Использование выведенных аналитически достаточных условий неустойчивости по Ляпунову и численного метода вихревых элементов позволяет проводить анализ при малых затратах вычислительных ресурсов.
E-mail: [email protected]; [email protected]; [email protected]
Ключевые слова: устойчивость по Ляпунову, условия устойчивости, профиль, галопирование, метод вихревых элементов.
Задача исследования устойчивости положений равновесия профиля в потоке является весьма актуальной при проектировании элементов строительных конструкций, подверженных ветровым нагрузкам, прямолинейных участков трубопроводов, находящихся в поперечном потоке жидкости или газа, трубчатых элементов систем охлаждения и теплообменников энергетических установок, отдельных проводов и расщепленных фаз линий электропередачи.
В современных нормативных документах по проектированию конструкций, испытывающих на себе действие потока среды (ветровых нагрузок, морских течений и пр.), указано [1], что в необходимых случаях в зависимости от условий возведения и эксплуатации сооружений следует учитывать воздействия, вызывающие аэродинамически неустойчивые колебания типа галопирования, т. е. колебания конструкции с одной степенью свободы поперек потока. Галопированию, как правило, подвержены протяженные конструкции, обтекание которых можно считать плоскопараллельным, поэтому можно ограничиться рассмотрением задачи обтекания профиля. Вопросы исследования устойчивости положения равновесия профиля в потоке (аэродинамической устойчивости) рассмотрены в специальной литературе [2, 3] и сводятся к проверке известного условия Глауэрта—Ден-Гартога: положение равновесия неустойчиво, если выполняется условие
G(a) = C'ya + Cxa < 0. (1)
Здесь и далее Cxa и Cya означают стационарные аэродинамические коэффициенты лобового сопротивления и подъемной силы соответственно, штрихом обозначена производная по углу атаки. Данное условие
получено как необходимое условие возбуждения больших колебаний профиля с одной степенью свободы [4, 5], а в работе [6] показано, что для случая галопирования оно также является и достаточным условием неустойчивости. При этом колебательная система предполагалась идеальной — демпфирование не рассматривалось.
Условие Глауэрта—Ден-Гартога удобно для практического использования, если известны зависимости стационарных аэродинамических характеристик профиля от угла атаки. При решении практических задач возникает проблема определения этих зависимостей, поскольку экспериментальные данные известны лишь для профилей, имеющих простую геометрическую форму. Кроме того, конструкция, взаимодействующая с потоком, не всегда может быть представлена в виде колебательной системы с одной степенью свободы. Также следует учитывать, что в реальных конструкциях присутствует демпфирование колебаний, что может повлиять на ее устойчивость. При оценке устойчивости положений равновесия задемпфированной колебательной системы с несколькими степенями свободы, в которой аэродинамические нагрузки воспринимает профиль произвольной формы, требуется применение как аналитических условий устойчивости, обобщающих условие (1) и расширяющих область его применения, так и методов вычислительной гидродинамики для определения аэродинамических характеристик профиля.
В настоящей работе предлагается подход к решению указанной задачи, называемый численно-аналитическим методом исследования устойчивости положений равновесия профиля в потоке. Данный метод может применяться для исследования систем с одной, двумя и тремя степенями свободы (движение по потоку, поперек потока и вращение) при наличии идеально упругих либо линейных вязкоупругих связей. Рассматривается обтекание профиля потоком несжимаемой среды.
Условия устойчивости профиля в потоке. В работе [6] рассмотрена модель движения профиля с тремя степенями свободы под действием аэродинамических сил и сил со стороны упругих либо вязко-упругих связей. Достаточное условие неустойчивости положения равновесия профиля с тремя степенями свободы имеет вид
W(а) = Оха(С'уа + Оха) + Суа(Суа - С'ха) < 0. (2)
В работе [7] проведена проверка адекватности условий (1) и (2) для различных профилей. В работах [8, 9] решена задача об исследовании устойчивости по Ляпунову положений равновесия профиля в потоке и получены необходимые и достаточные условия устойчивости и неустойчивости при более общих по сравнению с [6] ограничениях на характер движения профиля и параметры связей. Из полученных результатов для разных случаев движения профиля (первая, вторая или
третья степени свободы; одинаковые либо различные характеристики связей) следуют достаточные условия неустойчивости для систем с идеально упругими связями
в (а) < 0, W (а) < 0,
М(а) = С'уа + 3 ■ СХа < 0, (3) с линейными вязкоупругими связями
С» = С'уа + Сха + 2Ду < 0, (4)
МДа) = СУа + 3Сха + 2(ду + Дх) < 0, (5)
ГС» = (Сха + Дх)(СУа + Сха + 2^) + Суа(Суа - СХа) < 0, (6)
где дх и ду — безразмерные коэффициенты демпфирования связей. Отметим, что для неустойчивости положения равновесия профиля с тремя степенями свободы в случае одинаковых жесткостей связей достаточно выполнения любого из условий М(а) < 0, W(а) < 0; при различных жесткостях связей достаточным условием неустойчивости будет С(а) < 0. При этом условие (1) является более сильным, чем условие (3): учитывая положительность коэффициента лобового сопротивления Сха легко видеть, что при выполнении условия (3) условие (1) с необходимостью выполняется. Это же замечание относится к условиям (5) и (4).
Область неустойчивости для систем с вязкоупругими связями, описываемая условиями (4)—(6), уже, чем область неустойчивости для систем с идеально упругими связями, поэтому при исследовании устойчивости системы с малым демпфированием целесообразно пользоваться условиями С(а) < 0 и W(а) < 0. Удовлетворяющие этим условиям углы атаки профиля будут определять неустойчивые положения равновесия конструкции с некоторым запасом, что может быть весьма полезно в инженерных приложениях.
Особенностью приведенных достаточных условий неустойчивости профиля с идеально упругими связями является их инвариантность — они зависят лишь от стационарных аэродинамических коэффициентов профиля и их производных и не зависят от механических параметров системы — массы профиля, жесткостей связей, скорости набегающего потока.
Определение аэродинамических коэффициентов профиля и их производных, входящих в условия неустойчивости (1)—(6), возможно как экспериментальным, так и расчетным путем. На практике могут применяться, например, вихревые методы, позволяющие при сравнительно низких затратах вычислительных ресурсов моделировать обтекание профиля несжимаемым потоком и с достаточной точностью
определять силы, действующие на профиль со стороны потока. Метод вихревых элементов [10-12] — это бессеточный лагранжев метод численного моделирования нестационарного обтекания профилей, основанный на моделировании эволюции завихренности.
При этом решение двумерных уравнений динамики несжимаемой среды
V-V = 0, — + (V -V) ■ / = -V(^ Р) + vV2- W,
где / — скорость, p — давление, р = const — плотность, v — коэффициент кинематической вязкости среды, сводится к решению уравнения
д/ -— = Vx ((/ + W) х /) dt
относительно / — распределения завихренности в области течения. В последнем уравнении W — диффузионная скорость, учитывающая влияние взякости среды и вычисляемая по формуле
- (Vx /) х / W = v---.
|/12
Распределение завихренности, в свою очередь, моделируется набором вихревых элементов (ВЭ), для каждого из которых на плоскости задается положение / и интенсивность Г (i = 1,...,N). Скорость среды в любой точке течения вычисляется по закону Био—Савара по известным характеристикам ВЭ, а давление определяется c помощью аналога интеграла Коши—Лагранжа.
Метод вихревых элементов особенно эффективен при решении внешних задач обтекания, поскольку на каждом шаге расчета по времени выполнение уравнения неразрывности и граничного условия затухания возмущений на бесконечности происходит автоматически. Граничное условие на профиле обеспечивается генерацией ВЭ на поверхности обтекаемых профилей, а течение идеальной либо вязкой среды моделируется движением имеющихся ВЭ. Интенсивности рождаемых ВЭ находятся из решения соответствующей системы линейных алгебраических уравнений вида [А]{Г} = {b}, аппроксимирующей граничное условие на профиле, и в дальнейшем не меняются, а движение ВЭ описывается системой обыкновенных дифференциальных уравнений
dr * —* —*
-± = /(/) + W(/i), i = 1,...,N, dt
правые части которой представляют собой скорости движения ВЭ.
Высокая эффективность задачи метода вихревых элементов и его модификаций [9, 11, 12], выявленная при решении данной задачи, подтверждается результатами тестовых расчетов.
Численно-аналитический метод исследования устойчивости и примеры его применения. В основе предлагаемого метода лежит последовательное использование численного метода вихревых элементов и полученных аналитически условий неустойчивости профиля по Ляпунову (1)—(3). Алгоритм включает в себя 3 этапа.
1. Моделирование обтекания профиля методом вихревых элементов при различных углах атаки и определение стационарных коэффициентов лобового сопротивления и подъемной силы.
2. Аппроксимация зависимостей Сха(а) и Суа(а) гладкими кривыми.
3. Определение с помощью инвариантных достаточных условий углов атаки, соответствующих неустойчивым положениям профиля в потоке.
Для проверки работоспособности предложенного численно-аналитического метода и оценки его адекватности и эффективности были проведены тестовые расчеты по исследованию устойчивости положений равновесия профилей с сечениями в форме ромба и квадрата, для которых известны экспериментальные данные [7], и профиля с сечением обледенелого провода ЛЭП. Все параметры расчетных схем, а также особенности проведения вычислительных экспериментов и обработки результатов расчетов подробно описаны в работе [9].
Профиль с сечением в форме ромба. Исследуется устойчивость профиля, имеющего форму ромба с отношением длин диагоналей 1 : 0,75. На рис. 1 приведены характерные картины, на которой показаны вихревые следы и линии тока для углов атаки профиля а = 0° и а = 60°.
При определении аэродинамических коэффициентов профиля в каждом расчете выполнялось 5000 временных шагов, течение к этому моменту считалось установившимся, а стационарные значения аэродинамических коэффициентов получались путем усреднения нестационарных по последним 2500 шагам. Найденные значения отмечены
Рис. 1. Обтекание профиля в форме ромба
Рис. 2. Зависимости коэффициентов Сха и Суа от угла атаки (слева); графики функций С(а) и W(а) и амплитуды колебаний профиля в потоке (справа)
на левом графике (рис. 2) точками. Гладкие зависимости, показанные сплошными линиями, являются результатом аппроксимации полиномами Чебышева. На правом графике сплошными линями показаны зависимости С(а) и W(а), построенные по формулам (1) и (2).
В эксперименте по продувке упругозакрепленного профиля в аэродинамической трубе измерялись амплитуды колебаний профиля поперек потока. Полученные значения безразмерных амплитуд колебаний при различных углах атаки показаны точками рис. 2 (справа). Видно, что найденный при использовании численно-аналитического метода диапазон углов атаки 23° < а < 43°, в котором функции С(а) и W(а) принимают отрицательные значения, хорошо согласуется с областью углов, соответствующих возбуждению колебаний профиля с большой амплитудой (что является проявлением неустойчивости).
Профиль с сечением в форме квадрата. Аналогичное исследование было проведено для квадратного профиля. На рис. 3 показано обтекание профиля для углов атаки а = 0° и а = 20°.
На рис. 4 слева приведены вычисленные значения аэродинамических коэффициентов и их аппроксимация гладкими кривыми, справа — графики функций С (а) и W (а), наложенные на полученную в эксперименте в аэродинамической трубе зависимость амплитуды колебаний от угла атаки.
Хорошее совпадение диапазона углов атаки 0° < а < 15°, в котором С(а) < 0 и W(а) < 0, с областью углов, соответствующих неустойчивым положениям профиля в потоке, показывает эффективность и адекватность построенного численно-аналитического метода.
Рис. 3. Обтекание квадратного профиля
-I • - .•
Рис. 4. Зависимости коэффициентов Сха и Суа от угла атаки (слева); графики функций С(а) и W(а) и амплитуды колебаний профиля в потоке (справа)
Рис. 5. Обтекание профиля обледенелого провода ЛЭП
Профиль обледенелого провода ЛЭП. В литературе отмечается, что для проводов воздушных ЛЭП в условиях обледенения характерно явление пляски. Проведем исследование устойчивости положений равновесия с помощью разработанного численно-аналитического метода. На рис. 5 приведены характерная форма наледи, образующейся на круглом проводе, и картина обтекания при нулевом угле атаки (нулевой угол атаки соответствует естественному положению равновесия обледенелого провода на ЛЭП).
На рис. 6 представлены значения аэродинамических коэффициентов обледенелого провода при малых углах атаки, полученные в результате расчетов методом вихревых элементов, а также графики функций С(а) и W(а).
При углах атаки -6° < а < 9° выполняются условия С(а) < 0 и W(а) < 0, поэтому эти положения профиля неустойчивы. Данный результат хорошо согласуется с наблюдаемой при а = 0 неустойчивостью обледенелых проводов.
Рис. 6. Зависимости стационарных аэродинамических коэффициентов обледенелого провода от угла атаки (слева); графики функций О (а) и W(а) (справа)
Выводы. Предложенный численно-аналитический метод исследования устойчивости положений равновесия профиля в потоке дает возможность проведения эффективного анализа на этапе эскизного проектирования конструкций и организации экспериментов. Использование метода вихревых элементов позволяет с высокой эффективностью выполнять расчеты аэродинамических характеристик профилей, особенно при использовании параллельных вычислительных технологий.
Авторы благодарят Межведомственный суперкомпьютерный центр РАН за предоставленную возможность использования кластера МВС-100К.
СПИСОК ЛИТЕРАТУРЫ
1. Строительные нормы и правила: нагрузки и воздействия: СНиП 2.01.0785* / Госстрой СССР; введ. 01.01.87. — М.: ГП ЦПП, 2005. — 44 с.
2. Савицкий Г. А. Ветровая нагрузка на сооружения. — М.: Стройиздат, 1972.
— 110 с.
3. Симиу Э., Сканлан Р. Воздействия ветра на здания и сооружения. — М.: Стройиздат, 1984. — 360 с.
4. Den-Hartog J. P. Transmission LineYs Vibrations Due to Sleet // Transactions AIEE. — 1932. — Vol. 51. — P. 1074-1076.
5. Glauert H. The Rotation of an Aerofoil About a Fixed Axis // Reports & Memoranda / Great Britain Advisory Committee for Aeronautics (GBACA). — 1919.
— No. 595. — 8 p.
6. Ванько В. И. Математическая модель пляски провода ЛЭП // Изв. вузов. Энергетика. — 1991. — № 11. — С. 36-42.
7. Ванько В. И., Соловьева Е. В. Об условиях аэродинамической неустойчивости положений равновесия профилей // Прикладная механика и техническая физика. — 1996. — Т. 37, № 5. — С. 29-34.
8. Марчевский И. К. Об условиях устойчивости положения равновесия профиля в потоке // Вестник МГТУ им. Н.Э. Баумана. Естественные науки. — 2007. — № 4. — С. 29-36.
9. М а р ч е в с к и й И. К. Математическое моделирование обтекания профиля и исследование его устойчивости в потоке по Ляпунову: Дисс. ... канд. физ.-мат. наук. — М., 2008. — 119 с.
10. C o 11 e t G. -H., Koumoutsakos P. D. Vortex Methods: Theory and Practice.
— Cambridge University Press, 2008. — 328 p.
11. Андронов П. Р., Гувернюк С. В., Дынникова Г. Я. Вихревые методы расчета нестационарных гидродинамических нагрузок. — М.: Изд-во МГУ им. М.В. Ломоносова, 2006. — 184 с.
12. Дынникова Г. Я. Вихревые методы исследования нестационарных течений вязкой несжимаемой жидкости: Диса ...д-ра физ.-мат. наук. — М., 2011. — 269 c.
Статья поступила в редакцию 25.10.2011