Труды МАИ. Выпуск № 95
http://trudymai.ru/
УДК 519.62; 629.7.01, 623.465.7
Применение метода Ритца и метода конечных элементов к расчёту аэроупругих колебаний крылатой ракеты
Благодырёва О.В.
Корпорация « Тактическое ракетное вооружение», ул. Ильича, 7, Королёв, Московская область, 141080, Россия e-mail: [email protected]
Аннотация
Представлены математические модели аэроупругих колебаний крылатой ракеты, разработанные на основе метода Ритца и метода конечных элементов (МКЭ). Крыло ракеты рассматривается как упругая балка, работающая на изгиб с поперечным сдвигом и на кручение, а корпус и стабилизатор ракеты считаются абсолютно жёсткими. Аэродинамические нагрузки определяются на основе квазистационарной теории плоско-параллельного обтекания поперечных сечений крыла. Расчёты произведены в программной среде «WolframMathematica 8».
Ключевые слова: аэроупругие колебания, флаттер, метод Ритца, метод конечных элементов, моделирование колебаний летательного аппарата.
Введение
Исследование аэроупругих колебаний крылатой ракеты в дозвуковом потоке воздуха является одним из важнейших этапов её проектирования, поскольку
эффективность применения ракеты напрямую зависит от её аэродинамической устойчивости. Возникновение незатухающих упругих колебаний несущих поверхностей, органов управления, фюзеляжа или всей ракеты в целом в процессе полёта может привести к потере управляемости и даже разрушению.
Расчёт аэроупругих характеристик любого летательного аппарата представляет собой задачу, точное аналитическое решение которой получить в явном виде достаточно непросто, поскольку зачастую летательный аппарат представляет собой сложную конструкцию.
Для исследования аэроупругой устойчивости крылатой ракеты в рассматриваемой работе используются МКЭ [1] и метод Ритца [2]. При написании статьи были изучены дополнительные источники литературы. Пример решения задачи определения собственных частот по формуле Релея с помощью метода последовательных приближений приводится в работе [3].
Для статической прочности получены формулы оптимальных геометрических параметров панели, при которых не происходит потери упругой устойчивости [4].
Пример расчёта флаттера отдельных частей летательного аппарата в дозвуковом потоке на основе метода Ритца можно найти в статьях [5], [6].
Задаче об устойчивости крыла самолёта, находящегося в набегающем установившемся потоке воздуха с помощью МКЭ посвящены работы [7, 8]. В [9] рассматривается обтекание потоком газа прямоугольного крыла большого удлинения, моделируемого тонким стержнем. Дифференциальные уравнения изгиба, поперечного сдвига и кручения, а также изгиб-но-крутильных колебаний
крыла большого удлинения с учетом конусности, т.е. с учетом наклона образующих оболочки по отношению к продольной оси получены в [10].
В статье [11] представлена упрощенная методика расчета колебаний панелей фюзеляжа самолета, проведена верификация результатов расчетов на собственные формы колебаний прямоугольной пластины и отсека фюзеляжа.
Общая постановка и схема решения задачи анализа аэроупругой устойчивости складных рулей описана в [12].
Аналитические зависимости угла атаки, угла наклона вектора скорости и угла тангажа ракеты по времени на пассивном участке управляемого полета ракеты приведены в [13]. Аналитическое определение углового положения ракеты с учётом изгибных колебаний корпуса рассматривается в [14].
Нестационарные аэродинамические нагрузки (приращение давления и обобщенные силы), действующие на деформируемый тонкий профиль крыла при его произвольном движении в несжимаемом потоке с поперечными порывами ветра, определяются в [15].
Постановка задачи
Рассмотрим задачу об определении аэроупругой устойчивости упругой ракеты, совершающей продольное короткопериодическое движение в несжимаемом потоке воздуха. Ракету будем исследовать на основе описанной в [16] балочной схемы, включающей в себя фюзеляж, две прямые консоли крыла и отклоняемые органы управления - рули высоты. Приводы с проводкой управления будем считать отсоединёнными, и их действие на конструкцию заменим неизвестными управляющими силами.
Уравнения малых колебаний системы запишем в форме уравнений Лагранжа второго рода в обобщённых координатах. В соответствии с [17,18] в качестве обобщённых координат будем рассматривать нормальные координаты, представляющие движение по собственным формам колебаний свободной конструкции с закреплёнными органами управления.
Для решения задачи и анализа динамического поведения и устойчивости упругой системы [19,20,21] воспользуемся двумя методами: методом Ритца и методом конечных элементов, - и проведём сравнительный анализ полученных результатов расчёта на конкретном примере.
Ракета представляет собой систему, состоящую из фюзеляжа, двух прямых консолей крыла и двух консолей руля высоты, рассматриваемых здесь как недеформируемые цельноповоротные стабилизаторы (Рисунок 1).
г
Рисунок 1. Схема крылатой ракеты В возмущённом потоке воздуха на корпус и несущие поверхности ракеты действуют аэродинамические силы, влияющие на устойчивость движения ракеты.
Описание модели ракеты. Уравнение колебаний.
Кроме того, следует учесть продольное сжатие корпуса ракеты под влиянием силы тяги двигателя.
В общем случае довольно сложно описать возмущённое движение упругой ракеты, имеющей бесконечное число степеней свободы, с учётом геометрических связей. Ограничимся набором ^ обобщённых координат « (7), I = 1,...,^,
соответствующих ^ степеням свободы. Уравнения колебаний ракеты в обобщённых координатах составим в форме уравнений Лагранжа второго рода:
Ж
З— ^ —
=а,'=1, *, (1)
дд,
где Т - кинетическая энергия системы, — - потенциальная энергия системы, а - обобщённая сила, соответствующая обобщённой координате «
Обобщённые координаты « (7) являются скалярными функциями времени и
могут представлять собой любые перемещения и углы поворота системы, характеризующие все ^ степеней свободы системы. Симметричные относительно плоскости хОу колебания фюзеляжа описываются с помощью поперечных перемещения абсолютно жёсткой ракеты « (7) и угла поворота относительно центра
тяжести «2 (7). Горизонтальное оперение представляет собой недеформируемый
цельноповоротный стабилизатор - руль высоты. Угол отклонения руля относительно корпуса ракеты описывается углом поворота относительно оси вращения стабилизатора 8 (7):
8, (7 ) = 8: (7 ) + «3 (').
где 5пр (t) - заданное (программное) значение угла отклонения руля; (t) - отклонение руля за счёт упругости проводки управления.
Консоли крыла симметричным образом изгибаются и закручиваются. Предполагается, что упругие колебания консолей крыла происходят в соответствии с гипотезой плоских сечений без сдвига.
В случае симметричных колебаний достаточно рассматривать только одну консоль крыла. Полученные для неё выражения кинетической и потенциальной энергий и вариации работы аэродинамических сил удваиваются. Интерференцией между крылом и горизонтальным оперением пренебрегаем.
Симметричные колебания ракеты описываются функциями уф (х, t), V (t),
(р(2, t), представляющими собой поперечные перемещения оси фюзеляжа и оси крыла и угол закручивания крыла.
Кинематические условия абсолютно жёсткого сопряжения крыла с фюзеляжем ракеты имеют вид:
V (0, t) = V,, (0, t), ((0, t) = -^(0, t). (2)
Кинетическая энергия ракеты складывается из кинетических энергий фюзеляжа, крыла и горизонтального оперения:
Т = Тф + 2 Ткр + 2 Тго,
где
^ 'же ^ £
тф = т I Мф^фЛх, Ткр = -1(ту2 -1тауф + 3ф2)й2, (3)
2 -/„ 2 0
Тго = \ {т^ф (хс, 0 - (, О 4 + )>
] т
2
ас = Хц.ш. рл - ХС , = - ^^ (X > *) •
Здесь
т^ - погонная масса фюзеляжа, /н, 1хв - расстояния до носового и хвостового
сечений фюзеляжа, соответственно;
т - погонная масса крыла, а - расстояние между осью жёсткости крыла и линией центров масс крыла (положительное, если линия центров масс лежит за осью жёсткости), J - погонный массовый момент инерции крыла относительно оси жёсткости крыла, Ь - размах крыла;
т - масса руля высоты (стабилизатора), х - координата центра тяжести
руля, ас - расстояние между осью жёсткости руля высоты и линией центров масс
руля высоты (положительное, если линия центров масс лежит за осью жёсткости), хс - координата оси вращения руля, Jc - момент инерции консоли руля высоты,
Эс - угол поворота руля.
Потенциальную энергию ракеты также представим в виде суммы энергий всех её компонентов
П = Пф + 2 ПКр + 2 П го,
где
1 /хв 1 Ь 1
Пф =11 Е1фу;Чх, Пкр = - КБ1у"2 + ОТ^'2) dz, Пго = - ксд23 • (4)
2 - /„ 2 о 2
Здесь и Е1 - погонные изгибные жёсткости фюзеляжа и крыла, 01 -погонная крутильная жёсткость крыла, Пго - потенциальная энергия за счёт упругости проводки управления, - коэффициент жёсткости проводки управления.
В случае абсолютно жесткой проводки управления, отклонение руля д3 за счёт упругости проводки управления будет отсутствовать, и поворот на программный угол будет осуществляться без изменения потенциальной энергии, поэтому в рассматриваемой задаче Пго = 0.
Вариация работы аэродинамических нагрузок, действующих на корпус (фюзеляж), крылья и рули высоты будет иметь вид:
ЗАа/д = дАфа/д + 2дАкрд + 28Крд, (5)
где
1хв
дАф/д = | qyдvфdx, дА^ = AYcдvф (хе, 0 + ЛМ,,
дАкрд = | (AYдv + ЛМд() dz
(6)
Здесь ц (х, t) - подъемная сила, действующая в текущем сечении фюзеляжа:
Цу (х't) = ^ (х)"(X1:), где а - местный угол атаки; £0 - площадь миделя фюзеляжа; С" (х) - коэффициент подъёмной силы фюзеляжа.
- 1н
0
При обтекании воздушным потоком корпуса, установленного под малым положительным углом атаки а, возникает нормальная сила, пропорциональная углу атаки. Согласно теории тонких осесимметричных тел нормальная сила появляется только на участках корпуса с переменной площадью поперечного сечения ^. [22]
В рассматриваемой задаче поперечное сечение корпуса имеет постоянную площадь, а изменением площади поперечного сечения носовой части можно пренебречь, и, следовательно, нормальная сила равна нулю.
Аэродинамические нагрузки, действующие на крыло при использовании гипотезы плоского обтекания профиля, записываются в виде
2
т2
ДМ = ХЬ% (<р,у) + £¥-со*%Ъ2Вх (ф,у).
(7)
где
А [(, V ] = с> + ( с3;Ьр'-а') tgZ, Во [р,V] = тар + (т^Ьр'-т>')tgх,
А[ФА = ~су*+[с? +с;)ьф+{с;*ъф'- с**) ъ
Вх [ф,у\ = -тауу + {т^ + т"^Ьф + {т^Ьф'-та2у'^Ъ
(8)
Здесь х - угол стреловидности крыла.
В случае применения квазистационарной теории коэффициенты будут равны
- 7Г 7Г - - - - -
са=2я, С =-, <=-, <г =0, с"=с? =та=та/ =0. (9)
у У У z z V/
Аэродинамические нагрузки, действующие на руль высоты, записываются как
<
ЛУ = У аас + У :а2, ЛМ 2,с = ма ас + М?:;
где коэффициенты вычисляются по формулам:
Уа
рУ
рУ
с и : л ла рУ с •< г, „а л л: рУ о г.2—:
Бсса, У : БсЬссю;, ма
2 су г^ с с у ^ 2 2
■ЗА<, М? = БЫт:
с с 2
Здесь Б - площадь консоли руля, Ьс - средняя аэродинамическая хорда руля, ас - местный угол атаки руля высоты, а2 - частота вращения руля; Эс - угол
тангажа, ю2 - приведённая частота:
А, #с=*р-*фМ>
- : Ь
У
, а
V
ф(Хс>*) + ^ф(Хс>*)
Вариация работы аэродинамических нагрузок (5) в обобщённых координатах
приводится к виду
дАа/д = -^а/д,гдЯг •
(10)
¿=1
Сила тяги двигателя ракеты Т, согласно [23], представляет собой «следящую» силу, остающуюся в процессе поворота корпуса при колебаниях перпендикулярной поперечному сечению, в котором она приложена.
Эту силу можно разложить на две составляющие (рисунок 2): продольную Т и поперечную Т •
Тх - Т,
ту=-^а( с).
т.
т
Рисунок 2
-х
2
я
Продольная составляющая силы Тх в процессе деформации не меняет величину
и направление и обладает потенциалом ^ 0, равным по модулю работе А^ по
N 0
продольному сжатию корпуса и противоположным ей по знаку:
UN0 = АN0
^ 1хв
- | N0 (х ) V, dx,
- /„
где N° (х) - усилие, с которым корпус ракеты подвергается продольному сжатию.
Вариация работы продольного сжатия имеет вид:
1хв
дАN0 = | N0 (х) v'фдvаdx,
(11)
или в обобщённых координатах: дА 0 = 0 ддг.
¿=1
Эти усилия N0 (х) в каждом сечении находятся из уравнения равновесия сил, действующих на отсечённую часть (рисунок 3) при V = 0 в направлении оси х:
1хв
N0 ( х )+ | т (х ) adx = Т,
где а - ускорение движения ракеты.
В статье рассматривается случай, когда погонная масса фюзеляжа т (х) постоянна по
длине, и, следовательно,
mф(x)•a
dx
T
Т
N 0 ( х ) = — (+ х )
Рисунок 3
Вариация работы поперечной составляющей следящей силы имеет вид:
х
x
l
^ =-Ту'Ф{1в )Зуф (1яв). (12)
а
Или, если привести к обобщённым координатам: 8АТ = ^0 г 8чг.
у ¿=1 у
Уравнения колебаний ракеты составляются как уравнения Лагранжа (1), где 0 = йыы + ^о,г + йТу, .
Определение критической скорости флаттераракеты с закреплёнными рулями
высоты.
Метод Ритца
На основании метода Ритца неизвестные функции представим в виде разложений по обобщённым координатам:
уф (* ) = ЁЧ (* )ъ (х)'
¿=1
а
V (* *) = Х Ч (*) 1г (*)
¿=1
а
Н * * ) = Е Чг ( *) Фг ( * )
(13)
г=1
где - обобщённые координаты; ^, , ф - заданные функции, удовлетворяющие
условиям сопряжения, у = 1,.. ; л- - размерность задачи.
В силу условий сопряжения крыла с фюзеляжем функции по двум первым координатам будут иметь вид:
V
11 =1 А =1 ^ =0,
1 = Хц.т. - Х, ./2 = Хц.т., Ф2 = 1
Подставляя (13) в (3), (4), запишем кинетическую и потенциальную энергию следующим образом:
т 1хх<"'/'/ = , П = = 1чткч,
2 1=1 / =1 2 2 1=1 /=1 2
где
ky = J ЫфЛ' ' г ^'7 dx + 2-J| EI f' ' г f' ' 7 + ^ФФ}
L у
dz.
mj = J тфП{П1dx + 2 -J(m - m c (ЪФ + f jфi)+ Jфф])dz +
4 0
+ 2• ~{тЛ (xc)lj(xc)~(xc)(snp- v'j 00) +
+j c ((s; f - 2 ¿;фс )+фс) л] (xc))).
Здесь m - обобщённая масса, £y - обобщённая жёсткость, M - матрица
инерции К - матрица жёсткости, q = [g^,...,^]1 - вектор обобщённых координат,
дпр = const.
Вариация работы аэродинамических нагрузок (10) в матричной форме имеет
вид:
¿Aa/d=-6qT(Bq + Dq),
где B - матрица аэродинамической жёсткости, D - матрица аэродинамического демпфирования.
- 1н
о
L
Отклонение рулей высоты учитывается при записи вариации работы соответствующих аэродинамических сил:
где Ь - вектор-столбец аэродинамической жёсткости руля высоты, возникающей
из-за заданного (программного) угла отклонения руля; ^ - вектор-столбец
аэродинамического демпфирования руля высоты, возникающего из-за заданного (программного) угла отклонения руля.
Вариация работы продольного сжатия ЗА о (11) и вариация работы
поперечной составляющей следящей силы 5АТ (12) в обобщённых координатах в
матричном виде выглядят следующим образом:
5Амй = бдтСд, ЗАТу = -бд^.
Таким образом, вариацию работы силы тяги двигателя можно записать в виде:
ЗАТ = ЗА^ + 6АТ = бд^д,
где С = С - 8.
Метод конечных элементов
Для более точного исследования поведения изгибно-крутильных колебаний крыла дополнительно можно воспользоваться методом конечных элементов (МКЭ) и вычислить соответствующие матрицы инерции, жёсткости, аэродинамической жёсткости и аэродинамического демпфирования для крыла.
Разобьём крыло на достаточно большое число п отсеков - конечных элементов (КЭ) и будем считать в пределах каждого из них геометрические, жесткостные и массовые характеристики осреднёнными и постоянными.
Колебания крыла в пределах длины ^го КЭ описываются функциями перемещения V, угла поворота 3 и угла закручивания <:
V = с0 + с2<т+Сз<т, 3 = С + 2С2 4+ зс3 42, < = С4 + С5£.
V =
я
"к-1
(14)
где, индекс «к» относится к ^му КЭ, Я - некоторая характерная длина, общая для всех КЭ.
Вектор неизвестных коэффициентов С = [С0 С С2 С3 С4 С5 Г выражается
через вектор обобщённых д(к^ = 3к_х <рк_х ук 3к <]т координат ^го КЭ с
помощью матричного преобразования
С = Ак q
(к)
(15)
которое записывается на основании граничных условий на краях КЭ:
у = ^к-1, 3 = 3к-1, < = <-1 пРи 4 =0; V = ук, 3 = 3к, < = < пРи 4 =1.
Матрица А имеет вид:
А к =
1 0
-3 2 0 0
0
0 0
0 0 0 0 1
0 0 3 -2 0 0
0 0
0 0
у
<
I
к
Кинетическая энергия и потенциальная энергия изгибно-крутильных колебаний k-го КЭ имеют вид:
Л Zk Л Zk
ТГ = - J (mkv2-2mk<Tkv<p + Jk<p2)dz, П* =- J (ткЭ'2 + GJ Kpcp'2)dz, (16)
zk-1 zk-1
где m - погонная масса, Jk - погонный массовый момент инерции относительно оси жёсткости, <гк - расстояние от центра тяжести до центра жёсткости. С учётом (14) и (15) из (16) получаем
ТГ=^№М*Ч(1,> nf^qW^q« (17)
где M - матрица инерции k-го КЭ. K - матрица жёсткости k-го КЭ.
Погонная сила и погонный крутящий момент аэродинамических нагрузок, действующих на крыло при изгибно-крутильных колебаниях, определяются на основе гипотезы плоского обтекания поперечных сечений z = const по нестационарной теории для гармонических колебаний, описываемой формулами (7),
(8), (9).
Вариация работы аэродинамической нагрузки, действующей на k-ый КЭ
zk
SAKp = J (AYSv + AMS^) dz
zk-i
с учётом (2), (3), (4), (6) приводится к виду
SAf = -6q(A)T (BfqW +Dfq(*)). (18)
Здесь B¡р, DKP - матрицы аэродинамической жёсткости и аэродинамического демпфирования для k-го КЭ.
Выражения для кинетической и потенциальной энергий и вариации работы аэродинамических нагрузок для консоли получаются суммированием по всем КЭ:
Ткр = ^Ткр , Пкр = , 5Акр =
к , П = ^п , 5А = • (19)
Выражения (19) с учётом (17) и (18) принимают вид:
ТКР = ^4, Пкр = ^qтKкp Ч, ЗАкр = -бЧт (в*р д + Б*р д);
где д - вектор обобщённых координат крыла:
Ч = «90 % срх у2 32 ср2 ... V, Зп
Здесь матрицы инерции Мкр, жёсткости Кккр, аэродинамической жёсткости Вкр и аэродинамического демпфирования Бкр для крыла образуются из соответствующих матриц отдельных элементов.
Уравнение колебаний
Подставляя в (1) кинетическую энергию, потенциальную энергию и вариацию аэродинамических сил, полученные на основе метода Ритца или МКЭ, уравнение колебаний запишется в матричном виде следующим образом:
м<1+1)<1+(к-с+в^ = а])яр+ъ])яр, (20)
где матрицы инерции, жесткости, аэродинамического демпфирования и аэродинамической жёсткости включают в себя компоненты для каждой из составляющих частей ракеты:
М = Мф + 2Мкр + 2Мрл, К = Кф + 2 Ккр + 2 Крл, В = Вф + 2Вкр + 2Врл, Б = Бф + 2Бкр + 2Брл.
При закреплении недеформируемого руля с абсолютно жесткой проводкой управления 5 (7) = 0 уравнение (20) принимает вид:
Мс1 + Б4 + (К-С + В)д = 0. (21)
Для решения уравнения (21), его удобно преобразовать к уравнению первого порядка в фазовом пространстве координат:
Аг + Сг = 0, (22)
" М 0" " Б К - О + В" " Р"
где А = , с = , Г =
0 Е -Е 0
Для отыскания корней уравнения (22) решается задача о собственных значениях системы:
[ЛА + С]V = 0 ЛрУр] = \...^.
Нахождение корней Л. характеристического уравнения однородной системы
(22) позволяет судить об устойчивости движения ракеты на заданной высоте и при заданной скорости полёта. Решение характеристического уравнения представляет собой пары комплексно-сопряжённых Л . = а. + г ¿а., Л] = а. - г ^ и действительных
Л=ак собственных чисел, которые можно представить как частный случай комплексных.
При достаточно малых скоростях движение ракеты будет устойчивым, и действительные части всех корней будут отрицательными. При увеличении скорости полёта V может оказаться, что при некоторой критической скорости V = V действительная часть одного из корней обратится в нуль и при дальнейшем
увеличении V станет положительной. В этом случае система станет неустойчивой.
В зависимости от изменения корней при дальнейшем увеличении скорости может возникнуть динамическая неустойчивость (флаттер) или статическая неустойчивость (дивергенция). При а. > 0 и а ф 0 система будет динамически
неустойчива, при а. > 0 и а. = 0 - статически неустойчива. Граница устойчивости
определяется из условия а. = 0: на границе динамической устойчивости Л = +а.;
на границе статической устойчивости Л = 0. Величина а. < 0 характеризует запас
устойчивости по у -й форме.
Пример расчёта
В качестве примера рассмотрим ракету, геометрические, инерционно-массовые и жесткостные характеристики которой представлены в таблицах 1, 2, 3.
Фюзеляж Таблица 1.
1 [м] 1 хв [м] тФ [кг/м] Е1ф [Н-м]
2,54 3 178,68 2,95485х107
Крыло Таблица 2.
Ь [м] т [кг/м] а [м] 3 [кг-м] Е1 [Н-м] 01 [Н-м]
1 26,3543 -0,0103846 0,205251 84840 99586,5
Руль Таблица 3.
т [кг] [м] Л [кг-м2]
10,0844 0,00784153 0,0668728
Для ракеты в условиях свободного полёта первые две обобщённые координаты представляют собой поперечное перемещение абсолютно жёсткой ракеты и угол поворота ракеты относительно центра тяжести (приращение угла тангажа). Эти формы колебаний ортогональны как между собой, так и с формами упругих колебаний крыла. Формы перемещений твёрдого тела имеют нулевые собственные частоты (юх=ю2 = 0). Две низшие собственные частоты упругих колебаний ракеты, вычисленные по методу Ритца и с помощью МКЭ, приведены в таблице 4. При расчёте по методу Ритца для аппроксимации прогиба и угла закручивания консоли крыла ограничимся только двумя функциями:
1 ( 2 3 4 Л ( 2\
/ (—) = 1| б^2 - 4—- + —4 , Ф, (—) = [ 2 — - —.
ч Ь2 1 гу
1 1
2
V /
При этом ф (—) = 0, / (—) = 0.
Собственные частоты упругих колебаний крыла Таблица 4.
Метод Ритца МКЭ
п = 4 п = 5 п = 8 п = 9 п = 10
Изгиб крыла [Гц] 32.3905 32.2388 32.2382 32.2378 32.2377 32.2377
Кручение крыла [Гц] 176.449 176.122 175.72 175.286 175.228 175.186
Для ракеты с неподвижными, закреплёнными рулями построены графики изменения действительной и мнимой частей собственных частот:
• в зависимости от изменения скорости полёта V при нулевой силе тяги Т = 0 (рисунок 4);
• в зависимости от изменения силы тяги Т при нулевой скорости полёта V = 0 (рисунок 5).
Графики построены для расчётов, выполненных с помощью МКЭ (количество элементов п = 4,5,8,9,10) и метода Ритца. Из рисунков видно, что оба метода дают схожие результаты вычислений. При этом МКЭ даёт более точный результат, чем метод Ритца, и с увеличением числа конечных элементов значения собственных частот приближаются к точным значениям сверху. МКЭ показывает результаты более близкие к аналитическому решению задачи аэроупругости, так как при расчёте по методу Ритца для аппроксимации изгибно-крутильных колебаний крыла используются только две функции. С помощью увеличения числа аппроксимирующих функций точность результатов вычисления по методу Ритца можно повысить.
Рисунок 4. Графики изменения собственных частот в зависимости от изменения скорости полета V = 0 - 340 [ м/с]
при нулевой силе тяги Т = 0 [ Н ].
Рисунок 5. Графики изменения собственных частот в зависимости от изменения силы тяги Т = 0 —108 [Н ]
при нулевой скорости полета V = 0[ м/с ]. На рисунке 4 можно определить скорость полета, при которой возникает
динамическая неустойчивость (флаттер) - это область изменения знака
действительной части собственного значения, соответствующего крутильному тону
колебаний крыла, с отрицательного значения на положительное. Соответствующее
критическое значение скорости полета по методу Ритца составляет 254.88 [м/с] при
частоте 175.168 [Гц], а по МКЭ (при п = 10) 277.93 [м/с] при частоте 172.078 [Гц].
Аналогичным образом на рисунке 5 можно определить критическое значение силы
тяги.
Выводы
В статье представлен сравнительный расчет флаттера упругой ракеты с помощью метода Ритца и МКЭ: построены графики изменения собственных частот ракеты в зависимости от изменения скорости полета и силы тяги двигателя; определены критические скорости и области устойчивости полета ракеты соответственно для каждого из методов. МКЭ показывает результаты более близкие
к аналитическому решению задачи аэроупругости, так как при расчёте по методу Ритца используются только две аппроксимационные функции.
Библиографический список
1. Образцов И.Ф., Савельев Л.М., Хазанов Х.С. Метод конечных элементов в задачах строительной механики летательных аппаратов. - М.: Высшая школа, 1985. - 392 с.
2. Гроссман Е.П. Курс вибраций частей самолёта. - М.: Государственное издательство оборонной промышленности, 1940. - 312 с.
3. Замлинський О.С. Алгоритм автоматичного розрахунку частот i форм при коливаннях крила // 1нформацшш системи, мехашка та керування: науково-техшчний збiрник. 2011. №. 6. С. 89-95.
4. Кичеев В.Е. Проектирование сжатых панелей легких самолетов с использованием критерия подобия // Труды МАИ. 2007. № 27. URL: http://trudymai.ru/published.php?ID=34006
5. Благодырёва О.В. Расчёт на безопасность от флаттера крыла малого удлинения методом полиномов // Труды МАИ. 2013. № 68. URL: http: //trudymai.ru/published.php?ID=41717
6. Благодырёва О.В. Исследование флаттера композитного крыла // Труды МАИ. 2014. № 74. URL: http://trudymai.ru/published.php?ID=49345
7. Мартынова А.Д. Определение критической скорости флаттера и дивергенции крыла самолёта с использованием метода конечных элементов //
Молодёжный научно-технический вестник. 2014. № 10. http: //sntbul .bmstu.ru/archive. html
8. Войтышен В.С., Семенов В.Н. Решение задачи аэроупругости в переменных метода конечных элементов // Известия Коми научного центра УрО РАН. 2013. № 4(16). С.68-72.
9. Желтков В.И., Чыонг Ван Хуан. Определение критических скоростей прямого крыла большого удлинения // Известия Тульского государственного университета. Естественные науки. 2014. № 3. С. 71-80.
10. Гришанина Т.В. Расчет деформаций и колебаний крыльев большого удлинения с учетом конусности // Известия вузов. Авиационная техника. 2004. № 2. С. 10-13.
11. Разбегаева И.А. Расчет и оценка частотных характеристик обшивки фюзеляжа. Защита от шума // Труды МАИ. 2011. № 45. URL: http: //trudymai .ru/published.php?ID=25552
12. Неделин В.Г. Анализ конструктивно-технологических решений складных рулей с учетом требований аэроупругой устойчивости // Труды МАИ. 2012. №52. URL: http://trudymai.ru/published.php?ID=29424
13. Пархаев Е.С., Семенчиков Н.В. Некоторые вопросы оптимизации профиля крыла малоразмерного беспилотного летательного аппарата // Труды МАИ. 2015. № 80. URL: http ://trudymai .ru/published.php?ID=56884
14. Комиссаренко А.И. Аналитическое определение углового положения ракеты с учётом упругих колебаний // Вестник Московского авиационного института. 2015. Т. 22, № 2. С. 36-41.
15. Гришанина Т.В., Шклярчук Ф.Н. Неустановившиеся колебания деформируемого профиля крыла в несжимаемом потоке. // Известия вузов. Авиационная техника. 2009. № 2. С. 3-7.
16. Шклярчук Ф.Н. Аэроупругость самолета: Учебное пособие. М.: Издательство МАИ, 1985. - 77 с.
17. Гришанина Т.В., Шклярчук Ф.Н. Динамика упругих управляемых конструкций. - М.: Издательство МАИ, 2007. - 328 с.
18. Гришанина Т.В., Шклярчук Ф.Н. Избранные задачи аэроупругости. - М.: Издательство МАИ, 2007. - 48 с.
19. Рабинович Б.И. Прикладные задачи устойчивости стабилизированных объектов. - М.: Машиностроение, 1978. - 232 с.
20. Четаев Н.Г. Устойчивость движения. - М.: "НАУКА" Фармалит, 1990. -
176 с.
21. Микишев Г.Н., Рабинович Б.И. Динамика тонкостенных конструкций с отсеками, содержащими жидкость. - М.: Машиностроение, 1971. - 564 с.
22. Лебедев А.А., Чернобровкин Л.С. Динамика полета беспилотных летательных аппаратов. - М.: Машиностроение, 1973. - 616 с.
23. Гришанина Т.В., Шклярчук Ф.Н. Колебания упругих систем: Учебное пособие. - М.: Издательство МАИ, 2013. - 100 с.