Departamento de Ingenier´ Ingenier´ ıa Mec´ Mecanica anica ´ Facultad de Ciencias F´ısicas y Matem´ aticas aticas Universidad de Chile
Din´ Din´ amic amica a Estr Estruc uctu tura rall Apuntes para el curso ME706
Viviana Meruane
´Indice general Contenidos
II
1 Intr Introdu oducc cci´ i´ on on
2
1.1 Siste Sistema mass de un un grad grado o de liber libertad tad . . . . . . . . . . . . . . . . . . .
2
1.1. 1.1.11
Ecua Ecuaci ci´´on on de movimiento, funci´on on de transferencia . . . . . .
2
1.1.2
Polos olos del del sist sistem ema, a, frec frecue uenc ncia iass natu natura rale less y fact factor ores es de amortiguamiento . . . . . . . . . . . . . . . . . . . . . . . .
3
1.1. 1.1.33
Solu Soluci ci´´on on de la ecuaci´on on de movimiento . . . . . . . . . . . .
4
1.1. 1.1.44
Resi Residu duos os . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
1.1.5 Funci´on on de respuesta en frecuencia y funci´on on de respuesta a un impulso . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
1.1.6 1.1.6
Influenci Influencia a de camb cambios ios en la masa, masa, rigide rigidezz y amorti amortiguac guaci´ i´on on .
8
1.2 Sistemas Sistemas con multip multiples les grado gradoss de libertad libertad . . . . . . . . . . . . . .
10
1.2. 1.2.11
Ecua Ecuaci ci´´on on del sistema, funci´on on de transferencia . . . . . . . .
10
1.2.2 1.2.2
Polos, Polos, frecue frecuencia nciass naturale naturaless y factores factores de de amortigu amortiguamie amient ntoo
12
1.2.3 1.2 .3
Vectore ectoress modal modales, es, resid residuos uos . . . . . . . . . . . . . . . . . .
12
1.2.4
Funci´ on on de respuesta en frecuencia (FRF) y funci´ on on de respuesta a un impulso (IRF) . . . . . . . . . . . . . . . . .
14
Sistem Sistemas as sin amort amortigu iguami amien ento to y con amorti amortigu guami amien ento to proporcional . . . . . . . . . . . . . . . . . . . . . . . . . .
16
Ortogona Ortogonalida lidad d y coordenada coordenadass modales modales . . . . . . . . . . . .
18
1.2.5 1.2.6 1.2.6
ii
´ INDICE GENERAL
iii
2 El M´ etodo de Elementos Finitos 2.1
2.2
2.3
2.4
Elemento de barra . . . . . . . . . . . . . . . . . . . . . . . . . . .
23
2.1.1
Coordenadas locales y globales . . . . . . . . . . . . . . . .
24
Elemento de viga . . . . . . . . . . . . . . . . . . . . . . . . . . . .
25
2.2.1
Coordenadas locales y globales . . . . . . . . . . . . . . . .
27
Elemento de viga 3D . . . . . . . . . . . . . . . . . . . . . . . . . .
28
2.3.1
Rotaci´on . . . . . . . . . . . . . . . . . . . . . . . . . . . .
30
Ensamble . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
31
3 Integraci´ on num´ erica de las ecuaciones de movimiento 3.1 M´ etodos de integraci´on directa . . . . . . . . . . . . . . . . . . . .
3.2
21
35 35
3.1.1
El m´etodo de las diferencias centrales . . . . . . . . . . . .
36
3.1.2
El m´etodo de Wilson θ . . . . . . . . . . . . . . . . . . . . .
37
3.1.3
El m´etodo de Newmark . . . . . . . . . . . . . . . . . . . .
40
3.1.4
Ejemplo . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
42
Superposici´ on modal . . . . . . . . . . . . . . . . . . . . . . . . . .
43
3.2.1
43
Ejemplo . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4 Procesamiento de se˜ nales 4.1 La transformada de Fourier . . . . . . . . . . . . . . . . . . . . . . 4.1.1
45 46
Algunos par´ ametros importantes . . . . . . . . . . . . . . .
48
4.2 Errores y ventanas . . . . . . . . . . . . . . . . . . . . . . . . . . .
48
4.2.1
Aliasing . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
49
4.2.2
Leakage . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
49
4.2.3
Ventanas . . . . . . . . . . . . . . . . . . . . . . . . . . . .
50
4.3 Funciones de respuesta en frecuencia y funciones de coherencia . .
54
4.3.1
Efectos del ruido experimental . . . . . . . . . . . . . . . .
5 Estimaci´ on de par´ ametros modales
56
60
´ INDICE GENERAL
iv
5.1 M´ etodos de un grado de libertad . . . . . . . . . . . . . . . . . . .
5.2
61
5.1.1
Peak picking . . . . . . . . . . . . . . . . . . . . . . . . . .
61
5.1.2
Circle fitting . . . . . . . . . . . . . . . . . . . . . . . . . .
62
M´ etodos con multiples grados de libertad en el dominio de frecuencias 64
5.2.1
M´etodo de m´ınimos cuadrados no lineales, LSFD . . . . . .
64
5.2.2
Rational Fractional Polynomials . . . . . . . . . . . . . . .
65
5.3 M´etodos con multiples grados de libertad en el dominio del tiempo
68
5.3.1
El m´etodo Ibrahim . . . . . . . . . . . . . . . . . . . . . . .
68
5.3.2
El m´ etodo least-squares complex exponential (LSCE) . . .
70
5.4 Diagramas de estabilidad . . . . . . . . . . . . . . . . . . . . . . .
72
6 Medici´ on Experimental 6.1
74
An´alisis previo a las mediciones . . . . . . . . . . . . . . . . . . . .
74
6.1.1
Rango de frecuencias . . . . . . . . . . . . . . . . . . . . . .
75
6.1.2
Selecci´ on de la ubicaci´on de las respuestas . . . . . . . . . .
75
6.1.3
Selecci´ on de los puntos de excitaci´on . . . . . . . . . . . . .
76
6.1.4
Selecci´ on de los puntos de suspensi´on . . . . . . . . . . . .
77
6.2 El montaje experimental . . . . . . . . . . . . . . . . . . . . . . . .
78
6.3
6.4
6.2.1
Mecanismo de Excitaci´on . . . . . . . . . . . . . . . . . . .
78
6.2.2
Aceler´ ometro . . . . . . . . . . . . . . . . . . . . . . . . . .
79
6.2.3
Sensor de fuerzas . . . . . . . . . . . . . . . . . . . . . . . .
80
Selecci´on de la fuerza de excitaci´on . . . . . . . . . . . . . . . . . .
82
6.3.1
Excitaci´on sinusoidal . . . . . . . . . . . . . . . . . . . . . .
82
6.3.2
Excitaci´on aleatoria . . . . . . . . . . . . . . . . . . . . . .
83
6.3.3
Excitaci´on pseudo-aleatoria . . . . . . . . . . . . . . . . . .
84
6.3.4
Excitaci´on aleatoria en trenes (burst random) . . . . . . . .
84
6.3.5
Excitaci´on de impacto . . . . . . . . . . . . . . . . . . . . .
85
Evaluaci´on inicial de las FRFs medidas . . . . . . . . . . . . . . . .
86
´ INDICE GENERAL
v
6.4.1
Repetibilidad . . . . . . . . . . . . . . . . . . . . . . . . . .
87
6.4.2
Reciprocidad . . . . . . . . . . . . . . . . . . . . . . . . . .
87
6.4.3
Linealidad . . . . . . . . . . . . . . . . . . . . . . . . . . . .
87
6.4.4
Caracter´ısticas especiales de una FRF . . . . . . . . . . . .
87
7 Correlaci´ on Num´ erico-Experimental y Ajuste de Modelos 7.1 Parear modelos num´ ericos y experimentales . . . . . . . . . . . . .
89 92
7.1.1
T´ecnicas de reducci´ on . . . . . . . . . . . . . . . . . . . . .
92
7.1.2
T´ecnicas de expansi´ on . . . . . . . . . . . . . . . . . . . . .
94
7.2 T´ ecnicas de correlaci´ on . . . . . . . . . . . . . . . . . . . . . . . . .
95
7.2.1
Diferencia en las frecuencias de resonancia . . . . . . . . . .
96
7.2.2
Comparaci´ on visual de los modos . . . . . . . . . . . . . . .
96
7.2.3
Modal Assurance Criterion . . . . . . . . . . . . . . . . . .
97
7.2.4
Frequency Response Assurance Criterion . . . . . . . . . .
98
7.3 M´ etodos iterativos de ajuste de modelos . . . . . . . . . . . . . . .
98
7.3.1
Modos y frecuencias naturales . . . . . . . . . . . . . . . . . 101
7.3.2
Funciones de respuesta en frecuencia . . . . . . . . . . . . . 102
8 Ejemplos de An´ alisis Modal en Estructuras Reales
104
8.1 Estructura de barras . . . . . . . . . . . . . . . . . . . . . . . . . . 104 8.1.1 8.2
Modelo num´erico . . . . . . . . . . . . . . . . . . . . . . . . 108
Modelo a escala de un avion . . . . . . . . . . . . . . . . . . . . . . 110 8.2.1
Modelo num´erico . . . . . . . . . . . . . . . . . . . . . . . . 111
8.3 Viga de concreto reforzado . . . . . . . . . . . . . . . . . . . . . . . 113 8.3.1 8.4
Tubo de escape de un auto . . . . . . . . . . . . . . . . . . . . . . 116 8.4.1
8.5
Modelo num´erico . . . . . . . . . . . . . . . . . . . . . . . . 114
Modelo num´erico . . . . . . . . . . . . . . . . . . . . . . . . 117
El puente I-40 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119
´ INDICE GENERAL
1
8.5.1
Modelo num´erico . . . . . . . . . . . . . . . . . . . . . . . . 121
Bibliograf´ıa
124
Cap´ıtulo 1
Introducci´ on 1.1.
Sistemas de un grado de libertad
1.1.1.
Ecuaci´ on de movimiento, funci´ on de transferencia
Consideremos el sistema de un grado de libertad mostrado en la Figura 1.1. La ecuaci´ on de movimiento de este sistema viene dada por: mx ¨ + cx˙ + kx = f (t)
(1.1)
donde, m: masa k: coeficiente de rigidez c: coeficiente de amortiguaci´on x ¨, x˙ , x: aceleraci´on, velocidad y desplazamiento f : fuerza de excitaci´on externa t: tiempo
2
SISTEMAS DE UN GRADO DE LIBERTAD
3
� ���� �
����
�
Figura 1.1: Transformando la ecuaci´on 1.1 al dominio de Laplace y asumiendo que el desplazamiento y velocidad inicial son 0, se obtiene: (mp2 + cp + k )X ( p) = F ( p)
(1.2)
Lo que se puede escribir como: Z ( p)X ( p) = F ( p)
(1.3)
Z ( p) se denomina la rigidez din´amica.
Invirtiendo la ecuaci´on 1.2 ´ o 1.3 se obtiene la funci´ on de transferencia H ( p) = Z −1 ( p): X ( p) = H ( p)F ( p) H ( p) =
1.1.2.
(1.4)
1/m p2 + (c/m) p + (k/m )
(1.5)
Polos del sistema, frecuencias naturales y factores de amortiguamiento
El denominador de la ecuaci´on 1.5 representa la ecuaci´on caracter´ıstica del sistema. Las ra´ıces de esta ecuaci´on son: λ1,2 =
−
c 2m
± − c 2m
2
k m
(1.6)
Esta ecuaci´on permite la introducci´on de varios conceptos importantes. Si no hay amortiguamiento, el sistema en estudio es un sistema conservativo (c=0), y la frecuencias natural (rad/seg) del sistema sin amortiguamiento se define como: ωn =
k m
(1.7)
SISTEMAS DE UN GRADO DE LIBERTAD
4
Se define como amortiguamiento cr´ıtico, cc , al valor del amortiguamiento que hace que el termino dentro de la ra´ız en la ecuaci´on 1.6 sea cero: cc = 2m
k m
(1.8)
A partir de esta definici´on, se define la raz´on de amortiguamiento, ζ , como: c ζ = (1.9) cc Las ra´ıces presentadas en la ecuaci´on 1.6, se pueden escribir como dos ra´ıces complejas conjugadas: λ λ∗
= σ + jω d
(1.10)
= σ
(1.11)
− jω d
donde σ es el factor de amortiguamiento, ωd es la frecuencia natural del sistema amortiguado y ∗ representa el complejo conjugado
A partir de la definici´on anterior se pueden definir las siguientes relaciones: λ = ωn
=
ζ = σ ωd
1.1.3.
=
− − − − − ζ + j
1
ζ 2 ωn
ωd2 + σ 2
(1.13)
σ
(1.14)
ωd2 + σ 2
ζωn
= ωn
1
(1.12)
ζ 2
(1.15) (1.16)
Soluci´ on de la ecuaci´ on de movimiento
La soluci´on de la ecuaci´on 1.1 en el tiempo se puede determinar asumiendo una soluci´ on de la forma: x = Ae λt
(1.17)
donde A es una constante. Por lo tanto: x˙ = λAeλt
(1.18)
x ¨ = λ2 Aeλt
(1.19)
SISTEMAS DE UN GRADO DE LIBERTAD
5
Sustituyendo en la ecuaci´on 1.1 y asumiendo que las fuerzas externas son nulas, se obtiene: mλ2 Aeλt + cλAeλt + kAe λt = 0
(1.20)
Dividiendo por Aeλt se obtiene la ecuaci´on caracter´ıstica del sistema: mλ2 + cλ + k = 0
(1.21)
cuya soluci´on son las ra´ıces de la ecuaci´on 1.6. Por lo tanto, x = A 1 eλ t y x = A 2 eλ t , son soluciones del problema y su suma es la soluci´on general: 1
x = A 1 eλ t + A2 eλ 1
2
t
2
(1.22)
Las soluciones particulares dependen de las constantes, A 1 y A 2 , que se determinan a partir de las condiciones iniciales.
La ecuaci´on 1.22 se puede escribir como:
x = e −ζω n t A1 ejωd t + A2 e−jω d t
Utilizando las identidades: ejω d t = cosωd t + jsenωd t y e−jωd t = cosωd t la ecuaci´on 1.23 se puede escribir como: x = e −ζω n t (Bcosωd t + Dsenωd t)
(1.23)
− jsenωdt, (1.24)
ζ < 1 d u t i l p 0 m A
ζ = 1 d u t i l p m A
0
ζ > 1 d u t i l p m A
0 Tiempo
Figura 1.2:
SISTEMAS DE UN GRADO DE LIBERTAD
6
Dependiendo del valor de la raz´on de amortiguamiento el sistema se puede clasificar como sobreamortiguado (ξ > 1), con amortiguamiento cr´ıtico (ξ = 1) o con amortiguamiento d´ ebil (ξ < 1). Como se muestra en la Figura 1.2, la respuesta de un sistema sobreamortiguado es decreciente sin oscilaciones. En cambio, la respuesta de un sistema con amortiguamiento d´ebil es una respuesta oscilatoria decreciente. Amortiguamiento cr´ıtico es el caso de borde entre un sistema sobreamortiguado y con amortiguamiento d´ebil. Para el caso de una respuesta forzada, en donde f (t) de la ecuaci´on 1.1 es distinto de cero. La soluci´on viene dada por la soluci´on general cuando f = 0 m´as una soluci´ on particular X : x = A 1 eλ t + A2 eλ t + X 1
2
(1.25)
La forma de X depende de la forma de la fuerza de excitaci´on, la tabla a continuaci´on muestra varios casos particulares:
Forma de f (t) f (t) = A , constante f (t) = At f (t) = At 2 f (t) = Asen(ωt ) o Acos(ωt )
Soluci´on particular B Bt + D Bt 2 + D Bcos(ωt ) + Dsen(ωt )
Tabla 1.1:
1.1.4.
Residuos
Conociendo las ra´ıces del sistema, la ecuaci´on 1.5 se puede escribir como: H ( p) =
( p
−
1/m λ)( p λ∗ )
−
(1.26)
Ocupando fracciones parciales se tiene que: H ( p) =
A
( p
−
A∗ 1/m , A + con = λ) ( p λ∗ ) 2 jω d
−
(1.27)
En esta formulaci´on A y A∗ se denominan residuos .
1.1.5.
Funci´ on de respuesta en frecuencia y funci´ on de respuesta a un impulso
En las secciones anteriores se discute la relaci´ on entre la fuerza y la respuesta de un sistema de un grado de libertad en el dominio de Laplace. Esta relaci´ on
SISTEMAS DE UN GRADO DE LIBERTAD
7
tambi´en se puede expresar en el dominio de frecuencias o del tiempo. La funci´on de transferencia evaluada en el eje de frecuencias ( jω ) se denomina funci´on de respuesta en frecuencia (FRF): H ( p) p=jω = H (ω ) =
|
A∗
A
( jω
− λ) + ( jω − λ )
∗
(1.28)
En general, la contribuci´ on de la parte del complejo conjugado (o la parte de la frecuencia negativa) es despreciable cerca de la resonancia, ω = ω d . En consecuencia, la funci´on de respuesta en frecuencia para un grado de libertad frecuentemente se aproxima por: H (ω )
( jω A− λ)
(1.29)
La Figura 1.3 muestra un gr´afico de la amplitud y fase de una funci´on de respuesta en frecuencia. Notar que el peak en la FRF se obtiene cuando ω = ω d , a esta misma frecuencia la fase cambia en 180◦ .
−
d u t i l p m A
0
e s a F
-pi ωd
Figura 1.3: Aplicando una transformada inversa de Laplace a la expresi´on de la funci´on de transferencia se obtiene la funci´on de respuesta a un impulso. La respuesta a un impulso corresponde a la respuesta a un delta de Dirac en t = 0. Para un grado de libertad la funci´on de respuesta a un impulso es:
h(t) = e σt Aeiωd t + A∗ e−iωd t
(1.30)
SISTEMAS DE UN GRADO DE LIBERTAD
8
El residuo A define la amplitud inicial, σ la tasa de decrecimiento y ω d la frecuencia de oscilaci´ on. Notar que el gr´ afico de la funci´ on de respuesta a un impulso es equivalente al de la Figura 1.2.
1.1.6.
Influencia de cambios en la masa, rigidez y amortiguaci´ on
Las Figuras 1.4, 1.5 y 1.6 muestran como, la funci´ on de respuesta en frecuencia de un sistema de 1 grado de libertad, se ve afectada por cambios en la rigidez, amortiguamiento y masa. Un incremento en la rigidez resulta en una frecuencia de resonancia mayor y un valor menor de la FRF en el rango de frecuencias bajas. Debido al efecto dominante de la rigidez a bajas frecuencias, esta regi´ on de la FRF se denomina la linea de rigidez o “compliance line”.
) F R F ( g o l
frecuencia
Figura 1.4: Incremento de la rigidez Un incremento del amortiguamiento produce una peque˜na reducci´on de la frecuencia de resonancia. Sin embargo, la mayor contribuci´on es una disminuci´on de la amplitud de la funci´on de respuesta en frecuencia en las cercan´ıas a la resonancia. Tambi´en, la fase cambia m´ as suavemente. Si el amortiguamiento es cero, la magnitud en la resonancia se vuelve infinito y la fase disminuye repentinamente en 180 grados. Adicionalmente, los polos del sistema, λ, se vuelven puramente imaginarios y en amplitud iguales a la frecuencia natural sin amortiguamiento ωn = k/m .
SISTEMAS DE UN GRADO DE LIBERTAD
9
) F R F ( g o l
frecuencia
Figura 1.5: Incremento del amortiguamiento Un incremento de masa disminuye el valor de la frecuencia de resonancia. La amplitud de la FRF a altas frecuencias tambi´en disminuye. Debido al efecto dominante de la masa a altas frecuencias esta ´area de la FRF se denomina la “linea de masa”.
) F R F ( g o l
frecuencia
Figura 1.6: Incremento de la masa
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
1.2.
10
Sistemas con multiples grados de libertad
En esta secci´on vamos a extender los conceptos vistos en la secci´on anterior para el caso con multiples grados de libertad.
1.2.1.
Ecuaci´ on del sistema, funci´ on de transferencia
En la mayor´ıa de los casos el sistema no se puede describir como un sistema de un grado de libertad. La mayor´ıa de los sistemas consisten en un ensamble de un n´umero infinito de masas, resortes y amortiguadores. A continuaci´ on se mostrar´a como se relaciona la funci´on de transferencia de sistema con multiples grados de libertad con los par´ametros modales (frecuencias de resonancia y vectores modales). En el caso de un sistema con multiple grados de libertad la ecuaci´on de movimiento equivale a la ecuaci´on de un grado de libertad, pero matricial. Para ilustrar esto, se desarrollar´ a el caso de un sistema de dos grados de libertad mostrado en la Figura 1.7. � ���� �� ���
� ���� �� ��� ��
��
�� ��
�� ��
��
��
Figura 1.7: Ejemplo de un sistema de dos grados de libertad Las ecuaciones de movimiento para este sistema son: m1 x ¨1 + ( c1 + c2 )x˙ 1
− c2 ˙x2 + (k1 + k2 )x1 − k2x2 m2 x ¨2 + ( c2 + c3 )x˙ 2 − c2 ˙x1 + ( k2 + k3 )x2 − k2 x1
= f 1
(1.31)
= f 2
(1.32)
En notaci´on matricial:
m1 0 0 m2
x c + c2 c2 ¨1 + 1 x c2 c2 + c3 ¨2
−
−
x˙ 1 k + k2 k2 + 1 x˙ 2 k2 k2 + k3
−
−
x1 x2
=
f 1 (1.33) f 2
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
11
que se puede escribir como: M x ¨ + C x˙ + K x = f
{}
{}
{} {}
(1.34)
donde, M: matriz de masa K: matriz de rigidez C: matriz de amortiguaci´on
{x}: vector de respuesta {f }: vector de fuerzas Esta ecuaci´on tambi´ en describe el comportamiento de sistemas con m´as grados de libertad. La dimensi´ on de las matrices aumenta con el n´umero de grados de libertad. Transformando esta ecuaci´ on de movimiento en el dominio de Laplace, asumiendo que las velocidades y desplazamientos iniciales son cero, se obtiene: (M p2 + Cp + K )X ( p) = F ( p)
(1.35)
Lo que se puede escribir como: Z ( p)X ( p) = F ( p)
(1.36)
Z ( p) es la matriz de rigidez din´amica.
Invirtiendo la ecuaci´on 1.35 ´ o 1.36 se obtiene la matriz funci´on de transferencia H ( p): X ( p) = H ( p)F ( p) H ( p) = Z ( p)−1 =
(1.37) adj (Z ( p)) Z ( p)
|
|
(1.38)
donde: adj(Z ( p)) es la matriz adjunta de Z ( p) y Z ( p) es el determinante de Z ( p).
|
|
T La matriz adjunta viene dada por adj (Z ( p)) = [εij Z ij ] . Donde:
| |
|Z ij | es el determinante de Z ( p), eliminando la fila i y la columna j εij = 1, si i + j es par; = -1, si i + j es impar
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
1.2.2.
12
Polos, frecuencias naturales y factores de amortiguamiento
on caracter´ıstica del sistema viene dada por el denominador de la ecuaci´on La ecuaci´ 1.38, es decir, el determinante de Z ( p). Al igual que en el caso de un grado de libertad, las ra´ıces de esta ecuaci´ on o polos del sistema, definen las frecuencias naturales. Estas ra´ıces se pueden determinar resolviendo un problema de valores propios. Para transformar la ecuaci´on 1.35 en una ecuaci´on general de valores propios, se a˜nade la siguiente identidad:
( pM pM ) x = 0
−
(1.39)
{}
combinando las ecuaciones 1.35 y 1.39 se obtiene: ( pA + B ) y = f
{} { }
(1.40)
donde: A =
− { } 0 M
M , B= C
M 0 px , y = K x 0
y f =
{ }
0 f
(1.41)
Por lo tanto, los polos del sistema son los valores de p que satisfacen:
| pA + B| = 0
(1.42)
Notar que las ra´ıces de esta ecuaci´on, tambi´ en son las ra´ıces de la ecuaci´on Z = 0. Esta ecuaci´on genera 2n (n= n´umero de grados de libertad) valores propios complejos, que aparecen en pares de complejos conjugados:
| |
[Λ] =
λ1
..
.
0 λn
0
λ∗1
..
. λ∗n
=
σ1 + jω1
..
.
0 σn + jωn σ1
0
− jω1
..
. σn
− jωn
(1.43)
Como en el caso de un grado de libertad, la parte real del polo, σi , es el factor de amortiguamiento y la parte imaginaria, ωi , es la frecuencia natural amortiguada.
1.2.3.
Vectores modales, residuos
A cada valor propio del sistema le corresponde un vector propio. Para sistemas con multiples grados de libertad estos vectores propios introducen el concepto de
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
13
modos normales, formas modales o vectores modales, φi . Estos tambi´en aparecen en pares de complejos conjugados. Cada vector propio, se relaciona a un valor propio especifico.
Φ=
λ1 φ1 . . . φ1 . . .
λn φn φn
λ∗1 φ∗1 . . . λ∗n φ∗n φ∗1 . . . φ∗n
(1.44)
En general, los vectores propios contiene desplazamientos modales con valores complejos. Para el valor propio correspondiente, λi , se cumple que: (M λ2i + Cλi + K )φi = 0
(1.45)
De manera similar al caso de un grado de libertad, se puede introducir el concepto de residuos. Dado que λi , λ∗i son ra´ıces de la ecuaci´on caracter´ıstica del sistema, la ecuaci´ on 1.37 se puede escribir como: H ( p) =
adj (Z ( p)) n λi )( p i=1 E ( p
−
∗
− λi )
(1.46)
donde E es una constante. Utilizando fracciones parciales: n
H ( p) =
i=1
[A]i [A]∗i + ( p λi ) ( p λ∗i )
−
−
(1.47)
Los t´erminos [A]i y [A]∗i se denominan residuos , y vienen dados por: [A]i = P i adj (Z (λi ))
(1.48)
P i es una constante que depende del polo. Estudiando estos residuos se puede esclarecer la relaci´ on entre la matriz de funci´ on de transferencia y los vectores modales. Reescribiendo la ecuaci´on 1.38: Z ( p)adj (Z ( p)) = Z ( p) [I ]
|
|
(1.49)
Evaluando para p = λ i , se obtiene: Z (λi )adj (Z (λi )) = 0
(1.50)
Si consideramos una columna arbitraria de adj (Z (λi )) se cumple que: Z (λi ) adj (Z (λi ))
{
}k = 0
(1.51)
Esta ecuaci´on es equivalente a la ecuaci´on 1.45. Esto quiere decir la columnas de adj(Z (λi )) son proporcionales al i-´ esimo vector modal. Considerando que la matriz
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
14
adj(Z (λi )) es sim´etrica, entonces sus filas tambi´en son proporcionales al i-´esimo vector modal. Es decir, la matrix adj (Z (λi )) es de la forma:
adj (Z (λi )) = R i φi φT i = R i
φi,1 φi,1 φi,2 φi,1 .. .
φi,1 φi,2 . . . φi,1 φi,n φi,2 φi,2 . . . φi,2 φi,n .. .. .. . . . φi,n φi,2 . . . φi,n φi,n
φi,n φi,1
(1.52)
donde Ri es una constante asociada al vector φi . A partir de la ecuaci´on anterior, se tiene que cada residuo es de la forma: [A]i = Q i φi φT i
(1.53)
Y se obtiene finalmente que: n
H ( p) =
i=1
1.2.4.
Qi φi φT Q∗i φ∗i φ∗i T i + ( p λi ) ( p λ∗i )
−
−
(1.54)
Funci´ on de respuesta en frecuencia (FRF) y funci´ o n de respuesta a un impulso (IRF)
La funci´on de respuesta en frecuencia (FRF) es la funci´on de transferencia evaluada en el eje de frecuencias ( jω ): n
H ( jω ) =
i=1
Qi φi φT Q∗i φ∗i φ∗i T i + ( jω λi ) ( jω λ∗i )
−
−
(1.55)
La funci´on de respuesta a un impulso (IRF) viene dada por la transformada inversa de la funci´on de respuesta en frecuencia: N
h(t) =
∗
λr t (Qr φr φT + Q∗r φ∗r φ∗r T eλr t ) r e
(1.56)
r=1
Se debe notar que, en el caso de multiples grados de libertad, para cada valor de ω, H ( jω ) es una matriz. H ik ( jω ) corresponde a funci´ on de respuesta en frecuencia, cuando se excita la estructura en k y se mide la respuesta en i, o viceversa (H ( jω) es sim´etrica):
H ik ( jω ) = H ki ( jω ) =
X i ( jω ) X k ( jω ) = F k ( jω ) F i ( jω )
(1.57)
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
15
Las Figuras 1.8,1.9 y 1.10 muestran un ejemplo de las tres FRF para el sistema de dos grados de libertad de la Figura 1.7. Los “peaks” en las curvas son las frecuencias de resonancia, que tienen el mismo valor para todas las FRF, ya que son propiedades globales de la estructura. Por otro lado, las frecuencias donde las FRF tienen “ca´ıdas” se denominan antiresonancias. Los valores de las antiresonancias var´ıan para cada FRF, ya que son propiedades locales de la estructura.
)
1 1
H ( g o l
w1
w2
frecuencia
Figura 1.8: Funci´on de respuesta en frecuencia, excitaci´on en 1, respuesta en 1
)
2 1
H ( g o l
w1
w2
frecuencia
Figura 1.9: Funci´ on de respuesta en frecuencia, excitaci´ on en 1, respuesta en 2 ´o excitaci´on en 2, respuesta en 1
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
16
)
2 2
H ( g o l
w1
w2
frecuencia
Figura 1.10: Funci´on de respuesta en frecuencia, excitaci´on en 2, respuesta en 2 Se debe recalcar que la funci´on de respuesta en frecuencia posee valores complejos, y se puede representar de distintas maneras: parte real e imaginaria vs frecuencia, amplitud (usualmente en escala logar´ıtmica) y fase vs frecuencia, y por u ´ ltimo se puede graficar la parte real vs imaginaria con la frecuencia como un par´ ametro ´ variable dentro de la curva. Este u ´ ltimo se denomina diagrama de Nyquist, en la Figura 1.11 se muestra un ejemplo.
)
1 1
H ( g a m I
Real(H ) 11
Figura 1.11: Diagrama de Nyquist, H 11
1.2.5.
Sistemas sin amortiguamiento y con amortiguamiento proporcional
En las secciones anteriores se vio es caso m´as general, en donde, el amortiguamiento es aproximado como un amortiguamiento viscoso general. Estos sistemas entregan polos de valores complejos, vectores modales de valores complejos con fases distintas para cada vector y funciones de respuesta en frecuencia complejas. A continuaci´on
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
17
se estudiaran dos casos particulares: sin amortiguamiento y con amortiguamiento proporcional. Si la matriz de amortiguamiento, C, es cero, el sistema de ecuaciones en el dominio de Laplace es:
{ }
2
p M + K
x = f
{}
(1.58)
2
En este caso, las ra´ıces del polinomio caracter´ıstico p M + K , se encuentran f´acilmente resolviendo el siguiente problema de valores y vectores propios:
{ }
p2 M + K
x =0
(1.59)
La soluci´on a este problema entrega los valores propios, ωi2 , y los vectores propios correspondientes, φi . Se debe notar que en este caso, los vectores propios tienen valores reales. Es por esto, que en el caso sin amortiguamiento, los vectores modales se denominan vectores modales normales o reales .
−
La matriz de rigidez din´amica en el dominio de frecuencias ( p = jω ), Z ( jω ) =
−ω2 M + K
(1.60)
tiene valores reales. Por lo tanto, su inversa, la matriz de frecuencia en respuesta H ( jω), tambi´ en tiene valores reales. La expresi´on general 1.55 se puede simplificar a: n
H ( jω ) =
i=1
2 jω i Qi φi φT i 2 2 (ωi ω )
−
(1.61)
Notar que dado que H ( jω ) posee valores reales, entonces Qi es puramente imaginario. El caso de amortiguamiento proporcional es una forma hipot´etica de amortiguamiento. En su forma m´as sencilla, la matriz de amortiguamiento, C, toma la forma: C = αM + βK
(1.62)
donde α y β son constantes reales. En este caso el sistema de ecuaciones es:
{ } { }
p2 M + p(αM + βK ) + K
x
=
( p2 + pα)M + ( pβ + 1) K x
=
{f } {f }
(1.63)
(1.64)
Considerando la versi´on homog´enea de esta ecuaci´on y dividendo por ( pβ + 1):
p2 + pα pβ + 1
{ }
M + K
x =0
(1.65)
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
18
Esta ecuaci´ o n es similar al caso sin amortiguamiento. Los sistemas con amortiguamiento proporcional tienen polos complejos, λi , que cumplen: λ2i + λi α = λi β + 1
−ωi2
(1.66)
y modos normales, iguales a los del caso sin amortiguamiento. La complejidad de los c´ alculos para el caso de amortiguamiento proporcional es menor al caso de amortiguamiento general. Esta es la raz´on principal para la introducci´ on ´ del amortiguamiento proporcional. Este entrega un intermedio entre el caso sin amortiguamiento y el caso de amortiguamiento viscoso general. Con amortiguamiento proporcional, la funci´on de respuesta en frecuencia posee valores complejos y toma la forma: n
H ( jω ) =
i=1
1.2.6.
2 jω i Qi φi φT i 2 2 2 (σi + ωi ω ) 2 jω σi
−
(1.67)
−
Ortogonalidad y coordenadas modales
Los modos normales satisfacen un set de condiciones, conocidas como condiciones de ortogonalidad. Consideremos el problema de valores y vectores propios para el caso sin amortiguamiento o con amortiguamiento proporcional:
−
ωi2 M +
K φi = 0
Consideremos ahora dos modos arbitrarios r y s: Kφr
= ωr2 M φr
Kφs
= ωs2 M φs
pre-multiplicando la primera ecuaci´on por φT on por φT s y la segunda ecuaci´ r : φT s Kφr
= ωr2 φT s M φr
φT r K φs
= ωs2 φT r M φs
Transponiendo la segunda ecuaci´on y rest´andola a la primera, se tiene que:
− ωs2)φT s M φr = 0
(ωr2
(1.68)
Esta ecuaci´on nos da dos posibilidades diferentes: φT s M φr = 0
r = s
(1.69)
φT s M φr = 0
r = s
(1.70)
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
19
La primera ecuaci´on nos dice que dos modos distintos son ortogonales con la matriz de masa (tambi´ en con la matriz de rigidez) y que cuando se calcula el producto del mismo modo con la matriz de masa como en la segunda ecuaci´on, el producto no es cero. El valor de este producto se denomina masa generalizada mi mi = φ T i M φi
(1.71)
De manera similar se define la rigidez generalizada ki ki = φ T i Kφi
(1.72)
Dado lo anterior, se tiene que: ΦT M Φ = m
(1.73)
ΦT K Φ = k
(1.74)
Φ es la matriz de modos normales, m y k son matrices diagonales denominadas matriz de masa modal y matriz de rigidez modal respectivamente. Debido a que la matriz de modos normales esta sujeta a una normalizaci´on arbitraria, los valores de ki y mi no son ´unicos. Lo que si se puede demostrar es que la raz´ on ki /mi es u ´ nica y tiene un valor igual a ωi2 . Entre los distintos m´etodos de normalizaci´ on de los modos, ´el m´as relevante en an´alisis modal es la normalizaci´on con respecto a la matriz de masa. Los modos normalizados con respecto a la matriz de masa tiene la propiedad:
ΦT M Φ = I
(1.75)
ΦT K Φ = Ω
(1.76)
donde Ω es una matriz diagonal que contiene los valores caracter´ısticos ωi2 en la diagonal principal. Introduciendo la transformaci´ on modal x = Φq en la ecuaci´on 1.58 y preT multiplicando por Φ , se obtiene que: p2 ΦT M Φq + Φ T K Φq = ΦT f
(1.77)
( ω 2 m + k )q =
(1.78)
−
{} ΦT {f }
Para el caso de amortiguamiento proporcional, se define la matriz de amortiguamiento modal: c = ΦT C Φ
(1.79)
SISTEMAS CON MULTIPLES GRADOS DE LIBERTAD
20
Lo que da el siguiente sistema de ecuaciones desacopladas: ( p2 m + pc + k )q = ΦT f
{}
(1.80)
Cap´ıtulo 2
El M´ etodo de Elementos Finitos Fundamentalmente, existen dos tipos de m´ etodos de elementos finitos; el m´ etodo de las fuerzas, donde se asumen las fuerzas y se calculan los desplazamientos y el m´ etodo de los desplazamientos, donde se asumen los desplazamientos y se calculan las fuerzas. Este ´ ultimo es el m´as utilizado para an´alisis de vibraciones y ser´a el que se describir´a a continuaci´on.
El primer paso para construir un modelo en elementos finitos es discretizar la estructura en un n´umero de elementos. Se pueden definir tres familias de elementos. (1) Elementos unidimensionales (linea), (2) Elementos bidimensionales (planos) y (3) elementos tridimensionales (s´olidos). ��������� ���������������
��������� ����������������
�������� ������ ����� ����
��������� ������ ����
��������� ����������������
�������� �������� ������������ ����
Figura 2.1: Tipos de elementos finitos Para cada elemento se define u ˆ como los desplazamientos en los nodos y u(x), ε(x) y σ (x) como los desplazamientos, deformaciones y esfuerzos en el elemento, en 21
´ EL METODO DE ELEMENTOS FINITOS
22
funci´on de las coordenadas del elemento x. u ˆ y u se relacionan por la matriz de funciones de forma N (x), la que es caracter´ıstica para cada tipo de elemento. u(x) = N (x)ˆ u
(2.1)
la deformaci´on se puede expresar como: ε(x) = Du(x) = DN (x) u u ˆ = B (x)ˆ
(2.2)
≡B(x)
D es el operador diferencial espacial, que para el caso general viene dado por:
D =
∂ ∂x 1
0 ∂ ∂x 2
0 0
0
0 0 ∂ ∂x 3
∂ ∂x 2 ∂ ∂x 1
0
0 ∂ ∂x 3 ∂ ∂x 2
∂ ∂x 3
0 ∂ ∂x 1
T
(2.3)
Se puede demostrar que las matrices de rigidez y masa de un elemento viene dadas por: K =
B T HBdV
(2.4)
ρN T N dV
(2.5)
V
M =
V
Donde V es el volumen del elemento, ρ es la densidad del material y H es una matriz del material. La tabla a continuaci´ on entrega las expresiones de H para cada tipo de elemento. Problema
Matriz H
Barra
E
Viga
EI
Deformaciones planas
Esfuerzos planos
E 1−ν 2
1 ν ν 1 0 0
E (1−ν ) (1+ν )(1−2ν )
0 0 1−ν 2
1 ν (1−ν )
0
ν (1−ν )
1 0
0 0 1−2ν 2(1−ν )
ELEMENTO DE BARRA
Axis Ax isim´ im´etri etrico co
Placa en flexi´on on
Tridimensional
23
E (1 (1−ν ) (1+ν (1+ν )(1 )(1−2ν )
Eh 3 12(1−ν 2 )
ν (1−ν )
0
1 0
ν (1−ν )
ν (1−ν )
1 ν ν 1 0 0
(1−ν ) E (1 (1+ν (1+ν )(1 )(1−2ν )
ν (1−ν )
1
0 0 1−ν 2
1 ν (1−ν ) ν (1−ν )
0 0 0
0 0 1−2ν 2(1−ν )
0
ν (1−ν )
ν (1−ν ) ν (1−ν )
1 ν (1−ν )
1 0 0 0
0 0 0
Notaci´ on: E=Modulo de Young, ν = on: = Coeficiente de Poisson h=Espesor de la placa, I =Momento =Momento de inercia
2.1. 2.1.
Elem Elemen ento to de barr barra a
Consideremos el elemento de barra de Figura 2.2: Figura 2.2: ��
� � ���
�
�
� �
Figura 2.2: Elemento de barra L: largo del elemento A: area de la secci´on on E: Modulo de Young uˆi : desplazamiento en nodo i u(x): desplazamientos
ν (1−ν ) ν (1−ν )
0 1
0 0 0 1−2ν 2(1−ν )
0 0
0 0 0 0 1−2ν 2(1−ν )
0
0 0 0 0 0 1−2ν 2(1−ν )
ELEMENTO DE BARRA
24
Para el caso de una barra se definen dos funciones de forma: N i (x) = 1
− Lx
x L
N j (x) =
(2.6)
(2.7)
Los desplazamientos se pueden escribir como: u(x) =
Por lo tanto, N = B
=
N i
N j
1
− x/L
d N = dx
uˆi uˆj
= N u ˆ
x/L
(2.8)
−
(2.9)
1/L 1/L
(2.10)
Por ultimo las matrices de rigidez y masa viene dada por:
− − − − L
K l =
1/L E 1/L
x=0
L
M l =
ρ
x=0
2.1.1. 2.1.1.
1
x/L x/L
EA 1/L 1/L Adx = L
1
1 1
−
− 1 1
L/3 L/6 x/L x/L Adx = ρA L/6 L/3
(2.11)
(2.12)
Coorde Coordenad nadas as locales locales y global globales es
Consideremos una barra como la mostrada en la Figura 2.3. 2.3. Para determinar las matrices de rigidez y masa en funci´ on on de las coordenadas globales, es necesario definir las matrices de rotaci´on. on. �
� �
��
� �� �� α
��
Figura 2.3: Sistema de coordenadas locales y globales
ELEMENTO DE VIGA
25
De la Figura 2.3 2.3 se se tiene que:
u1 u2
=
cos α sin α 0 0 0 0 cos α sin α
ql
R
x1 y1 x2 y2
(2.13)
qg
Dado que la energ´ energ´ıa es invarian invariante te respecto del sistema de referencia: 1 T q K l q l 2 l
=
1 T q K g q g 2 g
(2.14)
1 T ˙q M l ˙q l 2 l
=
1 T ˙q M g ˙q g 2 g
(2.15)
De donde se obtiene que las matrices con respecto al sistema global vienen dadas por: K g
= RT K l R
(2.16)
M g
= RT M l R
(2.17)
2.2. 2.2.
Elem Elemen ento to de viga viga
Consideremos elemento de viga como el mostrado en la Figura 2.4 2.4.. � � �� θ �
��� θ� ��� �
� �
Figura 2.4: Elemento de viga L: Largo I: Momento de inercia de la secci´on on E: Modulo de Young
�
ELEMENTO DE VIGA
26
v(x): Desplazamiento vertical θ(x): Rotaci´o n en el eje z
Para el caso de un elemento de viga se definen las cuatro funciones de forma mostradas en la Figura 2.5.
� ��
�
−
�
� ��
� �
�
�
Figura 2.5: Funciones de forma para un elemento de viga Por lo tanto, se tiene que: N (x) = B (x) =
− 1
2 3 Lx 2
3 + 2 Lx 3
d N (x) = d2 x
−
x 6 L2
−
+
2 2 xL
12x L3
+
x3 L2
2 3 Lx 2
− L4 + L6x 2
− 6 L2
x3 L3
−
x2 L
+
x3 L2
(2.18)
− 12x − L2 + L6x (2.19) L 3
2
ELEMENTO DE VIGA
27
De donde se obtiene finalmente que: vi K l
=
M l
2.2.1. 2.2.1.
=
E I L3
ρA 420
θi
−
vj
θj
6L 4L2 6L 2L2
−12 6L2 −6L 2L − 12 −6L −6L 4L2 156 22L 54 −13L 22L 4L2 13L −3L2 54 13L 156 −22L −13L −3L2 −22L 4L2 12 6L 12 6L
(2.20)
(2.21)
Coorde Coordenad nadas as locales locales y global globales es
Consideremos una viga como la mostrada en la Figura 2.6. �
�
� θ
��
� �� �� α
��
Figura 2.6: Sistema de coordenadas locales y globales De la Figura 2.6 2.6 se se tiene que:
v1 θ1 v2 θ2
−
=
sin α cos α 0 0 0 0 0 0
0 1 0 0
−
R
0 0 0 0 sin α cos α 0 0
0 0 0 1
x1 y1 θ1 x2 y2 θ2
(2.22)
ELEMENTO DE VIGA 3D
28
Las matrices con respecto al sistema global vienen dadas por: K g
= RT K l R
(2.23)
M g
= RT M l R
(2.24)
2.3. 2.3.
Elem Elemen ento to de viga viga 3D
Un elemento de viga 3D combina un elemento de viga, barra e incorpora adem´as torsi´ on. on. La Figura 2.7 2.7 muestra muestra un esquema, cada nodo tiene 6 grados de libertad. � � � ��
�� �
�
��
�� �� ��� ��� ��
�
�� �
��
�
�
��
�
Figura 2.7: Elemento de viga 3D En este caso, las matrices de rigidez y masa del elemento de viga 3D en el sistema coordenadas local vienen dadas por:
ELEMENTO DE VIGA 3D
u1
K l =
EI L3
M l = ρAL
EA L
29
v1
0
0 0 0 0 0
12EI 12EI z L3
−EA
0
L
0 0 0 0 0
1 3
0 0 0 0 0 1 6
0 0 0 0 0
0 0 0 6EI y L2 −12EI 12EI z
L3
w1
0 0 12EI 12EI y L3
0
−6EI y
L2
0 0 0
0 0 0
−12EI 12EI y
6EI z L2
0
0 13 35
0 0 0
11L 11L 210
0
9 70
0 0 0
−13L 13L
210
L3
0
−6EI y
L2
0 0 13 35
0 −11L 11L
210
0 0 0 9 70
0 13L 13L 210
0
θu1
θv 1
0 0 0
0 0 2
L
0
4EI y L
0 0 0 0 0
0 0 0
−GI x
L
0 0 0 0 0 2
r 3
0 0 0 0 0 2
r 6
0 0
0
−6EI y
GI x L
θw 1 6EI y L2
0 0 0 4EI z L
0
−6EI z
6EI y L2
0
2EI y L2
0
0 0
L2
0 0 0 2EI z L2
0
−11L 11L
210
0 2
L 105
11L 11L 210
0 0 0
0 0 0
L2 105
−13
0 0 0
0
420
0 −L
2
140
0
13 420
−L2
140
u2 −EA
L
0 0 0 0 0 EA L
0 0 0 0 0
1 6
0 0 0 0 0 1 3
0 0 0 0 0
v2
0 −12EI 12EI z
L3
w2
θu2
0 0
0 0 0
−12EI 12EI y
−6EI z
0 0 0
L2
0
12EI 12EI z L3
3
L
0
0 0 0 −6EI z
0
0 9 70
0 0 0
13 420
0
13 35
0 0 0
−11L 11L
210
−GI x
L
6EI y L2
12EI 12EI y L3
L2
0 0 0
0
6EI y L2
0 0 9 70
0 −13
420
0 0 0 13 35
0 11L 11L 210
0
0 0 0 0 0 GI x L
0 0 0 0 0 2
r 6
0 0 0 0 0 2
r 3
0 0
θv2
0 0 −6EI y
L2
0
2EI y L2
0 0 0 6EI y L2
0
4EI y L
0
0 0 13L 13L 210
0 −L2
140
0 0 0 11L 11L 210
0 L2 105
0
θw2
0 6EI z L2
0 0 0 2EI z L2
0
−6EI z
L2
0 0 0 4EI z L
0 −13L 13L
210
0 0 0
−L2
140
0
−11L 11L
210
0 0 0
L2 105
ELEMENTO DE VIGA 3D
2.3.1.
30
Rotaci´ on
A continuaci´on se construir´a el operador de rotaci´on, que permite expresar las matrices con respecto a un sistema global de coordenadas. Sean,
p1
=
p2
=
x1 y1 z1 x2 y2 z2
las posiciones de los nodos del elemento en el sistema global de coordenadas.
−→
La direcci´on del eje neutral del elemento, ex se define por:
− p1 = p2 − p1 (2.25) L − p1 → → ey y − ez , es necesario definir un tercer punto p3 que junto a p1 y p2 Para definir − −→ ex =
p2 p2
definen el plano Oxz del elemento. Definiendo p3 se tiene que:
−→ ey
=
−→ ez
=
( p3 ( p3
− p1) × ( p2 − p1) − p1) × ( p2 − p1) → −ex × −→ ey
(2.26) (2.27)
× es el producto cruz entre dos vectores. Ahora se puede definir el operador de rotaci´ on que relaciona: xg = Rx l
(2.28)
donde xg es la posici´on de los ejes globales y xl es la posici´on de los ejes locales:
→− −→−→ R =
ex ey ey
T
(2.29)
Se tiene que: u v w
θu θv θw
= R
x y z
= R
θx θy θz
(2.30)
(2.31)
ENSAMBLE
31
De donde, se puede expresar el vector de desplazamientos del elemento en los ejes locales y globales como:
x1 y1 z1 θx1 θy1 θz1 x2 y2 z2 θx2 θy2 θz2
R 0 0 0
=
0 0 0 R 0 0 0 R 0 0 0 R T
u1 v1 w1 θu1 θv1 θw1 u2 v2 w2 θu2 θv2 θw2
(2.32)
Las matrices de rigidez y masa del elemento, expresadas en el sistema global de coordenadas quedan, K g
= T T K l T
(2.33)
M g
= T T M l T
(2.34)
2.4.
Ensamble
La matrices de masa y rigidez de un sistema es un ensamble de las matrices de cada elemento. Para realizar el ensamble, necesitamos saber que grados de libertad est´ an asociados que elemento en la estructura. Por ejemplo, consideremos el caso de una viga simple apoyada en ambos extremos que se muestra en la Figura 2.8. �
�
�
�
�
Figura 2.8: Esta viga esta compuesta por 5 elementos de viga, cada uno con 4 grados de libertad. La estructura en total tiene 5 elementos, 6 nodos y 12 grados de libertad (2 grados por nodo). La matriz de rigidez de la viga completa es un ensamble de las 5 matrices de rigidez, ubicadas en los grados de libertad correspondiente. Esto es, la matriz K 1 del primer elemento abarca los grados 1 a 4, la matriz K 2 los
ENSAMBLE
32
elementos 3 a 6, etc. La Figura 2.9 ilustra este ensamble. La matriz de masa se construye de manera equivalente.
��
� ��
��
�� �� �
��
Figura 2.9: Ensamble matriz de rigidez La condici´ on de borde, es que los grados de libertad de desplazamiento vertical en los nodos 1 y 6 est´an fijos. Para fijar grados de libertad, se eliminan las filas y columnas asociadas a estos grados en las matrices ensambladas. El siguiente c´odigo resuelve el problema en Matlab: %Propiedades de la viga E=3.2e10; %Pa modulo de elasticidad rho=2500; %densidad Kg/m3 A=0.05; %mˆ2 a ´rea Iz=1.66e-4; %mˆ4 inercia secci´ on %Definici´ o n de los elementos l=6; %largo total en metros n=5; %n´ umero de elementos le=l/n; %m largo de un elemento ne=(n+1)*2; %grados de libertad %Nodos con restricciones de libertad fix=[1 11]; %Nodos con restricciones %matriz de masa de un elemento Me=rho*A*le/420*[156 22*le 54 22*le 4*leˆ2 13*le 54 13*le 156 -13*le -3*leˆ2 -22*le
-13*le; -3*leˆ2; -22*le; 4*leˆ2];
ENSAMBLE
%matriz de rigidez de un elemento Ke=E*Iz/leˆ3*[12 6*le -12 6*le; 6*le 4*leˆ2 -6*le 2*leˆ2; -12 -6*le 12 -6*le; 6*le 2*leˆ2 -6*le 4*leˆ2]; %matrices iniciales K=zeros(ne,ne); M=zeros(ne,ne); %locel: grados de libertad asociados a cada elemento for i=1:n j=2*(i-1)+1; locel(i,1:4)=(j):(j+3); end %ensamble for i=1:n K(locel(i,:),locel(i,:))=K(locel(i,:),locel(i,:)) + Ke; M(locel(i,:),locel(i,:))=M(locel(i,:),locel(i,:)) + Me; end %fijar grados de libertad free=setdiff(1:ne,fix); K=K(free,free); M=M(free,free); %frecuencias naturales y modos [Phif,W]=eig(K,M); w=sqrt(diag(W))/2/pi; Phi=zeros(ne,ne-length(fix)); Phi(free,:)=Phif; %Gr´ a fica los primeros tres modos for i=1:3 figure set(gcf,’Position’,[100 100 500 400]) set(gca,’FontSize’,18) set(gca,’Units’,’centimeters’) set(gca,’Position’,[1.0 1.5 11 8]) set(gcf, ’PaperUnits’,’inches’); plot(0:ne/2-1,Phi(1:2:ne,i),’-b’,’LineWidth’,3)
33
ENSAMBLE
34
xlim([0 ne/2-1]) set(gca,’YTick’,[]); set(gca,’XTick’,[0:5]); set(gca,’XTicklabel’,[1:6]); xlabel(’Nodo’,’FontSize’,16) title([’Modo ’ num2str(i) ’, \omega= ’ num2str(w(i),3) ’Hz’],’FontSize’,16) grid on; end
La Figura 2.10 muestra los primeros tres modos de la viga. Modo 1,
1
2
= 9Hz
Modo 3,
ω
3
4 Nodo
5
6
1
2
= 81.6Hz
Modo 2,
ω
3
4
5
6
Nodo
Figura 2.10: Modos de la viga
1
2
= 36Hz
ω
3
4 Nodo
5
6
Cap´ıtulo 3
Integraci´ on num´ erica de las ecuaciones de movimiento Para un problema general, la ecuaci´on de movimiento viene dada por: M x ¨ + C x˙ + Kx = F
(3.1)
Donde M, C y K son las matrices de masa, amortiguaci´on y rigidez respectivamente. F es un vector de fuerzas externas, y x ¨, x, ˙ x son vectores de aceleraci´on, velocidad y desplazamiento respectivamente. Matem´aticamente, la ecuaci´on 3.1 representa un sistema de ecuaciones diferenciales lineales de segundo orden, las que en un principio, se pueden resolver por procedimientos est´ andar para ecuaciones diferenciales ordinarias. Sin embargo, estos procedimientos pueden ser muy costosos (en t´erminos de tiempo y cantidad de c´alculos necesarios) si las matrices son grandes. En el an´ alisis de elementos finitos, se utilizan algunos algoritmos eficientes para resolver las ecuaciones de movimiento, los que veremos en las secciones siguientes. Los m´ etodos que veremos se dividen en dos categor´ıas: integraci´on directa y superposici´on modal. Aunque ambas t´ecnicas puedan parecer muy distintas a primera vista, est´an fuertemente relacionadas, y la selecci´on de una u otra depende de su eficiencia para calcular la soluci´ on num´erica.
3.1.
M´ etodos de integraci´ on directa
En la integraci´ on directa las ecuaciones de movimiento son integradas por un procedimiento num´erico por pasos, el t´ ermino “directa” significa que no es necesario 35
´ ´ METODOS DE INTEGRACION DIRECTA
36
hacer transformaciones de las ecuaciones previo a la integraci´on num´erica. En los procedimientos de integraci´ on temporal, se asume que los vectores de aceleraci´on, velocidad y desplazamiento en el tiempo t = 0 ,x x(0), son ¨(0), x(0), ˙ conocidos. La soluci´on se divide en pasos de tiempo ∆ t.
3.1.1.
El m´ etodo de las diferencias centrales
El m´etodo de las diferencias centrales utiliza la aproximaci´on: x ¨(t ) =
1 (x(t ∆t2
− ∆t) − 2x(t) + x(t + ∆t))
(3.2)
El error en la expansi´on de 3.2 es del orden de ∆ t2 , para tener el mismo orden de error en la expansi´on de la velocidad se utiliza: x˙ (t) =
1 ( x(t + ∆t) 2∆t
− x(t − ∆t))
(3.3)
El desplazamiento en el tiempo t + ∆t se obtiene al considerar la ecuaci´ on de movimiento en t: M x ¨(t) + C x˙ (t) + Kx(t) = F (t)
(3.4)
Substituyendo, se obtiene:
1 1 M C x(t+∆t) = F (t) + ∆t2 2∆t
−
K
−
−
2 M x(t) ∆t2
1 M ∆t2
−
1 C x(t ∆t) 2∆t
−
(3.5) Se debe notar que para calcular x(t + ∆t), son necesarios x(t) y x(t ∆t). Por lo tanto, se debe utilizar un procedimiento especial de inicializaci´on. Dado que se x(0), las ecuaciones 3.2 y 3.3 se pueden utilizar para calcular conocen x ¨(0), x(0), ˙ x( ∆t):
−
−
x( ∆t) = x (0)
−
−
∆t2 x ∆tx˙ (0) + ¨(0) 2
(3.6)
En general, para llegar a soluciones estables se requiere un paso de tiempo relativamente peque˜no. De hecho, en el m´ etodo de las diferencias centrales es paso de tiempo debe ser menor que un valor cr´ıtico ∆tcr . El paso de tiempo cr´ıtico se determina como: ∆tcr =
T n π
(3.7)
´ ´ METODOS DE INTEGRACION DIRECTA
37
donde T n es el menor periodo en el sistema, entre el inverso de la mayor frecuencia natural y la frecuencia de las fuerzas externas. La Tabla 3.1 ilustra el procedimiento paso a paso para resolver las ecuaciones de movimiento mediante el m´ etodo de las diferencias centrales.
Tabla 3.1: Procedimiento paso a paso, m´ etodo de diferencias centrales C´ alculos iniciales:
1. Crear matrices de rigidez K, masa M y amortiguaci´on C. 2. Inicializar x ¨ (0), x˙ (0), x(0). 3. Seleccionar paso de tiempo ∆t, ∆t a0 =
1 ; ∆t2
a1 =
4. Calcular x( ∆t) = x (0)
−
≤ ∆tcr . Calcular constantes de integraci´on.
1 ; 2∆t
a2 = 2a0 ;
a3 =
1 a2
− ∆tx˙ (0) + a3 x¨(0)
ˆ = a 0 M + a1 C . 5. Formar matriz de masa efectiva M Para cada paso de tiempo:
1. Calcular fuerzas efectivas en t: ˆ (t) = F (t) F
− (K − a2M )x(t) − (a0 M − a1C )x(t − ∆t)
2. Determinar desplazamientos en t + ∆t: ˆ −1 ˆ x(t + ∆t) = M F (t) 3. Si se requiere, evaluar las aceleraciones y velocidades en t: x ¨(t) = a 0 (x(t ∆t) x˙ (t) = a 1 (x(t + ∆t)
−
3.1.2.
− 2x(t) + x(t + ∆t)) − x(t − ∆t))
El m´ etodo de Wilson
θ
En el m´etodo de Wilson θ se asume una variaci´on lineal de la aceleraci´on entre un paso de tiempo t y t + θ ∆t, donde θ 1, como se ilustra en la Figura 3.1.
≥
´ ´ METODOS DE INTEGRACION DIRECTA
38
Si θ = 1, el m´ etodo se reduce a un esquema de aceleraci´ on lineal. Pero se puede demostrar que para asegurar estabilidad incondicional es necesario utilizar θ 1,37, usualmente se emplea θ = 1,4
≥
�
� � θ ∆�
� � ∆�
��
��
∆�
��
θ∆�
Figura 3.1: Supuesto de aceleraci´on lineal, m´etodo de Wilson θ Denotemos τ como el incremento en el tiempo, donde 0 intervalo t a t + θ∆t se tiene que: x ¨(t + τ ) = x ¨(t) +
τ x(t + θ∆t) (¨ θ∆t
≤ τ ≤ θ∆t. Entonces en el
− x¨(t))
(3.8)
Integrando la ecuaci´on anterior se obtiene que, τ 2 x˙ (t + τ ) = x˙ (t) + τ x x(t + θ ∆t) ¨(t) + (¨ 2θ∆t
− x¨(t))
(3.9)
y τ 3 1 2 x(t + τ ) = x (t) + τ x˙ (t) + τ x x(t + θ∆t) ¨(t) + (¨ 2 6θ∆t
− x¨(t))
(3.10)
Usando 3.9 y 3.10, se tiene en el tiempo t + θ ∆t, x˙ (t + θ∆t) = x˙ (t) +
θ∆t x(t + θ∆t) + x (¨ ¨(t)) 2
θ2 ∆t2 x(t + θ∆t) = x(t) + θ∆tx˙ (t) + x(t + θ∆t) + 2¨ x(t)) (¨ 6
(3.11) (3.12)
De donde se puede expresar x ¨(t + θ ∆t) y x˙ (t + θ∆t) en t´erminos de x(t + θ ∆t): x ¨(t + θ∆t) =
6 (x(t + θ∆t) θ 2 ∆t2
x˙ (t + θ∆t) =
3 ( x(t + θ ∆t) θ ∆t
− x(t)) − θ∆6 t ˙x(t) − 2¨x(t)
− x(t)) − 2x˙ (t) − θ∆2 t x¨(t)
(3.13) (3.14)
´ ´ METODOS DE INTEGRACION DIRECTA
39
Para obtener la soluci´on de los desplazamientos, velocidades y aceleraciones en el tiempo t + ∆t, se considera la ecuaci´ on de movimiento en el tiempo t + θ∆t. Sin embargo, dado que la aceleraci´on varia linealmente, se utiliza un vector de fuerzas extrapolado linealmente, i e., la ecuaci´on utilizada es: ¯ (t + θ∆t) M x ¨(t + θ ∆t) + C x˙ (t + θ∆t) + Kx(t + θ∆t) = F
(3.15)
¯ (t + θ ∆t) = F (t) + θ (F (t + ∆t) F
(3.16)
− F (t))
Sustituyendo las ecuaciones 3.13 y 3.14 en 3.15, se obtiene una ecuaci´on en donde se puede despejar x(t + θ∆t). Luego, sustituyendo x(t + θ∆t) en 3.13, se obtiene x ¨(t + θ∆t) la que es usada en las ecuaciones 3.8 a 3.10, todas evaluadas en τ = ∆t para calcular x(t + ∆t), x(t ˙ + ∆t) y x ¨(t + ∆t). El algoritmo completo se da en la Tabla 3.2
Se debe notar que con este m´ etodo no se requiere inicializaciones especiales, dado que los desplazamientos, velocidades y aceleraciones en el tiempo t + ∆t se expresan en t´erminos de las mismas cantidades en el paso de tiempo anterior t. Tabla 3.2: Procedimiento paso a paso, m´etodo de Wilson θ C´ alculos iniciales:
1. Crear matrices de rigidez K, masa M y amortiguaci´on C. 2. Inicializar x ¨ (0), x˙ (0), x(0). 3. Seleccionar paso de tiempo ∆t y θ. Calcular constantes de integraci´on. a0 =
6 ; (θ∆t)2
−a2 ; a5 = θ
a1 =
a6 = 1
3 ; θ ∆t
−
3 ; θ
a2 = 2a1 ;
∆t a7 = ; 2
a3 =
θ∆t ; 2
a4 =
a0 θ
∆t2 a8 = 6
ˆ = K + a0 M + a1 C . 4. Formar matriz de rigidez efectiva K Para cada paso de tiempo:
1. Calcular fuerzas efectivas en t + θ ∆t: ˆ (t + θ∆t) = F (t) + θ(F (t + ∆t) F (t)) + M (a0 x(t) + a2 ˙x(t) + 2¨ F x(t)) +C (a1 x(t) + 2x˙ (t) + a3 x ¨(t))
−
´ ´ METODOS DE INTEGRACION DIRECTA
40
2. Determinar desplazamientos en t + θ∆t: ˆ −1 ˆ x(t + θ ∆t) = K F (t + θ∆t) 3. Evaluar aceleraciones, velocidades y desplazamientos en t + ∆t: x ¨(t + ∆t) = a 4 (x(t + θ∆t) x(t)) + a5 ˙x(t) + a6 x¨(t) x˙ (t + ∆t) = x˙ (t) + a7 (¨ x(t + ∆t) + x ¨(t)) x(t + ∆t) = x (t) + ∆tx˙ (t) + a8 (¨ x(t + ∆t) + 2¨ x(t))
−
3.1.3.
El m´ etodo de Newmark
El m´etodo de Newmark tambi´en se puede relacionar con un m´etodo de aceleraciones lineales. Utiliza los siguientes supuestos: x˙ (t + ∆t) = x˙ (t) + [(1
− δ )x¨(t) + δ x¨(t + ∆t)] ∆t
x(t + ∆t) = x(t) + x˙ (t)∆t +
− 1 2
(3.17)
α x ¨(t) + αx ¨ (t + ∆t) ∆t2
(3.18)
donde α y δ son par´ametros que determinan la precisi´on y estabilidad del m´etodo de integraci´ on. Cuando δ = 12 y α = 16 , las ecuaciones 3.17 y 3.18 corresponden al m´ etodo de aceleraciones lineales (que tambi´en se obtiene con θ = 1 e n el m´etodo de Wilson θ). Newmark propuso originalmente como condici´ on de estabilidad incondicional el esquema de promedio constante de la aceleraci´ on 1 1 (tambi´en denominado regla trapezoidal), en donde δ = 2 y α = 4 (ver Figura 3.2).
� � ∆� ��
� �
�
� � � ∆�
� � ∆�
Figura 3.2: Esquema de promedio de aceleraci´on constante de Newmark
´ ´ METODOS DE INTEGRACION DIRECTA
41
Adicionalmente a las ecuaciones 3.17 y 3.18, para la soluci´on de los desplazamientos, velocidades y aceleraciones en t + ∆t, se utilizan las ecuaciones de movimiento en t + ∆t: M x ¨(t + ∆t) + C x˙ (t + ∆t) + Kx(t + ∆t) = F (t + ∆t)
(3.19)
De la ecuaci´on 3.18 se puede despejar x ¨(t + ∆t) en funci´on de x(t + ∆t) y despu´es se puede reemplazar en la ecuaci´on 3.17. De esta manera se obtienen x ¨(t + ∆t) y x˙ (t + ∆t) en funci´ on del desplazamiento x(t + ∆t). Sustituyendo estas dos expresiones en la ecuaci´on 3.19, se puede despejar x(t + ∆t), de donde tambi´en se tienen x ¨(t + ∆t) y x˙ (t + ∆t). El algoritmo completo utilizando el algoritmo de Newmark se ilustra en la tabla 3.3. Se debe notar la cercana relaci´on entre el m´etodo de Wilson θ y el m´etodo de Newmark.
Tabla 3.3: Procedimiento paso a paso, m´etodo de Newmark C´ alculos iniciales:
1. Crear matrices de rigidez K, masa M y amortiguaci´on C. 2. Inicializar x ¨ (0), x˙ (0), x(0). 3. Seleccionar ∆t y los par´ametros α y δ . Calcular constantes de integraci´on. δ
≥ 0,50;
a0 =
1 ; α∆t2
∆t a5 = 2
a1 =
− δ α
2 ;
α
δ ; α∆t
≥ 0,25(0,5 + δ )2 a2 =
a6 = ∆t(1
1 ; α∆t
− δ );
a3 =
1 2α
a7 = δ ∆t
ˆ = K + a0 M + a1 C . 4. Formar matriz de rigidez efectiva K Para cada paso de tiempo:
1. Calcular fuerzas efectivas en t + θ ∆t: ˆ (t + ∆t) = F (t + ∆t) + M (a0 x(t) + a2 ˙x(t) + a3 x F ¨(t)) +C (a1 x(t) + a4 ˙x(t) + a5 x¨(t)) 2. Determinar desplazamientos en t + ∆t:
− 1;
a4 =
δ α
− 1
´ ´ METODOS DE INTEGRACION DIRECTA
42
ˆ −1 ˆ x(t + ∆t) = K F (t + ∆t) 3. Evaluar aceleraciones y velocidades en t + ∆t: x ¨(t + ∆t) = a 0 (x(t + ∆t) x(t)) a2 ˙x(t) x˙ (t + ∆t) = x˙ (t) + a6 x ¨(t) + a7 x ¨(t + ∆t)
−
3.1.4.
−
− a3x¨(t)
Ejemplo
Consideremos un sistema simple, cuyas ecuaciones de movimiento son las siguientes:
2 0 0 1
x ¨1 x ¨2
+
0,3 0,1
−
−0,1 0,2
x˙ 1 x˙ 2
+
6 2
−
− 2 4
x1 x2
=
0 10
(3.20)
Las frecuencias naturales del este sistema son 0.23Hz y 0.36Hz, por lo tanto el menor periodo del sistema es 2.78s. El paso de tiempo cr´ıtico es ∆tcr = 2,78 π = 0,88s. Utilizando un paso de tiempo 0.18s se obtienen, de acuerdo a los algoritmos descritos en las secciones anteriores, los resultados mostrados en la Figura 3.3. Los tres algoritmos de integraci´ on dan resultados similares, aunque el m´etodo de las diferencias centrales es que m´etodo que m´as se acerca a la soluci´on exacta. Se asumieron desplazamientos y velocidades iniciales nulas. Dif. Centrales
Wilson
Newmark
Exacta
3 2.5 2 1.5 ) t (
1
x
1 0.5 0 -0.5 -1 0
1
2 3 Tiempo (s)
4
Figura 3.3: Desplazamientos obtenidos
5
´ SUPERPOSICION MODAL
3.2.
43
Superposici´ on modal
En los m´ etodos de integraci´on directa se deben realizar muchas operaciones en cada paso de tiempo. Si se debe integrar en un numero alto de pasos de tiempo, puede ser m´as conveniente transformar las ecuaciones de movimiento en una forma que requiera de menos operaciones para cada paso de tiempo. Aprovechando las propiedades ortogonales de los modos propios, se puede definir la transformaci´on: x(t) = Φy (t)
(3.21)
De donde se obtiene la ecuaci´on de movimiento: y + P hiT C Φy˙ + P hiT K Φy = ΦT F ΦT M Φ¨
(3.22)
Las matrices ΦT M Φ, ΦT C Φ y ΦT K Φ son matrices diagonales, por lo que el sistema de ecuaciones 3.22 es un sistema de ecuaciones desacoplado. En consecuencia, se tienen n ecuaciones individuales de la forma:
mi ¨ yi + ci y˙i + ki yi = r i ri = φ T i F (t)
i = 1, 2, 3, . . . , n
(3.23)
La soluci´ on para cada una de las ecuaciones se puede obtener utilizando los algoritmos de integraci´on antes descritos o de acuerdo al procedimiento descrito en la secci´ on 1.1.3. Una vez calculada la soluci´on para cada una de las ecuaciones. La soluci´on general se obtiene al superponer las soluciones modales: n
x(t) =
φi yi (t)
(3.24)
i=1
3.2.1.
Ejemplo
Consideremos el sistema de ecuaciones 3.20. Para este sistema los modos normales son:
Φ=
−
0,5774 0,5774
−
−0,4082 0,8165
(3.25)
Reemplazando Φ en ecuaci´on 3.22, se obtiene:
1 0 0 1
y¨1 y¨2
+
0,1 0,0 0,0 0,25
− y˙1 y˙2
+
2 0 0 5
y1 y2
=
5,7735 8,1650
(3.26)
´ SUPERPOSICION MODAL
44
Este sistema de ecuaciones se puede resolver utilizando alg´ un algoritmo de integraci´ on. Utilizando diferencias centrales con el mismo paso de tiempo del ejemplo anterior, se obtiene la soluci´on para y(t). La soluci´on, x(t), viene dada por:
x(t) =
− 0,5774 0,5774
−
y1 (t) +
− 0,4082 0,8165
y2 (t)
(3.27)
La soluci´on obtenida se muestra en la Figura 3.4. Como es de esperarse, la soluci´on es la misma que si se utiliza diferencias modales con o sin transformaci´on modal de las ecuaciones. La ventaja de la transformaci´on modal es que se reduce el tiempo requerido para calcular la soluci´on a cada paso de tiempo. 3 2.5 2 1.5 ) t (
1
x
1 0.5 0 -0.5
Superposicion Modal Dif. Centrales
-1 0
1
2 3 Tiempo (s)
4
Figura 3.4: Desplazamientos obtenidos
5
Cap´ıtulo 4
Procesamiento de se˜ nales El procesamiento de digital se˜nales es una herramienta muy importante en el an´alisis de sistemas. En esta secci´on se dar´a un resumen de los conceptos b´asicos asociados al procesamiento de se˜nales. La se˜ nales, en general, se pueden clasificar de acuerdo a la tabla 4.1. Las se˜ nales estacionarias son aquellas cuyas propiedades promedio no var´ıan con el tiempo, pueden ser o deterministicas o aleatorias. El grupo m´ as importante de se˜nales deterministicas son las se˜ nales peri´odicas. Una funci´on pseudo-aleatoria es una se˜ nal aleatoria que se repite con un cierto periodo. Las se˜ nales no estacionarias se pueden dividir en continuas o transientes. Las se˜nales transientes se pueden definir como se˜nales que comienzan y terminar en cero en el periodo de observaci´on.
Se˜ nales Estacionarias Determin´ıstica Aleatoria Peri´ odica Cuasi-peri´ odica Pseudo-aleatoria
No estacionarias Continua Transiente
Tabla 4.1: Tipos de se˜nales Dado que el objetivo del procesamiento de se˜ nales es extraer el m´ aximo de informaci´ on de las se˜nales, es en general beneficioso estudiar las se˜nales en distintos dominios. Las se˜ nales medidas son, claramente, funciones en el dominio del tiempo. Para estudiar su contenido en frecuencias, es m´as f´acil examinar las se˜nales en el
45
LA TRANSFORMADA DE FOURIER
46
dominio de frecuencias. La transformada (inversa) de Fourier, permite transformar de manera sensilla una se˜nal en el tiempo a una se˜nal en frecuencia y viceversa. En consecuencia, la transformada de Fourier es uno de los temas m´as importante en procesamiento de se˜ nales.
4.1.
La transformada de Fourier
J.B. Fourier prob´o que una funci´on peri´odica en el tiempo se puede representar como una suma de componentes sinusoidales a frecuencias equiespaciadas: +∞
g (t) =
G(k ∆f )ej2πk∆ft
(4.1)
−∞
Los coeficientes de fourier viene dados por: 1 G(k ∆f ) = T con:
+T /2
g (t)e−j2πk∆f t dt
(4.2)
−T /2
t: tiempo k: entero que cuenta los pasos en frecuencia ∆f : espaciado de frecuencias o resoluci´on de frecuencia: ∆ f = 1/T j =
√ −1
T: periodo de tiempo: T = 1/∆f El set de valores G(k∆f ) se denomina espectro de la funci´on g(t). En general el espectro posee valores complejos.
Al utilizar computadores digitales, es necesario adquirir la se˜ nal continua en intervalos de tiempo. Esto significa que la se˜nal continua es representada por una se˜ nal discreta con valores a tiempos equidistantes. Considerando esto la transformada de Fourier queda como: 1 g (n∆t) = f s
+f s /2
G(f )ej2πf n∆t df
(4.3)
−f s /2
+∞
G(f ) =
−∞
con:
g (n∆t)e−j2πf n∆t
(4.4)
LA TRANSFORMADA DE FOURIER
47
n: entero contando el numero de pasos de tiempo
∆t: intervalo de muestreo: ∆t = 1/f s f s : frecuencia de muestreo: f s = 1/∆t En condiciones reales de medici´on experimental es imposible medir la se˜nal temporal hasta un tiempo infinito. Una parte de la se˜nal debe ser seleccionada. Se asume que la se˜ nal capturada se repite con un periodo T, entregando una funci´ on peri´odica. Combinando la hip´otesis de periodicidad con un muestreo temporal de la se˜nal, se obtiene la definici´on de la transformada discreta de Fourier:
1 g (n∆t) = f s
N s −1
G(k ∆f )ej2πnk/N s
(4.5)
k=0
1 G(k ∆f ) = N s
N s −1
g(n∆t)e−j2πnk/N s
(4.6)
n=0
con: N s : n´umero de datos: T = N s ∆t y f s = N s ∆f La evaluaci´on directa de la transformada discreta de Fourier requiere N s2 operaciones. Con la transformada r´apida de Fourier (FFT) se reduce el n´umero de operaciones a N s log2 N s . La transformada r´apida de Fourier es el n´ucleo de todos los procesadores de se˜ nal modernos. Consideremos como ejemplo una se˜nal sinusoidal, con periodo T=0.02s y amplitud A=1: x(t) = sen(2π 50t)
(4.7)
La se˜ nal es muestreada a una frecuencia de 10kHz y se considera el periodo de tiempo 0-0.2s. En la Figura 4.1 se ilustra el resultado de la transformada r´apida de Fourier de la se˜ nal.
ERRORES Y VENTANAS
48
1.5
1
0.5 d u t i l p m A
d u t i l p m A
0
1
0.5
-0.5
-1 0
0.05
0.1 Tiempo (s)
0.15
0.2
0 0
20
40 60 Frecuencia (Hz)
80
100
Figura 4.1: Espectro de una se˜nal peri´odica
4.1.1.
Algunos par´ ametros importantes
Aqu´ı se resumen algunos de los par´ametros mencionados en los p´arrafos anteriores: T : Periodo de tiempo en donde se adquieren N s muestras equiespaciadas de la se˜ nal a ser analizada. En la mayor´ıa de los analizadores de Fourier N s est´a restringidos a potencias de 2 (por ejemplo, 1024). Relaciones: T = N s ∆t = 1/∆f f s : Frecuencia de muestreo, frecuencia a la cual se adquieren y digitalizan los datos. Relaciones: f s = 1/∆t = N s ∆f ∆t : Intervalo de muestreo: intervalo de tiempo al cual la se˜nal es muestreada. Relaciones: ∆t = T /N s = 1/f s ∆f : Espaciado de frecuencias en el espectro. Relaciones: ∆ f = 1/T = f s /N s . En consecuencia, mediciones con un espaciado de frecuencias peque˜no (i.e. una resoluci´on en frecuencias alta) conllevan tiempos de medici´on mayores. f max : Frecuencia m´ axima: frecuencia mayor contenida o permitida en la se˜ nal temporal. Seg´un el teorema de Shannon: f max f s /2
≤
N s : N´ umero de muestras en el periodo de tiempo T: N s = T /∆t = f s /∆f
4.2.
Errores y ventanas
Durante el proceso de an´alisis digital de una se˜nal pueden ocurrir muchos errores. Errores t´ıpicos son sobrecargas, ruido digital, errores de cuantificaci´on, limitaciones del rango din´amico. Sin embargo, los dos errores principales son aliasing y leakage.
ERRORES Y VENTANAS
4.2.1.
49
Aliasing
El aliasing se produce debido al hecho que la se˜nal temporal debe ser muestreada. Componentes de alta frecuencia en la se˜nal pueden causar errores de amplitud y frecuencia en el espectro. Si la mayor frecuencia contenida en una se˜nal no cumple con el teorema de Shannon: f max f s /2, entonces las frecuencias por sobre f s /2 van a aparecer como frecuencias menores a f s /2. La Figura 4.2 muestra un ejemplo de Aliasing, se muestran tres senos con frecuencias de 1, 4 y 6 Hz muestreados a 5 Hz, a la frecuencia de muestreo los tres senos son id´ enticos.
≤
1 1 Hz 4 Hz 6 Hz
0.5 d u t i l p m A
0
-0.5
-1 0
0.2
0.4 0.6 Tiempo (s)
0.8
1
Figura 4.2: Aliasing: Seno a 1, 4 y 6 Hz, muestreados a 5 Hz (c´ırculos) Aliasing se puede evitar removiendo todos los componentes con frecuencias mayores a f s /2. Esto se puede lograr con una se˜nal de excitaci´on apropiada, aunque se logra generalmente utilizando un filtro pasa bajas. Dado que no existen filtros que remuevan todas las frecuencias altas a cero sin influenciar en las bajas, los filtros se fijan normalmente a un 40 % de f s .
4.2.2.
Leakage
El Leakage se origina debido a que los datos deben ser adquiridos en un periodo de observaci´on finito T . La transformada discreta de Fourier asume entonces que la se˜ nal es peri´odica con periodo T. Si esta condici´on no se cumple, se produce un error de “leakage”. La Figura 4.3 ilustra el espectro obtenido de una se˜nal tipo coseno, cuando la funci´on es peri´odica en T y cuando no lo es. En el segundo caso,
ERRORES Y VENTANAS
50
el espectro discreto no coincide con el real. El error en la hip´otesis de periodicidad produce errores importantes de amplitud y frecuencia. 1.5
d u t i l p m A
0
d u t i l p m A
1
0.5
1
2
0 0
3
2
4 6 Frecuencia (Hz)
8
10
2
4 6 Frecuencia (Hz)
8
10
Tiempo (s)
1.5
d u t i l p m A
0
d u t i l p m A
1
0.5
0.8 1.6 Tiempo (s)
2.4
0 0
Figura 4.3: Hip´otesis de periodicidad, Leakage La u ´nica soluci´on al problema de leakage, es asegurarse que la se˜nal es peri´odica o se observa completamente en el periodo de adquisici´on, lo que en general, es muy dif´ıcil de lograr. En sistemas perfectamente lineales, se puede lograr al excitarlos con una se˜nal peri´odica en el intervalo de tiempo considerado. Aumentar el tiempo de adquisici´ on, i.e. aumentar la resoluci´ on en frecuencias, ayuda a mejorar la periodicidad de la se˜ nal. El uso de ventanas de tiempo tambi´ en ofrecen una soluci´on parcial al problema de leakage.
4.2.3.
Ventanas
El uso de ventanas de tiempo no se puede evitar en el procesamiento digital de una se˜ nal. Al medir una se˜ nal temporal, solo una parte de la se˜nal total es considerada.
ERRORES Y VENTANAS
51
Esto equivale a multiplicar la se˜nal actual con una ventana de tiempo rectangular (Figuras 4.4(a) y 4.5(a)). Sin embargo, una mejor selecci´ on de la ventana puede reducir considerablemente el error debido a leakage. En general, se buscan ventanas que reduzcan las discontinuidades en los extremos de la se˜nal, dado que reducen el error por leakage al forzar la se˜nal a ser peri´odica. La selecci´on de una ventana de tiempo, es siempre un compromiso entre una buena estimaci´on de la amplitud y una buena resoluci´on espectral. La figura 4.4 muestra un conjunto de ventanas de tiempo que son usualmente empleadas. l
i
d u t i l p m A
i
d u t i l p m A
0
d u t i l p m A
0 Tiempo (s)
0 Tiempo (s)
(a) Rectangular
Tiempo (s)
(b) Hanning l
i
d u t i l p m A
(c) Hamming
d u t i l p m A
0
0 Tiempo (s)
(d) Gaussian
Tiempo (s)
(e) Flat-top
Figura 4.4: Diferentes ventanas de tiempo En la Figura 4.5 se muestra el espectro obtenido de una se˜nal peri´odica utilizando las distintas ventanas de tiempo. La se˜nal esta compuesta por un coseno de frecuencia f 1 y amplitud 1, y un coseno de frecuencia f 2 y amplitud 0.01. El coseno de frecuencia f 1 no es peri´odico en el intervalo de tiempo considerado. Se observa que solo algunas ventanas son capaces de separar bien ambas se˜nales. Es importante notar, que dado que la se˜nal no es peri´odica en el periodo considerado, la amplitud en el dominio de frecuencias depende de la ventana y de la ubicaci´on de la se˜nal no-peri´odica con respecto a los puntos de frecuencias muestreadas. En consecuencia, este error no puede ser corregido por un factor general. Tambi´en es importante destacar, que el uso de ventanas no rectangulares, reduce el total de energ´ıa en la se˜nal. Lo que conlleva a una reducci´on en la amplitud mostrada en el dominio de frecuencias. La cantidad de energ´ıa solo depende de la forma de la ventana y por lo tanto, este error puede ser compensado.
ERRORES Y VENTANAS
52
) d u t i l p m A ( g o l
d u t i l p m A
f1 f2
Tiempo (s)
Frecuencia (Hz)
(a) Rectangular
) d u t i l p m A ( g o l
d u t i l p m A
f1 f2
Tiempo (s)
Frecuencia (Hz)
(b) Hanning
) d u t i l p m A ( g o l
d u t i l p m A
f1 f2
Tiempo (s)
Frecuencia (Hz)
(c) Hamming
) d u t i l p m A ( g o l
d u t i l p m A
f1 f2
Tiempo (s)
Frecuencia (Hz)
(d) Gaussian
) d u t i l p m A ( g o l
d u t i l p m A
f1 f2
Tiempo (s)
Frecuencia (Hz)
(e) Flat-top
Figura 4.5: Diferentes ventanas, f (t) = cos (f 1 (2πt )) + 0,01cos(f 2 (2πt ))
ERRORES Y VENTANAS
53
Existen dos ventanas especiales que se utilizan test de impacto. En el caso de test de impacto la se˜ nal de entrada es una se˜nal tipo pulso y la respuesta es una combinaci´on de sinusoides que disminuye en el tiempo (Figura 4.6). Para la respuesta se usa una ventana exponencial w(t) = e−αt . En el caso de estructuras con baja amortiguaci´ on o si el tiempo de adquisici´ on es muy breve, la respuesta no llegar´ a a cero al final del bloque de tiempo, lo que causa discontinuidades y leakage. Al multiplicar esta respuesta con una ventana exponencial da como resultado una se˜nal que es casi zero al final del bloque de tiempo. En el caso de estructuras con alta amortiguaci´ on o con un tiempo de adquisici´ on alto, la se˜nal llega a cero antes del final del bloque de tiempo y el resto de lo medido es b´ asicamente ruido experimental. Al aplicar una ventana exponencial, en este caso, se reduce la contribuci´on del ruido. Dado que la fuerza es de corta duraci´ on, cualquier ruido durante el resto de la se˜ nal no es deseable. Aplicar una ventana que sea igual a la ventana exponencial durante el pulso, suavemente desciende a cero justo despu´es del pulso y permanece igual a cero durante el resto del bloque de tiempo, ofrece una soluci´on al problema del ruido (Figura 4.7).
Figura 4.6: Se˜nal t´ıpica de impacto y respuesta al impacto
Figura 4.7: Ventanas de fuerza y exponencial
FUNCIONES DE RESPUESTA EN FRECUENCIA Y FUNCIONES DE COHERENCIA
54
La combinaci´ on de una venta de fuerza y una ventana exponencial, a˜nade amortiguamiento al sistema. La funciones de respuesta en frecuencia resultantes van a contener el efecto de este amortiguamiento extra. No es f´ acil remover este amortiguamiento de las funciones de respuesta en frecuencia, sin embargo, puede ser considerado en la etapa de estimaci´on de par´ametros.
4.3.
Funciones de respuesta en frecuencia y funciones de coherencia
Sea F (f ) el espectro en frecuencia de una se˜nal de entrada f (t) y X (f ) es espectro en frecuencia de una se˜nal de salida x(t), la funci´ on de respuesta en frecuencia (FRF), H (f ), entre ambas se˜ nales se define como: H (f ) =
X (f ) F (f )
(4.8)
Al calcular H (f ) con la expresi´on anterior se corre el riesgo que existan t´erminos donde F (f ) sea cero. Por lo tanto, en la practica se utilizan maneras alternativas de calcular H (f ), utilizando las potencias espectrales: H 1 (f ) =
X (f ) F ∗ (f ) GXF = F (f ) F ∗ (f ) GF F
(4.9)
H 2 (f ) =
X (f ) X ∗ (f ) GXX = F (f ) X ∗ (f ) GF X
(4.10)
El principal motivo para estimar las funciones de respuesta en frecuencia con las ecuaciones anteriores es la reducci´on del ruido no correlacionado en las se˜nales de entrada o salida al promediar.
FUNCIONES DE RESPUESTA EN FRECUENCIA Y FUNCIONES DE COHERENCIA
55
En la pr´ actica, la funci´ on de respuesta en frecuencia es estimada con valores promedio de las potencias espectrales, ˆ F F G
ˆ XX G
ˆF X G
ˆ XF G
=
=
=
=
1 N a 1 N a 1 N a 1 N a
N a
(GF F )n
(4.11)
(4.12)
(4.13)
(4.14)
n=1 N a
(GXX )n
n=1 N a
(GF X )n
n=1 N a
(GXF )n
n=1
donde N a es el numero de promedios (el ensayo se repite N a veces), lo que entrega una aproximaci´on de m´ınimos cuadrados de H (f ). Dado que las funciones de respuesta en frecuencia se obtiene de una aproximaci´on de m´ınimos cuadrados, se puede definir un coeficiente de correlaci´on. En este caso, la correlaci´ on se denomina funci´ on de coherencia y es una medida del error de m´ınimos cuadrados. La coherencia se define por: γ 2 =
ˆF X G
2
ˆ F F G ˆ XX G
=
H 1 (f ) H 2 (f )
(4.15)
La coherencia varia entre 0 y 1. Un valor de 1, indica una relaci´on perfectamente lineal entre las se˜nales de entrada y salida por sobre todos los promedios. Una coherencia menor a uno, se puede deber a uno de los siguientes motivos:
Ruido no correlacionado en las mediciones de f (t) y/o x(t) No-linealidades del sistema en investigaci´on Leakage en el an´alisis Desfase en las mediciones no compensado en el an´alisis.
FUNCIONES DE RESPUESTA EN FRECUENCIA Y FUNCIONES DE COHERENCIA
4.3.1.
56
Efectos del ruido experimental
Consideremos que las se˜ nales de excitaci´ on y respuesta “verdaderas” son e(t) y r (t), respectivamente, las se˜ nales medidas son: f (t) = e(t) + m(t)
(4.16)
x(t) = r(t) + n(t)
(4.17)
donde m(t) y n(t) representa representan ruido no correlacionado. La ecuaci´ on anterior implica que las potencias espectrales de la se˜nal son: GF F
= GEE + GMM + GEM + GME
GXX
= GRR + GN N + GRN + GN R
GF X
= GER + GEN + GMR + GMN
GXF
= GRE + GRM + GN E + GN M
(4.18)
La funci´on de respuesta en frecuencia verdadera viene dada por: H (f ) =
|H (f )|2
=
GRE GRR = GEE GER
(4.19)
GRR GEE
De la ecuaci´on anterior se tiene: GER GRE = G RR GEE
(4.20)
Bajo el supuesto razonable que las se˜nales de ruido no est´an correlacionadas entre ellas ni con las se˜ nales de excitaci´ on y de respuesta, entonces GEM = GME = GRN = GN R = GEN = GN E = GMR = GRM = GMN = GN M = 0, se pueden estudiar tres casos de inter´es.
FUNCIONES DE RESPUESTA EN FRECUENCIA Y FUNCIONES DE COHERENCIA
57
Caso 1. Ruido en la se˜ nal de excitaci´ on, sin ruido en la respuesta Dado que n(t) = 0, ecuaci´on 4.18 se convierte: GF F
= GEE + GMM
GXX
= GRR
GF X
= GER
GXF
= GRE
Evaluando H 1 de acuerdo a la ecuaci´on 4.9: GRE H H 1 = = GEE + GMM 1 + GMM /GEE
(4.21)
(4.22)
Se ve la verdadera funci´on de respuesta en frecuencia es subestimada ya que el denominador es mayor a uno. Evaluando ahora H 2 de acuerdo a la ecuaci´on 4.10: GRR GRE = = H (4.23) GER GEE lo que demuestra que, bajo la hip´otesis de ruido no correlacionado, H 2 es insensible a ruido en la se˜nal de excitaci´on. H 2 =
De lo anterior se deduce que si la se˜nal de la fuerza esta contaminada por ruido no correlacionado, un estimador preferente para H es H 2 .
La funci´on de coherencia en este caso viene dada por, 1 γ 2 = 1 + GMM /GEE
(4.24)
como es esperado, la coherencia depende de la raz´ on entre la se˜ nal y el ruido: mientras menor sea la raz´ on entre el ruido y la se˜nal, m´ as cercana a uno es la coherencia.
Caso 2. Sin ruido en la se˜ nal de excitaci´ on, ruido en la respuesta Ecuaci´ on 4.18 es en este caso: GF F
= GEE
GXX
= GRR + GN N
GF X
= GER
GXF
= GRE
(4.25)
FUNCIONES DE RESPUESTA EN FRECUENCIA Y FUNCIONES DE COHERENCIA
58
Evaluando H 1 de acuerdo a la ecuaci´on 4.9: H 1 =
GRE = H GEE
(4.26)
Por lo tanto, H 1 es insensible a ruido no correlacionado en la respuesta. Evaluando ahora H 2 de acuerdo a la ecuaci´on 4.10:
GRR + GN N GRE (1 + GN N /GRR ) GN N H 2 = = = H 1 + GER GEE GRR
(4.27)
La ecuaci´on anterior muestra que H 2 sobreestima a la verdadera FRF, dado que el valor en el par´entesis es mayor a 1.
De lo anterior se deduce que si la respuesta esta contaminada por ruido no correlacionado, un estimador preferente para H es H 1 . La coherencia viene dada por, γ 2 =
1
(4.28)
1 + GN N GRR
Caso 3. Ruido en ambas se˜ nales En este caso, GF F
= GEE + GMM
GXX
= GRR + GN N
GF X
= GER
GXF
= GRE
(4.29)
De donde se tiene: H 1 =
GRE H = GEE + GMM 1 + GMM /GEE
GRR + GN N GN N H 2 = = H 1 + GER GRR
(4.30)
(4.31)
y la funci´on de coherencia viene dada por, γ 2 =
1 (1 + GN N GRR )(1 + GMM /GEE )
(4.32)
FUNCIONES DE RESPUESTA EN FRECUENCIA Y FUNCIONES DE COHERENCIA
59
Es evidente que: H 1 < H < H 2
(4.33)
por lo tanto una buena estimaci´ on para H es la media geom´ etrica de estas dos cantidades. Esta estimaci´on se denomina H v H v =
H 1 H 2 = H
1 + GMM /GRR 1 + GMM /GEE
(4.34)
A partir de los resultados anteriores se pueden sacar algunas conclusiones. En las resonancias la se˜ nal de excitaci´ on es particularmente susceptible al ruido, ya que se requiere una peque˜na fuerza para generar desplazamientos significativos, por otro lado, la respuesta tiene una raz´on se˜ nal-ruido baja. Por lo tanto, en las zonas cercanas a la resonancia se puede esperar una mejor estimaci´ on si se utiliza el estimador H 2 . En contraste, en las zonas lejanas a la resonancia (y espec´ıficamente cerca de las antiresonancias) la respuesta es m´as sensible a contaminaci´on por ruido. Entonces en estas zonas se obtiene una mejor estimaci´on con H 1 .
Cap´ıtulo 5
Estimaci´ on de par´ ametros modales Los m´etodos de identificaci´on de par´ametros buscan extraer la informaci´on modal de una estructura a partir de mediciones experimentales. Estos m´etodos se clasifican principalmente en m´etodos en el dominio del tiempo o m´ etodos en el dominio de frecuencias. Los m´ etodos en el dominio del tiempo siempre se pueden utilizar, ya sea para respuesta libre o forzada (con o sin conocimiento de las fuerzas). Por otro lado, los m´ etodos en el dominio de frecuencias se pueden utilizar s´olo en casos de vibraciones forzadas y cuando las fuerzas son conocidas.
Para cada dominio hay m´etodos que utilizan informaci´on de un s´olo punto de medici´on y otros que utilizan la informaci´on de varios puntos simult´aneamente. En cada uno de estos casos pueden haber una o varias fuerzas de excitaci´on externas. Lo que lleva a la siguiente clasificaci´on: 1. Una respuesta debido a una fuerza (SISO) single-input single-output 2. Varias respuestas debido a una fuerza (SIMO) single-input multiple-output 3. Varias respuestas debido a varias fuerzas (MIMO) multiple-input multipleoutput 4. Una respuesta debido a varias fuerzas (MISO) multiple-input multiple-output En el dominio de tiempo, las respuestas contienen naturalmente informaci´on acerca del contenido en frecuencias, aunque est´a “escondida”, por lo que no es posible definir a priori cuantas resonancias hay presentes en un cierto periodo de tiempo. En
60
´ METODOS DE UN GRADO DE LIBERTAD
61
consecuencias, los m´etodos en el dominio del tiempo deben estimar simult´aneamente varias resonancias de la estructura y para los m´etodos SIMO y MIMO varios modos de vibraci´on. Estos m´etodos son conocidos como m´etodos de multiples grados de libertad (MDOF). En el dominio de frecuencia, dado que los “peaks” de las resonancias son visibles, es posible realizar una identificaci´on modo por modo. Estos m´etodos se denominan m´etodos de un grado de libertad (SDOF). A continuaci´on se detallar´an algunos de los m´ etodos de identificaci´on de par´ametros m´ as utilizados.
5.1.
M´ etodos de un grado de libertad
En general, la respuesta din´amica de un sistema es una superposici´on de muchos de sus modos. Sin embargo, si en una cierta banda de frecuencia se asume que un solo modo es importante, lo par´ ametros de este modo se pueden determinar separadamente. Los m´etodos basados en este supuesto se denominan “m´etodos de un grado de libertad”, estos m´ etodos son una manera r´ apida de obtener los par´ ametros modales ya que requieren poco esfuerzo y memoria computacional.
5.1.1.
Peak picking
El m´ etodo “peak picking” es talvez el m´ etodo m´as sencillo dentro de los de un grado de libertad. Este m´ etodo utiliza la curva de la FRF en la cercan´ıa a una resonancia como si fuese la curva de un sistema de un grado de libertad. El procedimiento del m´etodo peak picking es el siguiente:
Estimar la frecuencia natural La frecuencia natural se selecciona como la frecuencia del m´aximo en la curva de FRF, ωr = ω m´ax . Estimar el amortiguamiento Para estimar el amortiguamiento, se ubican las frecuencia ubicadas a cada lado del m´aximo identificado y que corresponden a una amplitud de la FRF αm´ax igual a (ver Figura 5.1). La raz´on de amortiguamiento se puede estimar 2 como,
√
ωb2 ωa2 ζ r = 4ωr2
−
≈ ωb2−ωrωa
(5.1)
Estimar la constante modal La funci´ on de respuesta en frecuencia a un sistema de un grado de libertad
´ METODOS DE UN GRADO DE LIBERTAD
62
cerca de la resonancia viene dada por (ecuaci´on 1.29), H (ωr )
( jωrA− λr ) = ( jωr − σAr − jω r ) = − σAr = ζ rAωr
(5.2)
Conocido el valor del m´ aximo de la FRF, la constante A se puede estimar como: αm´ax ζ r ωr
FRF (dB)
| α (ω r ) | | α (ω ) | r
2
ω b
Frecuencia
ω r
ω a
Figura 5.1: Peak picking
5.1.2.
Circle fitting
El m´ etodo de “circle fitting” es un de m´etodo de un grado de libertad, puede estimar los polos del sistema y los modos normales (reales o complejos). Este m´ etodo se basa en el hecho que la funci´ on de respuesta en frecuencia de un sistema de un grado de libertad describe un circulo en el diagrama de Nyquist (ver Figura 1.11). Si la influencia de otros modos se aproxima por una constante compleja R + jI , la funci´ on de respuesta en frecuencia cerca de la resonancia, ωr , se puede escribir como, H (ω )
+ R + jI −σr +U j+( jV ω − ωr )
(5.3)
Para determinar los par´ametros de H, primero se selecciona un conjunto de puntos en las cercan´ıas de la resonancia seleccionada. Luego se ajusta un circulo a estos puntos, como se ilustra en la Figura 5.2.
´ METODOS DE UN GRADO DE LIBERTAD
63
Figura 5.2: Circle fitting La frecuencia natural amortiguada, ωr , corresponde al punto de m´axima raz´on de cambio de ´angulo entre puntos (m´axima distancia angular entre puntos) o con un ´ angulo de fase cercano al ´ angulo de fase del centro del circulo, para modos bien separados la diferencia entre estos ´angulos ser´a peque˜ na.
La raz´ on de amortiguamiento, ζ r , se puede estimar como, ζ r =
ω2 ω1 ωr (tan(θ1 /2) + tan(θ2 /2))
−
(5.4)
donde, ω1 y ω2 son dos frecuencias a ambos lados de ωr θ1 y θ2 son los ´angulos entre el radio de ω1 y ω2 y el radio de ωr . El di´ ametro del circulo y la posici´on angular de la frecuencia natural contiene la informaci´ on del residuo complejo U + jV : φ =
tan(α) =
√ U 2 + V 2
(5.5)
σr U V
(5.6)
donde φ es el di´ametro del circulo y α es el ´ angulo entre la linea que conecta el centro del circulo con la frecuencia natural y el eje imaginario.
´ METODOS CON MULTIPLES GRADOS DE LIBERTAD EN EL DOMINIO DE FRECUENCIAS
64
5.2.
M´ etodos con multiples grados de libertad en el dominio de frecuencias
5.2.1.
M´ etodo de m´ınimos cuadrados no lineales, LSFD
El m´ etodo de m´ınimos cuadrados no lineales en el dominio de frecuencias, es un m´ etodo para estimar los polos y modos normales (si se utilizan multiples respuestas). Se basa en el modelo modal en el dominio de frecuencias. La funci´on de respuesta en frecuencia entre una respuesta en el punto i y una excitaci´on en el punto k se puede aproximar por: N m
H ik (ω ) =
r=1
φir Lrk φ∗ir L∗rk + jω λr jω λ∗r
−
−
+ U Rik
ik − LR 2 ω
(5.7)
Gi k(ω)
donde U Rik y LRik son los residuos superiores e inferiores respectivamente. Los residuos aproximan el efecto de modos bajo y sobre el rango de frecuencias en estudio. H ik (ω) es la funci´on de respuesta en frecuencia medida experimentalmente, mientras que el lado derecho de la ecuaci´ on es el modelo modal considerando N u par´ ametros desconocidos λr , φir , Lrk , U Rik , LRik , indicado en la funci´on: Gik (ω ) = Gik (ω, λr , φir , Lrk , U Rik , LRik ) r=1,...,N m
|
(5.8)
La diferencia entre la funci´on de respuesta en frecuencia medida y estimada viene dada por: eik (ω ) = H ik (ω )
− Gik (ω)
(5.9)
El error total en el rango de frecuencias de inter´es es: N f
E ik =
eik (ωf )e∗ik (ωf )
(5.10)
f =0
Considerando todas las funciones de respuesta en frecuencia entre N i entradas y N o respuestas, el error total es: N 0 N i
E =
i=1 k=1
E ik
(5.11)
´ METODOS CON MULTIPLES GRADOS DE LIBERTAD EN EL DOMINIO DE FRECUENCIAS
65
Los par´ametros desconocidos se obtienen al imponer que ´estos minimicen el error total: ∂E ∂λr
∂E ∂LRik
= 0
(5.12)
.. .
(5.13)
= 0
(5.14)
El sistema de ecuaciones anterior es altamente no-lineal, pero puede resolverse iterativamente como un problema linealizado (con una expansi´on de primer orden). Debido a las desventajas conocidas de estos algoritmos, como la necesidad de tener buenos valores de partida, la velocidad de convergencia limitada, el riesgo de divergencia, etc., estos m´ etodos nunca se hicieron populares. Sin embargo, pueden ser una herramienta u ´ til para mejorar la precisi´on de un modelo modal existente. Si los polos del sistema y los factores de participaci´on modal ya fueron estimados, la ecuaci´ on 5.7 se convierte en un set de ecuaciones lineales en funci´ on de los par´ametros desconocidos: φir , U Rik y LRik . El set de ecuaciones resultantes es relativamente f´acil de resolver.
5.2.2.
Rational Fractional Polynomials
La idea detr´as de este m´ etodo es expresar la funci´on de respuesta en frecuencia como el cuociente entre dos polinomios. Estableciendo la relaci´on entre los coeficientes de los polinomios y los par´ametros modales, se puede llegar a la identificaci´on de los par´ametros.
La definici´ on de la funci´ on de respuesta en frecuencia como el cuociente de dos polinomios viene de la definici´ on de la funci´ on de transferencia en el dominio de Laplace, N (s) a0 + a1 s + a2 s2 + . . . + am sm H (s) = = D (s) b0 + b1 s + b2 s2 + . . . + sn
(5.15)
En este caso el orden del denominador, n, es mayor que el del numerador, m, por 2. Para simplificar se define: p0 (s) = 1,
p1 (s) = s,
p2 (s) = s 2 ,
...,
pm (s) = s m
(5.16)
q 0 (s) = 1,
q 1 (s) = s,
q 2 (s) = s 2 ,
...,
qm (s) = s m
(5.17)
´ METODOS CON MULTIPLES GRADOS DE LIBERTAD EN EL DOMINIO DE FRECUENCIAS
66
Entonces la funci´on de respuesta en frecuencia se puede expresar como: H (ω ) =
m k=0 ak pk ( jω ) n k=0 bk q k ( jω )
(5.18)
Para el an´alisis subsiguiente se asume que hay mediciones en p frecuencias ω1 , ω2 , . . ., ω p . Tambi´en se asume que hay p frecuencias negativas ω p , ω p−1 , . . ., ω−1 . Se puede demostrar para las frecuencias negativas se cumple que, H (ω−i ) = H ( jωi ) = H ∗ ( jω i )
(5.19)
−
El prop´osito de las frecuencias negativas quedar´a claro m´as adelante. Para comenzar la identificaci´on de los coeficientes se define la funci´on de error, e(ω ) = H (ω )
− H ˜ (ω)
(5.20)
˜ (ω) es la donde H (ω) es la FRF estimada como el cuociente de dos polinomios y H FRF experimental. Lo que lleva a: e(ω ) =
m k=0 ak pk ( jω ) n k=0 bk q k ( jω )
− H ˜ (ω)
(5.21)
Este error es una funci´ on no-lineal de los coeficientes ak y bk . Para facilitar el an´ alisis el error se re-define como: n
eˆ(ω ) = e (ω )
m
bk q k ( jω ) =
k=0
k=0
n−1
˜ (ω ) ak pk ( jω ) H
−
bk q k ( jω ) + q n ( jω ) (5.22)
k=0
El error total para todas las frecuencias se puede escribir en forma matricial: E = E =
{e(ω p), . . . , e(ω 1 ), e(ω1 ), . . . , e(ω p)}T [P ]{A} − [Q]{B } − {W } −
(5.23)
−
(5.24)
´ METODOS CON MULTIPLES GRADOS DE LIBERTAD EN EL DOMINIO DE FRECUENCIAS
67
donde,
[P ] =
[Q] =
A
=
B
=
p0 ( jω − p ) p1 ( jω− p ) . . . pm ( jω− p ) .. . ... ... ... p0 ( jω −1 ) p1 ( jω−1 ) . . . pm ( jω−1 ) p0 ( jω 1 ) p1 ( jω1 ) . . . pm ( jω1 ) .. . ... ... ... p0 ( jω p ) p1 ( jω p ) . . . pm ( jω p )
˜ ( jω − p )q 0 ( jω − p ) H ˜ ( jω − p )q 1 ( jω − p ) . . . H .. . ... ... ˜ ( jω −1 )q 0 ( jω −1 ) H ˜ ( jω −1 )q 1 ( jω −1 ) . . . H ˜ ( jω 1 )q 0 ( jω 1 ) ˜ ( jω 1 )q 1 ( jω 1 ) H H ... .. . ... ... ˜ ( jω p )q 0 ( jω p ) ˜ ( jω p )q 1 ( jω p ) H H ...
{a0 , a1, . . . , am}T {b0 , b1 , . . . , bn 1}T {H ˜ ( jω p)q n( jω p), . . . , H ˜ ( jω ˜ ( jω p )q n ( jω p )}T , H
˜ ( jω− p )q m ( jω− p ) H ... ˜ ( jω−1 )q m ( jω −1 ) H ˜ ( jω1 )q m ( jω 1 ) H ... ˜ ( jω p )q m ( jω p ) H
−
W =
−
−
˜
−1 )q n ( jω −1 ), H ( jω1 )q n ( jω 1 ), . . .
La magnitud total del error se puede definir como:
{ }H {E }
J = E
(5.25)
donde H denomina la transpuesta conjugada. Para determinar los coeficientes A y B que minimizan la magnitud del error, las siguientes derivadas parciales deben ser iguales a cero:
{ }
{ }
∂J ∂J = =0 ∂ A ∂ B
{ }
(5.26)
{ }
Lo que lleva al siguiente sistema de ecuaciones para estimar los coeficientes:
1 ∗ H T 2 ([P ] [P ] + [ P ] [P ] ) H T
( Re([P ] [Q]))
−
−H Re([P ]H [QT ])
1 2 ([Q]
[Q] + [ Q] [Q]∗ )
{ } A B
{ }
Re([P ]H W ) = Re([Q]H W )
{ } { }
(5.27)
´ METODOS CON MULTIPLES GRADOS DE LIBE RTA D EN EL DOMINIO DEL TIEMPO
68
Los par´ ametros modales se puedes determinar a partir de los coeficientes A y B , al representar el cuociente entre polinomios como fracciones parciales:
{ }
{ }
n/1
H (ω ) =
rk
jω
k=1
rk∗
− pk + jω − pk ∗
(5.28)
donde pk es el k-esimo polo y rk el k -esimo residuo.
5.3.
M´ etodos con multiples grados de libertad en el dominio del tiempo
Los m´ etodos en el dominio del tiempo, a diferencias de los algoritmos en el dominio de frecuencias, utilizan informaci´on de la respuesta en el tiempo. Estos algoritmos ajustan los datos a la funci´on de respuesta a un impulso. La funci´on de respuesta al impulso para un sistema con multiples grados de libertad viene dada por:
N
h(t) =
∗
λr t (Qr φr φT + Q∗r φ∗r φ∗r T eλr t ) r e
(5.29)
r=1 N
=
∗
(Ar eλr t + A∗r eλr t )
(5.30)
r=1 2N
=
Ar eλr t , para r > N : A r = A ∗r ; λr = λ ∗r
(5.31)
r=1
5.3.1.
El m´ etodo Ibrahim
El m´etodo de Ibrahim, tambi´ en conocido como ITD, construye un problema de valores y vectores propios a partir de la respuesta en el tiempo del sistema. La soluci´ on de este problema permite derivar las frecuencias naturales, factores de amortiguamiento y los modos normales. La respuesta en el tiempo en un punto i, para un instante de tiempo tj se puede expresar como la suma de respuestas individuales de cada modo: 2N
xi (tj ) =
r=1
φir eλr tj
(5.32)
´ METODOS CON MULTIPLES GRADOS DE LIBE RTA D EN EL DOMINIO DEL TIEMPO
69
Considerando q puntos de medici´on y L instantes de tiempo, se tiene que,
x1 (t1 ) . . . x1 (tL ) .. .. .. . . . xq (t1 ) . . . xq (tL )
q×L
o bien,
=
φ11 .. . φq1
. . . φ1L .. .. . . . . . φq2N
q ×2N
eλ t .. . 1
eλ
1
t
2N 1
. . . eλ tL .. .. . . . . . eλ N tL 1
2
2N ×L
X = ΦΛ
(5.33)
(5.34)
donde Λ es una matriz constituida por los elementos e λi tj . Esta ecuaci´on por si sola es insuficiente para determinar los par´ametros modales. Consideremos entonces un segundo set de L puntos, desfasados en ∆t con respecto al set anterior: 2N
xi (tj + ∆t) =
2N
λr (tj +∆t)
φir e
r=1
=
φir eλr ∆t eλr tj
(5.35)
r=1
Definiendo: xi (tj )∗ φ∗ir
= xi (tj + ∆t)
(5.36)
= φir eλr ∆t
(5.37)
Siguiendo el mismo procedimiento anterior, se obtiene un segundo set de ecuaciones: X ∗ = Φ∗ Λ
Definiendo una matriz As de q
(5.38)
× q con la siguiente propiedad:
As Φ = Φ∗
(5.39)
Premultiplicando la ecuaci´on 5.34 por As : As X = A s ΦΛ = Φ∗ Λ
(5.40)
Sustituyendo 5.38: As X = X ∗
(5.41)
Si X y X ∗ son conocidos, As se puede despejar como: As = X ∗ X +
(5.42)
donde X + es la pseudoinversa de X y viene dada por: X + = X T (X T X )−1 asumiendo que L > q .
´ METODOS CON MULTIPLES GRADOS DE LIBE RTA D EN EL DOMINIO DEL TIEMPO
70
Recordando que As Φr = Φ∗r , se tiene que: As Φr
As
λr ∆t
−e
I Φr
= Φr eλr ∆t =
(5.43)
{0}
(5.44)
La ecuaci´ on anterior es una ecuaci´ on de valores propios. Los modos normales corresponden a los vectores propios y las frecuencias naturales y factores de amortiguamiento se deducen de los valores propios: eλr ∆t
= e(σr +jω dr )∆t = α r + jγ r
(5.45)
αr
= eσr ∆t cos(ωr ∆t)
(5.46)
γ r
= eσr ∆t sin(ωr ∆t)
(5.47)
σr
ln(αr2 + γ r2 ) 2∆t
=
arctan ωdr
=
γ r αr
∆t
(5.48)
(5.49)
Finalmente, los par´ametros din´amicos se obtienen como: ωr
=
ζ r
=
5.3.2.
−
2 ωdr + σr2
σr
2 + σ 2 ωdr r
(5.50) (5.51)
El m´ etodo least-squares complex exponential (LSCE)
Tenemos que la respuesta a un impulso IRF, se puede expresar como:
2N
h(t) =
r=1
Ar eλr t
(5.52)
´ METODOS CON MULTIPLES GRADOS DE LIBE RTA D EN EL DOMINIO DEL TIEMPO
71
Si consideramos que los datos son adquiridos en intervalos de tiempo k ∆ (k = 0, 1, . . . , 2N ), se tiene la siguiente serie 2N
h(k ∆) =
Ar eλr k∆
(k = 0, 1, . . . , 2N )
(5.53)
r=1 2N
hk
=
Ar zrk
zrk = e λr k∆
(k = 0, 1, . . . , 2N ),
(5.54)
r=1
Cada termino de la suma posee valores reales, aunque los residuos A r y los polos λ r tienen valores complejos. Se puede demostrar que existe un polinomio de coeficientes reales de manera que: β 0 + β 1 zr + β 2 zr2 + . . . + β 2N zr2N = 0
(5.55)
Esta ecuaci´on se conoce como la ecuaci´on de Prony. Dado que hay 2N +1 ecuaciones en la expresi´on 5.58, se pueden multiplicar las ecuaciones por un coeficiente β correspondiente y sumarlas para formar siguiente ecuaci´on: 2N
2N
β k hk
2N
β k
=
r=1
r=1
r=1
2N
2N
Ar
=
r=1
Ar zrk
(5.56)
β k zrk
(5.57)
r=1
Esta u ´ ltima expresi´ on sabemos que es igual a 0, si zr es una ra´ız de la ecuaci´on 5.55. Lo que nos lleva a una relaci´on simple entre los coeficientes β y los datos de la respuesta a un impulso: 2N
β k hk = 0
(5.58)
r=1
A partir de la ecuaci´on anterior se pueden determinar los valores de β k . Conocidos los coeficientes β k se puede resolver la ecuaci´on 5.55 y extraer las ra´ıces z r . A partir de zr se calculan las frecuencias naturales y razones de amortiguamiento: ωr ζ r
=
1 ∆
=
− ln(zr zr )
ln zr ln zr∗
2ωr ∆
(5.59)
∗
(5.60)
DIAGRAMAS DE ESTABILIDAD
72
Para identificar los modos propios del sistema, podemos reescribir la ecuaci´on 5.54 de la siguiente forma:
1 z1 .. .
1 z2 .. .
z12N −1
z22N −1
... ... .. .
1 z2N .. .
2N −1 . . . z2N
A1 A2 .. . A2N
=
h0 h1 .. . h2N −1
(5.61)
La soluci´on de este sistema de ecuaciones entrega los residuos y por lo tanto los modos normales.
5.4.
Diagramas de estabilidad
En todos los algoritmos anteriores se debe definir a priori cual es el numero de polos que se desea estimar. Lo que muchas veces no es una decisi´on simple. Si se utilizan menos polos que lo actuales, el ajuste no ser´a adecuado. En cambio, si se definen m´as polos que los reales, los algoritmos van a entregar polos “computacionales” que no corresponden a polos reales del sistema, sino que son polos que tratan de modelar el ruido en los datos. Una metodolog´ıa muy ´util para determinar el numero de polos reales en un sistema son los diagramas de estabilidad.
Los diagramas de estabilidad son una herramienta muchas veces imprescindible en un an´ alisis modal experimental, ellos ayudan a separar los polos reales de los polos computacionales. Estos diagramas se obtienen, al repetir el an´alisis modal incrementando el orden del sistema (n´umero de polos asumidos). Para cada orden se estiman los polos, los resultados son presentados gr´aficamente en el diagrama de estabilidad (Figura 5.3). En el eje vertical se encuentra el orden y el eje horizontal representa la frecuencia natural del polo estimado. En general, los polos reales aparecen a la misma frecuencia en el diagrama, independiente del orden del sistema. En cambio, la frecuencia de los polos computacionales, var´ıa al aumentar el orden del sistema. El diagrama de la Figura 5.3 se representan los polos estables (en frecuencia y amortiguamiento) con un circulo y con una cruz los polos inestables. Normalmente, se define que un polo es estable si la frecuencia no var´ıa m´ as del 1 % de su magnitud y la raz´ on de amortiguamiento no var´ıa m´a s del 5%.
DIAGRAMAS DE ESTABILIDAD
73
40
30 n e d r 20 O
10
0
20
40 60 80 Frecuencia (Hz)
100
Figura 5.3: Diagrama de estabilidad
120
Cap´ıtulo 6
Medici´ on Experimental El primer paso en un an´ alisis modal experimental es adquirir las funciones de respuesta frecuencia experimentales de una estructura. El m´etodo m´as usual es excitar la estructura con una fuerza conocida y medir de forma simultanea la fuerza y las respuestas en la estructura. Como resultado, se obtiene un grupo de FRFs que pueden ser utilizadas posteriormente para derivar los par´ametros modales de la estructura utilizando algunos de los algoritmos descritos en el capitulo anterior.
6.1.
An´ alisis previo a las mediciones
Al preparar un montaje experimental para un an´alisis modal se debe considerar el prop´osito del experimento, los datos requeridos (FRFs o par´ametros modales), la precisi´ on requerida, etc. Para ello se necesita la mayor cantidad de informaci´ on posible de la estructura, la que se puede obtener de experimentos anteriores en estructuras similares o de modelos num´ericos de la estructura. A continuaci´on se describir´an algunas herramientas u ´ tiles para definir el montaje experimental de acuerdo a los requerimientos. Desde un punto de vista pr´actico, los criterios siguientes se deben cumplir en un buen dise˜ no de un montaje experimental: Correspondencia: Los modos medidos experimentalmente deben corres-
ponder a los modos reales, los que desafortunadamente son desconocidos. Sin embargo, experimentos previos en estructuras similares o un modelo en elementos finitos pueden ayudar a estimar los modos. Adicionalmente, el test debe producir modos claramente distinguibles. La independencia del los
74
´ ANALISIS PREVIO A LAS MEDICIONES
75
modos esta directamente relacionada con el rango de la matriz de vectores propios Φ. Excitaci´ on: El montaje debe incluir un sistema de excitaci´on que garantice
que todos los modos de inter´es son excitados. Identificaci´ on: Los datos medidos deben contener la informaci´on necesaria
para identificar los par´ametros de inter´ es. Por lo tanto, el dise˜no del montaje depende del prop´osito del experimento.
actica, se requiere visualizar los modos obtenidos, Visualizaci´ on: En la pr´ de manera de evaluar de precisi´ on de ´estos y para compararlos con modos estimados. La visualizaci´on tambi´ en es importante para detectar posibles discrepancias. Robustez: Dado que el montaje est´a basado en experimentos previos o en
modelos num´ ericos, donde ambos contiene errores, ´este debe ser robusto: No debe ser muy sensible a estos errores. Por lo tanto, es preferible alg´un grado de redundancia. Accesibilidad: Las ubicaciones seleccionadas para medir la respuesta y para
excitar la estructura deben ser accesibles.
6.1.1.
Rango de frecuencias
Es importante definir bien el rango de frecuencias para identificar todos los modos de inter´es. Si se sabe cuantos modos se quieren identificar, el rango de frecuencias se puede definir con un an´alisis previo de un modelo en elementos finitos. Desde luego se debe considerar un margen de error. En caso que no se cuente con un modelo en elementos finitos, el rango de frecuencias se puede definir in situ. Al intentar en una primera instancia con una rango de frecuencias amplio, luego contar los peaks de una de las FRFs para definir la frecuencia m´axima que se debe medir para identificar n modos.
6.1.2.
Selecci´ on de la ubicaci´ on de las respuestas
En un an´alisis modal se deben definir suficientes puntos de medici´on de manera que los modos identificados sean independientes. A continuaci´on se describir´ an algunos de los algoritmos disponibles para determinar la ubicaci´on ´optima de los puntos de medici´on. En los algoritmos a continuaci´on se asume que se cuenta con los modos normales obtenidos a trav´es de alg´un m´etodo num´erico. La independencia de los modos se puede determinar a partir del rango de la matriz Φ. Este rango equivale al rango de la matriz Φ T Φ (la matriz de Fisher).
´ ANALISIS PREVIO A LAS MEDICIONES
76
La contribuci´on de cada grado de libertad i al rango de la matriz de Fisher, viene dado por el siguiente indicador (independencia efectiva):
EfIi = diag Φ ΦT Φ
−1
ΦT
(6.1)
Al ir eliminando iterativamente los grados de libertad con menor EfI i y re-calculando los nuevos EfI, se llega a un set ´ optimo de puntos de medici´ on con respecto a la independencia de los modos de inter´es.
Otro indicador de la independencia entre modos es el MAC (Modal Assurance Criterion), el que se define como:
2
φT φj M AC ij = T i φi φi φT j φ j
(6.2)
Donde φi es el iesimo modo normal. MAC es un factor que mide la correlaci´ on entre dos modos. Un valor de 0 indica que no hay correlaci´ on, mientras que un valor de 1 indica dos modos perfectamente correlacionados. Si se calcula la correlaci´ on para todos los pares de modos, se puede construir la matriz MAC. Si los modos son 100 % independientes, los valores fuera de la diagonal de la matriz son todos cero y en la diagonal son todos 1. Los puntos de medici´on se determinan al encontrar al menor set de puntos de medici´on que minimize los valores fuera de la diagonal de la matriz MAC. En general, se aceptan valores de hasta un 20 % de correlaci´on cruzada. Este es el m´etodo m´as utilizado en la selecci´ on de los puntos de medici´on.
6.1.3.
Selecci´ on de los puntos de excitaci´ on
Los puntos de excitaci´on de la estructura se deben seleccionar de manera garantizar que todos los modos de inter´es sean excitados adecuadamente. Un modo en particular va a estar bien excitado si la fuerza se aplica en un punto de alto desplazamiento. El m´etodo m´as usual para seleccionar los puntos de excitaci´on es a trav´es de los “driving point residues” (DPR). El residuo A ikr esta definido por la expresi´on de la funci´on de respuesta en frecuencia en t´erminos de los par´ ametros modales:
´ ANALISIS PREVIO A LAS MEDICIONES
N
H ik (ω ) =
r=1 N
=
r=1
77
Qr φir φkr Q∗r φ∗ir φ∗kr + jω λr jω λr
−
−
Aikr A∗ikr + jω λr jω λr
−
−
(6.3)
(6.4)
De donde se define el DPR del punto i como: Aiir
φ2ir = 2 jωr
(6.5)
Al estudiar los DPR para todas puntos candidatos de excitaci´on y para todos los modos de inter´es, se obtiene informaci´on ´ util para la selecci´ on de los puntos de excitaci´ on. En general, los grados de libertad con DPR altos para el mayor numero posible de modos van a ser buenos puntos de excitaci´on.
6.1.4.
Selecci´ on de los puntos de suspensi´ on
Dado que es muy dif´ıcil imitar condiciones de borde reales en un laboratorio, en la mayor´ıa de las mediciones experimentales la estructura se monta de manera de simular una condici´ on libre (sin condiciones de borde). Esta condici´ on se logra al suspender la estructura con materiales blandos como resortes o el´asticos. Con este arreglo se producen modos de cuerpo r´ıgido. Si las frecuencias naturales de los modos r´ıgidos son lejanas a la primera frecuencia natural de la estructura, la medici´ on de la FRF no se ver´a afectada por esta condici´on de borde. En la Figura 6.1 se muestra un montaje de una placa suspendida por cuatro resortes.
Figura 6.1: Una placa suspendida por cuatro resortes.
EL MONTAJE EXPERIMENTAL
78
Al contrario de los puntos de excitaci´ on, los puntos de suspensi´ on se deben seleccionar en puntos con muy bajas amplitud de movimiento. La t´ecnica, por lo tanto es similar a la de la secci´ on anterior, pero ahora se deben seleccionar los grados de libertad con DPR bajos.
6.2.
El montaje experimental
6.2.1.
Mecanismo de Excitaci´ on
La primera parte de un montaje experimental es el mecanismo de excitaci´ on que aplica una fuerza de suficiente amplitud y contenido en frecuencia a la estructura. Existen diferentes mecanismo capaces de excitar la estructura. Los dos m´as comunes son el shaker (Figura 6.2) y el martillo (Figura 6.3).
Un shaker electromagn´etico, tambi´ en conocido como shaker electrodin´amico, es el tipo de shaker m´as utilizado en an´alisis modal. Consiste en un im´an, un bloque libre y una bobina. Cuando una se˜nal el´ ectrica pasa por la bobina, se genera una fuerza proporcional a la corriente y a la densidad de flujo magn´ etico generado, la que mueve al bloque. Un shaker electromagn´etico tiene un rango amplio de frecuencias y de amplitudes. Para frecuencias bajas y altas amplitudes de excitaci´on, se puede utilizar un shaker electro-hidr´aulico. Generador de la señal
Analisador
Amplificador de poder
Shaker Sensor de fuerzas
Signal flow
Acondicionador de la señal
Aceler ómetro
Figura 6.2: Set-up experimental con excitaci´on por shaker Un martillo es un aparato que produce una fuerza de excitaci´ on de la forma de un pulso. Consiste en la punta del martillo, un sensor de fuerzas, una masa y el mango. La punta del martillo se puede cambiar para alterar su dureza. Materiales t´ıpicos para puntas de martillo son; goma, pl´astico y acero. La dureza de la punta, en conjunto con la dureza de la superficie de la estructura, est´ a directamente
EL MONTAJE EXPERIMENTAL
79
relacionada con el rango de frecuencias contenida en la se˜nal de la fuerza. Mientras mayor es la dureza, mayor es el rango de frecuencias contenido.
Figura 6.3: Set-up experimental con excitaci´on por martillo
6.2.2.
Aceler´ ometro
El aceler´ometro es el sensor m´as com´un utilizado en an´ alisis modal. Mide la aceleraci´ on de un punto en la estructura, la se˜nal de salida viene en la forma de voltaje. Esta se˜nal es transformada por un acondicionador de la se˜nal antes de ser procesada por otro hardware o software. El tipo m´as com´ un de aceler´ometros son los piezoel´ectricos, ilustrados en la Figura 6.4. Estos sensores contienen un cristal piezoel´ectrico en su interior, ´este cristal produce una carga el´ ectrica al ser deformado. Al seleccionar un aceler´ ometro hay distintos factores importantes que se deben considerar, estos son; el rango de frecuencias, sensibilidad, masa y estabilidad bajo cambios de temperatura. La sensibilidad de un aceler´ ometro determina la raz´ on entre la se˜nal medida y el ruido. Mientras m´as alta es la sensibilidad del aceler´ometro m´as precisas son las mediciones. La masa del aceler´ometro tambi´en es muy importante, ya que puede modificar las caracter´ısticas de la estructura. Mientras menor es la masa mejor, aunque esto significa muchas veces una menor sensibilidad.
EL MONTAJE EXPERIMENTAL
80
Resorte
pre-comprimido Masa Cr istal
+ q –
Cuerpo
Figura 6.4: Esquema de un aceler´ometro piezoel´ectrico
Sensor
Aceler ómetro
Estructur a Montaje Estructur a
Figura 6.5: Montaje de un aceler´ometro La precisi´on de un aceler´ometro tiene una alta dependencia con la manera en que es montado sobre la estructura. En la Figura 6.5 se muestra el montaje t´ıpico de un aceler´ ometro. La flexibilidad del montaje afecta las mediciones del sensor, mientras m´as r´ıgido es el montaje mayor es la precisi´ on. Existen diferentes mecanismos de montaje; tornillo, adhesivo, magn´etico, capa delgada de cera, entre otros. Al atornillar el sensor a la estructura se obtienen los mejores resultados, aunque esto significa un montaje experimental m´as complejo. En la pr´actica, s´olo se utilizan sensores atornillados en estructuras que deben ser monitoreadas de forma continua.
6.2.3.
Sensor de fuerzas
El sensor de fuerzas es otro tipo de sensor utilizado en an´alisis modal. Al igual que un aceler´ometro, un sensor de fuerzas piezoel´ectrico produce una carga o voltaje proporcional a la fuerza aplicada. A diferencia de un aceler´ometro, un sensor de fuerzas no tiene una masa inercial en su interior, si no que debe ser comprimido o estirado f´ısicamente para generar la carga (Figura 6.6). Para un montaje con
EL MONTAJE EXPERIMENTAL
81
un shaker, el sensor de fuerzas se conecta entre la superficie de la estructura en el shaker. En el caso de un martillo, el sensor se ubica en la punta del martillo y es comprimido cuando el martillo impacta la superficie. Las consideraciones para seleccionar un sensor de fuerzas son similares a las de un aceler´ometro: rango de frecuencias, sensibilidad, masa y estabilidad bajo cambios de temperatura.
F
Cr istales
+
q –
F
Figura 6.6: Esquema de un sensor de fuerzas piezoel´ectrico
´ ´ SELECCION DE LA FUERZA DE EXCITACION
6.3.
82
Selecci´ on de la fuerza de excitaci´ on
La precisi´on de las mediciones experimentales depende en gran parte de la forma de la fuerza de excitaci´ on utilizada. Aunque te´ oricamente una FRF no debiese depender de la excitaci´ on utilizada, en la pr´actica la precisi´ on y calidad de una FRF depende, entre otros factores, en la elecci´ on de la fuerza de excitaci´ on. Es por lo tanto, muy importante en un an´alisis modal experimental evaluar todos los m´etodos posibles de excitaci´on de manera de elegir el m´as adecuado.
6.3.1.
Excitaci´ on sinusoidal
La excitaci´ on sinusoidal es la forma m´as tradicional en an´alisis modal. La fuerza contiene una sola frecuencia en un tiempo y la excitaci´on cambia de una frecuencia a otra con un paso dado, permitiendo a la estructura excitar un modo a la vez. Esta excitaci´ on es efectiva para excitar una estructura con niveles de vibraci´on altos, para caracterizar no-linealidades de una estructura y para excitar modos normales de una estructura amortiguada.
La Figura 6.7 presenta una se˜ nal sinusoidal t´ıpica y su contenido en frecuencias.
Figura 6.7: Se˜nal sinusoidal
´ ´ SELECCION DE LA FUERZA DE EXCITACION
83
En este caso la din´amica de una estructura se descompone f´ısicamente por frecuencia y posici´on. Al sintonizar cerca de una frecuencia natural, la respuesta de la estructura es dominada por ese modo de vibraci´on. Esta segregaci´on natural proporciona una via para la identificaci´on directa de par´ametros, la que usualmente tiene una buena raz´on se˜ nal/ruido. Esta es una posibilidad que otros m´ etodos de excitaci´ on no ofrecen.
6.3.2.
Excitaci´ on aleatoria
Una se˜ nal de excitaci´on aleatoria es una se˜nal aleatoria estacionaria que sigue una distribuci´ on Gaussiana. Contiene todas las frecuencias en el rango seleccionado. Para una estructura que tiene un comportamiento no-lineal, una excitaci´on aleatoria tiende a linealizar este comportamiento. Por lo tanto la funci´ on de respuesta en frecuencia derivada de una excitaci´ on puramente aleatoria va a ser una versi´on linealizada de la FRF. Esta FRF, aunque no contiene informaci´ on sobre las nolinealidades, es en realidad una funci´ on muy ´util y es percibida como la “mejor” estimaci´ on para las FRFs. Sin embargo, el hecho que para una excitaci´on aleatoria ni la fuerza ni las respuestas son peri´odicas, puede llevar a errores de leakage.
La Figura 6.8 presenta una se˜ nal aleatoria t´ıpica y su contenido en frecuencias.
Figura 6.8: Se˜nal aleatoria
´ ´ SELECCION DE LA FUERZA DE EXCITACION
6.3.3.
84
Excitaci´ on pseudo-aleatoria
Una se˜ nal de excitaci´on pseudo-aleatoria es una se˜nal aleatoria estacionaria que consiste en frecuencias discretas formadas m´ultiplos de la resoluci´on en frecuencia utilizada en la transformada de Fourier. Es una se˜nal peri´odica con amplitud fija y fase aleatoria. Esta excitaci´on elimina el problema de leakage de la se˜nal aleatoria y usualmente entrega FRFs precisas. La Figura 6.9 presenta una se˜ nal pseudo-aleatoria t´ıpica y su contenido en frecuencias.
Figura 6.9: Se˜ nal pseudo-aleatoria
6.3.4.
Excitaci´ on aleatoria en trenes (burst random)
Una mejor se˜ nal que la pseudo-aleatoria es la denominada se˜nal aleatoria en trenes. Una se˜ nal aleatoria en trenes se crea al encender y apagar una se˜nal aleatoria, de esta manera la medici´ on comienza con una excitaci´ on aleatoria y continua luego que la excitaci´ on se ha apagado. El espectro de la se˜nal aleatoria en trenes tiene amplitud y fase aleatorias y contiene suficiente energ´ıa en todo el rango de frecuencias. Al seleccionar adecuadamente el tiempo de apagado de la se˜nal, esta asemeja una se˜nal pseudo-aleatoria pero sin la necesidad de esperar el decaimiento de la parte transiente de la respuesta. La proporci´on de tiempo que esta apagada la se˜ nal de excitaci´on se debe seleccionar de manera que la respuesta de la estructura
´ ´ SELECCION DE LA FUERZA DE EXCITACION
85
sea cero al final del periodo de medici´on, este tiempo depende principalmente del nivel de amortiguamiento en la estructura. La Figura 6.10 presenta una se˜ nal de excitaci´on aleatoria en trenes, su contenido en frecuencias y la respuesta t´ıpica.
Figura 6.10: Se˜nal aleatoria en trenes
6.3.5.
Excitaci´ on de impacto
La se˜ nal en el tiempo de una se˜nal de fuerza debido a un impacto, es un pulso con un contenido en frecuencias no controlable. En t´ erminos de hardware involucrado, la excitaci´on por impacto es relativamente simple comparado con la excitaci´on con un shaker. Es conveniente de usar y muy portable para mediciones en terreno y laboratorio. Debido a que no existe conexi´on f´ısica entre la fuerza de excitaci´on y la estructura, el test de impacto evita el problema de su interacci´on. La principal desventaja del test de impacto es la dificultad de controlar el nivel de la fuerza y su contenido en frecuencias. Esto puede afectar la raz´on se˜ nal/ruido en las mediciones, resultando en datos de baja calidad.
En la Figura 6.11 se muestra un t´ıpica se˜nal de impacto.
´ EVALUACION INICIAL DE LAS FRFS MEDIDAS
86
Figura 6.11: Se˜nal de impacto
6.4.
Evaluaci´ on inicial de las FRFs medidas
La calidad de un an´ alisis modal recae principalmente en la calidad de las FRFs medidas. Aunque algunos m´etodos de an´alisis modal pueden minimizar ciertos errores en la informaci´ on obtenida experimentalmente, ning´un m´etodo puede rectificar errores fundamentales o mediciones err´ oneas. Las propiedades modales obtenidas de FRFs err´oneas son susceptibles a errores inaceptables. Debido a esto se vuelve esencial verificar la calidad de las funciones de respuesta en frecuencia obtenidas experimentalmente. Las FRFs medidas deben cumplir b´ asicamente con dos requerimientos: (1) la estructura satisface los supuestos requeridos en an´alisis modal y (2) los errores de sistema y humanos son minimizados o eliminados. B´asicamente en an´alisis modal la estructura debe cumplir con reciprocidad, invariabilidad en el tiempo y linealidad. Si estos supuestos no son verificados, las propiedades modales obtenidas no ser´an confiables. Aunque existen algunos m´ etodos para asistir en la b´usqueda de posibles errores en las FRFs medidas, algunos de estos errores no pueden ser identificados. Por ejemplo, no es posible detectar a partir de las FRFs si el sistema de medici´on esta calibrado correctamente.
´ EVALUACION INICIAL DE LAS FRFS MEDIDAS
6.4.1.
87
Repetibilidad
La verificaci´on m´as simple, pero no menos ´util, es la repetibilidad de las mediciones. Esto es, asegurar que el comportamiento din´amico de la estructura y de todo el montaje experimental no depende del tiempo. Para una fuerza de excitaci´ on y un punto de medici´on seleccionados, una estructura lineal deber´ıa dar resultados id´enticos en cada medici´ on. Para un par de puntos de entrada y salida (fuerza y respuesta), se pueden realizar mediciones a intervalos de tiempo. T´ıpicamente, una FRF seleccionada se mide antes y despu´es que se termine de medir toda la estructura. Este proceso que parece trivial, puede ser muy ´util para verificar no s´olo que el comportamiento de la estructura es constante, sino tambi´en que las condiciones de medici´on no han cambiado durante todo el experimento.
6.4.2.
Reciprocidad
Una estructura lineal e invariante en el tiempo debe tener la propiedad de reciprocidad. Esto significa, que una FRF deber´ıa ser id´entica si se intercambian la ubicaci´on de la fuerza y la respuesta. Te´oricamente, esta propiedad se puede derivar a partir de la simetr´ıa de las matrices de masa, rigidez y amortiguaci´on. Debido a esta simetr´ıa, la matriz FRF, que es la inversa de la matriz din´amica, tambi´en va a ser sim´ etrica. La propiedad de reciprocidad de la FRF se puede utilizar para evaluar la fiabilidad y exactitud de los datos medidos.
6.4.3.
Linealidad
El supuesto m´as com´ un en an´ alisis modal es que la estructura se comporta de manera lineal. Sin este supuesto, el an´alisis modal no tiene sentido. Una manera de verificar la linealidad de una estructura es asegurar que la medici´on de las FRFs es independiente de la amplitud de las fuerzas de excitaci´ on. Para verificar esto, se pueden repetir mediciones en una misma ubicaci´on, pero variando la amplitud de la fuerza de excitaci´on.
6.4.4.
Caracter´ısticas especiales de una FRF
A partir de la teor´ıa de an´alisis modal, se pueden derivar algunas caracter´ısticas de las FRFs, las que se pueden utilizar para evaluar la calidad de los datos medidos. Esta evaluaci´on puede detectar errores en las mediciones. La primera caracter´ıstica es la de una FRF directa, en donde el punto de excitaci´on es el mismo punto de medici´on. En este caso se espera siempre ver una antiresonancia
´ EVALUACION INICIAL DE LAS FRFS MEDIDAS
88
entre dos resonancias. Por lo tanto, si esta caracter´ıstica no se observa en la medici´on de una FRF directa, es muy probable que los sensores de fuerza y de respuesta no est´en en realidad en la misma coordinada. Para una estructura empotrada, a bajas frecuencias, la caracter´ıstica predominante de la estructura es su rigidez est´ atica. Por lo tanto, al comienzo de una FRF, se deber´ıa observar una linea de rigidez antes de la primera resonancia, como se muestra en la Figura 6.13. Por otro lado, para una estructura libre, la caracter´ıstica predominante a bajas frecuencias es la masa y la inercia. Esto significa que se debe observar una linea de masa al inicio de la FRF, como se muestra en la Figura ??.
) B d ( F R F
0
Frecuencia (log)
Figura 6.12: FRF directa de una estructura empotrada
–150 –200 –250 d(
B
)
–300 F F
R
–350 –400 –450
0
1
2
3
4
Frecuencia (log)
Figura 6.13: FRF directa de una estructura libre
5
Cap´ıtulo 7
Correlaci´ on Num´ erico-Experimental y Ajuste de Modelos En el dise˜ no de estructuras mec´anicas, el comportamiento din´amico es un tema de vital importancia. La vida ´util bajo cargas c´ıclicas, niveles de vibraci´on o ruido, interacci´on entre sistemas de control y la vibraci´ on de la estructura, . . . son restricciones relevantes en el dise˜no. Sin embargo, el an´alisis din´amico de una estructura no es directo. Los par´ametros modales se pueden determinar de manera experimental o por m´etodos num´ ericos. En general, los resultados experimentales entregan informaci´on de la estructura s´olo para la configuraci´on experimental. Un modelo en elementos finitos, en cambio, permite predecir el comportamiento din´ amico de la estructura bajo distintas condiciones de borde y carga, pero la confiabilidad de un modelo en elementos finitos muchas veces no est´a garantizada. Los m´etodos de ajuste de modelos permiten verificar y corregir estos modelos en elementos finitos por medio de los datos experimentales. El resultado de un ajuste de modelo es un modelo en elementos finitos que es m´as confiable para predicciones futuras. El ajuste de un modelo num´erico busca desarrollar un modelo en elementos finitos que entregue predicciones confiables de la din´amica de una estructura mec´anica. En la Figura 7.1 se ilustra el esquema general de un m´etodo de ajuste de modelos.
89
´ ´ CORRELACION NUMERICO-EXPERIMENTAL Y AJUSTE DE MODELOS
90
���������� ����
������ ��������
��������
���� ��������
�����
���� ��������� ��
�����������
��������� ��
��������� � �����
����������
�� �����������
�� ��� ����������� ��������� � ����� �� ���
�����������
�����
���
������ ��� ������ ��������
Figura 7.1: Esquema general de un ajuste de modelo El procedimiento comienza con la construcci´on de un modelo en elementos finitos. La estructura se divide en elementos conectados por nodos. Cada nodo tiene uno o m´ as grados de libertad. Los grados de libertad representan desplazamientos y deformaciones de la estructura en forma discretizada. Este modelo est´a representado por tres matrices: de rigidez K, masa M y amortiguaci´ on C. A partir de estas matrices se pueden determinar las propiedades modales de la estructura (ver secci´on 1.2). Por otro lado, se realiza un an´ alisis modal experimental de la estructura. En donde se miden las funciones de respuesta en frecuencia en distintos puntos de la estructura. Luego los par´ametros modales se pueden identificar a partir de un m´etodo de identificaci´on de par´ametros.
´ ´ CORRELACION NUMERICO-EXPERIMENTAL Y AJUSTE DE MODELOS
91
El ajuste de modelos comienza primero al parear los puntos de medici´ on experimentales con nodos del modelo num´erico. Generalmente, los puntos experimentales no coinciden completamente con los nodos num´ericos. Primero, los nodos num´ericos y los puntos experimental pueden diferir en su ubicaci´on f´ısica. Esto se puede solucionar con una buena comunicaci´on entre la persona que realiza el modelo num´erico y el que realiza las mediciones experimentales. Segundo, y m´as importante, un modelo en elementos finitos tiene muchos mas grados de libertad que los que se miden experimentalmente. Para algunas t´ ecnicas de correlaci´ on es necesario que los grados medidos experimentalmente y los grados del modelo num´erico coincidan. Para resolver esta incompatibilidad con las mallas, se puede o reducir el modelo num´erico a los grados de libertad medidos experimentalmente, o expandir los datos experimentales al n´umero de grados de libertad num´ericos. La reducci´on de matrices o expansi´on datos experimentales es un obst´aculo importante en ajuste de modelos. Luego de coincidir ambos modelos, el procedimiento continua con la correlaci´ on entre los datos num´ ericos y experimentales. Existen distintas t´ecnicas de correlaci´on que se pueden utilizar. Si la correlaci´on es buena el algoritmo termina, y el modelo en elementos finitos se considera suficientemente bueno. En caso contrario (y m´as probable), si la correlaci´on es mala, el modelo en elementos finitos debe ser corregido por medio de un procedimiento de ajuste. Los datos num´ericos con los experimentales muchas veces no son compatibles por alguna de las razones siguientes: Los grados libertad experimentales no coinciden con los grados de libertad del modelo en elementos finitos: Esto se puede solucionar
parcialmente con algoritmos de reduci´on/expansi´on. olo el n´ umero El conjunto de datos experimentales es incompleto: No s´ de grados experimentales esta limitado, sino tambi´ en el n´umero de modos que se pueden determinar experimentalmente. La medici´on de las funciones de respuesta en frecuencia se realiza en rango de frecuencias limitado. Por lo tanto, los modos fuera de este rango no pueden ser identificados. Dependiendo de la cantidad de informaci´on experimental, los resultados del ajuste de modelos pueden no ser u ´ nicos. Los datos experimentales est´ an contaminados con ruido: En general,
se esperan errores del orden del 3 % en las frecuencias naturales y del 10 % en las formas modales. El uso de informaci´on experimental contaminada con ruido puede llevar a una ajuste de modelos inadecuado. El amortiguamiento no se puede incluir de manera precisa en el on de amortiguamiento est´a inherentemente modelo num´ erico: La informaci´
presente en los datos experimentales, pero usualmente no se considera en el
´ PAREAR MODELOS NUMERICOS Y EXPERIMENTALES
92
modelo en elementos finitos. Esta discrepancia puede causar errores en el modelo ajustado.
7.1.
Parear modelos num´ ericos y experimentales
En la mayor´ıa de los casos de n´ umero de grados de libertad del modelo num´erico excede por mucho el n´umero de grados de libertad medidos. Los modelos en elementos finitos requieren mallas finas para poder entregar resultados precisos. No es practico , y muchas veces no es posible, medir todos los grados de libertad en la estructura real: Muchos de los grados de libertad son internos y no se puede acceder a ellos para realizar mediciones.
Los grados de libertad de rotaci´on son dif´ıciles de medir. Para el prop´ osito de an´alisis modal, no es necesaria una malla muy fina de puntos de medici´on. Sin embargo, muchos algoritmos de ajuste de modelos requieren una correspondencia uno-a-uno entre los grados de libertad anal´ıticos y experimentales. Para lograr esto se utilizan algoritmos de reducci´on de modelos num´ ericos y de expansi´on de los datos experimentales.
7.1.1.
T´ ecnicas de reducci´ on
Las t´ecnicas de reducci´on expresan las ecuaciones de movimiento en t´erminos de los grados de libertad medidos experimentalmente. En una primera fase, se identifican los grados de libertad anal´ıticos para cada grado de libertad experimental. Esto define a los grados de libertad “activos”. El resto de los grados de libertad anal´ıticos se denominan grados de libertad “eliminados”:
X f =
X a X d
donde: X a = grados de libertad activos X d = grados de libertad eliminados X f = set completo de grados de libertad
(7.1)
´ PAREAR MODELOS NUMERICOS Y EXPERIMENTALES
93
Las t´ecnicas de reducci´on establecen una relaci´on entre los grados de libertad activos y los grados de libertad eliminados por medio de una matriz de transformaci´on T d o T f : X d
= T d X a
(7.2)
X f = T f X a
(7.3)
I T d
T f =
(7.4)
T f se utiliza para construir matrices de masa y rigidez reducidas M r y K r respectivamente. Para un sistema sin amortiguamiento se cumple la siguiente relaci´ on: X f T M X f = X aT M r X a
(7.5)
X f T KX f = X aT K r X a
(7.6)
Sustituyendo X f por T f X a se obtiene: M r
= T f T M T f
(7.7)
K r
= T f T KT f
(7.8)
Las matrices del sistema reducido entregan s´olo una descripci´on aproximada del comportamiento din´amico “exacto” del sistema con las matrices completas.
Se sabe que todos los grados de libertad deben cumplir con la ecuaci´ on de movimiento:
K
2
− ω M
X f = Z (ω )X f = F f
(7.9)
Para el sistema reducido no deben haber fuerzas externas en los grados de libertad eliminados. Dividiendo el sistema de ecuaciones en los grados de libertad activos y los eliminados:
Z aa Z da
Z ad Z dd
X a X d
=
F a 0
(7.10)
despejando X d se obtiene que: X d =
−Z dd1 Z daX a −
(7.11)
´ PAREAR MODELOS NUMERICOS Y EXPERIMENTALES
94
De donde se deduce que: T d =
−1
−Z dd(ωr )
Z da (ωr )
(7.12)
Se debe notar que Z est´a evaluada a una frecuencia ω r , en consecuencia las matrices M r y K r dar´an resultados exactos para esta frecuencia. Para las otras frecuencias M r y K r dan s´olo resultados aproximados, siendo peores al alejarse de ωr . Un caso particular es la bien conocida reducci´on de Guyan, en donde se utiliza ωr = 0.
7.1.2.
T´ ecnicas de expansi´ on
las t´ecnicas de expansi´on expanden los modos experimentales desde el set de grados de libertad activos al set de grados de libertad completo del modelo anal´ıtico. La evaluaci´ on de los grados de libertad no medidos en los modos experimentales se lleva a cabo la mayor´ıa de las veces utilizando las relaciones entre los grados activos y eliminados del modelo num´ erico. El m´etodo de reducci´on presentado en la secci´on anterior se puede utilizar tambi´ en como un m´etodo de expansi´on ya que define la relaci´on entre los grados de libertad activos y eliminados en el modelo anal´ıtico: X d = T d X a
(7.13)
La misma matriz de transformaci´on se puede utilizar para los modos experimentales: φde = T d φae
(7.14)
Donde φde representa los grados de libertad eliminados de un modo (o matriz de modos) experimental y φae los grados de libertad activos.
A continuaci´on se presentan m´etodos de expansi´on alternativos.
Mezcla de vectores modales Este m´ etodo simplemente completa los grados de libertad faltantes en los modos experimental con los valores correspondientes de los modos anal´ıticos: φf e
=
φae φda
(7.15)
φf e representa un modo experimental (o matriz de modos) con el set completo de grados de libertad, φda es el modo anal´ıtico en los grados de libertad no medidos y φae es el modo experimental en los grados de libertad medidos.
´ ´ TECNICAS DE CORRELACION
95
Antes de completar el modos experimental con los valores del modo anal´ıtico, ´estos deben tener la misma escala de modo que:
φae = φaa
(7.16)
M´ etodo de coordenadas modales Este m´ etodo define los modos experimentales como una combinaci´on lineal de los modos anal´ıticos. Los coeficientes de la combinaci´on lineal son calculados con los grados de libertad activos: φae
= φaa q
a q = φa+ a φe
(7.17)
(7.18)
La matriz q se utiliza luego para evaluar los grados de libertad experimental faltantes: φde = φ da q
(7.19)
Interpolaci´ on Los grados de libertad faltantes tambi´en se pueden aproximar por medio de m´etodos de interpolaci´on. Una ventaja de este m´etodo es que no se necesitan datos num´ ericos para expandir los datos experimentales. Este m´etodo depende fuertemente de la conexi´on entre los distintos grados de libertad. En consecuencia, cada interpolaci´ on depende mucho del caso en estudio, lo que dificulta su uso pr´actico. En la pr´actica, estos m´ etodos se aplican en estructuras que consisten en componentes unidimensionales. Tambi´en se pueden utilizar para estimar grados de libertad de rotaci´ on a partir de los grados de libertad de traslaci´ on, pero este procedimiento es num´ ericamente poco estable.
7.2.
T´ ecnicas de correlaci´ on
El prop´osito de los m´ etodos de ajuste de modelo es lograr que el modelo num´erico se acerque lo m´ as posible a los datos experimentales. Para evaluar el nivel de correlaci´on entre los datos num´ ericos y experimentales es necesario definir que par´ametros se van a comparar y en que forma. A continuaci´on, se describir´ an algunas t´ecnicas usuales de correlaci´on num´erica-experimental.
´ ´ TECNICAS DE CORRELACION
7.2.1.
96
Diferencia en las frecuencias de resonancia
La manera m´as sencilla de verificar la correlaci´on es comparando las frecuencias naturales anal´ıticas y experimentales. La diferencia m´axima permitida depende en la precisi´on de las frecuencias naturales.
7.2.2.
Comparaci´ on visual de los modos
Una manera de visualizar la correlaci´ on entre dos modos son los denominados ◦ gr´ aficos de 45 . En estos gr´aficos se gr´afica un modo experimental versus el modo anal´ıtico correspondiente en un gr´afico x-y. Si ambos modos son id´enticos y tienen la misma escala, entonces todos los puntos se encuentran sobre una linea de 45 ◦ que parte desde el origen. La distancia entre los puntos y la l´ınea de 45◦ da una indicaci´on de la correlaci´on entre los modos. En la Figura 7.2 se ilustra un gr´afico de 45◦ t´ıpico.
� � � � � � � � � � � � � � � � ���� ���������
Figura 7.2: Gr´afico de 45◦ Una manera de asegurar que los modos experimentales tengan la misma escala que los modos anal´ıticos es utilizando el “Factor de Escala Modal (MSF)”. El MSF mide el factor de escala entre dos modos, los modos experimentales se pueden escalar a los modos anal´ıticos multiplicados por el MSF correspondiente: φ∗e,i M SF i
= φe,i M SF i
·
=
φT a,i φe,i φT e,i φe,i
(7.20) (7.21)
´ ´ TECNICAS DE CORRELACION
97
Donde φ∗e,i y φe,i es el i´esimo modo experimental escalado y original, φa,i es el i´esimo modo anal´ıtico, y M SF i es el factor de escala modal para el modo i. Al multiplicar el modo experimental por el MSF correspondiente, tambi´ en se soluciona el problema que los modos anal´ıtico y experimental pueden estar con un desfase de 180◦ .
7.2.3.
Modal Assurance Criterion
El “Modal assurance criterion (MAC)” expresa la correlaci´on que es visualizada en un gr´afico de 45◦ en un s´olo n´ umero. Se define como:
M AC ij =
φT a,i φe,j
φT a,i φa,i
2
φT e,j φe,j
(7.22)
Donde φa,i es el iesimo modo anal´ıtico y φe,j es el jesimo modo experimental. Un valor de 0 indica que no hay correlaci´ on, mientras que un valor de 1 indica dos modos perfectamente correlacionados. Al ordenar todos los valores M AC ij en una matriz, la diagonal deber´ıa tener valores altos (en general mayores a 0.8) para una buena correlaci´on. Una ventaja del MAC es que la correlaci´ on no depende de la escala de los modos, sino que s´olo de forma de ´estos. En la Figura 7.3 se representa gr´aficamente una matriz de valores MAC t´ıpica. 1
116 0.9
58.9 ) z 56.6 H (
0.8
s 44.7 e l a t 39 n e m31.2 i r e p28.4 x E s 26.3 o d24.2 o M 19
0.7 0.6 0.5 0.4 0.3 0.2
15.4
0.1
13 12.7 15.2 19.4 24.2 26.3 28.3 31.4 38.4 43.8
57 58.7 116
Modos numéricos (Hz)
Figura 7.3: Matriz de valores MAC
´ METODOS ITERATIVOS DE AJUSTE DE MODELOS
7.2.4.
98
Frequency Response Assurance Criterion
Tal como el MAC mide la correlaci´ on de dos vectores modales, el “Frequency Response Assurance Criterion (FRAC)” mide la correlaci´on entre dos funciones de respuesta en frecuencia. Se define de manera equivalente al MAC como:
FRAC pq =
2 ω2 a e H ω H ω ( ) ( ) pq ω=ω1 pq ω2 ω 2 a a e e ω=ω1 H pq (ω )H pq (ω ) ω=ω1 H pq (ω )H pq (ω )
(7.23)
Donde H pq es la funci´ on de respuesta en frecuencia cuando la excitaci´ on es en el punto p y se mide en el punto q. Los indices a y e se refieren a anal´ıtico y experimental respectivamente. Por ´ultimo, ω1 y ω2 es el rango de frecuencias en que se quiere correlacionar las FRFs. Como en el caso del MAC, un valor del FRAC igual a 1 equivale a dos FRFs perfectamente correlacionadas, mientras que un valor cercano a 0 indica una mala correlaci´ on.
7.3.
M´ etodos iterativos de ajuste de modelos
Los m´etodos m´as utilizados en ajuste de modelo son m´etodos iterativos basado en gradientes. Estos m´etodos minimizan una funci´on error ε que relaciona los par´ ametros experimentales con los entregados por el modelo num´erico. Esta funci´on de error puede estar compuesta de distintos residuos: La diferencia entre las frecuencias naturales anal´ıticas y experimentales. La diferencia entre los modos anal´ıticos y experimentales. El balance de fuerzas a una frecuencia especifica. La diferencia entre las funciones de respuesta en frecuencia anal´ıticas y experimentales. El problema se resuelve con una aproximaci´on de Taylor de primer orden
ε(β ) = ε (β 0 ) +
i
∂ε ∆β i ∂β i
(7.24)
´ METODOS ITERATIVOS DE AJUSTE DE MODELOS
99
Donde β = β 1 , β 2 , . . . βm es un vector de par´ ametros a ser actualizados. ∆β i = j+1 j j β i β i , β i es el valor actual del par´ametro y β ij+1 es el valor estimado.
{
−
}
Las derivadas parciales de primer orden son evaluadas num´ericamente. Definiendo una matriz S que contenga las derivadas parciales del residuo con respecto a los par´ ametros β :
S =
∂ε1 ∂β 1 ∂ε2 ∂β 1 .. . ∂εn ∂β 1
∂ε1 ∂β 2 ∂ε2 ∂β 2 .. . ∂εn ∂β 2
... ...
..
.
...
∂ε1 ∂β n ∂ε2 ∂β n ... ∂εn ∂β n
(7.25)
la ecuaci´on 7.24 se puede escribir como: ∆ε = S ∆β,
∆ε = z E
− zA(β )
(7.26)
Donde z E and z A (β ) son los datos experimentales y anal´ıticos respectivamente, ∆ ε es el vector de residuos. El vector de par´ametros β es estimado al minimizar la siguiente funci´on objetivo, J (β ) = (S ∆β
− ∆ε)T (S ∆β − ∆ε)
(7.27)
Este problema de minimizaci´on se resuelve por m´ınimos cuadrados. La soluci´on viene dada por, ∆β = S + ∆ε
(7.28)
Donde S + es la pseudoinversa de la matriz S : S n+×m
=
S −1 (S T S )−1 S T S T (S T S )−1
n = m n>m n
(7.29)
donde n es el n´ umero de mediciones y m el n´ umero de par´ametros. La ecuaci´on 7.28 se puede escribir de manera equivalente como, β j+1 = β j + S (β j )+ ∆ε(β j )
(7.30)
Donde ∆β i = β ij+1 β ij , β ij es el valor actual del par´ ametro y β ij+1 es el valor estimado. El valor final de β se obtiene por un procedimiento iterativo; la ecuaci´on 7.30 es iterada hasta que se alcanza convergencia.
−
´ METODOS ITERATIVOS DE AJUSTE DE MODELOS
100
La ecuaci´ on 7.27 pondera de igual forma a todos los datos experimentales. Sin embargo, ser´ıa preferible dar mayor peso a los datos que se saben son m´as precisos. Por ejemplo, se sabe que las frecuencias naturales se pueden determinar de manera m´ as precisa que los modos, o que los primeros modos contienen menos contaminaci´on por ruido que los modos de frecuencias m´as altas. Para compensar estas diferencias se pueden utilizar matrices de ponderaci´on. La funci´ on a minimizar considerando una matriz de ponderaci´on viene dada por, J (β ) = (S ∆β
− ∆ε)T W εε (S ∆β − ∆ε)
(7.31)
Donde W εε es una matriz definida positiva. En general, se define W εε como una matriz diagonal cuyos elementos son la inversa de la varianza de los datos experimentales. La soluci´on del problema de minimizaci´on viene dada por, j+1
β
j
j T
j
−1
= β + S (β ) W εε S (β )
S (β j )T W εε ∆ε(β j )
(7.32)
El problema iterativo definido en la ecuaci´on 7.32 puede estar mal condicionado. En un problema mal condicionado no se puede asegurar la convergencia. Una soluci´on al problema de mal condicionamiento es la regularizaci´ on de las ecuaciones. La regularizaci´ on de Tikonov es el m´ etodo mas utilizado. En este caso, se a˜nade un termino adicional a la funci´on objetivo para penalizar por el cambio relativo de los par´ametros. De esta manera se incentivan cambios peque˜nos de los par´ametros en cada paso, favoreciendo la convergencia. La funci´on objetivo queda como, J (β ) = (S ∆β
− ∆ε)T W εε (S ∆β − ∆ε) + ∆β T W ββ ∆β
(7.33)
Donde W ββ = I α es una matriz diagonal cuyos elementos son iguales al par´ametro de regularizaci´on α. La soluci´on a la ecuaci´on 7.33 es,
β j+1 = β j + S (β j )T W εε S (β j ) + W ββ
−1
S (β j )T W εε ∆ε(β j )
(7.34)
El valor del par´ametro de regularizaci´ on α se selecciona utilizando la curva ‘L’ de Hansen. La curva de ‘L’ es un gr´ afico de la norma de ∆β T W ββ ∆β versus la T norma del residuo (S ∆β ∆ε) W εε (S ∆β ∆ε). El valor ´optimo del par´ametro de regularizaci´on corresponde al punto de m´axima curvatura en la curva del gr´afico log-log de la curva de ‘L’.
−
−
Existen algunos m´etodos est´andar para construir la matriz de gradientes S . A continuaci´ on se describir´an estos m´ etodos cuando los siguientes datos son utilizados: (1) modos y frecuencias naturales y (2) funciones de respuesta en frecuencia.
´ METODOS ITERATIVOS DE AJUSTE DE MODELOS
7.3.1.
101
Modos y frecuencias naturales
Un m´ etodo bien establecido es el “inverse eigen-sensitivity method” (IESM). El que correlaciona valores y modos propios. La funci´on de error se define como, λE,k λA,k λE,k
−
ελ,k
=
εφ,k
= φE,k
(7.35)
− φA,k
(7.36)
Donde λk , φk son los kesimo valor y modo propio respectivamente. Los sub´ındices A, E se refieren a anal´ıtico y experimental respectivamente.
Las derivadas de los valores y modos propios, se pueden calcular a partir de las relaciones dadas por Fox y Kapoor (1968),
∂λA,k ∂K = φ T A,k ∂β i ∂β i
∂φA,k = ∂β i
∂M λA,k φA,k ∂β i
−
N
φT A,j
=r j=1,
∂K ∂M λA,k φA,k ∂β i ∂β i (λA,k λA,j )
−
−
φA,j
(7.37)
(7.38)
K y M son las matrices de rigidez y masa del modelo num´erico. La ecuaci´on 7.38 requiere el calculo de todos los modos propios. Aunque, las derivadas se pueden aproximar con los primeros t´erminos de la serie en el ecuaci´ on 7.38.
Las derivadas de las matrices de rigidez y masa, utilizando diferencias finitas:
∂K = ∂β i
K (β i + δ ) δ
∂M ∂β i
M (β i + δ ) δ
=
− K (β i) − M (β i)
donde δ es un n´ umero peque˜ no.
∂K ∂M y , se pueden estimar ∂β i ∂β i
(7.39)
(7.40)
´ METODOS ITERATIVOS DE AJUSTE DE MODELOS
102
La matriz de gradientes S y el vector de residuos viene dados por,
S =
∂λA,1 /λA,1 , ∂β 1 ∂φA,1 ∂β 1 .. . ∂λA,m /λA,m , ∂β 1 ∂φA,m ∂β 1
...
..
.
... ...
(λE,1 λA,1 )/λA,1 φE,1 φA,1 .. .
−
∆ε =
7.3.2.
...
∂ λA,1 , /λA,1 ∂β n ∂φA,1 ∂β n .. . ∂ λA,m , /λA,m ∂β n ∂φA,m ∂β n
−
(λE,m λA,m )/λA,m φE,m φA,m
−
−
(7.41)
(7.42)
Funciones de respuesta en frecuencia
El m´ etodo denominado “response function method (RFM)” propuesto por Lim y Edwins (1990) utiliza directamente los funciones de respuesta en frecuencia. La funci´on de error est´a definida como, εα = α jA
− αE j
(7.43)
Donde αj representa la jesima columna de la matriz de FRFs α. Los sub´ındices A, E se refieren a anal´ıtico y experimental respecticamente.
Lim y Edwins demostraron la siguiente relaci´on, ∂
αjA
−
j αE
∂β i
= α A
∂Z j α ∂β i E
(7.44)
Z es la matriz de rigidez din´amica. Sus gradientes a una frecuencia dada, ω , se pueden determinar como,
∂ K ω 2 M ∂Z ∂K = = ∂β i ∂β i ∂β i
−
− ω2 ∂M ∂β i
(7.45)
K y M son las matrices de rigidez y masa del modelo num´erico. Los gradientes de ∂K ∂M las matrices de rigidez y masa, and , se pueden determinar por diferencia ∂β i ∂β i finitas, como se describe en la ecuaciones 7.39 y 7.40.
´ METODOS ITERATIVOS DE AJUSTE DE MODELOS
103
La matriz de gradientes S y el vector de residuos viene dados por,
S =
− −
∆ε =
∂Z j αA (ω1 ) α (ω1 ), ∂β 1 E .. . ∂Z j αA (ωnf ) α (ωnf ), ∂β 1 E
αjE (ω1 ) j αE (ωnf )
− .. .
j αA (ω1 )
− αAj (ωnf )
∂ Z j ... , αA (ω1 ) α (ω1 ) ∂β n E .. .. . . ∂ Z j . . . , αA (ωnf ) α (ωnf ) ∂β n E
−
−
(7.46)
(7.47)
Donde ω1 , ω2 , . . . ωnf representan las frecuencias seleccionadas para definir la funci´on de error.
Cap´ıtulo 8
Ejemplos de An´ alisis Modal en Estructuras Reales 8.1.
Estructura de barras
La estructura consiste en un ensamble tridimensional de barras est´ aticamente indeterminado, compuesto por 43 barras y 20 uniones. Las dimensiones de la estructura son: largo 3m, ancho 0.5m y alto 0.5m. Las barras est´an formadas por tubos de aluminio con un di´ametro de 22mm y un espesor de 1mm. En la Figura 8.1 se muestra el montaje experimental, la estructura esta suspendidas por resortes blandos para simular una condici´on de borde libre. Los puntos de medici´ on, excitaci´ on y suspension se definieron con la ayuda de un modelo num´erico. Este modelo utiliza elementos de viga 3D para modelar las barras y masas concentradas en las uniones. La configuraci´on o´ptima de puntos de medici´on se defini´o al minimizar los valores fuera de la diagonal de la matriz MAC. Este procedimiento asegura que los modos sean independientes entre s´ı. Para efectos de visualizaci´on, se utilizaron aceler´ometros tri-axiales en los puntos seleccionados. Como resultado se midieron 20 3 grados de libertad, un aceler´ometro triaxial por uni´on.
×
La ubicaci´ on de los puntos de excitaci´ on y suspensi´ on se seleccion´o de acuerdo a un estudio de los driving point residues (DPR). Los puntos de excitaci´ on son los grados de libertad con DPR altos para el mayor numero de modos. En la Figura 8.2 se muestra los DPR promedio (en los primeros 12 modos) para los 60 grados de libertad medidos. Los primeros cuatro m´aximos corresponden a las cuatro esquinas inferiores de la estructura cuando son excitadas en la direcci´ on
104
ESTRUCTURA DE BARRAS
105
x. Los cuatro m´aximos siguientes corresponden a las esquinas superiores cuando son excitadas en la direcci´on z. De acuerdo a estos resultados se seleccionaron dos puntos de excitaci´on (ver Figura 8.3(a)), el primero en una esquina inferior en la direcci´ on x y el segundo en una esquina superior en la direcci´on z .
z y
x
Figura 8.1: Montaje experimental de la estructura de barras
4.5
x 10
-3
↓ ↓ ↓ ↓
4 3.5 3 R P D2.5 e g a r 2 e v A 1.5
↓
↓
↓
↓
1 0.5 0 0
10
20
30
40
50
60
Degree of Freedom
Figura 8.2: Driving point residues promedio para los posibles puntos de excitaci´on
ESTRUCTURA DE BARRAS
106
z y x
(a) Puntos de excitaci´ on
(b) Puntos de suspensi´ on
Figura 8.3: Puntos de excitaci´on y suspensi´on seleccionados Los puntos de suspensi´ on se seleccionaron buscando puntos donde fuese factible suspender la estructura y que tuviesen un DPR promedio bajo. En la Figura 8.4 se muestran los puntos candidatos para la suspensi´on y sus valores DPR. La ubicaci´on f´ısica de los puntos seleccionados se puede visualizar en la Figura 8.3(b). El rango de frecuencias de medici´on es de 0 a 256Hz con una resoluci´on de 0.25Hz. Este rango contiene todos los modos globales y unos pocos modos locales. En la Figura 8.5 se muestran los 12 modos globales identificados. 4.5
x 10
-3
4 3.5 3 R P D2.5 e g a r 2 e v A 1.5 1 0.5 0 0
↓ 10
↓ 20
↓ 30
40
↓ 50
60
Degree of Freedom
Figura 8.4: Driving point residues promedio para los posibles puntos de suspensi´on
ESTRUCTURA DE BARRAS
107
(a) Modo 1, ω = 13,03Hz
(b) Modo 2, ω = 15,37Hz
(c) Modo 3, ω = 19,05Hz
(d) Modo 4, ω = 24,23Hz
(e) Modo 5, ω = 26,31Hz
(f) Modo 6, ω = 28,42Hz
(g) Modo 7, ω = 31,17Hz
(h) Modo 8, ω = 39,04Hz
(i) Modo 9, ω = 44,71Hz
(j) Modo 10, ω = 56,64Hz
(k) Modo 11, ω = 58,87Hz
(l) Modo 12, ω = 115,6Hz
Figura 8.5: Doce modos globales identificados
ESTRUCTURA DE BARRAS
8.1.1.
108
Modelo num´ erico
El modelo num´ erico se construy´o en Matlab utilizan elementos de viga 3D para las barras e inercias concentradas para las uniones. Cada barra se model´o con cuatro elementos de viga, como se muestra en la Figura 8.6. Los grados de libertad de rotaci´ on fueron eliminados por medio de una metodolog´ıa de reducci´on de Guyan. En total, el modelo tiene 447 grados de libertad y 192 elementos. Se identifican 12 modos globales con frecuencias entre 12Hz y 125Hz, los modos siguientes son modos locales de las barras.
1
2
3
4
Figura 8.6: Modelo en elementos finitos de cada barra
En la Figura 8.7 se muestra la correlaci´ on entre los modos num´ ericos y experimentales, donde ωE y ωA son las frecuencias naturales experimentales y anal´ıticas respectivamente. El valor MAC indica la correlaci´ on entre las formas modales num´ ericas y experimentales. La m´axima diferencia en frecuencia es del 2,3 % y el MAC m´ınimo es 0,86. Es importante destacar que la correlaci´on num´ erico-experimental mostrada, es una vez realizado el ajuste del modelo num´erico.
ESTRUCTURA DE BARRAS
ω
109
= 12.72 Hz
ω
= 13.03 Hz
ω
A
ω
E
MAC = 0.97
ω
= 15.37 Hz
ω
MAC = 0.98
ω
ω
= 24.23 Hz
ω
E
MAC = 0.98
= 31.43 Hz
A
ω
E
= 26.31 Hz
ω
E
MAC = 0.92
ω
ω
= 56.64 Hz
MAC = 0.88
= 43.84 Hz
= 39.04 Hz
ω
A
= 44.71 Hz
E
MAC = 0.97
(i) Modo 9 = 115.8 Hz
= 58.8 Hz
ω
= 58.87 Hz
ω
A
E
MAC = 0.86
Numérico
(j) Modo 10
(f) Modo 6
(h) Modo 8
E
MAC = 0.93
ω
MAC = 0.98
= 56.93 Hz
A
= 28.42 Hz
E
= 38.43 Hz
(g) Modo 7
= 28.33 Hz
A
(e) Modo 5
E
MAC = 0.92
(c) Modo 3 ω
A
= 31.17 Hz
MAC = 0.99
= 26.29 Hz
A
(d) Modo 4 ω
= 19.05 Hz
E
(b) Modo 2
= 24.25 Hz
A
A
E
(a) Modo 1 ω
= 19.38 Hz
= 15.25 Hz
A
(k) Modo 11
A
= 115.6 Hz
E
MAC = 0.97
Experimental
(l) Modo 12
Figura 8.7: Comparaci´on entre los modos num´ ericos y experimentales
MODELO A ESCALA DE UN AVION
8.2.
110
Modelo a escala de un avion
La estructura simula el comportamiento din´ amico de un avion. Est´ a formada por vigas de acero. El fuselaje consiste en una viga recta de secci´on rectangular, mientras que vigas planas forman las alas y cola. Las turbinas est´an representadas por dos masas conectadas a las alas. El montaje experimental se ilustra en la Figura 8.8. La estructura est´a suspendida por tres resortes y un shaker electrodin´amico la excita con una se˜nal aleatoria en al punta del ala. La respuesta es capturada por 11 aceler´ometros.
Figura 8.8: Montaje experimental del avion a escala
El test es realizado en una rango de frecuencias de 0 130Hz. Se obtuvieron modos operacionales mediante en m´ etodo en el dominio del tiempo (stochastic subspace identification method). Los modos identificados se muestran en la Figura 8.9.
−
MODELO A ESCALA DE UN AVION
111
(a) Modo 1, ω = 18,44Hz
(b) Modo 2, ω = 39,84Hz
(c) Modo 3, ω = 81,34Hz
(d) Modo 4, ω = 83,44Hz
(e) Modo 5, ω = 91,18Hz
(f) Modo 6, ω = 99,32Hz
Figura 8.9: Primeros seis modos experimentales del avion a escala
8.2.1.
Modelo num´ erico
El modelo num´erico se construy´ o en Matlab con elementos de viga 3D para el fuselaje, alas y cola e inercias concentradas para las turbinas. El modelo num´erico tiene 43 elementos y 252 DOF, como se ilustra en la Figura 8.10. Las propiedades iniciales del modelo num´erico fueron ajustadas para coincidir con los modos y frecuencias naturales. La Figura 8.11 presenta la correlaci´on entre los modos num´ ericos y experimentales, donde ωE y ωA son las frecuencias naturales experimentales y anal´ıticas respectivamente. El MAC mide la correlaci´on entre los modos num´ ericos y experimentales. El MAC m´ınimo tiene un valor de 0,82 y la diferencia m´axima de frecuencias es del 2 ,0 %.
MODELO A ESCALA DE UN AVION
112
4 0
3 7
1 6
1 5
1 4
1 3
1 2
1 1
1 0
9
42
1 7
1 8
1 9
2 1
2 0
2 2
2 3
2 4
2 5
2 7
2 6
2 8
2 9
3 9
3 6
3 0
3 5
3 8
3 1
3 4
3 3
3 2
21
9 41
8
7
20
6
5
4
3
2
8
1
43
Figura 8.10: Modelo en elementos finitos
MAC = 0.98 ω
= 18.43 Hz
A
,
ω
= 18.44 Hz
E
ω
= 39.8 Hz
,
A
ω
= 39.84 Hz
E
(a) Modo 1
(b) Modo 2
MAC = 0.82
MAC = 0.91
= 82.75 Hz, = 83.44 Hz
A
MAC = 0.99
E
= 89.83 Hz, = 91.18 Hz
A
E
Numérico (d) Modo 4
(e) Modo 5
MAC = 0.93 ω
= 82.17 Hz
A
,
ω
= 81.34 Hz
E
(c) Modo 3
(f) Modo 6
Figura 8.11: Comparaci´on entre los modos num´ ericos y experimentales
VIGA DE CONCRETO REFORZADO
8.3.
113
Viga de concreto reforzado
La estructura consiste en una viga de concreto de largo 6m y una secci´on de 0,2 0,25m2 . Seis vigas de acero con un di´ametro de 16mm refuerzan la viga. 30 vigas de acero con un di´ametro de 8mm refuerzan transversalmente la viga, como se muestra en la Figura 8.12. El peso total de la viga es de 750kg.
×
Figura 8.12: Secci´on de la viga y refuerzos
Figura 8.13: Montaje experimental para la viga reforzada
VIGA DE CONCRETO REFORZADO
114
La Figura 8.13 muestra el montaje experimental. La viga est´ a suspendida por resortes blandos, simulando una condici´ on de borde libre. La aceleraci´ on verticalmente en 62 puntos distribuidos equidistantes en la cara superior a lo largo de la viga. Un martillo de impacto excita a la viga en uno de las esquinas, excitando tambi´ en los modos torsionales. Las mediciones son realizadas en un rango de frecuencias de 0 625Hz . Se identifican los modos operacionales por medio de un algoritmo en el dominio del tiempo (stochastic subspace identification method). La Figura 8.14 muestra los primeros seis modos identificados.
−
(a) Modo 1, ω=22.76Hz
(b) Modo 2, ω=64.85Hz
(c) Modo 3, ω =126.9Hz
(d) Modo 4, ω=186.4Hz
(e) Modo 5, ω=208.4Hz
(f) Modo 6, ω=302.8Hz
Figura 8.14: Primeros seis modos experimentales
8.3.1.
Modelo num´ erico
El modelo num´erico se construy´ o en Matlab utilizando elementos de viga 2D, como se muestra en la Figura 8.15. Tiene 30 elementos de viga y 62 grados de libertad. Luego de ajustar el modelo num´ erico las propiedades de la viga son: Modulo de Young 3,60 1010 P a, densidad 2500kg/m3 , area 0,05m2 e inercia de la secci´on 1,93 10−4 m4 . El modelo num´erico no incluye grados de libertad torsionales. Por lo tanto, solo se consideran los modos de flexi´on. La Figura 8.16 presenta la
·
·
VIGA DE CONCRETO REFORZADO
115
correlaci´ correlac i´on on num´ erica-experimental erica-exp erimental de los primeros cuatro modos de flexi´on. on. La m´ axima diferencia en frecuencias es del 2 ,46 % y el MAC m´ axima m´ınimo es 0,996. x
1
2
x
x
3
x
4
x
5
x
6
7
x
x
8
x
9
x
10 11 12 13 14 15 16
x
x
x
x
x
x
x
17 18 19 20
x
x
x
x
21 22
x
x
23 24 25 26
x
x
x
x
27 28 29 30
x
x
x
x
Figura 8.15: Modelo en elementos finitos de la viga reforzada
ω
= 23.32 Hz,
A
ω
E
= 22.76 Hz, MAC= 0.99 8
= 64.29 Hz,
ω
(a) Modo 1
= 126 Hz,
A
E
E
= 64.85 Hz, MAC= 0.997
(b) Modo 2
= 126.9 Hz, MAC= 0.997
Numérico (c) Modo 3
ω
A
= 208.4 Hz,
A
E
= 208.4 Hz, MAC= 0.996
Experimental (d) Modo 4
Figura 8.16: Comparaci´on on entre los modos num´ ericos ericos y experimentales exper imentales
TUBO DE ESCAPE DE UN AUTO
8.4. 8.4.
116
Tubo de esca escape pe de de un un aut auto o
En la Figura 8.17 Figura 8.17 se se muestra el montaje experimental. La estructure corresponde a un tubo de escape de un auto. Las dimensiones son: largo 2.3m 2.3 m, ancho 0.45m 0.45m. El tubo de escape esta hecho de acero y tiene un di´ametro de 38mm 38mm.. Resortes blandos suspenden a la estructura y un shaker electrodin´amica la excita. La respuesta es capturada por 20 aceler´ ometros ometros distribuidos en la estructura. Las mediciones se realizan en un rango de frecuencias de 0-1024 H z con una resoluci´on on en frecuencias de 0.25H z . La Figura 8.18 Figura 8.18 muestra muestra los primeros cuatro modos identificados.
Figura 8.17: Montaje experimental del tubo de escape
(a) Modo 1, ω=27.2Hz
(b) Modo 2, ω=63.4Hz
(c) Modo 3, ω=169Hz
(d) Mo do 4, ω=291Hz
Figura 8.18: Primeros cuatro modos experimentales
TUBO DE ESCAPE DE UN AUTO
8.4.1.
117
Modelo num´ erico erico
El modelo mo delo num´erico erico se construyo const ruyo en Matlab, Ma tlab, el tubo tub o de escape se aproxim´o como la uni´ on de varios tubos de diferentes di´ametros, on ametros, el espesor de los tubos fue definido inicialmente como 0.8mm 0.8 mm.. En la Figura 8.19 Figura 8.19 se se muestra la geometr´ıa ıa aproximada, aproximad a, las dimensiones est´an an en mil´ mil´ımetros y los ´angulos angulos en grados.
2 . 5 4
300
5 0 2
5 9 1
10
2 3 5
45 57
75
81
75
60 47
38
4 4
38
0 3
4 4
597
140
20 40 45
190
60 45 25
360
15
4 4
3 0
Figura 8.19: Geometr´ Geometr´ıa aproximada aproximada utilizada en el modelo n´umerico umerico
El modelo num´ num´erico erico mostrado en la Figura 8.20 Figura 8.20 fue construido en Matlab con elem elemen ento toss de viga viga 2D e iner inerci cias as conc concen entr trad adas as para para las las masa masas. s. Se util utiliz izar aron on 11 propieda propiedades des gen´ gen´ericas ericas para el acero; acero; Modulo de Young 2,1 10 P a, densidad 3 7800kg/m 7800kg/m , coeficiente de Poisson 0, 0 ,3. El modelo tiene 47 elementos de viga y 3 inercias, con 144 grados de libertad.
·
Se ajust´o el espesor del tubo, el peso de las inercias y las propiedades del material para coincidir con los modos experimentales. En la Figura 8.21 se presenta la correlaci´ correlac i´on on num´ erica-experimental erica-exp erimental de los primeros cuatro modos de flexi´on. on. La m´ axima axima diferencia en frecuencias es del 2.78 % y el MAC m´ m´ınimo es 0.98.
TUBO DE ESCAPE DE UN AUTO
118
Figura 8.20: Modelo en elementos finitos
ω =27.2 E
Hz,
ω =26.5 A
Hz, MAC= 0.988
ω
=63.4 Hz,
(a) Modo 1
=169 Hz,
E
A
=63.5 Hz, MAC= 0.986
A
(b) Modo 2
=171 Hz, MAC= 0.994
Numérico (c) Modo 3
ω
E
=291 Hz,
E
A
=287 Hz, MAC= 0.98
Experimental (d) Modo 4
Figura 8.21: Correlaci´on num´ erico-experimental de los modos
EL PUENTE I-40
8.5.
119
El puente I-40
La estructura corresponde al antiguo puente I-40 sobre Rio Grande en New Mexico. En la Figura 8.22 se ilustra una vista en elevaci´on de la porci´on del puente que fue sometida a un an´ alisis modal. La secci´ on del puente consiste en tres tramos. Los dos tramos en los extremos son de igual largo, 39m, el tramo central tiene 49.7m de largo. En la Figura 8.23 se muestra una vista en secci´ on del puente, el puente esta conformado por una secci´on de concreto soportado por dos vigas laterales de acero y tres vigas intermedias. Las cargas de las vigas intermedias se transmiten a las vigas laterales por medio de vigas cruzadas ubicadas en intervalos de 6.1m. 39.9 m
49.7 m
39.9 m
y
Pier 3
Pier 2
z
Pier 1
Figura 8.22: Vista en elevaci´on de la secci´on del puente estudiada
0.178m
Stringers (21 WF 62)
L 5x5x5/16 Bracing m 5 0 . 3
r e d r i G e t a l P
2.06m
36 WF 182 or 36 WF 150
Floor Beam Plate Girder
2.29m
2.29m
2.29m
2.29m
Figura 8.23: Secci´on del puente I-40
2.06m
EL PUENTE I-40
120
Durante las mediciones el puente fue excitado por medio de un shaker hidr´aulico y se midi´o la respuesta por medio de 26 aceler´ometros como se muestra en la Figura 8.24. La estructura fue excitada por una se˜nal aleatoria de frecuencias entre 2Hz y 12Hz. Los modos se obtienen a partir de un m´ etodo de identificaci´on en el dominio del tiempo. En la Figura 8.25 se muestran los primeros seis modos experimentales.
Shaker
Figura 8.24: Mediciones experimentales
(a) Modo 1, ω=2.48Hz
(b) Modo 2, ω =2.95Hz
(c) Modo 3, ω=3.49Hz
(d) Modo 4, ω=4.08Hz
(e) Modo 5, ω=4.17Hz
(f) Modo 6, ω =4.64Hz
Figura 8.25: Primeros seis modos experimentales
EL PUENTE I-40
8.5.1.
121
Modelo num´ erico
El modelo num´ erico se construy´o utilizando una secci´on aproximada del puente, la que se muestra en la Figura 8.26, la que consiste en los siguientes componentes:
1. Una secci´on de concreto de espesor constante, con un area equivalente a la mostrada en la Figura 8.23. 2. Dos vigas de acero laterales. 3. Tres vigas de acero intermedias. 4. Vigas de acero cruzadas. El modelo num´ erico se construy´o en Matlab. La secci´on de concreto y las vigas laterales se modelaron con elementos de placa. Las vigas de acero intermedias y cruzadas se modelaron con elementos de viga 3D. Se utilizaron propiedades de material gen´ericas:
· 1010P a, densidad: ρ =
Concreto: Modulo de Young: E = 3,0
coeficiente de Poisson: ν = 0,2.
2300kg/m3 ,
Acero: Modulo de Young: E = 2,1 1011 P a, densidad: ρ = 7800kg/m3 , coeficiente
·
de Poisson: ν = 0,3.
261’’
21 WF 62
8.7’’
36 WF 182 or 36 WF 150 ’ ’ 0 2 1
90’’
90’’
0.375’’
Y 21’’ or 24’’ 1.5” or 2.625” (top & bottom)
Figura 8.26: Secci´on aproximada del puente
X
EL PUENTE I-40
122
La Figura 8.27 muestra el modelo num´ erico. Tres resortes conectados en la parte superior del ´ ultimo pilar simulan la rigidez a˜nadida por la secci´on siguiente del puente que comparte este pilar. Las constantes de los resortes vienen dadas por: kx = 3,17 106 N/m
· ky = 1,26 · 106 N/m kz = 4,29 · 106 N/m
Donde kx , ky y kz son las rigideces de los resortes en las direcciones x, y y z. En la parte inferior de los pilares se restringen los 6 grados de libertad. En el extremo derecho del puente se restringen todos los grados de libertad de traslaci´ on y las rotaciones con respecto a los ejes Y y Z. Se eliminaron los grados d libertad del modelo num´erico por medio del m´etodo de Guyan, como resultado se tienen 711 grados de libertad. Las propiedades iniciales del modelo num´erico fueron ajustadas para coincidir con los resultados experimentales. En la Figura 8.28 se presenta la correlaci´on num´ericaexperimental de los primeros seis modos. La m´axima diferencia en frecuencias es del 2.45 % y el MAC m´ınimo es 0.979.
Figura 8.27: Modelo en elementos finitos
EL PUENTE I-40
123
= 2.48 Hz,
ω
A
= 2.48 Hz, MAC= 0.99 7
ω
= 3.02 Hz,
ω
A
E
(a) Modo 1
= 3.58 Hz,
ω
A
E
= 4.18 Hz,
ω
A
(c) Modo 3
= 4.14 Hz,
A
= 4.08 Hz, MAC= 0.979
ω
E
(d) Modo 4
= 4.17 Hz, MAC= 0.982
E
éo
(e) Modo 5
E
(b) Modo 2
= 3.50 Hz, MAC= 0.994
ω
= 2.96 Hz, MAC= 0.992
ω
= 4.70 Hz,
A
= 4.63 Hz, MAC= 0.981
E
(f) Modo 6
Figura 8.28: Comparaci´on entre los modos num´ ericos y experimentales