УЧЁНЫЕ ЗАПИСКИ Ц АГИ
Т о м XV 198 4 №2
УДК 629.735.33.015.4:533.6.013.422 629.7.018.7:53.087
ИССЛЕДОВАНИЕ ФЛАТТЕРА НА ОСНОВЕ ЧАСТОТНЫХ ИСПЫТАНИЙ ПРИ ДОКРИТИЧЕСКИХ РЕЖИМАХ
Б. Д. Брянцев
В работе кратко излагаются основные алгоритмы построения уравнений колебаний упругого летательного аппарата в потоке воздуха и оценки флаттерных характеристик по результатам частотных испытаний при скоростях потока, значительно меньших критической скорости флаттера (Ккр). Приводятся расчетные и экспериментальные данные, показывающие возможности алгоритмов и их преимущества по сравнению с алгоритмами других методов оценки Ккр по результатам испытаний при скоростях потока, меньших Укр.
На рис. 1 дана упрощенная схема испытаний. Испытываемый объект снабжен системой из п датчиков, преобразующих колебания объекта в электрические сигналы, параметры которых в виде комплексных амплитуд уъ . . . , уп измеряются с помощью информационно-измерительной системы (ИИС), управляемой от ЦВМ, и запоминаются на внешних накопителях (магнитофонных лентах, дисках и т. д.). Объект возбуждается с помощью силовозбудителей (СВз, . . . , СВ,.), управляемых блоками подбора сил (БПС^ . . . , БПСГ) и задающим генератором. Величины входных колебаний контролируются с помощью датчиков, преобразующих возбуждающие силы
Рис. 1
в электрические сигналы fu . . /,.. В ходе эксперимента информа-
ция о параметрах установившихся вынужденных колебаний объекта обычно собирается в виде векторов комплексных амплитуд сигналов датчиков Ут (со,) = \_уи . . . , уп_\ и векторов возбуждающих гармонических сил /7Т (ш;[= /Л_|, при частотах возбуждения шг,
г= 1,..
Приближенно уравнения установившихся вынужденных колебаний упругого ЛА в потоке воздуха можно записать в виде:
[- С + (В, + Р V2 В,) + ] ш (А, + ? К 0,) ] Х = Ф, (1)
где С, В0, А, — к X ¿-мерные матрицы масс, жесткости и демпфирования конструкции, Вх и А — А X ¿-мерные матрицы аэродинамических коэффициентов жесткости и демпфирования соответственно, а А" и Ф — ¿-мерные векторы обобщенных координат и обобщенных возбуждающих сил соответственно.
Связь между векторами X и У, ^ и Ф может быть выражена формулами:
Х = Ы, Ф = МР,
где Ь — ¿Хя-мерная матрица, а М — £ X г-мерная матрица.
Основная задача обработки данных эксперимента состоит в том, чтобы определить матрицы Ь, М, С, В0, Ви />0, которые составляют математическую модель колебаний упругого ЛА в потоке воздуха.
Эту задачу (задачу идентификации) можно решать различными методами, например, с помощью алгоритмов, изложенных в [1-4]. В данной работе рассматривается двухэтапная процедура идентификации. На первом этапе обрабатываются результаты испытаний базового варианта объекта, для которого справедливы упрощающие идентификацию предположения, например, предположение о том, что испытываемый объект по его динамическим характеристикам можно представить системой линейных независимых осцилляторов [5—7].
На втором этапе, используя полученную во время первого этапа информацию, оцениваются те коэффициенты уравнений, которые отличаются от соответствующих коэффициентов базового варианта. Опыт использования изложенных ниже алгоритмов показал, что по сравнению с [3, 4] такой подход позволяет в несколько раз уменьшить влияние случайных ошибок измерений на оценки критической скорости флаттера.
При идентификации уравнений движения базового варианта по результатам наземных частотных испытаний известными методами (см. [5—7]) определяются: и1 — «-мерные векторы показаний датчиков при собственных колебаниях конструкции без потока (собственные формы), обобщенные массы — ¡¡.¡, собственные частоты—2, и коэффициенты демпфирования а1 (1=\ ,...,£) для к тонов собственных колебаний объекта. Далее, приняв в качестве координатных векторов «-мерные векторы Е1 (¿—1 можно определить
матрицы С, В0, по формулам:
с = фт [\ц\]Ф, я0 = Фт 1\2гц\] ф, я0 = Фч\2^\]Ф;
матрицы Ь и М — по формулам:
Ь = (Ет Е)-! Ет, М = Е\,
где £=[еи ей], е, —«-мерный вектор показаний датчиков при единичном отклонении конструкции по г-й обобщенной координате, Е9 — г X ¿-мерпая матрица форм обобщенных координат в точках возбуждения, Ф = Ч*-1, Ч* = [1Г1, .. . ¥й], ¥г — ¿-мерный вектор коэффициентов разложения собственного вектора и1 по векторам ег, . . . , гк. В принципе векторы е,, . . . ,ек могут задаваться произвольно.
При таком подходе существует произвол в выборе координатных векторов и его можно избежать, если в качестве координатных векторов принять, например, пронормированные на единичную обобщенную Массу собственные векторы ег = 11^\/~
Легко показать (см. [8]), что матрицы В — ?У2В1 и £> = р V можно определить по результатам измерений вынужденных колебаний объекта в потоке с помощью метода наименьших квадратов по формуле:
[ШЬ] = Ие (с?^т) [ие (гг7) Г1, (2)
где = • • у ?4], ^=ф,+(®? с-у«*1О0 — в0)х1,г=[г1...........zй],
г] = [Х], ]<а1 Х]\, Х^¿-мерный вектор комплексных амплитуд вынужденных колебаний и Ф,- — ¿-мерный вектор возбуждающих сил при возбуждении с частотой (ог, ¿=1, . . . . , 5; 5 > ¿. Решение (2) соответствует минимальному значению величины:
У = 2 II Ф/ — [— С + (Но + В) + /“г (А + &) ] Х^2.
г=1
Если экспериментальные данные имеют при разных сог одинаковые относительные ошибки векторов Х1 и Т7,, то рационально ввести весовые коэффициенты 5 (юг) = 1/| Ф( |, где |Фг| — модуль вектора Фг. В этом случае вместо Х1 и Фг в выражении (2) надо принять векторы
ф;=5(ш/).ф/.
Получив несколько матриц В1 и А при скоростях потока У1 (г'= 1, . . . , /га; /«> 1), можно найти среднеквадратичные значения матриц и А, например, по формулам:
/ т \—1 т'® / т \ — 1 ш
5.= 2 5,; 0! = 2! Р, У,
\/=1 / г = 1 \г = 1 / ¡=1
Таким образом, основная задача эксперимента об определении матриц Ь, М, С, В0, О0, В1} может быть решена с помощью приведенных выше алгоритмов.
Необходимо отметить, что выражение для £ соответствует наилучшей вероятностной оценке вектора обобщенных координат X по вектору показаний датчиков У, если случайные ошибки по-
казаний датчиков не коррелированы и дисперсии показаний датчиков одинаковы. В более общем случае наилучшую оценку вектора X можно получить (см. [9]), если в качестве Ь использовать выражение:
А = (£т Е)~1 ЁКт, (3)
где V — «Хя-мерная матрица собственных векторов (УКТ = /), а а2 = [\з2\] — п X я-мерная диагональная матрица собственных значений корреляционной матрицы ошибок показаний датчиков Уъ • • • > У„! а Е — а-1 V7 Е. Поскольку получение экспериментальных данных для вычисления корреляционной матрицы ошибок показаний датчиков перед каждым измерением требует больших затрат времени, то приближенно в качестве матрицы V принимается единичная матрица /, а в качестве матрицы а2 принимается диагональная матрица дисперсий показаний датчиков, которая определяется по результатам специальных испытаний и принимается неизменной в течение эксперимента. При таком подходе оценка X будет наилучшей, если ошибки показаний датчиков не коррелированы и отношение среднеквадратичных ошибок показаний датчиков сохраняется в ходе эксперимента.
Основные требования к эксперименту при наземных частотных испытаниях вытекают из методов определения иг, Ц, а(-, ¡¿г, г=1, .. .,к. Число датчиков п должно быть не меньше числа учитываемых при обработке тонов колебаний А, а расположение датчиков по конструкции должно быть таким, чтобы сЫ (Ч*т Ет Е Ч1) Ф 0 и имел наибольшую величину. При испытаниях в потоке достаточно провести измерения векторов У и Т7 при £ значениях частот ш, например, при резонансных частотах, соответствующих каждому тону колебаний. После того, как из экспериментальных данных получены уравнения (1), задача об определении границ флаттера в пространстве параметров р и V свелась к стандартной задаче об определении границ флаттера по известным уравнениям движения. Такой подход может быть использован в качестве экстраполяционного метода определения Укр. Результаты применения этого метода приведены
Г таблице, где Уисп= УИсп!Укр— относительная скорость потока, по результатам испытаний при которой производилась оценка относительного значения критической скорости флагтера—Укр. р = = УКР.Р1УКР и оценка частоты ¡флаттера — Рфл.р, а Рфл — частота флаттера, полученная в прямом эксперименте. Из таблицы видно, что: в среднем точность оценки Укр растет с приближением Уисп к Укр; метод дает удовлетворительные результаты оценок Укр как 11 при Уша, в 1,5 — 2 раза меньших Укр, так и при 1/исп « Укр; наибольшая точность оценок Укр получена в „простейшем“ случае (№ 1,. .., № 6), а относительно большие ошибки в случае с полной ДПМ (№ 7, . . . , № 10) можно частично объяснить влиянием тонов, не учтенных при обработке результатов эксперимента; метод работает как при учете минимального числа степеней свободы (два тона), так и при учете значительного числа степеней свободы (до шести тонов). Для сравнения с другими экстраполяционными методами на рис. 2 приведены интегральные распределения вероятности ошибки определения 1/кр, полученные путем обработки результатов
ЮЗ
№ ^исп у •'кр.р Р фл. р Краткая информация об объекте и „расчете“
1 0.51 0,973 8,66 Полужесткая, двухстепенная модель:
2 0,63 1,02 8,60 Укр = 31 м/с; Яфл = 8,4 Гц
3 0,73 1.012 8,49
4 0,81 1,012 8,39 Расчет с учетом двух тонов
5 0,86 1,00 8,37
6 0,9 0,996 8,32
7 0,36 0,84 6,88 Полная ДПМ тяжелого самолета:
8 0.44 0,86 6,15 Укр = 68,7 м/с; Рфл = 7,0 Гц
9 0,58 0,96 6,89 Расчет с учетом двух симметричных тонов
10 0,73 0,94 6,82 изгиба крыльев
11 0,58 0,93 6,02 ДПМ консольно защемленного крыла:
12 0,71 0,95 6,04 Укр = 62 м/с; Рфл = 6,2 Гц
13 0.82 1,002 5,54 Расчет с учетом четырех тонов
14 1,00 1,002 5,49
15 1,08 1,013 5,65
16 0,46 0,98 5,3 ДПМ консольно защемленного крыла:
17 0,66 1,04 5,0 ^кр = 54>5 м/с; Рфл = 6,2 Гц
18 0,81 1,035 5,5 Расчет с учетом трех тонов
19 0,935 0,998 6,3
20 0,4 0,935 9,97 ДПМ г. о. тяжелого самолета:
21 0,56 0,912 10,2 ^кр = 41 м/с; РфЛ =10 Гц
22 0,794 0,99 9,98 Расчет с учетом трех тонов
23 0,79 1,018 10,13 ДПМ г. о. тяжелого самолета: Ккр = 38 м/с; Яфл = 10 Гц Расчет с учетом шести тонов
испытаний моделей в дозвуковой аэродинамической трубе. Сравнивались три метода: метод № 1 (см. [10]), основанный на экстраполяции функции
+!/*£,+ Vі В2, /=■ (\^кр) = 0,
вычисляемой по собственным частотам <о1( ш2 и коэффициентам затухания а2, двух собственных тонов, определяющих флаттер; метод № 2 (см. [11]), основанный на экстраполяции коэффициентов характеристического полинома, и рассмотренный в данной работе метод № 3. Сравнение проводилось по условной ошибке экстраполяции 5 1/кр.р = (Икр.р—Укр)1(Укр — 1/исп). При использовании второго метода учитывалось от двух до четырех собственных тонов, а при использовании третьего метода учитывалось от двух до шести тонов колебаний испытываемого объекта. Из графиков видно, что третий метод дает в несколько раз меньшие ошибки определения 1/кр, чем другие методы.
Существенное преимущество третьего метода состоит в том, что на основании результатов испытаний одного варианта объекта можно производить некоторые виды параметрических исследований, например, оценивать влияние дополнительных сосредоточенных грузов, жесткостей отдельных элементов конструкции, или учитывать наличие дополнительных демпферов и дополнительных кинематических связей в конструкции объекта. Пусть £* —матрица форм обобщенных координат в точках ри . . . , р(. Легко показать (см. [12]), что введение системы дополнительных пружин, связывающих ТОЧКИ Рь имеющей Ь X ¿-мерную матрицу коэффи-
циентов влияния К, приведет к изменению матрицы жесткости конструкции В0 на величину Д В0 — Е*т КЕ*. Аналогичным образом можно вычислить эквивалентные изменения матрицы С при введении системы грузов и матрицы Д> при введении системы дополнительных демпферов, связывающих точки ри . . . , р(. Учет кинематических связей, задаваемых уравнением 8Ур—0, где 5—ЛХ^-мерная матрица кинематических связей, а Ур—/-мерный вектор перемещений в точках ри . . . , рг, приводит к необходимости исключения „лишних“ координат с помощью £ X (£ — $)-мерной матрицы линейного преобразования /?, где я = Иё (5).
Приведенные выше алгоритмы для параметрических исследований точны для систем с £ степенями свободы. Однако реальный ЛА можно представить такой системой только приближенно и поэтому возможны ошибки при оценке влияния перечисленных выше модификаций объекта. Для иллюстрации этих эффектов на
Исходный Вариант №тах=0,027
'-х==
ЬЯтах=0,02
1=Т=
5.
0 N V тах 0,15 Ь
N V N г 1
Т^Г
Рис. 3
^тах = 0,ОШ
й*пахт№
рис. 3 приведены расчетные примеры с балками единичной длины. В расчетах учитывалось различное число собственных тонов исходной балки. После введения связей получены значения собственных частот 2, балок с различными условиями заделок. Поскольку наибольшие ошибки возможны при введении кинематических связей, то только эти связи и рассматривались. Расчеты показали, что в случаях с безмоментными связями результаты оценки собственных частот хорошо стыкуются с точными решениями, а в случаях с моментными связями (примеры №№ 4, 5, 7, 8) относительные ошибки оценок собственных частот Д2 могут быть весьма значительны. Характер зависимости Д2 от числа учитываемых тонов виден на рис. 4, где ^— число упругих тонов балки после введения кинематических связей.
Для иллюстрации эффективности методики в условиях реального эксперимента на рис. 5 приведены результаты оценки влияния на величину критического скоростного напора ^кр жесткостей пружин К и С для жесткого отсека крыла, установленного между двумя стенками. Все три кривые <7кр (С) получены „расчетом“ по результатам испытаний без потока и в потоке при <7ИСП ~ 0,55 при значении жесткостей пружин К = 1 и С= 1. Из графика видно, что
с
К
1,5 ~
>
Эксперимент
0,25
0,5
Точка испытании.
Точки испытании [К=1,0)
О
ОМ 0,6 о,в
о 11 с
о
г
3 к
Рис. 5
Рис. 6
полученные результаты удовлетворительно согласуются с результатами прямого эксперимента.
На рис. 6 приведены результаты оценки влияния жесткости проводки К на <7кр для ДПМ горизонтального оперения тяжелого самолета.
В этом эксперименте „расчет“ проводился с учетом трех симметричных тонов упругих колебаний. Из графика видно, что результаты „расчета“ и прямого эксперимента удовлетворительно согласуются в диапазоне реальных величин жесткости проводки управления.
Опыт использования рассмотренного метода и приведенные выше результаты показывают, что метод может быть применен в практике исследований флаттера и с его помощью можно получать оценки критической скорости флаттера и оценивать влияние некоторых параметров конструкции и потока по результатам испытаний при скоростях потока, в 1,5 — 2 раза меньших Укр.
1. GaukrogerD. R. SkingleC. W., Heron К. H. An application of system identification to flutter testing. — Journal of Sound and Vibration. 1980, 72 (2).
2. Heron K. H., Gaukroger D, R., Skingle C. W. The derivation of equation of motion from response data and its application in flutter testing. 1973, RAE Technical Report 73051.
3. И л ь и ч е в В. Д. Линейная идентификация в аэроупругости.— Ученые записки ЦАГИ, 1972, т. III, № 4.
4. Назаров В. В. Идентификация неконсервативной упругой конструкции.—Ученые записки ЦАГИ, 1972, т. III, № 4.
5. Gauzy Н. Vibration testing by harmonic exitation. — Manual on Aeroelasticity, vol. 1. Structural analysis, 1964 (AGARD).
6 Beatrix Ch. La determination experimental des caractéristiques vibration des structures ONERA NOTE Technique № 212 (1973).
7. Шибанов P. А. Метод анализа результатов частотных испытаний.—Труды ЦАГИ, вып. 1188, 1970.
8. Брянцев Б. Д. Определение аэродинамических коэффициентов уравнений колебаний упругого летательного аппарата по результатам частотных испытаний в потоке. — Труды ЦАГИ, 1978, вып. 1926.
ЛИТЕРАТУРА
9. Андерсен Т. Статистический анализ временных рядов.— М.: Мир, 1976,
10. Zimmerman N. H., Weissenburger J. T. Prediction о! flutter onset speed based on flight testing at subcritical speed. — Journal of Aircraft, at july — august, 1964.
11. Брянцев Б. Д., Карклэ П. Г. Некоторые результаты определения критической скорости флаттера экстраполяционными методами.—Труды ЦАГИ, 1976, вып. 1772.
12. Breit bach E. Bestimmung der Kopplungsglieder und der Bewegungsgleichungen einer aus mehreren elastomechenischen Teilsystemen Zusammengesetzten Flugkonstruktion. — Zeitschrift für Flugwissenschaften № 9/10, 1970.
Рукопись поступала 24/IX 1982