В этом примере показано, как оценить частотно-ответные функции (FRF) и модальные параметры на основе экспериментальных данных. В первом разделе описан смоделированный эксперимент, который возбуждает систему с тремя степенями свободы (3DOF) с последовательностью ударов молотком и регистрирует результирующее смещение. Для трех режимов структуры оценивают частотно-ответные функции, собственные частоты, коэффициенты демпфирования и векторы формы моды. Вторая секция оценивает векторы формы режима из оценок частотно-ответной функции в эксперименте с лопаткой ветровой турбины. Конфигурация измерения лопаток турбины и результирующие формы режимов визуализируются. В этом примере требуется панель инструментов идентификации системы (TM).
Возбуждение молотка с одним входом/одним выходом
Серия ударов молотком возбуждает систему 3DOF, и датчики регистрируют результирующие смещения. Система пропорционально демпфируется, так что демпфирующая матрица представляет собой линейную комбинацию матриц массы и жесткости.

Импортируйте данные для двух наборов измерений, включая сигналы возбуждения, ответные сигналы, временные сигналы и функции частоты-отклика истинности земли. Первый набор ответных сигналов, Y1, измеряет смещение первой массы, и второй, Y2, измеряет вторую массу. Каждый сигнал возбуждения состоит из десяти конкатенированных ударов молотка, и каждый ответный сигнал содержит соответствующее смещение. Продолжительность каждого сигнала удара составляет 2,53 секунды. Аддитивный шум присутствует в сигналах возбуждения и отклика. Визуализируют первый канал возбуждения и отклика первого измерения.
[t,fs,X1,X2,Y1,Y2,f0,H0] = helperImportModalData(); X0 = X1(:,1); Y0 = Y1(:,1); helperPlotModalAnalysisExample([t' X0 Y0]);

Вычислите и постройте график FRF для первого канала возбуждения и отклика с точки зрения динамической гибкости, которая является мерой смещения над силой [1]. По умолчанию FRF вычисляется путем усреднения спектров оконных сегментов. Поскольку каждое возбуждение молотка затухает по существу до следующего возбуждения, можно использовать прямоугольное окно. Укажите датчик как смещение.
winLen = 2.5275*fs; % window length in samples modalfrf(X0,Y0,fs,winLen,'Sensor','dis')

FRF, оцененный по умолчанию 'H1' оценщик содержит три выдающихся пика в измеряемой полосе частот, соответствующих трем гибким модам вибрации. Когерентность близка к единице вблизи этих пиков и низкая в антирезонансных областях, где отношение сигнал/шум измерения отклика низкое. Когерентность вблизи одного указывает на высокую оценку качества. 'H1' оценка является оптимальной, когда шум существует только при выходном измерении, тогда как 'H2' оценщик оптимален, когда аддитивный шум имеется только на входе [2]. Вычислите 'H1' и 'H2' оценки для этого FRF.
[FRF1,f1] = modalfrf(X0,Y0,fs,winLen,'Sensor','dis'); % Calculate FRF (H1) [FRF2,f2] = modalfrf(X0,Y0,fs,winLen,'Sensor','dis','Estimator','H2');
При значительном шуме измерения или слабом возбуждении параметрические методы могут предложить дополнительные возможности для точного извлечения FRF из данных. Метод «подпространства» сначала соответствует модели «состояние-пространство» данным [3], а затем вычисляет свою функцию «частота-отклик». Порядок модели state-space (равный количеству полюсов) и наличие или отсутствие прохождения могут быть определены для конфигурирования оценки state-space.
[FRF3,f3] = modalfrf(X0,Y0,fs,winLen,'Sensor','dis','Estimator','subspace','Feedthrough',true);
Здесь FRF3 оценивается путем подгонки модели состояния-пространства, содержащей проходной член и оптимального порядка в диапазоне 1:10. Сравнение расчетных FRF с использованием 'H1', 'H2' и 'subspace' методы к теоретической FRF.
helperPlotModalAnalysisExample(f1,FRF1,f2,FRF2,f3,FRF3,f0,H0);

Оценщики выполняют сравнительно близкие пики ответа, в то время как 'H2' оценщик переоценивает отклик на антирезонансы. На согласованность не влияет выбор оценщика.
Затем оцените собственную частоту каждого режима с помощью алгоритма выбора пика. Алгоритм подбора пиков представляет собой простую и быструю процедуру идентификации пиков в FRF. Это локальный метод, поскольку каждая оценка генерируется из одной частотно-ответной функции. Это также метод одной степени свободы (SDOF), так как пик для каждого режима рассматривается независимо. В результате для каждого FRF генерируется набор модальных параметров. На основе предыдущего графика задайте диапазон частот от 200 до 1600 Гц, который содержит три пика.
fn = modalfit(FRF1,f1,fs,3,'FitMethod','PP','FreqRange',[200 1600])
fn =
1.0e+03 *
0.3727
0.8525
1.3707
Собственные частоты составляют приблизительно 373, 853 и 1371 Гц. Постройте график реконструированного FRF и сравните его с измеренными данными с помощью modalfit. FRF реконструируется с использованием модальных параметров, оцененных из матрицы частотно-ответной функции, FRF1. Звонить modalfit снова без выходных аргументов для создания графика, содержащего реконструированную FRF.
modalfit(FRF1,f1,fs,3,'FitMethod','PP','FreqRange',[200 1600])

Реконструированная FRF согласуется с измеренной FRF первого канала возбуждения и отклика. В следующем разделе рассматриваются два дополнительных места возбуждения.
Возбуждение Roving Hammer
Вычислите и постройте график FRF для ответов всех трех датчиков по умолчанию 'H1' оценщик. Укажите тип измерения как 'rovinginput' так как у нас есть возбуждение бродячего молота.
modalfrf(X1,Y1,fs,winLen,'Sensor','dis','Measurement','rovinginput')

В предыдущем разделе один набор модальных параметров был вычислен из одного FRF. Теперь оцените модальные параметры, используя алгоритм наименьших квадратов (LSCE). Алгоритмы LSCE и LSRF генерируют единый набор модальных параметров путем одновременного анализа множества ответных сигналов. Это глобальные методы множественной степени свободы (MDOF), так как параметры для всех режимов оцениваются одновременно из множества функций частотно-отклика.
Алгоритм LSCE генерирует вычислительные режимы, которые физически не присутствуют в структуре. Используйте схему стабилизации для определения физических режимов путем анализа устойчивости полюсов по мере увеличения числа режимов. Собственные частоты и коэффициенты демпфирования физических мод имеют тенденцию оставаться на том же месте или являются «стабильными». Создайте схему стабилизации и выведите собственные частоты тех полюсов, которые стабильны по частоте.
[FRF,f] = modalfrf(X1,Y1,fs,winLen,'Sensor','dis','Measurement','rovinginput'); fn = modalsd(FRF,f,fs,'MaxModes',20, 'FitMethod', 'lsce'); % Identify physical modes

По умолчанию полюса классифицируются как стабильные по частоте, если собственная частота полюсов изменяется менее чем на один процент. Полюса, которые стабильны по частоте, далее классифицируются как устойчивые по демпфированию для изменения коэффициента демпфирования менее чем на пять процентов. Оба критерия могут быть скорректированы на различные значения. Исходя из расположения устойчивых полюсов, выбирайте собственные частоты 373, 852,5 и 1371 Гц. Эти частоты содержатся в выходных данных fn из modalsd, наряду с собственными частотами других стабильных по частоте полюсов. Для получения хороших оценок модальных параметров с использованием алгоритма LSCE обычно требуется более высокий порядок моделей, чем количество физически присутствующих режимов. В этом случае порядок моделей из четырех режимов указывает на три устойчивых полюса. Интересующие частоты встречаются в первых трех столбцах в 4-м ряду fn.
physFreq = fn(4,[1 2 3]);
Оценить собственные частоты и демпфирование и построить график восстановленных и измеренных FRF. Укажите четыре режима и физические частоты, определенные на диаграмме устойчивости. 'PhysFreq'. modalfit возвращает модальные параметры только для указанных режимов.
modalfit(FRF,f,fs,4,'PhysFreq',physFreq)

[fn1,dr1] = modalfit(FRF,f,fs,4,'PhysFreq',physFreq)
fn1 =
1.0e+03 *
0.3727
0.8525
1.3706
dr1 =
0.0008
0.0018
0.0029
Затем вычислите FRF и постройте график стабилизации для второго набора ударов молотка с датчиком в другом месте. Измените критерий устойчивости на 0,1 процента для частоты и 2,5 процента для демпфирования.
[FRF,f] = modalfrf(X2,Y2,fs,winLen,'Sensor','dis','Measurement','rovinginput'); fn = modalsd(FRF,f,fs,'MaxModes',20,'SCriteria',[0.001 0.025]);

При более строгих критериях большинство полюсов классифицируются как не стабильные по частоте. Полюса, которые стабильны по частоте и демпфированию, близко совпадают со усредненными FRF, что позволяет предположить, что они присутствуют в измеренных данных.
physFreq = fn(4,[1 2 3]);
Извлекают модальные параметры для этого набора измерений и сравнивают с модальными параметрами для первого набора измерений. Укажите индексы ведущей точки FRF, соответствующие месту совпадения измерений возбуждения и отклика. Собственные частоты совпадают в пределах доли процента, а коэффициенты демпфирования - в пределах менее четырех процентов, что указывает на то, что модальные параметры согласуются от измерения к измерению.
[fn2,dr2] = modalfit(FRF,f,fs,4,'PhysFreq',physFreq,'DriveIndex',[1 2])
fn2 =
1.0e+03 *
0.3727
0.8525
1.3705
dr2 =
0.0008
0.0018
0.0029
Tdiff2 = table((fn1-fn2)./fn1,(dr1-dr2)./dr1,'VariableNames',{'diffFrequency','diffDamping'})
Tdiff2 =
3x2 table
diffFrequency diffDamping
_____________ ___________
2.9972e-06 -0.031648
-5.9335e-06 -0.0099076
1.965e-05 0.0001186
Параметрический способ для оценки параметров режима может обеспечить полезную альтернативу подбору пиков и способу LSCE, когда имеется шум измерения в FRF или FRF показывает высокую плотность режима. Подход рациональной функции наименьших квадратов (LSRF) соответствует общей передаточной функции знаменателя множественному входу, множественному выходу FRF и, таким образом, получает единую глобальную оценку модальных параметров [4]. Процедура использования подхода LSRF аналогична процедуре для LSCE. Схему стабилизации можно использовать для идентификации стационарных режимов и извлечения модальных параметров, соответствующих выявленным физическим частотам.
[FRF,f] = modalfrf(X1,Y1,fs,winLen,'Sensor','dis','Measurement','rovinginput'); fn = modalsd(FRF,f,fs,'MaxModes',20, 'FitMethod', 'lsrf'); % Identify physical modes using lsfr physFreq = fn(4,[1 2 3]); [fn3,dr3] = modalfit(FRF,f,fs,4,'PhysFreq',physFreq,'DriveIndex',[1 2],'FitMethod','lsrf')
fn3 =
372.6832
372.9275
852.4986
dr3 =
0.0008
0.0003
0.0018

Tdiff3 = table((fn1-fn3)./fn1,(dr1-dr3)./dr1,'VariableNames',{'diffFrequency','diffDamping'})
Tdiff3 =
3x2 table
diffFrequency diffDamping
_____________ ___________
-7.8599e-06 0.007982
0.56255 0.83086
0.37799 0.37626
Заключительное примечание о параметрических методах: метод оценки FRF («подпространство») и метод оценки модальных параметров ('lsrf') аналогичны тем, которые используются в панели инструментов идентификации системы для подгонки динамических моделей к сигналам временной области или к функциям частотного отклика. Если эта панель инструментов доступна, можно определить модели, соответствующие данным, с помощью таких команд, как tfest и ssest. Оценить качество идентифицированных моделей можно с помощью compare и resid команды. После проверки качества модели их можно использовать для извлечения модальных параметров. Это кратко показано с использованием оценщика состояния-пространства.
Ts = 1/fs; % sample time % Create a data object to be used for model estimation. EstimationData = iddata(Y0(1:1000), X0(1:1000), 1/fs); % Create a data object to be used for model validation ValidationData = iddata(Y0(1001:2000), X0(1001:2000), 1/fs);
Определите непрерывную модель состояния-пространства 6-го порядка, содержащую член прохождения.
sys = ssest(EstimationData, 6, 'Feedthrough', true)
sys =
Continuous-time identified state-space model:
dx/dt = A x(t) + B u(t) + K e(t)
y(t) = C x(t) + D u(t) + e(t)
A =
x1 x2 x3 x4 x5 x6
x1 4.05 -1765 149.8 -1880 -49.64 -358
x2 1764 -0.3332 2197 -232.5 -438.3 -128.4
x3 -152.4 -2198 2.85 4715 255.9 547.5
x4 1879 228.2 -4713 -15.9 -1216 -28.79
x5 59.42 440.9 -275.5 1217 35.05 -8508
x6 363.7 120.2 -545.4 -44.02 8508 -92.45
B =
u1
x1 -0.1513
x2 -1.911
x3 4.439
x4 -3.118
x5 -0.9416
x6 -8.039
C =
x1 x2 x3 x4 x5 x6
y1 3.135e-05 2.511e-06 8.634e-06 -1.416e-05 2.218e-06 -6.271e-06
D =
u1
y1 7.564e-09
K =
y1
x1 3.513e+07
x2 -3.244e+06
x3 -3.598e+07
x4 -1.059e+07
x5 1.724e+08
x6 7.521e+06
Parameterization:
FREE form (all coefficients in A, B, C free).
Feedthrough: yes
Disturbance component: estimate
Number of free coefficients: 55
Use "idssdata", "getpvec", "getcov" for parameters and their uncertainties.
Status:
Estimated using SSEST on time domain data "EstimationData".
Fit to estimation data: 99.3% (prediction focus)
FPE: 1.235e-16, MSE: 1.189e-16
Оцените качество модели, проверив ее соответствие данным проверки.
clf
compare(ValidationData, sys) % plot shows good fit

Использовать модель sys для вычисления модальных параметров.
[fn4, dr4] = modalfit(sys, f, 3);
Понимание динамического поведения лопаток ветровых турбин важно для оптимизации эффективности работы и прогнозирования отказа лопаток. В этом разделе анализируются экспериментальные данные модального анализа для лопатки ветровой турбины и визуализируются формы режимов лопатки. Ударник возбуждает лопатку турбины в 20 местах, а эталонный акселерометр измеряет отклики в 18 месте. В основании лопасти установлен алюминиевый блок, и лопасть возбуждается в наклонной ориентации, перпендикулярной плоской части лопасти. FRF собирается для каждого местоположения. Данные FRF любезно предоставлены Лабораторией систем структурной динамики и акустики Массачусетского университета в Лоуэлле. Во-первых, визуализировать пространственное расположение мест измерения.

Загрузить и составить график оценок FRF лопаток ветротурбины для местоположений 18 и 20. Увеличьте изображение первых нескольких пиков.
[FRF,f,fs] = helperImportModalData(); helperPlotModalAnalysisExample(FRF,f,[18 20]);

Первые два режима выглядят как пики около 37 Гц и 111 Гц. Постройте график стабилизации для оценки собственных частот. Первые два значения, возвращаемые для порядка модели 14, стабильны по частоте и коэффициенту демпфирования.
fn = modalsd(FRF,f,fs,'MaxModes',20);
physFreq = fn(14,[1 2]);

Затем извлеките формы режима для первых двух режимов с помощью modalfit. Ограничьте посадку диапазоном частот от 0 до 250 Гц на основе предыдущего графика.
[~,~,ms] = modalfit(FRF,f,fs,14,'PhysFreq',physFreq,'FreqRange',[0 250]);
Формы моды количественно определяют амплитуду движения для каждой моды структуры в каждом местоположении. Для оценки вектора формы моды необходима одна строка или столбец матрицы частотно-ответной функции. На практике это означает, что либо необходимо возбуждение в каждом месте измерения конструкции (в данном случае - бродячий молоток), либо необходимо измерение отклика в каждом месте. Формы режимов могут быть визуализированы путем изучения мнимой части FRF. Постройте график водопада воображаемой части матрицы FRF для расположений на одной стороне лопатки. Ограничьте диапазон частот максимум 150 Гц для проверки первых двух режимов. Пики графика представляют формы режимов.
measlocs = [3 6 9 11 13 15 17 19 20]; % Measurement locations on blade edge
helperPlotModalAnalysisExample(FRF,f,measlocs,150);

Формы, обозначенные на графике контуром вершин, представляют первый и второй изгибающие моменты лезвия. Затем постройте график величины векторов формы моды для тех же мест измерения.
helperPlotModalAnalysisExample(ms,measlocs);

В то время как амплитуды масштабируются по-разному (векторы формы моды масштабируются до единицы моды А), контуры формы моды сходятся по форме. Форма первого режима имеет большое смещение наконечника и два узла, где амплитуда колебаний равна нулю. Второй режим также имеет большое смещение наконечника и имеет три узла.
В этом примере анализировались и сравнивались моделируемые наборы данных модального анализа для 3DOF системы, возбуждаемой бродящим молотком. Оценивали собственную частоту и демпфирование с помощью схемы стабилизации и алгоритмов LSCE и LSRF. Модальные параметры были согласованы для двух наборов измерений. В отдельном случае использования модовые формы лопатки ветровой турбины визуализировали с использованием мнимой части матрицы FRF и векторов модовой формы.
Спасибо доктору Питеру Авитабиле из Лаборатории систем структурной динамики и акустики Массачусетского университета в Лоуэлле за содействие сбору экспериментальных данных по лопаткам ветровых турбин.
[1] Брандт, Андерс. Анализ шума и вибрации: анализ сигналов и экспериментальные процедуры. Чичестер, Великобритания: Джон Уайли и сыновья, 2011.
[2] Вольд, Ховард, Джон Кроули и Г. Томас Роклин. «Новые способы оценки функций частотного отклика». Звук и вибрация. Том 18, ноябрь 1984, стр. 34-38.
[3] Питер Ван Куражи и Барт Де Моор. «N4SID: Подпространственные алгоритмы идентификации комбинированных детерминированно-стохастических систем». Автоматически. том 30, январь 1994 года, стр. 75-93.
[4] Оздемир, А. А., С. Гумуссой. «Оценка передаточной функции в панели инструментов идентификации системы с помощью векторного фитинга». Материалы 20-го Всемирного конгресса Международной федерации автоматического управления. Тулуза, Франция, июль 2017 года.