Геометрический анализ при исследовании колебаний тела сложной конфигурации в потоке среды
В.А. Самсонов, Д. В. Беляков
Аннотация - Работа посвящена построению и исследованию математической модели автоколебаний аэродинамического маятника в потоке среды. В качестве модели воздействия среды на тело принята модель квазистатического обтекания пластинки средой. Согласно этой гипотезе, аэродинамические силы, действующие на тело, прикладываются в центре давления. В рассматриваемой задаче центр давления является подвижным относительно пластинки. Получены уравнения движения для рассматриваемого тела. Проведен переход к новым безразмерным переменным. Показано нарушение единственности при определении угла атаки. Проведен параметрический анализ областей неоднозначности. Найдены все стационарные точки, являющиеся решениями уравнений равновесия. Показано, что в наиболее характерном положении равновесия, соответствующем состоянию покоя областей неоднозначности нет. Проведено исследование
устойчивости положения равновесия, соответствующего состоянию покоя, в котором реализован критерий Гурвица и изображена область устойчивости. Показано, что силы аэродинамического воздействия для тел с одними формами могут способствовать развитию автоколебаний, а для других затуханию. В математическом пакете MATLAB 18 написан комплекс программ, позволяющий строить области устойчивости и проводить численное интегрирование уравнений, описывающих колебания тела, для того, чтобы подтвердить адекватность построенной модели.
Ключевые слова - Автоколебания, флаттер, маятник
I. Введение
Поиски, сопровождающие развитие техники, всегда сопряжены с риском принятия неверных решений, приводящих к авариям, а иногда и к катастрофам. Середина 30-х годов прошлого века стала периодом бурного развития авиации. В это время по миру прокатилась волна аварий при испытаниях скоростных самолётов. В то время, наиболее острой проблемой был флаттер - один из видов автоколебаний, ведущих к разрушению самолета в воздухе, который возникал при увеличении скорости полета. С точки зрения теории колебаний летательный аппарат представляет из себя
Статья получена 14 июля 2019
Д. В. Беляков, к.т.н., доцент, Московского Авиационного Института (Национального исследовательского университета), г. Москва (e-mail: [email protected])
В. А. Самсонов, д.ф.м.н., профессор МГУ им. М. В. Ломоносова, г. Москва (e-mail: [email protected])
автоколебательную систему, источником энергии в которой служит набегающий поток, на который реагирует конструкция самолета. При интенсивных колебаниях аэродинамические напряжения могут достигнуть разрушающих значений. М.В. Келдыш показал, что флаттер это форма потери устойчивости, он имеет резонансную природу. Был разработан метод расчета критической скорости флаттера. Чтобы избежать резонанса при движении крыла в воздушном потоке, он предложил соответствующим образом перераспределить массы вдоль крыла и так расположить упругие элементы, чтобы избежать совпадения собственных частот колебаний крыла с частотами вынуждающих внешних сил.
При детальном изучении флаттера почти во всех случаях обнаруживаются нелинейные
аэродинамические эффекты. Однако в ряде ситуаций оказалось возможным успешно решить задачу на основе линейных аналитических подходов. Рассмотрим одну модельную задачу, которая иногда обсуждается в литературе, в которой можно использовать такие подходы.
II. Постановка задачи Тело состоит из пластинки АВ и прикрепленной к ней массивной державки ВД. Центр масс тела Д может двигаться только по прямой ГЕ, державка может поворачиваться вокруг точки Д (см. рисунок 1). Введем неподвижную систему координат. Будем считать, что силы деформации элементов крепления зависят от отклонений линейным образом и сводятся к восстанавливающей силе и возвращающему моменту, и что в положении покоя маятник ориентирован по потоку. Аэродинамические силы, приложенные к телу, приняты в соответствии с эмпирической теорией
стационарного обтекания плоской пластины. [3],[4], [6],[7].
f < Д / < -V г S 3 в
>
\ Pa \= p(a)V2 = 0.5 pac y (a)Vj
где a - угол атаки между вектором
V
и пластинкой
(2)
Рис. 1 Рассматриваемое тело
Аэродинамические силы, действующие на пластинку, разложим на две составляющие: сила
сопротивления SA, направленная против воздушной
скорости ^ точки А относительно потока среды, и
подъемная сила РА, направленные ей ортогонально. При этом величины аэродинамических сил равны:
5(аУА2 = 0.5 p<Jcx (а)ГА2
р, 5 - аэродинамические функции углов атаки, сх , су -безразмерные аэродинамические функции, р -плотность воздуха, < - площадь одной пластинки.
В рассматриваемой модели предполагается, что центр давления пластинки точку А можно считать подвижной относительно пластинки. Сдвиг центра давления описывается функцией I(а), которая описывает расстояние между центром давления А и геометрическим центром пластинки В. Зависимость определена из продувок прямоугольных пластинок с заданным удлинением в аэродинамической трубе [5]. Аналогичная задача, только с неподвижным центром давления рассматривалась в работе [1].
Составим уравнения движения рассматриваемого тела. В качестве обобщенных координат, определяющих положение тела, введем координату у центра масс Д, совпадающего с серединой пластинки и угол отклонения стержня от горизонтали. Тогда теорема о движении центра масс в проекции на ось ОУ и теорема об изменении кинетического момента будут иметь вид:
my = s(a)VA (l(а) 3 sin 3 - r3cos 3 - y) -
- p(a)VA (r 3 sin 3 + l (а) 3 cos 3 + V) - ky (1)
J3 = r Vjz(a) - l(a)VA п(а) - с 3 Где
т(а) = p(а) sin а - s(а) cos а п(а) = p(а) cos а + s(а) sin а
- это аэродинамические функции нормальной и тангенциальной сил. Кинематические соотношения,
связывающие VA ,а с y, y, 3, 3', имеют вид:
VA sin а = l(а) 3 - y sin 3 + V cos 3
VA cos а = r3 + y cos 3 + V sin 3 Построенные уравнения содержат параметры:
m, J, s^), р(а), l (а), Va , V, k, с, r
Введем новые безразмерные переменные: y
Y = — безразмерная координата центра масс, b
V
т = — t безразмерное время,
^ bS
Q = безразмерная угловая скорость,
V
U = ~A безразмерная скорость центра давления,
При переходе к новым безразмерным переменным уравнения (1) преобразуются к виду:
MY = сх (а)и (RQcos 3 + s^Qsrn 3 -Y) + +cy (а)и(RQsin 3 + s^)Q cos 3 + 1) - KY
13 = U2 (RcT (а) - £(а)сп (а)) - C3
Здесь M, K, I, е(а), R это безразмерные коэффициенты :
(3)
M=
m
kb
bc
I=
, K =-- , C = 2
0.5pa 0.5 paV 0.5paV2
J l (a) r
-,s(a) = ,R = —
0.5 pab b b
(4)
После простых преобразований кинематические соотношения (2) будут иметь вид:
U sin а = s(а)Q - Y sin 3 + cos 3
U cos а = RQ + Y cos 3 + sin 3 Построенные уравнения содержат параметры:
M, I, c», Cy^s^U, K, C
Таким образом, построена математическая модель колебаний пластинки, которая представляется системой безразмерных уравнений (3)-(4), которая имеет меньшее
число параметров. Можно исследовать эту систему и потом вернуться к старым переменным.
III. НЕОДНОЗНАЧНОСТЬ ОПРЕДЕЛЕНИЯ УГЛА АТАКИ При численном интегрировании уравнений движения (3), нужно сначала определить угол атаки из соотношений (4). Разделим первое уравнение (4) на второе и, избавляясь от U , получим: s(a)Q -Y sin 3 + cos 3
(5)
b)
Q<0,RQ + Y cos S + sin S < 0
tga =
a)
Q > 0, RQ + Y cos S + sin S > 0
с)
Q<0,RQ + Y cos S + sin S > 0
RQ + Y cos3 + sin3
Умножим правую и левую части уравнения (5) на знаменатель
RQ +Ycos3 + sin3 ф 0
Если знаменатель все-таки обращается в нуль, то тогда
п
может быть потеряно решение a = —.
Получим нелинейное уравнение (6), равносильное соотношению (5): (RQ + Y cos 3 + sin 3)tga -
-Y cos 3-sin 3 = s(a)Q (6)
Мы можем решить уравнение (6) и найти угол атаки a при различных значениях фазовых переменных Y, Q,3. Поверхность a = a (Y, Q,3) можно изобразить только в четырехмерном пространстве. Мы будем изображать поверхность a = a(Q, Y) только при
фиксированных значениях угла 3. Приведем графическую интерпретацию результатов решений уравнения (6) для пластинки шириной один метр с удлинением Я = 8, из которой следует нарушение единственности определения угла атаки при некоторых значениях Y, Q, 3.
а) Q<0,RQ + Y ^3 + зш3<0
Рис. 2. Неоднозначность углов атаки. Таким образом, существуют значения У,Q, 3., при
которых поверхность а=а(0., У) имеет складку, внутри которой кинематические соотношения (4) имеют ровно по три решения. Участок области
неоднозначности для поверхности а = а(0, У) при
п
фиксированном значении рисунке 3 для рассматриваемых пластинок.
S = — изображен на 4
ш il ш
à
0.5 0 -0.5
п
Рис.3. Решения уравнения (6) в случае S = — .
4
При уменьшении значения угла ориентации тела 3 до нуля на рассматриваемой поверхности, эта тенденция сохраняется, и все равно остаются области неоднозначности. Построенная для случая 3 = 0 поверхность изображена на рисунке 4. Она
симметрична относительно точки (0,0, —) и имеет по краям две складки.
Рис. 4 Решения уравнения (6) при 3 = 0. Области неоднозначности.
Возникает вопрос: какое решение при прохождении области неоднозначности является истинным, и по какому правилу можно выбирать нужную ветвь при моделировании? Ведь неоднозначность определения угла атаки означает, что неоднозначно определятся аэродинамические функции и сдвиг центра давления. Учитывая то, что характер обтекания маятника и силового воздействия со стороны среды носит непрерывный характер, то значения угла атаки а нужно выбирать, сохраняя его непрерывность. Таким образом, математическая модель колебаний пластинки, (3)-(4) стала действительно замкнутой и определенной. Более подробный алгоритм выбора нужного значения угла атаки а рассмотрен в статье (2).
Приступим к проведению анализа положений равновесий для нашего тела.
IV. РЕШЕНИЕ УРАВНЕНИЙ РАВНОВЕСИЯ
При Y = 0,0 = 0 из кинематических соотношений (4), учитывая, что 0* = const, Y* = const следует, что
tga* = ctg0*, откуда а* +0* = П. Запишем уравнения равновесия для системы (3)
с (a*) -kY* = 0
(7)
Проведем анализ решений системы (7) для прямоугольных пластинок, имеющих удлинение Х = 8 при уменьшении коэффициента жесткости С от больших значений до нуля.
Второе уравнение системы (7) всегда имеет п
решение а = —, которое соответствует стационарной
точке в = 0, У = 0. Правая часть второго уравнения системы (8) симметрична относительно этой
точки и имеет период 2п. Положения равновесия, зависящие от взаимодействия линейной
восстанавливающей силы с аэродинамическими силами, выведем графически.
При С >>1 имеем только одну стационарную точку в = 0, У = 0.
При уменьшении значения С, прямая в правой части второго уравнения системы (7) становится более пологой, добавляются новые положения равновесия, и число стационарных точек увеличивается до трех.
Каждая из этих точек разветвляется еще на две (см. рис. 5) и мы имеем уже пять положений равновесия.
* O -5
-1O
-1.5
-O.5
O
teta
O.5
1.5
RC'(a-)-S(a')Cn (a) = С -a)
Рис. 5. 5 положений равновесия
По мере уменьшения значения С этот процесс ветвления
положений равновесия продолжается. Число положений равновесия увеличивается. Все они симметричны относительно стационарной точки в = 0, У = 0.
В случае С=0 имеем счетное множество стационарных решений, как в [4,7,8].
Продолжим анализ стационарных точек в случае г = 0 (пластинка ориентирована поперек потока). Если
2
С>>1 имеем одно тривиальное положение равновесия в На этом анализ простейших положений равновесия состоянии покоя в = 0, У = 0. При уменьшении С можно закончить. Далее будет проведено исследование
число симметричных стационарных точки растет (см. рис. 6, 7).
устойчивости наиболее характерного положения равновесия: в = 0, Y = 0.
1O б
* O
-1O
-2 -1.б -1 -O. б O Oi 1 1.б 2 teta
Рис. б. Решение уравнения равновесия в случае r = 0
V. ИССЛЕДОВАНИЕ УСТОЙЧИВОСТИ СОСТОЯНИЯ ПОКОЯ
Будем считать, что маятник совершает малые колебания около положения равновесия У = 0,3 = 0. Исследуем его устойчивость.
Для исследования устойчивости воспользуемся критерием Гурвица для системы четвертого порядка. Линеаризуем уравнения движения (3) - (4), в случае,
п
ко гда
Y ^ 0 и S ^ 0 a ^ ■
2
1O-
> O
-1O
(8)
-2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5 2 2. teta
Рис. 7. Решение уравнения равновесия в случае r = 0
В случае, когда C=0, имеем счетное множество стационарных точек.
Таким образом, у рассматриваемого тела есть интересные и неочевидные изолированные положения равновесия. Все найденные стационарные точки симметричны относительно положения в = 0, Y = 0, что подтверждает симметрию рассматриваемой задачи. Для всех найденных стационарных точек уравнение (6) преобразуется к виду: tga = ctge . Это уравнение
при фиксированном в имеет только одно решение, которое не попадает в область неоднозначности. Таким образом, целесообразно проводить исследование устойчивости в рассматриваемых стационарных точках.
^ MY + ( с;; + ()) Y + Ky + +R(cy ' (') + с( j))S + cy' (') S = 0 ^ (RcJ (|)'(§)c„ (f))Y +S +
+(R2cT (|) -Rs'(-|)c„ (-|))S + +(RcT' (- '( -си (- + СS = 0
Характеристиче ское уравнение системы (8) имеет вид:
a0Л4 + а1Л3 + a2l2 + a3Л + a4 = 0
Сформулируем условия устойчивости, учитывая, что
<7) = ^ (7)+
а0 = MI > 0 (Выполнено всегда)
a! = M^'-^-s' (-)сп (-Ж)+с/ -)I>0
a2 = М(сг ' (-)R-s'(-)сп (-)+C)+KI>0 a3 = К(ст ' -)R2-s'(-X (-)R)+ +сТ (-Х^ (-XC+с/ (-)R-s' (-с -))>0
а4 = К(с/ (-)R-s' (-)сп (-)+C)>0
A3 =
0 aj a3
aj a3 0
> 0
Области устойчивости изобразим безразмерных параметров (R С ).
(9)
(10)
(11)
(12) (13)
на плоскости
б
При выполнении условия С < с^ (—)Я-е'( — )сп (П)
всегда нарушается условия (10),(12) и нулевое решение системы (3) неустойчиво.
Во всех остальных случаях условия (9)-(13) равносильны системе неравенств:
С > с/ (|)R-,'((П) А 3 > 0
(14)
Далее на рисунке 8 зеленым цветом изображена область устойчивости для прямоугольных пластинок, имеющих удлинение Я = 8.
Рис. 8. Область устойчивости
Граница области устойчивости на рисунке 8 является прямой и соответствует появлению нулевого корня в характеристическом уравнении. С ростом скорости потока условие (12) нарушается. Такая форма потери устойчивости сопровождается ветвлением положений равновесия.
Исследование показало, что силы
аэродинамического воздействия для тел с одними формами могут способствовать развитию автоколебаний, а для других затуханию.
Перейдем к следующему этапу данной работы -построению программы, которая реализует численное интегрирование системы уравнений движения (3).
уравнений (3)-(4), описывающих колебания пластинки. При поиске численного решения используется процедура ode45, реализующая методы Рунге-Кутта четвертого и пятого порядка с переменным шагом. При поиске численного решения аэродинамические функции интерполируются кубическим сплайном.
Программа для решения задачи о колебании пластинки имеет следующий интерфейс. С помощью ползунка задается интервал интегрирования, значение удлинения пластинки можно выбрать из четырех значений в выпадающем меню. Значения М, I вводятся в соответствующие поля. Значения Я и С вводятся путем выбора точки на графике области устойчивости. Начальные условия также выбираются графически. Полученные путем интегрирования точки У = изображаются как траектории в конфигурационном пространстве системы.
Для моделирования колебаний пластинки возьмем удлинение, равное пяти. Зададим значения М = 10 , / = 10 (Я, С) введем графически:
Рис. 9. Графический ввод Я,С
начальные условия введем графически в окрестности состояния покоя.
VI. ЧИСЛЕННОЕ РЕШЕНИЕ УРАВНЕНИЙ ДВИЖЕНИЯ В математическом пакете МАТЬАВ 18 написана программа, реализующая численное интегрирование
4. Предложен комплекс программ, позволяющий строить области устойчивости и проводить численное интегрирование уравнений движения для исследуемой модели.
2 -1-1-1--1-i-
1.5 -1 ■
0.5 ■
"3=-D~:-■--="
-0.5 --1 -1.5 -
.2 -'-'-1--'-'-
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2 TETA
Рис 10 Ввод начальных условий Видно, что изображающая точка покидает окрестность начала координат, что иллюстрирует неустойчивость основного равновесия и соответствует переходу к колебаниям.
РЕШЕНИЕ ЗАДАЧИ О КОЛЕБАНИЯХ МАЯТНИКА
8 6 -4 . 2 -
)- 0 --2 _ -4 --6 --8 _
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2
teta
Рис. 11. Неустойчивость основного равновесия VII. ЗАКЛЮЧЕНИЕ
Таким образом, в работе:
1. Создана математическая модель колебаний аэродинамического маятника.
2. Показано нарушение единственности при определении угла атаки. Проведен параметрический анализ областей неоднозначности.
3. Проведены численные исследования устойчивости, в которых реализован критерий Гурвица и изображены области устойчивости при различных значениях параметров С, R. Аналогичные задачи по математическому моделированию проводились также в работах [4,6,7].
Библиография
[1] Самсонов В.А., Цыпцын С. В. Обобщение задачи о классическом флаттере упругой конструкции в потоке среды. Сборник: Механика в авиации и космонавтике, 1995, издательство Машиностроение г. Москва стр. 135
[2] Локшин Б.Я., Самсонов В.А. Авторотационные и автоколебательные режимы движения аэродинамического маятника 2013 г. Прикладная математика и механика, том 77, № 4, с. 501-513
[3] Локшин Б.Я. , Привалов В.А., Самсонов В.А. "Введение в задачу о движении точки и тела в сопротивляющейся среде". Издательство Московского университета. 1986.
[4] Паршин Д.Е. Самсонов В.А. «Качественный анализ в задаче о движении аэродинамического маятника.»1992г. (МГУ НИИмеханики, отчет 419).
[5] Табачников В.Г. «Стационарные характеристики крыльев на малых скоростях во всем диапазоне углов атаки.» Труды ЦАГИ 1974 г. выпуск 1621
[6] Strickland J.N. , Smith T. and Sun K. ,"A vortex Model of the Darrious Turbine: An analitycal and experimental stady", Final Report submitted to Sandia Laboratories on contract # 13-5602, 1982
[7] Paraschivoiu J., "Double Multiple Stremeamtube model with Recent Improvements", Journal of Energy, vol.7 no.3
Geometrical analysis in the study of body
oscillations of complex configuration in the
medium flow
V.A. Samsonov, D.V. Belyakov
Abstract - The work is devoted to the construction and study of a mathematical model of self-oscillations of an aerodynamic pendulum in a medium flow. The model of quasi-static flow of the plate by the medium was adopted as a model of the impact of the medium on the body. According to this hypothesis, the aerodynamic forces acting on the body are applied at the center of pressure. In the problem under consideration, the center of pressure is movable relative to the plate. The equations of motion for the body in question are obtained. The transition to new dimensionless variables. The violation of uniqueness in determining the angle of attack is shown A parametric analysis of ambiguity domains is carried out. All stationary points that are solutions of the equilibrium equations are found. It is shown that in the most characteristic position of equilibrium, corresponding to the state of rest, there are no ambiguity regions. The stability of the equilibrium position corresponding to the state of rest was studied, in which the Hurwitz criterion was implemented and the stability region was depicted. It is shown that aerodynamic forces for bodies with some forms can contribute to the development of self-oscillations, and for others, attenuation. In the mathematical package MATLAB 18, a set of programs is written, which makes it possible to build areas of stability and to carry out numerical integration of the equations describing body oscillations in order to confirm the adequacy of the constructed model.
Keywords - Auto-oscillations, flutter, pendulum
REFERENSES
[1] Samsonov V.A., Tsiptsyn S.V. Generalization of the problem of the classical flutter of an elastic structure in a flow of a medium. Collection: Mechanics in Aviation and Cosmonautics, 1995, publishing house Mechanical Engineering, Moscow p. 135
[2] Lokshin B.Ya., Samsonov V.A. Autorotation and self-oscillation modes of motion of an aerodynamic pendulum 2013. Applied Mathematics and Mechanics, Vol. 77, No. 4, p. 501-513
[3] Lokshin B.Ya. , Privalov V.A., Samsonov V.A. "Introduction to the problem of the motion of a point and body in a resisting medium". Moscow University Press. 1986
[4] Parshin D.E. Samsonov V.A. "Qualitative analysis in the problem of the motion of an aerodynamic pendulum." 1992. (MSU Research Institute of Mechanics, report 419).
[5] Tabachnikov V.G. "Stationary characteristics of the wings at low speeds in the entire range of angles of attack." Proceedings of TsAGI 1974 Issue 1621
[6] Strickland J.N. , Smith T. and Sun K., "A d'Arrious Turbine Turbine: An anesthesia and experimental stady", # 13-5602, 1982
[7] Paraschivoiu J., "Double Multiple Stremeamtube model with Recent Improvements", Journal of Energy, vol.7 no.3