Анализ пространственно спектра колебаний

[Schwungrad]

Требуется провести качественный анализ пространственных мод колебаний корпуса электродвигателя. «Качественный анализ» обозначает: определить какие пространственные моды вносят наибольший вклад в колебания при критической скорости вращения двигателя, которая известна.
Пространственные моды окружности имеют вид (изображение взято из Интернета):



Измеряются x/y/z компоненты ускорений 11 датчиками с круговым шагом 22.5° (см. картинку) в сечении:



Проблема (вопрос) состит в том, что ускорения на части окружности незвестны, т.к. мотор пристегнут ремнем к испытательному стенду и, поэтому, на части окружности было невозможно закрепить датчики. Графически это можно пояснить так: известна синяя линяя на картинке. Зная только синюю линию требуется произвести анализ пространственно спектра. Понятно, что данные неполные и придется иметь дело с так называемым aliasing. Но других данных на данный момент нет и непонятно как корректно сделать оценку пространственного спектра в Матлабе.




В данном случае измерение проводились по времени в 11 точках окружности. Спектр временных сигналов анализируется в программе Testlab и экспортируется в Matlab для последующего анализа пространственного спектра.


Результат спектрального анализа по времени (waterfall диаграмма) в одной из точек:



Во всех точках выделяется наиболее критическая 48-я временная гармоника (относительно оборота вала эл. двигателя) и экспортируется в Матлаб для анализа пространственных мод (только по окружности). И вот тут у меня возникает вопрос: как этот анализ проводить ибо данные по окружности неполные?

OBaMA (Order Based Modal Analysis) проведенный в Testlab не дает ответ на вопрос, какая пространственная мода дает наибольший вклад:



Тут, скорее, понятно, что ничего не понятно.

[Schwungrad]

Прошу высказать соображения о корректности довольно примитивного метода оценки вклада пространственных мод в данном случае. Это первое, что приходит в голову.

1. Предположим, что возбуждается максимум первые 6 пространственных мод. С учетом того, что измерение проводилось всего в 11 точках окружности, а также, что практически более высокие пространственные моды слабо возбуждаются из физических соображений, это разумное предположение.

2. Решение (разложение в ряд Фурье) будем рассматривать относительно угловой пространственной координаты и длины радиус-вектора перемещений в плоскости окружности.

3. С учетом 1, решение на окружности будет иметь вид (здесь x - угловая координата, R - длина радиус-вектора из центра окружности):

$$\[R(x)=\sum_{k=1}^{N_m} a_k\cos((k-1)x+\varphi_k), \quad N_m=6,\, \varphi_1=0\]$$

4. При фиксированной частоте колебаний входными данными (результатами измерений) являются значения радиус-вектора и ассоциированной угловой координаты в 11 точках окружности $x_p, \, R_p, \quad p=1\ldots 11$. Это легко вычисляется по результату гармонического анализа (преобразование Фурье по времени), проведенного в Testlab, где выделяется доминирующая 48-я гармоника.

5. Подставив входные данные амплитуд и фаз в 11 точках $x_p, \, R_p, \quad p=1\ldots 11$ в 3., получим нелинейную алгебраическую систему 11 уравнений относительно 11 неизвестных. Решив эту систему в Matlab, найдем неизвестные амплитуды и фазы пространственных мод: $a_k, \, \varphi_k$.

6. Повторив решение для каждой частоты, получим полный пространственно-временной спектр.

[wrest]

Schwungrad в сообщении #1527235 писал(а):

Тут, скорее, понятно, что ничего не понятно.

Глазом хорошо видна 3-я мода: пучности на 0, 120 и 240 градусов. Т.е. решение должно показать, что она есть и вклад большой.

[realeugene]

Schwungrad в сообщении #1527268 писал(а):

5. Подставив входные данные амплитуд и фаз в 11 точках $x_p, \, R_p, \quad p=1\ldots 11$ в 3., получим нелинейную алгебраическую систему 11 уравнений относительно 11 неизвестных. Решив эту систему в Matlab, найдем неизвестные амплитуды и фазы пространственных мод: $a_k, \, \varphi_k$.

Нет, не так. Каждая нарисованная вами мода кроме нулевой на самом деле состоит из двух ортогональных мод, повёрнутых друг относительно друга на угол $\varphi_m=\frac{\pi}{2 m}$. У каждой из этих мод своя собственная комплексная аимлитуда колебаний, т. е. размерность задачи всего 22 действительных числа. Выполнив преобразование Фурье во временной области для 11 датчиков вы тоже получите 11 комплексных чисел, но вам ни в коем случае нельзя терять информацию про фазу после преобразования Фурье во временной области, иначе задача окажется недоопределённая. И вам лучше всего простраственные моды описывать именно комплексными амплитудами ввиду линейности получаемой аппроксимирующей функции над полем комплексных чисел. Фазы во времени вам могут быть не интересны, но вот фазовые сдвиги между парами мод с одним номером существенны и могут приводить к вращающимся по кругу поляризациям колебаний, аналогично линейным и круговым поляризациям электромагнитной волны. Вы же сейчас рассматриваете только линейные поляризации колебаний.

[Schwungrad]

realeugene в сообщении #1527307 писал(а):

Schwungrad в сообщении #1527268 писал(а):

5. Подставив входные данные амплитуд и фаз в 11 точках $x_p, \, R_p, \quad p=1\ldots 11$ в 3., получим нелинейную алгебраическую систему 11 уравнений относительно 11 неизвестных. Решив эту систему в Matlab, найдем неизвестные амплитуды и фазы пространственных мод: $a_k, \, \varphi_k$.

Нет, не так. Каждая нарисованная вами мода кроме нулевой на самом деле состоит из двух ортогональных мод, повёрнутых друг относительно друга на угол $\varphi_m=\frac{\pi}{2 m}$. У каждой из этих мод своя собственная комплексная аимлитуда колебаний, т. е. размерность задачи всего 22 действительных числа.

Спасибо за ответ. Дайте пожалуйста ссылку где можно почитать о суперпозиции двух ортогональных мод, повёрнутых друг относительно друга на угол $\varphi_m=\frac{\pi}{2 m}$? Желательно какой-нибудь учебник попроще, ведь тема довольно классическая.

Перед тем как написать свое сообщение, я потрудился нарисовать нормальные моды в Matlab. Я исходил из того, что каждая мода определяется волновым числом, а ее вклад - двумя действительными числами (одним комплексным): амплитудой и фазой.

По этим 3-м числам нарисовал каждую из пространственных мод. Вот, например, если k=5

Код Матлаба:


Результат:


Вроде строится и без ортогональных мод. Хотя, видимо, я тут чего-то недопонимаю. Например, мне непонятно про

Цитата:

Фазы во времени вам могут быть не интересны

Ведь времени больше нет после FFT преобразования по времени. При выделенной временной гармонике (48-я дает наибольший вклад) и фиксированной частоте в каждой точке измерения после временного FFT остается амплитуда и фаза, которые необходимо корректно учесть для спектрального анализа в пространственной области.

[realeugene]

Schwungrad в сообщении #1527356 писал(а):

Дайте пожалуйста ссылку где можно почитать о суперпозиции двух ортогональных мод, повёрнутых друг относительно друга на угол $\varphi_m=\frac{\pi}{2 m}$?

Не дам. Не помню конкретные учебники. Но именно таков гармонический базис на окружности. Поверьте. Причём механикам, наверное, это должно быть непривычно и сложно понять, что функции, которые по этому базису вами раскладываются, изначально существенно комплексные. Так как исходные точки у вас уже Фурье-компоненты на какой-то определённой частоте действительного сигнала. В электромагнетизме и обработке сигналов это давно привычнее.

Schwungrad в сообщении #1527356 писал(а):

Вроде строится и без ортогональных мод.

Вот смотрите. У той моды, которую вы построили на графике, есть нули, т. е. места, где голубая линия пересекает красную окружность. В этих точках амплитуда колебаний нулевая. Вторая ортогональная ей мода с $m=5$ - это колебания с максимумами в тех местах, гду у вашей первой моды нули При фазовом сдвиге между этими модами $\pi/2$ нарисованная вами звёздочка будет вращаться с постоянной скоростью по окружности.

Удобнее при этом использовать комплексный гармонический на окружности базис $e^{i m \varphi}$, где у вас $\left|m\right|\le5$. Обратите внимание, что $m$ - целые, в том числе, отрицательные числа.

-- 27.07.2021, 22:26 --

Schwungrad в сообщении #1527356 писал(а):

При выделенной временной гармонике (48-я дает наибольший вклад) и фиксированной частоте в каждой точке измерения после временного FFT остается амплитуда и фаза, которые необходимо корректно учесть для спектрального анализа в пространственной области.

Верно, на фазы после ПФ зависят от начального момента отсчёта времени. Который выбирается произвольно, и от выбора которого никакие выводы зависеть не могут.

Преобразование Фурье - ПФ. А FFT - "Fast Fourier Transform", алгоритм быстрого вычисления преобразования Фурье. Точнее, целый класс таких алгоритмов.

[ozheredov]

realeugene в сообщении #1527367 писал(а):

Но именно таков гармонический базис на окружности.

Верно ли, что под разлагаемыми по базису функциями подразумеваются те, кто: 1) $\displaystyle 2\pi$-периодичные и 2) $f(x) \geqslant -r$ (чтобы радиус-вектор отрицательным не стал)?

[realeugene]

ozheredov в сообщении #1527377 писал(а):

Верно ли, что под разлагаемыми по базису функциями подразумеваются те, кто: 1) $\displaystyle 2\pi$-периодичные и 2) $f(x) \geqslant -r$ (чтобы радиус-вектор отрицательным не стал)?

1 - да, 2 - нет. На окружности вообще нет радиус-вектора. Разлагать такие функции можно в классический ряд Фурье по синусам и косинусам, но удобнее по базису из комплексных экспонент, особенно, если функция комплекснозначная.

Недостаток базиса ТС легко видеть для $m=1$. При помощи такой базисной функции можно представить биения вдоль некоторой оси. Например, по горизонтали. Или по вертикали. Или по диагонали. Но невозможно одной такой базисной функцией представить вращение с эксцентриситетом. Для которого тоже $m=1$.

[Schwungrad]

realeugene в сообщении #1527367 писал(а):

Schwungrad в сообщении #1527356 писал(а):

Дайте пожалуйста ссылку где можно почитать о суперпозиции двух ортогональных мод, повёрнутых друг относительно друга на угол $\varphi_m=\frac{\pi}{2 m}$?

Не дам. Не помню конкретные учебники. Но именно таков гармонический базис на окружности.

Совершенно верно и ссылки, конечно же, не нужно, т.к. это набор "синусов" и "косинусов" без учета фаз и один из другого получается поворотом на $\varphi_k=\frac{\pi}{2 k}$. Но формула разложения по ортогональным модам, очевидно, эквивалентна приведенной мной по косинусам, но с фазой. Сумма ортогональных мод представима

$$A_k \cos(kx) + B_k\sin(kx)= a_k\cos(kx + \varphi_k),\quad a_k=\sqrt{A_k^2 + B_k^2},\, \varphi_k=-\arccos{\frac{A_k}{a_k}}$$

Ортогональный базис, да, именно такой.





Вы советуете: искать решение, считая $A_k,\, B_k$ комплексными числами, т.е. разлагать по ортогональному базису контур, определяемый 11 комплексными числами (результат FFT по времени в 11 точках).

В коммерческих профильных пакетах, например Simcenter Testlab, для спектрального анализа используется только FFT заточенное на степень двойки. И, насколько мне известно, только один алгоритм - Butterfly algorithm.

[realeugene]

Schwungrad в сообщении #1527380 писал(а):

Сумма ортогональных мод представима

$$A_k \cos(kx) + B_k\sin(kx)= a_k\cos(kx + \varphi_k),\quad a_k=\sqrt{A_k^2 + B_k^2},\, \varphi_k=-\arccos{\frac{A_k}{a_k}}$$

$A_k$ и $B_k$ - комплексные числа, так что, нет. Разность их фаз имеет смысл фазового сдвига между косинусной и синусной модами колебаний, и без него вам не обойтись в любом случае.

-- 28.07.2021, 02:59 --

Schwungrad в сообщении #1527380 писал(а):

В коммерческих профильных пакетах, например Simcenter Testlab, для спектрального анализа используется только FFT заточенное на степень двойки. И, насколько мне известно, только один алгоритм - Butterfly algorithm.

Оценка спектра - это отдельная наука. Не забудьте про оконные функции. Окно должно быть одним на все каналы. Надеюсь, у вас сбор сигналов всех датчиков синхронный.

Матлаб умеет делать быстрое преобразование Фурье с длиной не равной степени двойки. Но не рекомендую, когда-то вылазили какие-то глюки. Алгоритм Кули-Тьюки привычнее и лучше отлажен скорее всего.

[Schwungrad]

Благодарю realeugene за то, что помогли мне разобраться с волновым решением.

Пространственный спектр в данном случае будет конечно же комплексный, т.к. преобразование Фурье применяется к вектору комплексных чисел, состоящим из амплитуд и фаз (действительной и мнимой части), полученных в результате применения FFT по времени к исходному сигналу (ускорение, измеренное в 11 точках окружности). Поэтому, пространственный спектр содержит положительные и отрицательныые волновые числа, а искомое разложение имеет вид:

$$\sum_{k}\left( a_k^{1}e^{i(n\omega t + kx +\varphi_k^{1})}+a_k^{2}e^{i(n\omega t - kx +\varphi_k^{2})} \right)$$

то есть решение состоит из пар волн с амплитудами и фазами $a_k^{1},\, \varphi_k^{1}$ и $a_k^{2},\, \varphi_k^{2}$, распространяющихся по и против часовой стрелке.

Оставив действительную часть, получим радиальные перемещения в виде:

$$R(t,x) = R_0 + \sum_{k}\left[ a_k^{1}\cos{(n\omega t + kx +\varphi_k^{1})}+a_k^{2}\cos{(n\omega t - kx +\varphi_k^{2})} \right]$$

Таким образом, исходное разложение (если учитывать только первые 6 пространственных мод) содердит 24 неизвестных (амплитуды и фазы распространяющихся направо и налево волн). Или с учетом того что при $k=1$ фазы можно считать равными 0, остается 22 неизвестных.

Распространение пары волн легко изображается в Matlab:



В результате получаются две волны:



Их можно представить в виде характеристик 1 и 2:

[realeugene]

Schwungrad в сообщении #1527553 писал(а):

$$\sum_{k}\left( a_k^{1}e^{i(n\omega t + kx +\varphi_k^{1})}+a_k^{2}e^{i(n\omega t - kx +\varphi_k^{2})} \right)$$

А что такое $x$?

Нет, по окружности у вас одна циклическая координата - сам угол. Соответственно, спектр по углу дискретный.

И попробуйте думать сразу в терминологии комплексных амплитуд. Решать систему линейных уравнений проще, чем заниматься нелинейной оптимизацией. После нахождения комплексных амплдитуд перейдёте к энергиям, просуммировав квадраты модулей амплитуд гармоник для одинаковых $m$. Эти энергии вам и интересны как результат анализа.

Впрочем, анализ амплитуд и фаз для одного $m$ даст дополнительную информацию о том, как именно эти гармоники вибрируют: вдоль некоторой оси или вращаясь. То есть для каждого $\left| m\right| > 0$ можно посчитать отношение энергий колебаний по главным осям и направление максимума колебаний. Ещё можно двигать окно анализа по времени и анализировать дрейф этих энергий и углов со сременем.

[Schwungrad]

$x$ это угловая координата. Спектральный анализ (и сами измерения ускорений) проводится в трех сечениях раздельно. Здесь Сечение 1: mp_01, ... mp_11; Сечение 2: mp_12, ... mp_22 и Сечение 3: mp_23 ... mp_33. Сначала измерение проводилось в Сечении 1, потом датчики переклеивались и измерялось Сечение 2, и затем - Сечение 3.

Точка mp_34 (сензор физически) присутствует во всех измерениях и используется в качестве phase reference для "склейки" измерений во всех трех сечениях. Это необходимо для определения собственных форм колебаний (модального анализа).



$x$ это то, что в Testlab обозначает Theta:

[Schwungrad]

Используя разложение, рекомендованное realeugene, получены первые результаты спектрального анализа.

В "терминологии комплексных амплитуд" пространственное представление по $x$ имеет вид:



Подставив в это разложение точки, $x_p$, в которых проводилось измерение, получим $N_p=11$ уравнений:



Записав уравнения относительно действительных и мнимых частей неизвестных:



Окончательно получим систему 22-х линейных уравнений относительно неизвестных $y$:



Здесь в правой части - результат преобразования Фурье по времени (при фиксированной частоте), то есть правая часть вычисляется и импортируется из Testlab.

Реализовав решение этой системы уравнений в Matlab (пока рассмотрена только нормальная составляющая ускорений) получились вот такие спектры.

40-я временная гармоника в Сечении 1 (середина статора):

Оси: частота в [Hz], Wave number, Amplitude в [g]. Датчики сертифицированы только до 10 kHz, поэтому частота отображается по 10500 Hz.



Оранжевый спектр соответствует малошумному мотору, который был отобран как "golden sample" по результатам тестирования вибрации непосредственно на сборочной линии, а синий спектр - соответствует шумному мотору, который рекламировали и прислали назад с просьбой установить причину шумности. Причина пока непонятна, но четко видно что амплитуды вибрации критических гармоник (40-й и 48-й) у рекламированного мотора заметно выше.

40-я временная гармоника в Сечении 2 (Non drive end):



40-я временная гармоника в Сечении 3 (Drive end):



48-я временная гармоника в Сечении 1 (середина статора):



48-я временная гармоника в Сечении 2 (Non drive end):



48-я временная гармоника в Сечении 3 (Drive end):

[realeugene]

Schwungrad в сообщении #1528092 писал(а):

Записав уравнения относительно действительных и мнимых частей неизвестных

Зачем? Матлаб умеет решать сразу комплексные системы линейных уравнений.