Tarea 4
Dinámica Estructural Experimental 2021-1
Integrantes: Esteban Balboa, Marilyn Ramos, Marcelo Sanhueza
Profesor: Claudio Oyarzo, PhD
Fecha: 04 de Junio 2021
Notación:
E : Modulo elástico( N /m2)
f n : Frecuencia natural , en Hz
G xx ( ω ) : Densidad auto−espectral de la respuesta
I y : Inercia , eje debil
K c :Rigidez de la columna
k v : Rigidez local del sistema
k i : Rigidez local de columna i
K u : Matriz de rigidez del sistema
L : Largode la columna
L : Matriz cinemática
M a : Masaacelerometros
M c : Masa columna
mi : Masa pisoi
M placa : Masa placa
M u : Masa del Sistema
S x ( ω ) : FFT de la respuesta
¿
S x (ω) :Conjugado de la FFT de la respuesta
T n : Periodo natural , en s
T d : Periodo Amortiguado , en s
W n =Frecuencia natural , en rad /s
2
1. Descripción del Modelo Físico.
La estructura ensayada es la que se muestra esquemáticamente en la Figura 1, compuesta
por 4 columnas metálicas rectangulares (divididas en tres partes iguales por las conexiones
columna-placa), tres placas cuadradas de madera MDF (pisos), y tres masas adicionadas
(golillas de acero) en la parte superior de cada placa. Las masas adicionadas se ubicaron en
el centro de masa y área (en este caso coinciden) de la placa de madera. Por simplicidad, las
conexiones se consideraron empotradas y la base de la estructura no se toma en cuenta ya
que el modelo se analiza desde la conexión que se encuentra sobre ella.
Figura 1: Modelo físico de estructura ensayada.
En la Tabla 1 se presentan las dimensiones y masas de los elementos de la estructura.
Tabla 1:
Dimensiones y masa de elementos.
Columnas Placa 1 Placa 2 Placa 3
Ancho (cm) 2.9 12.5 12.5 12.5
Largo (cm) 20 12.5 12.5 12.5
Espesor (cm) 0.9 0.9 0.9 0.9
Masa (g) 40.67 350 349 350
3
Nota: La masa de la placa considera la masa adicional y conexiones.
Respecto de la distribución de masas del modelo, las placas y las masas adicionadas se
consideran como masas puntuales en cada grado de libertad (gdl), tributando 100% al gdl
correspondiente y las columnas tributando la mitad de su masa al gdl respectivo, según se
observa en la Figura 2.
Figura 2: Distribución de masas de columnas (izquierda), placas y masa adicionada
(derecha).
Los elementos que proporcionan rigidez lateral a la estructura son las columnas y la placa,
pero para los cálculos, se considera la placa como infinitamente rígida.
Para el cálculo teórico del periodo natural se considera el modelo matemático de la Figura 3
con las siguientes simplificaciones de modelo físico:
Largo de columnas entre placas
Apoyo perfectamente empotrado
Columnas axialmente indeformables
Placa infinitamente rígida.
4
Figura 3: Modelo matemático
La Tabla 2 muestra las propiedades de la columna para calcular la rigidez del modelo.
Tabla 2:
Propiedades de las columnas
Propiedades columnas
−12
Inercia (m4) 1 ,76 x 10
11
Modulo elástico (N/m2) 2 , 1 x 10
2. Modelo teórico
Teóricamente, para encontrar las frecuencias y modos de vibrar que presenta la estructura,
se deben calcular las matrices tanto de masa y rigidez. La matriz de masa es la que se
presenta a continuación:
[ ]
m1 0 0
M u= 0 m 2 0
0 0 m3
5
Donde m 1, m 2 y m 3 se obtuvieron mediante las siguientes ecuaciones (2.1-2.3):
MC
m1= M Placa + ∗2∗4 + M a 2.1
3∗2
MC
m2= M Placa + ∗2∗4 + M a 2.2
3∗2
MC
m 3= M Placa + ∗4+ M a 2.3
3∗2
Luego, la matriz de rigidez se obtiene a partir de la ecuación 2.4:
T
K u=L ∗k v∗L 2.4
Donde k v se obtiene mediante la matriz que se presenta a continuación.
[ ]
k1 0 0 0 0 0
0 k2 0 0 0 0
0 0 k3 0 0 0
k v= ∗2
0 0 0 k4 0 0
0 0 0 0 k5 0
0 0 0 0 0 k6
En la matriz k v , k i son las matrices de rigidez de cada columna, son todas iguales y se
pueden obtener mediante la siguiente ecuación.
k 1=
12∗E∗I
h
3 [
∗
1 −1
−1 1
2.5
]
Luego la matriz cinemática L, se construye a partir de la Figura 4.
6
[]
Figura 4: Representación de gdls locales de columnas.
0 0 1
0 1 0
0 0 1
0 1 0
0 1 0
1 0 0
L=
0 1 0
1 0 0
1 0 0
0 0 0
1 0 0
0 0 0
Aplicando la ecuación 2.4, se obtiene la matriz de rigidez del sistema K u .
[ ]
1.0324 −0.5162 0
K u= −0.2178 0.4355 −0.2178 *104
0 −0.2178 0.2178
Las frecuencias naturales y formas modales se obtuvieron a través del software Matlab,
utilizando el comando “eig ( K u, M u)”, el cual arrojó las frecuencias naturales presentadas
en la Tabla 3.
Tabla: 3
Frecuencias naturales del modelo teórico
Modo 1 Modo 2 Modo 3
Frecuencia
natural (Hz) 4,827 13,343 18,893
7
A continuación, se presenta la matriz normalizada de formas modales entregadas por
Matlab.
[ ]
−0.4883 1.0492 −0.7819
∅ = −0.8708 0.3624 1.0300
−1.0648 −0.9241 −0.5751
En la Figura 5 se observa una representación esquemática de las formas modales de cada
modo.
Figura 5: Representación gráfica de formas modales.
Si se normaliza la matriz ∅ al gdl 1, se obtiene lo siguiente:
[ ]
1 1 1
∅ = 1,801 0,442 −1.255
2,246 −0,801 0,560
Anteriormente, se realizó un ensayo experimental de la estructura, en el que se
determinaron las frecuencias modales del sistema con la aplicación de un golpe de un
matillo. Los resultados son los de la Tabla 4.
Tabla 4:
Frecuencias naturales determinada en ensayo del martillo
Modo 1 Modo 2 Modo 3
Frecuencia
natural (Hz) 4,11 13,91 20,93
3. Descripción del código Matlab
De acuerdo con Brandt (2010), existen diversas excitaciones a las cuales se puede ver
expuesta la estructura. Algunas de estas excitaciones, es fácil medirlas y obtener su espectro
de frecuencia, como es el caso de las excitaciones inducidas por un martillo o un vibrador.
Por otro lado, hay excitaciones de las cuales no es posible tener un registro, un ejemplo
serían las vibraciones ambientales. Una de las desventajas de las vibraciones ambientales es
que se desconoce su naturaleza, y por lo tanto son señales aleatorias, provocando que el
efecto de leakage se presente en los espectros de frecuencias.
Para disminuir el efecto leakage se divide la medición en bloques con el fin de obtener
promedios representativos de la muestra. A cada bloque se le aplica una ventana para
8
asegurar la periodicidad de la señal y así poder calcular las FFT. Sin embargo, al aplicar la
ventana, se pierden datos ya que la ventana transforma la señal en una que empieza y
finaliza en cero. Para solucionar lo anterior y no perder información importante, los bloques
se traslapan un 50%. Finalmente, para cada bloque se calcula la densidad auto-espectral de
la respuesta y se promedian con el objetivo de obtener el espectro completo en el dominio
de la frecuencia.
La densidad auto-espectral según Avitabile (2001) se calcula como la multiplicación de la
FFT de la respuesta por su conjugado, ecuación 4.1.
¿
G xx ( ω )=S x ( ω ) S x ( ω ) ( 4.1 )
Se implementó en Matlab un código que realizara el procedimiento descrito anteriormente,
conocido como el Método de Welch para estimar el espectro de la señal de respuesta
medido por los acelerómetros 1, 2 y 3. Pero antes, se debe definir el tamaño de los bloques
con los que se va a trabajar. Para esto, se realizó un análisis de sensibilidad del parámetro
N: “número de datos en cada bloque”, con el objetivo de identificar su efecto en la gráfica
del espectro y en los valores de las frecuencias naturales. Se graficó (Figura 6) sólo para el
acelerómetro 1 para el rango de frecuenticas de interés.
Figura 6: Representación gráfica de análisis de sensibilidad de N.
En la Figura 6, se observa que, a mayor valor de N, la curva del espectro es más limpia en
comparación al menor valor de N, es decir, el espectro se define mejor (peaks más claros) si
9
se utilizan más datos en cada bloque, pero se obtienen menos promedios representativos, lo
cual implica que no se eliminan los posibles errores del registro.
El valor de las frecuencias naturales cambia para cada valor de N. Los resultados se
encuentran en la Tabla 5.
Tabla 5:
Frecuencias naturales para los primeros 3 modos al variar N
Acelerómetro 1
2^1 2^1 2^1 2^1
2^14
0 1 2 3
fn1 (Hz) 4 4 4,5 4,25 4,375
14,7
fn2 (Hz) 14 15 15 5 14,88
fn3 (Hz) 22 22 22 22 22
De acuerdo a la Tabla 5, se decidió trabajar con un valor N de 8192 ya que se definen bien
lo peaks del espectro y además los valores de las frecuencias naturales se acercan a los
calculados en el ensayo con el martillo. Se podría haber elegido el valor N de 1024, ya que
es el que más se asemeja a los valores experimentales previos, pero la gráfica del espectro
se define pobremente.
Una vez identificado el valor de N, para desarrollar el código, es necesario cargar el registro
y separar los datos en vectores como aparece en la Figura 7. La frecuencia de muestreo fue
de 2048 Hz y se cuenta con 7.300.000 de datos.
Figura 7: Parte Nº1 de código Matlab.
Luego, se separan los bloques considerando un N de 8192 datos de la forma en que se
observa en la Figura 8. Se define el vector t_i que contiene las posiciones de los datos para
luego extraerlos desde los vectores A1, A2 y A3. Por ejemplo, en el primer bloque,
t_i=[1:8192], para el segundo bloque, t_i=[4097:12288]. Puede ser que el último bloque no
contenga los N datos, por esta razón va desde t_i=[mitad bloque anterior : total de datos
(7.300.000)].
10
Figura 8: Parte Nº2 de código Matlab.
Se elimina el último elemento de cada celda Blot, BloA1, BloA2 y BloA3 para poder
convertir los arreglos en matrices, donde las filas corresponden al dato “n” y las columnas
al bloque “k”, según Figura 9.
11
Figura 9: Parte Nº3 de código Matlab.
Para aplicar de manera correcta la ventana de Hanning, es necesario restar a cada uno de los
datos de cada bloque, la respectiva media aritmética de cada bloque para centrar los datos al
0, como se muestra en la Figura 10.
Figura 10: Parte Nº4 de código Matlab.
Luego, es posible aplicar la ventana de Hanning a cada bloque (Figura 11), utilizando el
comando hann(N) con el objetivo de tener señales periódicas que empiecen y finalicen en 0.
Figura 11: Parte Nº5 de código Matlab.
12
En la Figura 12, se observa el efecto de la aplicación de la ventana en un segmento
cualquiera del acelerómetro 1.
Figura 12: Ventana de Hanning en el segmento 1500 del acelerómetro 1.
Una vez verificado el efecto de la ventana de Hanning, es posible calcular las FFT y el
respectivo conjugado de cada uno de los datos, para finalmente obtener Gxx, definido en la
ecuación 4.1. Las FFT y sus conjugados, se obtienen con el comando “fft” y “conj” como
se observa en la Figura 13.
Figura 13: Parte Nº6 de código Matlab.
13
Una vez obtenidos todos los Gxx, para acelerómetros 1, 2 y 3, se promedian para cada N,
como se observa en la Figura 14. Lo anterior, significa que el valor de un punto del espectro
final se obtiene mediante un promedio simple de cada uno de los puntos de los espectros
individuales.
Figura 14: Parte Nº7 de código Matlab.
Finalmente, se calcula el espectro de frecuencias para cada acelerómetro mediante el
comando “pwelch” en función de los vectores A1, A2 y A3, y el número de datos N (Figura
15).
Figura 15: Parte Nº8 de código Matlab.
Para poder comparar la forma de los espectros obtenidos mediando el código implementado
y el comando pwelch, se normalizaron los espectros Gxx1, Gxx2 y Gxx3 respecto a los
valores pw1, pw2 y pw3, así se obtienen las dos curvas superpuestas para observar si la
forma coincide. En la Figura 16, se observa que la forma de las dos curvas coincide, aunque
no es tan relevante dicho hecho porque lo importante es que los valores en el eje X
coincidan y no necesariamente las amplitudes. En un rectángulo negro se encierran los
peaks correspondientes a las frecuencias de los modos excitados.
14
Figura 16: Representación gráfica de espectros de cada acelerómetro normalizados.
En la Figura 17, se muestra un zoom al recuadro negro de la Figura 16. Esta vez, se
graficaron los valores Gxx1, Gxx2 y Gxx3 sin normalizar para comprobar que lo que
importa es el valor de los peaks en el eje X y no las amplitudes.
15
Figura 17: Zoom de Figura 15.
4. Determinación de propiedades modales mediante datos experimentales
Los resultados de las frecuencias naturales obtenidas para los primeros 3 modos de la estructura son
los que se encuentran en la Tabla 6.
Tabla 6:
Resultados de frecuencias naturales para los primeros 3 modos según código implementado y
comando pwelch
Código implementado Comando pwelch
A1 A2 A3 A1 A2 A3
fn1 (Hz) 4.25 4.25 4.25 4.25 4.25 4.25
fn2 (Hz) 14.75 14.75 14.75 14.75 14.75 14.75
fn3 (Hz) 22 22 22 22 22 22
16
Como se observa en la tabla anterior, los resultados obtenidos del código implementado en
MATLAB, son iguales a los resultados obtenidos de función pwelch de MATLAB.
Además, se obtuvo el error relativo (ecuación 4.1) entre el cálculo teórico de las frecuencias (Tabla
3) y las calculadas en el ensayo de vibraciones ambientales. Los resultados se muestran en la Tabla
7. De manera similar, se comparan los resultados experimentales calculados en los ensayos del
matillo (Tabla 4) con los del ensayo de vibraciones ambientales, donde el error relativo se obtuvo a
partir de la ecuación 4.2. Los resultados se muestran en la Tabla 8.
f n ,teo−f n , código
Error relativo= ∙ 100 %(4.1)
f n ,teo
f martillo−f n ,código
Error relativo= ∙100 % (4.2)
f martillo
Tabla 7:
Error relativo entre frecuencias obtenidas mediante ensayo de vibraciones ambientales y
frecuencias teóricas de la estructura
Error relativo
A1 A2 A3
fn1 (Hz) 12% 12% 12%
fn2 (Hz) 11% 11% 11%
fn3 (Hz) 16% 16% 16%
Tabla 8:
Error relativo entre frecuencias obtenidas mediante ensayo de vibraciones ambientales y
frecuencias experimentales mediante ensayo de martillo
Error relativo
A1 A2 A3
fn1 (Hz) 3,4% 3,4% 3,4%
fn2 (Hz) 6,0% 6,0% 6,0%
fn3 (Hz) 5,1% 5,1% 5,1%
Se observa que los errores respecto al modelo teórico, son mayores al 10% por lo que se considera
poco preciso. En cambio, al comparar los resultados entre ambos experimentos, se puede concluir
que se semejan bastante (errores menores al 6%), lo cual implica que es necesario ajustar el modelo
teórico para sea más representativo del sistema.
17
5. Comentarios y Conclusiones
De acuerdo al análisis anterior, se logró implementar en MATLAB el Método de Welch,
determinando la densidad auto-espectral de la respuesta y las frecuencias naturales del
sistema de 3 grados de libertad. Los resultados obtenidos por la implementación del método
y los arrojados por el comando “pwelch” de MATLAB resultaron ser idénticos.
Según el análisis de sensibilidad realizado, es necesario considerar no sólo que las
frecuencias naturales coincidan según experimentos realizados anteriormente, sino también
tomar en cuenta la definición del espectro para poder tener peaks bien definidos y así
disminuir el leakage. Podría darse el caso de que los resultados obtenidos en el experimento
del martillo hayan tenido una distorsión producto de ese fenómeno y como los datos no
fueron tratados, se obtienen frecuencias erróneas producto de la no periodicidad de la señal
medida, requisito crucial para calcular la FFT.
Adicionalmente, se compraron las frecuencias modales de este ensayo, con el modelo
teórico del sistema y el ensayo anterior del matillo. Se concluye que los resultados de los 2
ensayos experimentales se asemejan entre sí y que al modelo teórico hay que realizarle
ajustes para que sea más representativo del sistema, pues los errores relativos entre este
ensayo y el modelo teórico son mayores a 10%, en cambio, respecto al ensayo del matillo,
son menores que 6%.
6. Referencias
[1]. Avitabile, P. (2018). Modal Testing. USA: john wiley & sons Ltd.
[2]. Brandt, A. (2010). Noise and vibration analysis. John Wiley & Sons Ltd.
18