UNIVERSIDADE TECNOLÓGICA FEDERAL DO PARANÁ
CAMPUS DE CURITIBA
DEPARTAMENTO DE PESQUISA E PÓS-GRADUAÇÃO
PROGRAMA DE PÓS-GRADUAÇÃO EM ENGENHARIA MECÂNICA
E DE MATERIAIS - PPGEM
Francisco José Doubrawa Filho
CONTROLE DE VIBRAÇÃO FLEXIONAL EM SISTEMAS
GIRANTES UTILIZANDO NEUTRALIZADORES DINÂMICOS
VISCOELÁSTICOS
Curitiba
JUNHO - 2008
Francisco José Doubrawa Filho
CONTROLE DE VIBRAÇÃO FLEXIONAL EM SISTEMAS
GIRANTES UTILIZANDO NEUTRALIZADORES DINÂMICOS
VISCOELÁSTICOS
Dissertação apresentada como requisito parcial
à obtenção do título de Mestre em Engenharia,
do Programa de Pós-Graduação em Engenharia Mecânica e de Materiais, Área de concentração em Mecânica dos Sólidos e Vibrações,
do Departamento de Pesquisa e Pós Graduação,
do Campus de Curitiba da UTFPR.
Orientador:
Prof. Carlos Alberto Bavastri, Dr. Eng.
Curitiba
JUNHO - 2008
TERMO DE APROVAÇÃO
Francisco José Doubrawa Filho
CONTROLE DE VIBRAÇÃO FLEXIONAL EM SISTEMAS
GIRANTES UTILIZANDO NEUTRALIZADORES DINÂMICOS
VISCOELÁSTICOS
Esta Dissertação foi julgada para a obtenção do título de mestre em engenharia, Área de
concentração em Mecânica dos Sólidos e Vibrações, e aprovada em sua forma final pelo Programa de Pós-Graduação em Engenharia Mecânica e de Materiais.
Prof. Neri Volpato, Dr. Eng.
Coordenador de Curso
BANCA E XAMINADORA
Prof. Carlos Alberto Bavastri, Dr. Eng.
Universidade Tecnológica Federal do Paraná
Prof. Jucélio Tomás Pereira, Dr. Eng.
Universidade Tecnológica Federal do Paraná
Prof. Hans Ingo Weber, Dr. Eng.
Pontifícia Universidade Católica do Rio de
Janeiro
Prof. Marco Antônio Luersen, Dr. Eng.
Universidade Tecnológica Federal do Paraná
Curitiba, 30 de Junho de 2008
iii
In Memoriam:
Dedico este trabalho ao meu pai Francisco José Doubrawa,
Mecânico autodidata que esteve à frente do seu tempo e sempre me incentivou no gosto pela matemática e física.
iv
AGRADECIMENTOS
A minha família pela compreensão e apoio nos momentos de ausência.
Ao meu orientador Prof. Carlos Alberto Bavastri que aceitando um desafio, conseguiu com paciência, amizade e esforço fazer com que este trabalho se concretizasse.
Ao Prof. Jucélio Tomás Pereira outro gigante incansável mestre, que com seu infalível “kit multimídia” conseguiu transmitir relevantes informações nas suas aulas.
Ao Prof. Marco Antonio Luersen, pelos conhecimentos transmitidos em otimização
e pela presença nas defesas e apresentações.
Especialmente ao meu amigo Hideraldo Luis Vasconcelos do Santos pelas medições, suporte técnico e companheirismo, sem os quais as dificuldades seriam muito
maiores.
Ao pessoal do LAVIB, principalmente Aleksander Kokot e Eduardo Afonso Ribeiro
pela boa vontade nas montagens e medições.
Ao Programa de Pós-Graduação em Engenharia Mecânica e de Materiais, pela
infra-estrutura tornada disponível para o desenvolvimento deste e de outros trabalhos, e ao seu pessoal técnico e administrativo, que de uma forma ou de outra,
acabou por participar deste trabalho.
A WEG Equipamentos Elétricos S.A. - Máquinas, na pessoa do Sr. Fredemar Rüncos, por autorizar minha ausência do trabalho nos períodos de aula.
Ao MCT/FINEP/FNDCT - Chamada PROMOVE - Laboratórios de Inovação Convênio 4931/06 e à empresa WEG Equipamentos Elétricos S.A. pelo apoio financeiro no desenvolvimento deste trabalho.
A todos que contribuíram para esta dissertação.
v
DOUBRAWA FILHO, Francisco José, Controle de Vibração Flexional em Sistemas Girantes Utilizando Neutralizadores Dinâmicos Viscoelásticos, 2008, Dissertação (Mestrado em
Engenharia) - Programa de Pós-Graduação em Engenharia Mecânica e de Materiais, Universidade Tecnológica Federal do Paraná, Curitiba, 128p.
RESUMO
Sistemas rotativos estão sujeitos à forças geradas pelo desbalanceamento residual a qual é
proporcional ao quadrado da rotação. Não é incomum máquinas modernas de alta velocidade
operarem acima da rotação crítica. Próximo a esta rotação, em sistemas com baixo amortecimento, o fator de amplificação pode levar o rotor a operar em níveis de vibração excessiva ou
mesmo à destruição. Neutralizadores dinâmicos de vibração são dispositivos amplamente utilizados na atenuação de ruído e vibração em estruturas não girantes. Normalmente neutralizadores dinâmicos são construídos utilizando materiais viscoelásticos cujas propriedades dependem
da temperatura e da freqüência. Um dos modelos que melhor descreve este comportamento é
o que faz uso de derivada fracionária e quatro parâmetros. O sistema rotativo, chamado primário, pode ser modelado utilizando parâmetros modais obtidos no domínio da freqüência do
espaço de estado para uma dada temperatura de trabalho. Utilizando uma metodologia similar
à desenvolvida pelo grupo PISA o sistema composto (sistema girante + neutralizadores) pode
ser modelado em um subespaço modal do sistema primário. O modelo à parâmetros equivalentes generalizados utilizado nos neutralizadores permite que o sistema composto seja resolvido
utilizando as coordenadas generalizadas do sistema primário apenas, apesar dos mesmos introduzirem novos graus de liberdade. O projeto ótimo dos neutralizadores é implementado num
subespaço modal do espaço de estado do sistema primário utilizando um algoritmo de otimização não linear. A função objetivo é definida pela norma Euclidiana do vetor de máximos
absolutos das assim chamadas coordenadas principais no subespaço definido. Os neutralizadores devem ser fixados ao sistema primário utilizando um mancal flutuante localizado em um
ponto modal ativo, cujo deslocamento é elevado para o modo a ser controlado. Uma metodologia de projeto ótimo de neutralizadores para redução da resposta vibratória flexional em uma
faixa larga de freqüências é apresentada e resultados numéricos e experimentais são produzidos
e discutidos.
Palavras-chave: Máquinas Rotativas, Controle de Vibração, Neutralizadores Dinâmicos
vi
DOUBRAWA FILHO, Francisco José, Controle de Vibração Flexional em Sistemas Girantes Utilizando Neutralizadores Dinâmicos Viscoelásticos, 2008, Dissertação (Mestrado em
Engenharia) - Programa de Pós-Graduação em Engenharia Mecânica e de Materiais, Universidade Tecnológica Federal do Paraná, Curitiba, 128p.
ABSTRACT
All rotating systems are subjected to residual unbalance forces proportional to speed squared. In modern high speed machines is usual operate above the first critical speed. Near to this
speed in system with low damping, the amplification factor can lead the rotor to a high vibration
level condition or even to the destruction. Vibration neutralizers are well known devices and
successfully applied to reduce vibrations and noise on several dynamic non rotating structures. Usually the neutralizers are constructed using viscoelastic material which has frequencytemperature dependent properties. One of the models that best describes this behavior is the
four parameters fractional derivative. The rotating system here called primary system can be
modeled using modal parameters obtained in the frequency domain space state model for a given temperature. In a similar way to the general methodology developed by PISA group, the
compound system (rotating system + dynamic absorbers) can be modeled in a state space modal subspace of the primary system. The dynamic absorbers itself modeled using generalized
equivalent parameters allows the compound system to be solved using only the generalized coordinates of the primary system, even if they introduce new degrees of freedom. In this modal
subspace of the primary system the optimal design of the dynamic viscoelastic absorbers will be
performed using a non linear optimization algorithm. The cost function is defined by the Euclidean norm of the so called principal coordinates in the defined modal subspace. The absorbers
should be attached to a floating bearing located in a modal active point where the displacement
is highest for the modes to be controlled. A methodology to optimal design of dynamic vibration neutralizers to reduce a flexural unbalance response in a simple rotor, in a wide frequency
band, is presented and numerical-experimental results are produced and discussed.
Keywords: Rotating machinery, Vibration control, Dynamic neutralizers
vii
SUMÁRIO
RESUMO
v
ABSTRACT
vi
LISTA DE FIGURAS
xii
LISTA DE TABELAS
xvi
LISTA DOS PRINCIPAIS SÍMBOLOS
xvii
LISTA DE ABREVIATURAS E SIGLAS
xviii
1 INTRODUÇÃO
1
2 REVISÃO TEÓRICA
7
2.1
2.2
Equações de Lagrange . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
7
2.1.1
Princípio do Trabalho Virtual . . . . . . . . . . . . . . . . . . . . . .
8
2.1.2
Massa e Rigidez Generalizadas para Sistemas Lineares . . . . . . . . .
9
Elementos do Rotor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
2.2.1
Disco . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
13
2.2.2
Eixo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
15
2.2.3
Mancais . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
19
2.2.4
Excitação Tipo Desbalanceamento . . . . . . . . . . . . . . . . . . . .
20
2.2.5
Montagem das Matrizes Globais do Sistema . . . . . . . . . . . . . . .
22
Sumário
2.3
viii
O Material Viscoelástico . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
24
2.3.1
Relação Constitutiva Generalizada . . . . . . . . . . . . . . . . . . . .
24
2.3.2
Efeito da Temperatura . . . . . . . . . . . . . . . . . . . . . . . . . .
25
2.3.3
Efeito da Freqüência . . . . . . . . . . . . . . . . . . . . . . . . . . .
27
2.3.4
Outros Efeitos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
28
2.3.5
Efeito Combinado de Temperatura e Freqüência . . . . . . . . . . . . .
28
Sistemas com Um Grau de Liberdade . . . . . . . . . . . . . . . . . . . . . .
29
2.4.1
Modelo Viscoso . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
29
2.4.2
Modelo Viscoelástico . . . . . . . . . . . . . . . . . . . . . . . . . . .
32
Controle de Vibração . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
34
2.5.1
Rigidez Dinâmica na Base . . . . . . . . . . . . . . . . . . . . . . . .
35
2.5.2
Sistemas com Dois Graus de Liberdade . . . . . . . . . . . . . . . . .
35
2.6
Parâmetros Equivalentes Generalizados . . . . . . . . . . . . . . . . . . . . .
39
2.7
Múltiplos Graus de Liberdade, Sistemas Girantes . . . . . . . . . . . . . . . .
42
2.7.1
Ortonormalização . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
44
2.7.2
Diagrama de Campbell . . . . . . . . . . . . . . . . . . . . . . . . . .
46
2.7.3
Solução simplificada . . . . . . . . . . . . . . . . . . . . . . . . . . .
47
2.7.4
Solução no Domínio da Freqüência . . . . . . . . . . . . . . . . . . .
49
2.4
2.5
3 MODELO DO SISTEMA COMPOSTO
53
3.1
Implementação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
53
3.2
Resposta em Freqüência . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
55
3.3
Inclusão dos Neutralizadores Diretamente no Modelo do Eixo . . . . . . . . .
56
3.3.1
57
Diagrama de Campbell . . . . . . . . . . . . . . . . . . . . . . . . . .
4 PROJETO DO NEUTRALIZADOR PARA DESBALANCEAMENTO 62
4.1
Massa do Neutralizador . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
62
Sumário
ix
4.2
Função Objetivo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
63
4.3
Implementação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
64
4.4
Detalhamento do processo de projeto . . . . . . . . . . . . . . . . . . . . . . .
65
4.5
Sistema Primário - Identificação . . . . . . . . . . . . . . . . . . . . . . . . .
68
4.6
Sistema Auxiliar - Montagem . . . . . . . . . . . . . . . . . . . . . . . . . . .
71
4.7
Simulação Numérica Sem Neutralizadores . . . . . . . . . . . . . . . . . . . .
72
4.8
Exportação das Matrizes . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
73
4.9
Parâmetros do Modelo de Derivada Fracionária do Material Viscoelástico . . .
73
4.10 Simulação Numérica Com Neutralizadores . . . . . . . . . . . . . . . . . . . .
74
4.11 Construção Física do Neutralizador . . . . . . . . . . . . . . . . . . . . . . . .
76
4.11.1 Outros Materiais . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
78
4.11.2 Resposta em Freqüência . . . . . . . . . . . . . . . . . . . . . . . . .
80
4.11.3 Variação dos Parâmetros do Material Viscoelástico com a Temperatura
81
5 VALIDAÇÃO EXPERIMENTAL
83
5.1
Resposta em Freqüência do Sistema em Repouso . . . . . . . . . . . . . . . .
83
5.2
Resposta ao Desbalanceamento . . . . . . . . . . . . . . . . . . . . . . . . . .
86
5.3
Tabela de Resultados . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
89
5.4
Discussão dos Resultados . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
89
5.4.1
Teste de Resposta . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
89
5.4.2
Órbita . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
90
6 CONCLUSÕES
91
7 SUGESTÕES PARA TRABALHOS FUTUROS
94
Referências Bibliográficas
95
Sumário
x
ANEXO A -- SISTEMA NÃO GIRANTE
99
A.1 Introdução dos Neutralizadores no Modelo . . . . . . . . . . . . . . . . . . . . 101
A.1.1 Massa dos Neutralizadores . . . . . . . . . . . . . . . . . . . . . . . . 103
A.2 Otimização Multi-Modal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106
A.3 Funções Objetivo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107
ANEXO B -- ALGORITMOS DE OTIMIZAÇÃO
109
B.1 “Quasi-Newton” . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109
B.2 Algoritmo Genético . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110
ANEXO C -- MANCAIS DE ROLAMENTO
112
C.1 Rigidez . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112
ANEXO D -- FREQÜÊNCIA NATURAL DO NEUTRALIZADOR
113
D.1 Medição . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113
ANEXO E -- TRANSDUTORES DE DESLOCAMENTO
115
E.1 “Proximitor” . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115
ANEXO F -- MOVIMENTO ELÍPTICO
F.1
116
Definição . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116
ANEXO G -- PARÂMETROS MODAIS
117
G.1 Modos de Vibrar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117
G.1.1 Matriz [Φ] Complexa e Truncada 10 × 4 . . . . . . . . . . . . . . . . . 118
G.1.2 Matriz [Ψ] Complexa Truncada 10 × 4 . . . . . . . . . . . . . . . . . . 118
G.1.3 Matriz [Λ] Complexa Truncada 4 × 4 . . . . . . . . . . . . . . . . . . 118
ANEXO H -- PARÂMETROS DO MANCAL FLUTUANTE
119
Sumário
xi
H.1 Obtenção . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119
ANEXO I -- EQUIPAMENTOS DE MEDIÇÃO
120
I.1
Analisador de Sinais . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120
I.2
Martelo Instrumentado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121
I.3
Acelerômetro . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121
I.4
Sistema de Aquisição de Sinal . . . . . . . . . . . . . . . . . . . . . . . . . . 122
I.5
Proximitor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124
ANEXO J -- ROTORDIN
125
J.1
Introdução: . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 125
J.2
Capacidades . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126
xii
LISTA DE FIGURAS
Figura 1.1
Máquina de alta rotação, cortesia ASI Robicon. . . . . . . . . . . . .
Figura 2.1
Disco com velocidade Ω sobre AB. (fonte LALANNE e FERRARIS
1
(2001)) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
14
Figura 2.2
Modelo de viga discretizada. . . . . . . . . . . . . . . . . . . . . . .
16
Figura 2.3
Modelo para rigidez e amortecimento do mancal. . . . . . . . . . . .
19
Figura 2.4
Massa residual mu . . . . . . . . . . . . . . . . . . . . . . . . . . . .
21
Figura 2.5
Condições de equilíbrio para 2 graus de liberdade por nó. . . . . . .
22
Figura 2.6
Montagem da matriz global, para 6 elementos 8 × 8. . . . . . . . . .
23
Figura 2.7
Variação com temperatura (freqüência constante). . . . . . . . . . .
26
Figura 2.8
Nomograma material Isodamp C-1002. . . . . . . . . . . . . . . . .
27
Figura 2.9
Sistema massa mola 1gl, modelo viscoso. . . . . . . . . . . . . . . .
30
Figura 2.10
Função resposta em freqüência adimensional. . . . . . . . . . . . . .
31
Figura 2.11
Rigidez dinâmica adimensional. . . . . . . . . . . . . . . . . . . . .
31
Figura 2.12
Sistema massa-elemento resilente 1gl. . . . . . . . . . . . . . . . . .
32
Figura 2.13
Função resposta em freqüência adimensional. . . . . . . . . . . . . .
33
Figura 2.14
Área de carregamento: a) axial, b) cisalhamento. . . . . . . . . . . .
34
Figura 2.15
Modelo viscoelástico no domínio da freqüência. . . . . . . . . . . .
35
Figura 2.16
Gráfico da rigidez dinâmica Kb (Ω) na base. . . . . . . . . . . . . . .
36
Figura 2.17
Pontos fixos para vários amortecimentos. . . . . . . . . . . . . . . .
37
Figura 2.18
Desempenho dos modelos viscoelástico e viscoso. . . . . . . . . . .
38
Figura 2.19
Modelo com material viscoelástico no domínio da freqüência. . . . .
39
Figura 2.20
Equivalência de um sistema primário no domínio da freqüência. . . .
40
LISTA DE FIGURAS
xiii
Figura 2.21
Sistema dois graus de liberdade no domínio da freqüência. . . . . . .
41
Figura 2.22
Diagrama de Campbell. . . . . . . . . . . . . . . . . . . . . . . . .
47
Figura 2.23
Resposta ao desbalanceamento, RotorDin. . . . . . . . . . . . . . .
52
Figura 2.24
Resposta ao desbalanceamento, LALANNE e FERRARIS (2001). . .
52
Figura 3.1
Diagrama de Campbell interno. . . . . . . . . . . . . . . . . . . . .
58
Figura 3.2
Detalhe do diagrama de Campbell interno. . . . . . . . . . . . . . .
59
Figura 3.3
Detalhe do diagrama de Campbell interno a 12500 (rpm). . . . . . .
59
Figura 3.4
Campbell final com neutralizadores. . . . . . . . . . . . . . . . . . .
60
Figura 3.5
Resposta com neutralizadores de neoprene e Ωrpm = 0. . . . . . . . .
61
Figura 3.6
Resposta com neutralizadores de borracha butílica. . . . . . . . . . .
61
Figura 4.1
Diagrama de otimização.
. . . . . . . . . . . . . . . . . . . . . . .
65
Figura 4.2
Diagrama detalhado de otimização. . . . . . . . . . . . . . . . . . .
67
Figura 4.3
Bancada de teste do LAVIB. . . . . . . . . . . . . . . . . . . . . . .
68
Figura 4.4
Detalhe do dispositivo de segurança. . . . . . . . . . . . . . . . . .
70
Figura 4.5
Geometria do sistema primário do LAVIB. . . . . . . . . . . . . . .
71
Figura 4.6
Detalhe da fixação dos neutralizadores. . . . . . . . . . . . . . . . .
71
Figura 4.7
Diagrama de Campbell do sistema primário. . . . . . . . . . . . . .
72
Figura 4.8
Diagrama de Campbell do sistema primário (detalhe). . . . . . . . .
72
Figura 4.9
Resposta do sistema primário ao desbalanceamento. . . . . . . . . .
73
Figura 4.10
Nomograma do material viscoelástico. . . . . . . . . . . . . . . . .
74
Figura 4.11
Resposta ao desbalanceamento com e sem controle. . . . . . . . . .
75
Figura 4.12
Resposta linear com e sem controle. . . . . . . . . . . . . . . . . . .
75
Figura 4.13
Geometria do projeto físico do neutralizador de 150 (g). . . . . . . .
77
Figura 4.14
Modelo sólido do neutralizador. . . . . . . . . . . . . . . . . . . . .
77
Figura 4.15
Desempenho do neoprene. . . . . . . . . . . . . . . . . . . . . . . .
79
Figura 4.16
Desempenho da borracha butílica 20%. . . . . . . . . . . . . . . . .
79
LISTA DE FIGURAS
xiv
Figura 4.17
desempenho da borracha natural. . . . . . . . . . . . . . . . . . . .
80
Figura 4.18
desempenho do material EAR isodamp C-1002. . . . . . . . . . . .
80
Figura 4.19
Resposta em freqüência com neutralizadores e Ωrpm = 0. . . . . . .
81
Figura 4.20
Resposta na temperatura de 318 (K). . . . . . . . . . . . . . . . . .
82
Figura 5.1
Medição da resposta em repouso. . . . . . . . . . . . . . . . . . . .
83
Figura 5.2
Esquema de medição FRF em repouso. . . . . . . . . . . . . . . . .
84
Figura 5.3
Resultado em repouso (inertância). . . . . . . . . . . . . . . . . . .
84
Figura 5.4
Resultado completo em repouso. . . . . . . . . . . . . . . . . . . .
85
Figura 5.5
Comparação de resultados em repouso. . . . . . . . . . . . . . . . .
85
Figura 5.6
Detalhe da fixação dos transdutores, em destaque eixo z. . . . . . . .
86
Figura 5.7
Órbita de referência a 1800 (rpm). . . . . . . . . . . . . . . . . . . .
87
Figura 5.8
Órbita sem neutralizadores na rotação crítica a 3600 (rpm). . . . . .
88
Figura 5.9
Órbita com neutralizadores na rotação crítica a 3600 (rpm). . . . . .
88
Figura A.1
Viga e modelo equivalente generalizado. . . . . . . . . . . . . . . . 101
Figura A.2
Acoplamento das coordenadas principais {P(Ω)}. . . . . . . . . . . 103
Figura A.3
Otimização multi-modal. . . . . . . . . . . . . . . . . . . . . . . . . 107
Figura B.1
Comparação de técnicas de otimização. . . . . . . . . . . . . . . . . 111
Figura C.1
Deflexão devido compressão dos elementos. . . . . . . . . . . . . . 112
Figura D.1
Montagem para medição da massa aparente. . . . . . . . . . . . . . 113
Figura D.2
Gráfico massa aparente obtido para um neutralizador de 150 (g). . . . 114
Figura E.1
Transdutores em configuração ortogonal XY. . . . . . . . . . . . . . 115
Figura F.1
Movimento elíptico (órbita). . . . . . . . . . . . . . . . . . . . . . . 116
Figura G.1
Forma de vibrar do primeiro modo forward. . . . . . . . . . . . . . . 117
Figura H.1
Função resposta medida experimentalmente. . . . . . . . . . . . . . 119
Figura I.1
Analisador HP 3560A. . . . . . . . . . . . . . . . . . . . . . . . . . 120
Figura I.2
Martelo PCB 086C04. . . . . . . . . . . . . . . . . . . . . . . . . . 121
LISTA DE FIGURAS
xv
Figura I.3
Acelerômetro PCB 352C68. . . . . . . . . . . . . . . . . . . . . . . 122
Figura I.4
Acelerômetro PCB 352C68. . . . . . . . . . . . . . . . . . . . . . . 122
Figura I.5
Placa de entrada SCXI 1531. . . . . . . . . . . . . . . . . . . . . . . 123
Figura I.6
Placa de aquisição DAQCArd 6062E. . . . . . . . . . . . . . . . . . 124
Figura I.7
Sensor PS 1002. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124
Figura I.8
Condicionador GS 501. . . . . . . . . . . . . . . . . . . . . . . . . 124
Figura J.1
Tela principal. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126
Figura J.2
Geometria e modelo. . . . . . . . . . . . . . . . . . . . . . . . . . . 127
Figura J.3
Diagrama de Campbell. . . . . . . . . . . . . . . . . . . . . . . . . 127
Figura J.4
Resposta ao desbalanceamento. . . . . . . . . . . . . . . . . . . . . 128
Figura J.5
Modo 1, direção “forward”. . . . . . . . . . . . . . . . . . . . . . . 128
xvi
LISTA DE TABELAS
Tabela 2.1
Tabela de modelos equivalentes. . . . . . . . . . . . . . . . . . . . .
41
Tabela 4.1
Freqüências ótimas para diversos materiais. . . . . . . . . . . . . . .
78
Tabela 5.1
Resultados obtidos. . . . . . . . . . . . . . . . . . . . . . . . . . . .
89
xvii
LISTA DOS PRINCIPAIS SÍMBOLOS
δW - Trabalho virtual.
γr - Deslocamento virtual.
F - Função dissipativa de Rayleigh.
ω - Velocidade angular.
Ω - Rotação em torno do eixo de coordenadas y, freqüência angular.
ψ, θ e φ - Ângulos de giro sobre Z, X eY , respectivamente.
ρ - Massa específica.
σ - Tensão.
ε - Deformação.
ν - Coeficiente de Poisson.
α(T ) - Fator de deslocamento, dependente da temperatura.
β - Expoente fracionário característico do material viscoelástico.
η - Fator de perda material viscoelástico.
α(Ω) - Resposta em freqüência, receptância.
ε - Razão de freqüências.
βe - Fator de forma.
, - Por definição.
ζ - Razão de amortecimento.
ℑ - Função transformada de Fourier.
[λ ] - Matriz de autovalores.
δ + iυ - Autovalor complexo.
{φ }, {ϕ} - Vetores solução no domínio do tempo.
[θ ], [ψ] - Matrizes de autovetores.
[Θ], [Ψ] - Matrizes de autovetores ortonormalizados.
R - Domínio dos números reais.
µ - Razão de massas.
δm - Direção de descida.
δc - Folga radial.
δh - Compressão.
δo - Ovalização.
xviii
LISTA DE ABREVIATURAS E SIGLAS
BFGS - Algoritmo de otimização de Broyden - Fletcher - Goldfarb - Shanno.
CEFET-PR - Centro Federal de Educação Tecnológica do Paraná.
CNPq - Conselho Nacional de Desenvolvimento científico e tecnológico.
EPDM - Etileno-Propileno-Dieno (borracha).
FRF - Função Resposta em Freqüência.
UFPR - Universidade Federal do Paraná.
LAVIB - Laboratório de Vibrações, UTFPR.
LVA - Laboratório de Vibrações e Acústica - UFSC.
PISA - Grupo de Pesquisa Integrada em Sistemas Vibrantes e Acústicos.
RMS - “Root Mean Square”, valor médio quadrático.
TTL/CMOS - “Transistor-Transistor-Logic”/“Complementary Metal Oxide Semicondutor”
UFSC- Universidade Federal de Santa Catarina.
UTFPR - Universidade Tecnológica Federal do Paraná.
WLF - Fator de deslocamento de Williams-Landel-Ferry.
CPU - “Central Processing Unit”, unidade central de processamento.
1
1
INTRODUÇÃO
Máquinas de alta rotação tem um vasto campo de aplicação em compressores, bombas,
turbo-geradores, imposto por um mercado que exige potências maiores e tamanhos cada vez
menores. Como exemplo, a Figura 1.1, apresenta um motor de indução moderno de alta rotação
20.000 (rpm), potência 5 (MW). Vibrações em sistemas rotativos sempre estão presentes em
maior ou menor grau, devido à virtual impossibilidade de se atingir um balanceamento perfeito
de suas partes rotativas. A energia gerada por uma massa desbalanceada residual em rotação
pode ser suficiente para excitar modos de flexão com pouco amortecimento e provocar respostas
significativas em seus elementos de restrição à translação (mancais).
Figura 1.1: Máquina de alta rotação, cortesia ASI Robicon.
Algumas abordagens para solucionar problemas relacionados com altos níveis de vibração
podem ser adotadas:
1. Agir sobre a excitação:
Diminuindo a excitação diminui-se a resposta. Em certas aplicações uma melhoria no
balanceamento exige a utilização de equipamentos de alta sensibilidade e mesmo assim,
imprecisões mecânicas e não linearidades do sistema estabelecem um limite prático para
a qualidade do balanceamento.
2. Modificação estrutural:
A intervenção estrutural pode ser utilizada de forma a alterar características do sistema
Capítulo 1 Introdução
2
como rigidez e amortecimento, e conseqüentemente, suas rotações características. Porém
em aplicações atuais, nem sempre é possível alterar suficientemente o sistema por força
de dimensões padronizadas e outros tipos de restrições mecânicas.
3. Operar o sistema longe de suas rotações características:
Sistemas modernos exigem uma capacidade de operação em múltiplas rotações e muitas
vezes acima de uma ou mais rotações críticas. Somente em aplicações específicas este
tipo de solução é viável, e mesmo assim, para atingir a rotação de regime é necessário
passar pela rotação crítica, o que sempre é uma tarefa difícil.
4. Controle de vibração:
(a) O isolamento de vibração permite reduzir a propagação e amplificação de vibração entre a máquina e o suporte, porém entre eixo e mancal ainda podem ocorrer
deslocamentos suficientes elevados para danificar selos e o próprio mancal.
(b) O uso de neutralizadores dinâmicos, uma outra forma de modificação estrutural,
permite reduzir a resposta em uma região da freqüência na qual o sistema possui
uma ou várias freqüências naturais. Estes dispositivos auxiliares podem ser classificados em dois grupos, ativos e passivos. Os dispositivos ativos requerem sensores
e processadores, apresentando por isso um custo adicional ou maior se comparado
aos passivos. Os dispositivos passivos, se devidamente projetados, deverão atenuar
a resposta vibratória, permitindo a operação precisa do sistema, com a desvantagem de não se dispor da informação do resultado do controle em operação. Um
dos tipos de dispositivos que podem ser empregados, os neutralizadores dinâmicos,
são o objeto de estudo desta proposta. Para estruturas não girantes ESPÍNDOLA e
SILVA (1992) desenvolveu uma teoria geral para projetar de forma ótima este tipo
de dispositivo, utilizando materiais viscoelásticos como elemento resilente.
O projeto de neutralizadores dinâmicos para sistemas não girantes vem evoluindo na sua implementação e projeto. Segundo KRENEV e REZNIKOV (1993), para reduzir o movimento
oscilatório em navios, o almirante Makarov em 1897 desenvolveu o conceito de neutralizador
de vibrações utilizando a transferência de água entre tanques, mas foi FRAHM (1909), que definitivamente teve seu nome associado ao desenvolvimento deste tipo de neutralizador. Outro
pioneiro, HARTOG (1956), estudou modelos simples de um grau de liberdade massa-molaamortecedor viscoso e propôs o método dos pontos fixos, visando a obtenção de parâmetros
ótimos através da sintonização do neutralizador.
Capítulo 1 Introdução
3
Inicialmente os neutralizadores eram projetados sem amortecimento, tinham uma faixa de
operação muito estreita e não raro apresentavam problemas de fadiga. Por outro lado, o amortecimento viscoso, modelo utilizado como elemento dissipador de energia vibratória, era de
difícil implementação prática.
Modelos de neutralizadores utilizando material viscoelástico foram apresentados por SNOWDON (1959), aplicados a sistemas de um e dois graus de liberdade, onde o material viscoelástico
substitui a mola e o amortecedor viscoso. A grande vantagem deste tipo de material é o elevado
amortecimento que pode ser obtido e a versatilidade para projetar neutralizadores dinâmicos de
tamanhos variados. Existem varias abordagens para o modelamento do comportamento deste
tipo de material, dentre as quais não se pode deixar de citar os modelos simples propostos
por Kelvin e Maxwell LAKES (1999). Estes modelos não conseguem descrever o comportamento do material numa escala larga de freqüência. Modelos mais precisos, como estudado por
ROGERS (1983), utilizando derivadas fracionárias descrevem com a aproximação necessária
à aplicações em banda larga, as variações do módulo de elasticidade e o fator de perda com
a freqüência e temperatura. PRITZ (1996) estudou o modelo utilizando derivada fracionária
com quatro parâmetros. Materiais viscoelásticos são amplamente utilizados na engenharia em
controle de vibrações e ruído irradiado. Para tal fim, o conhecimento preciso das características
dinâmicas deste material é fundamental. LOPES et al. (2004) desenvolve e utiliza uma metodologia própria para a determinação das características do módulo de cisalhamento G e fator de
perda η, em função da freqüência e da temperatura, apresentando-as em forma de nomograma
conforme LOPES et al. (2004) e ESPÍNDOLA et al. (2006). Aplicações bem sucedidas deste
modelo, utilizado em neutralizadores dinâmicos para controle de vibração, podem ser encontrados nos trabalhos de ESPíNDOLA et al. (2005), CRUZ (2004), BAVASTRI et al. (2006) e
FERREIRA (2005) entre outros.
Um dos problemas decorrentes da utilização do material viscoelástico em neutralizadores
dinâmicos é a grande variação de seus parâmetros de rigidez e amortecimento com a freqüência
de operação e a temperatura, o que praticamente inviabiliza métodos analíticos de otimização
de parâmetros, tendo em vista a complexidade inerente à própria formulação. Visando resolver
este problema, as técnicas de otimização não linear tem sido aplicadas com sucesso na obtenção de parâmetros ótimos para neutralizadores, como realizado por KITIS (1983) ao minimizar
a resposta vibratória de uma viga engastada simples utilizando neutralizadores dinâmicos. O
método que foi empregado naquele trabalho resolve o sistema completo (no espaço de configurações), cujo problema pode se tornar extremamente pesado, do ponto de vista computacional,
quando o sistema apresenta um número de graus de liberdade elevado, situação comum em
sistemas rotativos modelados com elementos de viga.
Capítulo 1 Introdução
4
Para este tipo de problema, ESPÍNDOLA (1992) e SILVA (1991) apresentaram uma generalização para o projeto ótimo de neutralizadores dinâmicos, propondo o conceito de parâmetros equivalentes generalizados no qual a dinâmica do sistema composto (sistema primário
com neutralizadores dinâmicos) pode ser descrito apenas em função das coordenadas generalizadas do sistema primário, sendo possível então realizar as análises em um sub-espaço modal
do sistema primário, com um número reduzido de equações. O sistema composto foi suposto
como completamente desacoplado, e o controle foi realizado modo a modo, onde um ou mais
neutralizadores eram projetados para controlar um modo especificado.
O truncamento das matrizes modais, como proposto por ESPÍNDOLA e SILVA (1992) e
apresentado em BAVASTRI (1997), permite sem perda de precisão na faixa de freqüência proposta, uma redução no tamanho do sistema e conseqüente ganho de performance na solução.
Quando se propõe o controle de vários modos simultaneamente, utilizando um ou mais neutralizadores, é essencial que se estabeleça a função objetivo a qual é submetida a um algoritmo
de otimização determinístico, evolucionário ou combinação de ambos. Em função do tipo de
problema, que pode apresentar mínimos locais, a combinação de métodos se mostra bastante
eficiente. Enquanto o algoritmo evolucionário busca uma aproximação para o mínimo global,
o determinístico se encarrega de garantir a localização exata do ponto de mínimo da função
objetivo.
É possível aplicar os conceitos de controle passivo ótimo de vibração flexional em sistemas
girantes. Para tal fim, é necessário adaptar a teoria de estruturas não girantes, já que em geral
é necessário modelar o sistema composto em um sub-espaço modal do espaço de estado do
sistema primário.
O estudo de sistemas girantes, partiu de modelos muito simples, como o trabalho publicado por RANKINE (1869) no qual eram previstas precessão e deflexão sem limites. Conforme
GENTA (2005), De Laval1 , posteriormente em 1889, construiu com sucesso máquinas centrífugas e turbinas que operavam acima da primeira rotação crítica. DUNKERLEY (1894) publicou
um estudo sobre vibrações em eixos e utilizou pela primeira vez o termo “rotação crítica”.
Outro pioneiro JEFFCOTT (1919), modelou um rotor de forma simples mas consistente até a
rotação crítica. Grandes avanços foram obtidos por STODOLA (1927) na teoria da operação
super-crítica (acima da rotação crítica). SMITH (1933) publicou um artigo sobre estabilidade
introduzida pelo amortecimento. Posteriormente, MYKLESTAD (1944), PROHL (1945) e outros desenvolveram o método de cálculo da rotação crítica através de matrizes de transferência,
1 Carl
Gustaf Patrik de Laval (Orsa, Dalarna, 9 de Maio de 1845 - Estocolmo, 2 de Fevereiro de 1913) foi
um engenheiro e inventor sueco, que fez grandes contribuições no projeto da turbina a vapor e na maquinaria de
ordenha do leite.
Capítulo 1 Introdução
5
que ainda encontra aplicação atualmente.
O modelo numérico para sistemas rotativos está muito bem desenvolvido, podendo-se citar
vários autores como LALANNE e FERRARIS (2001), GENTA (2005) e VANCE (1988), entre
outros. Estes autores descrevem um sistema girante através de equações diferenciais gerais do
movimento, utilizando matrizes consistentes de massa, rigidez e amortecimento considerando o
efeito giroscópico. Estas matrizes, por sua vez, podem ser obtidas através de modelos discretos
finitos, em geral para elementos finitos de viga, disco e mancais, como proposto por LALANNE
e FERRARIS (2001) e que possibilita simulações com resultados precisos.
A atenuação de vibrações em sistemas girantes pode ser ativa, tipo balanceamento, como
publicado por ALAUZE et al. (2001) e ZHOU e SHI (2001), onde um sistema composto por
sensores e uma “CPU”, adquire sinais de velocidade, fase e amplitude de vibração e calcula
em tempo real, o deslocamento de massas distribuídas em diferentes raios através de servomecanismos. Mancais de levitação magnética, também tem se mostrado outro método ativo
de redução de vibrações, como apresentado por HOPE et al. (1998). Porém, devido ao alto
custo da implementação, sua utilização ainda é bastante restrita. Outro tipo de controle ativo
adaptativo foi apresentado por TAMMI (2007), utilizando atuador e software. Neste caso o
atuador é composto pelas bobinas eletromagnéticas que compõem o mancal.
Outro tipo de abordagem, o controle passivo de vibrações, se faz, por exemplo, introduzindo
amortecimento nos mancais do sistema, através de mantas de material viscoelástico, como proposto em FERREIRA (2005) e em CHÁVEZ (2003), ou ainda, introduzindo um fino filme de
óleo na capa externa dos rolamentos que compõem os mancais do sistema (“squeeze-film”),
como mostrado por ZEIDAN (1995). Nesses casos a resposta vibratória também é atenuada,
mas a aplicação requer a modificação do projeto dos mancais.
Com a evolução das técnicas de projeto do eixo de máquinas rotativas, é possível predizer suas rotações críticas com precisão. Em diversas situações, pode ser necessário operar a
máquina em rotações acima da primeira e até mesmo da segunda ou terceira rotações críticas,
como no caso de turbo-geradores, ou através de inversores de freqüência. Isto significa que em
um dado instante, o eixo poderá estar sujeito a uma amplificação da vibração, provocada pela
coincidência da rotação com sua freqüência característica de flexão. Esta rotação é comumente
denominada de rotação crítica. Dependendo da situação, se a passagem pela rotação crítica não
for rápida o suficiente ou se o amortecimento no mancal não for elevado , o nível de vibração
resultante pode ocasionar danos permanentes em mancais e selos mecânicos, com a necessidade
de parada e reparo da máquina. Em geral, máquinas com mancal de rolamento, devido ao baixo
amortecimento inerente, são sempre projetadas para operar abaixo da primeira rotação crítica.
Capítulo 1 Introdução
6
Neste trabalho, desenvolve-se um projeto ótimo de neutralizadores dinâmicos, utilizando
material viscoelástico como elemento resilente. Os mesmos são acoplados ao eixo do rotor da
máquina através de um falso mancal (mancal flutuante), possibilitando a operação ou a transição
suave ao passar por uma ou mais rotações críticas. Desta forma, pretende-se minimizar a amplificação da resposta vibratória à flexão ampliando a faixa de trabalho da máquina e garantindo
a durabilidade dos elementos que compõem o rotor.
Para que se possa projetar tal neutralizador é necessário obter o modelo modal do sistema
girante para o qual é utilizado um código numérico próprio. A partir destes parâmetros e utilizando uma teoria equivalente à apresentada por BAVASTRI (1997), propõe-se uma metodologia geral (aplicável em sistemas girantes modeláveis através de tais parâmetros) para o projeto
ótimo de neutralizador dinâmico de vibração viscoelástico para sistemas rotativos, controlando
uma ou várias rotações críticas. Obtidos os parâmetros modais do sistema primário, o sistema
composto (sistema primário + neutralizadores dinâmicos) é modelado em um sub-espaço modal no espaço de estado do sistema primário. Utilizando técnicas de otimização não linear, são
variadas as freqüências naturais dos neutralizadores, buscando a minimização da resposta do
sistema composto, dentro da faixa de freqüências proposta.
Para validação do procedimento, são construídos e ensaiados, neutralizadores dinâmicos
viscoelásticos a serem montados em um mancal flutuante, sobre um rotor simples. O material
viscoelástico a ser utilizando deve ter suas características dinâmicas perfeitamente conhecidas.
O modelo matemático proposto, para predizer o comportamento dinâmico do sistema composto (primário + neutralizadores), foi implementado em uma linguagem de alto nível. Os
resultados numéricos obtidos são apresentados, comparados com as medições experimentais
realizadas e discutidos. São também apresentadas as conclusões e sugestões para trabalhos
futuros.
No corrente projeto, o sistema girante se limitará a eixos e discos cilíndricos e simétricos,
modelados através de expressões para energia cinética e potencial, aplicadas à equação de Lagrange. As forças consideradas são devido ao desbalanceamento e os mancais do tipo rolamento
de alta rigidez, são considerados pouco amortecidos por hipótese. Os elementos do sistema são
discretizado através do método de elementos finitos. Através de matrizes elementares obtidas
da formulação de Lagrange, são montadas as matrizes globais de massa, rigidez, amortecimento
e giroscópica. Estas matrizes globais são aplicadas à equação geral do movimento que, quando
solucionada através de um conjunto de problemas de autovalor, permite obter os parâmetros
modais do sistema girante primário. A obtenção dos parâmetros, como descrita, está implementada em um código numérico próprio denominado “RotorDin” (ver J.1).
7
2
REVISÃO TEÓRICA
Neste capítulo são revisados os modelos matemáticos necessários para o claro entendimento
do propósito e dos meios utilizados no desenvolvimento deste trabalho, com ênfase em sistemas
girantes e neutralizadores dinâmicos.
Os principais tópicos revisados são:
• Obtenção dos elementos da equação do movimento para sistemas rotativos modelados
através do método dos elementos finitos utilizando as equações de Lagrange.
• Um modelo para o comportamento do material viscoelástico em função da freqüência e
da temperatura utilizando derivada fracionária.
• Os modelos de um e dois graus de liberdade para um sistema simples massa-molaamortecedor e massa-material resilente, bem como os princípios básicos de controle de
vibração.
• A teoria geral de parâmetros equivalentes generalizados.
• O modelo de múltiplos graus de liberdade do sistema girante, o problema de autovalores
e sua solução. O diagrama de Campbell e a resposta ao desbalanceamento.
2.1
Equações de Lagrange
Através de um equacionamento baseado em grandezas escalares (trabalho e energia), Lagrange desenvolveu uma teoria para descrever o movimento de sistemas dinâmicos. Para descrever completamente o movimento de um sistema, são necessárias tantas coordenadas independentes quanto os graus de liberdade considerados.
A um conjunto de coordenadas físicas pode ser imposta uma ou mais transformações lineares de maneira a representar as mesmas para um sistema referencial arbitrário. Tais coordenadas
são denominadas “generalizadas”.
Capítulo 2 Revisão Teórica
2.1.1
8
Princípio do Trabalho Virtual
Para um sistema de N partículas qualquer em equilíbrio, podem ser impostos deslocamentos virtuais infinitesimais δ r j , fisicamente realizáveis, de tal modo que o trabalho virtual δW
realizado pelas forças Fj responsáveis pelos deslocamentos seja nulo. Ou seja
N
δW = ∑ Fj δ r j = 0
Eq. 2.1
j
Trabalho Virtual em Coordenadas Generalizadas
Para um sistema qualquer com n graus de liberdade, é possível expressar um deslocamento
r em função de n coordenadas generalizadas q, como mostrado por MEIROVITCH (1990),
como
r j = r j (q1 , q2 , ..., qn ) j = 1, N.
Eq. 2.2
Um deslocamento virtual δ r j pode ser obtido da expressão
n
∂rj
δ qi .
i=1 ∂ qi
δrj = ∑
Eq. 2.3
O trabalho virtual δW realizado, associado à força Fj , através da força generalizada
n
Qi =
∂rj
∑ F j ∂ qi ,
Eq. 2.4
j=1
pode ser expresso através da relação
n
δW = ∑ Qi δ qi .
Eq. 2.5
i=1
Este conceito foi estendido por D’Alembert para sistemas dinâmicos, a partir da segunda lei de
Newton, aplicada à i-ésima partícula de massa mi , com fi representando forças internas ou de
restrição
Fi − fi = mi r̈.
Eq. 2.6
Aplicando o princípio do trabalho virtual em um sistema de n coordenadas generalizadas, a
energia cinética T de um sistema de partículas é dada pela expressão
N
1
T = ∑ mi ṙi2 ,
i 2
Eq. 2.7
Capítulo 2 Revisão Teórica
9
obtém-se
n
δWT =
∑
k=1
d ∂T
∂T
−
− Qk δ qk .
dt ∂ q̇k ∂ qk
Eq. 2.8
Em um sistema conservativo, o trabalho realizado só depende dos estados inicial e final de sua
energia potencial U. Não há dissipação no caminho. Em termos de coordenadas generalizadas
(arbitrando-se o estado inicial como energia potencial nula), tem-se
W = −U(q1 , q2 , ..., qn ).
Eq. 2.9
Portanto, o trabalho virtual em termos de energia potencial, é dado pela expressão
δWU = − ∑
k
∂U
δ qk .
∂ qk
Eq. 2.10
0
Forças não conservativas fi também efetuam trabalho. Definindo-as em termos de coordenadas
generalizadas tem-se
n
Qk = ∑ f i
0
i=1
∂ ri
,
∂ qk
Eq. 2.11
logo
δWQ = ∑ Q̄k δ qk .
Eq. 2.12
∂T
∂U
d ∂T
−
+
= Q̄k ,
dt ∂ q̇k ∂ qk ∂ qk
Eq. 2.13
k
A equação de Lagrange completa fica
ou na sua forma reduzida, com o Lagrangeano definido como L = T −U e
∂U
∂ q̇
=0
d ∂L
∂L
−
= Q̄k .
dt ∂ q̇k ∂ qk
2.1.2
Eq. 2.14
Massa e Rigidez Generalizadas para Sistemas Lineares
Quando o deslocamento das coordenadas generalizadas fica restrito à vizinhança do seu
ponto de equilíbrio as equações do movimento da estrutura podem ser consideradas lineares,
como mostrado em ESPÍNDOLA (1992).
Energia Potencial
A equação da energia potencial U (Equação 2.10) para um sistema de n graus de liberdade
em função de suas coordenadas generalizadas, pode ser representada em torno de um ponto de
Capítulo 2 Revisão Teórica
10
equilíbrio através da série de Taylor
n ∂U
1
U = U0 + ∑
qj +
2
j=1 ∂ q j 0
n
n
∑∑
j=1 l=1
∂ 2U
∂ q j ∂ ql
q j ql + · · ·
Eq. 2.15
0
As derivadas são calculadas na posição de equilíbrio, por conveniência denominado ponto “0”
(zero). Considerando que a posição de equilíbrio corresponde à mínima energia potencial do
sistema, sua derivada primeira pode ser considerada nula. O termo constante U0 também pode
ser considerado nulo no ponto de referência. Desprezando os termos superiores da série de
Taylor, a expressão da energia potencial é
1
U=
2
n
n
∑∑
j=1 l=1
∂ 2U
∂ q j ∂ ql
q j ql .
Eq. 2.16
0
Através da comparação da expressão acima com a da energia associada a um sistema massamola de um grau de liberdade U = 12 kx2 , a rigidez generalizada k jl fica definida
k jl ,
∂ 2U
∂ q j ∂ ql
.
Eq. 2.17
0
Assim, pode-se escrever genericamente a energia potencial elástica do sistema através da equação
1
U = {q}T [K]{q}.
2
Eq. 2.18
Energia Cinética
Em um sistema de n graus de liberdade, em função de suas coordenadas generalizadas q, a
derivada temporal da posição r, da i-ésima partícula é
n
ṙi =
∂ ri
∂ ri
∑ ∂ q j q̇ j + ∂t .
Eq. 2.19
j=1
A expressão geral para a energia cinética é
T=
1 N
∑ miṙi2.
2 i=1
Eq. 2.20
Para um sistema de N partículas em termos de coordenadas generalizadas é dada pela relação
"
#
n n
1 N
∂ ri ∂ ri
∂ ri n ∂ ri
∂ ri ∂ ri
T = ∑ mi ∑ ∑
q̇ j q̇l + 2
Eq. 2.21
∑ ∂ q j q̇ j + ∂t ∂t ,
2 i=1
∂t j=1
j=1 l=1 ∂ q j ∂ ql
ou, separando em três parcelas
T = T2 + T1 + T0 .
Eq. 2.22
Capítulo 2 Revisão Teórica
11
Reescrevendo a expressão geral (Equação 2.22) para um sistema cujas restrições são independentes do tempo. Neste caso então T1 = 0 e T0 = 0 e tem-se
n
1
T = T2 =
2
n
∑ ∑ q̇ j q̇l
j=1 l=1
N
∂ ri ∂ ri
∑ mi ∂ q j ∂ ql
i=1
!
.
Eq. 2.23
Assim a massa generalizada m jl pode então ser definida como
N
m jl = ∑ mi
i=1
∂ ri ∂ ri
∂ 2 T2
=
,
∂ q j ∂ ql
∂ q̇ j ∂ q̇l
Eq. 2.24
e a energia cinética em forma matricial é
1
T = {q̇}T [M]{q̇}.
2
Eq. 2.25
O segundo termo (T1 ) da Equação 2.21, introduz o efeito giroscópico. Definindo a força
N
f r = ∑ mi
i=1
∂ ri ∂ ri
,
∂t ∂ qr
Eq. 2.26
e desenvolvendo sua derivada temporal em série de Taylor no ponto de equilíbrio, obtém-se
n
∂ fr
1 n n ∂ fr
ql + ∑ ∑
qi ql + ...,
2 i=1 l=1 ∂ qi ∂ ql
l=1 ∂ ql
fr = fr (0) + ∑
Eq. 2.27
aproximando para o primeiro termo
n
fr w
∂ fr
∑ ∂ ql ql
Eq. 2.28
∂ fr
,
∂ ql
Eq. 2.29
∑ frl ql .
Eq. 2.30
l=1
definindo-se
frl =
obtém-se
n
fr =
l=1
Desta forma, a parcela componente T1 da energia cinética é
n
T1 =
∑ fr q̇ j ,
Eq. 2.31
j=1
ou em forma matricial
T1 = {q̇}T [F]{q}.
Eq. 2.32
Capítulo 2 Revisão Teórica
12
Da expressão matricial, pode-se deduzir que
f jl =
∂ 2 T1
∂ 2 T1
e fl j =
.
∂ ql ∂ q̇ j
∂ q j ∂ q̇l
O terceiro termo, T0 , componente da Equação 2.21, segundo ESPÍNDOLA (1992), refere-se a
uma alteração na rigidez devido ao efeito de campos longitudinais (força centrífuga) e não depende da velocidade generalizada, mas somente da coordenada generalizada, pode ser reescrito
como
T0 =
1 N
∂ ri ∂ ri
mi
.
∑
2 i=1 ∂t ∂t
Eq. 2.33
Forças Dissipativas
A função dissipativa de Rayleigh modela o trabalho realizado por forças devido ao amortecimento viscoso, coeficientes ci j e circulatórias (1 ), hi j , proporcionais à velocidade generalizada, como apresentado por MEIROVITCH (1990)
F=
n n
1 n n
c
q̇
q̇
+
i
j
i
j
∑∑
∑ ∑ hi j q̇iq j .
2 i=1
j=1
i=1 j=1
Eq. 2.34
Equação de Lagrange na Forma Final
Para um sistema qualquer no qual atua uma força generalizada F̄ na k-ésima coordenada
generalizada, incorporando-se todos os termos desenvolvidos e utilizando a definição do Lagrangeano L = T −U, obtém-se
∂L ∂F
d ∂L
−
+
= F̄k .
dt ∂ q̇k ∂ qk ∂ q̇k
Eq. 2.35
Substituindo os termos obtidos para a energia cinética T , a energia potencial U e de dissipação
F na 2.35, obtém-se o conjunto de equações que descrevem o movimento de um sistema linear
com n graus de liberdade
n
m
q̈
+
(c
+
g
)
q̇
+
(k
+
h
)q
= F̄i i = 1, 2, . . . , n,
i
j
j
i
j
i
j
j
i
j
i
j
j
∑
Eq. 2.36
j=1
onde gi j , fi j − f ji , ou na forma matricial
[M]{q̈} + ([C] + [G]){q̇} + ([K] + [H]){q} = {F̄},
1A
Eq. 2.37
matriz de circulação H está relacionada ao efeito da circulação de fluxo, como o que ocorre em um perfil
de asa aeronáutico, por exemplo.
Capítulo 2 Revisão Teórica
13
onde
• [M] - Matriz de inércia (simétrica).
• [C] - Matriz de amortecimento viscoso (simétrica).
• [G] - Matriz giroscópica (anti-simétrica).
• [K] - Matriz de rigidez (simétrica).
• [H] - Matriz de circulação (anti-simétrica).
• {F̄} - Vetor de excitação (forças generalizadas).
• {q} - Vetor de coordenadas generalizadas.
• {q̇} =
d
dt {q}
• {q̈} =
d2
{q}dt 2
2.2
2.2.1
- Vetor de velocidades generalizadas.
Vetor de acelerações generalizadas.
Elementos do Rotor
Disco
O disco é caracterizado como um segmento de cilindro rígido, de massa Md , cuja energia
cinética pode ser expressa através da relação
Td =
1
Md ṙ2 + Iω 2 ,
2
Eq. 2.38
onde ṙ representa a velocidade de translação do sistema de referência, móvel em relação ao sistema inercial, ω a velocidade angular e I a inércia, referentes ao centro de massa C. Considerase a rotação Ω livre em torno do eixo de coordenadas y e nulo o deslocamento nesta direção.
Capítulo 2 Revisão Teórica
14
Figura 2.1: Disco com velocidade Ω sobre AB. (fonte LALANNE e FERRARIS (2001))
Para o cálculo do vetor de velocidade angular ω é necessário transportar as coordenadas
do centro de massa C para um sistema referencial global através de sucessivos giros ψ, θ e φ ,
sobre os eixos Z, X eY , respectivamente (Figura 2.1). Através de uma matriz de transformação,
segundo LALANNE e FERRARIS (2001), é possível obter o vetor velocidade angular ω


−ψ̇ cos θ sin φ + θ̇ cos φ


.
Eq. 2.39
ω =
φ̇
+
ψ̇
sin
θ


ψ̇ cos θ sin φ + θ̇ sin φ
Logo, considerando que o disco é axi-simétrico (Idz = Idx ) e que a velocidade angular é constante, obtém-se
Td =
1
Md (u̇2 + ẇ2 ) + Idx (ψ̇ 2 cos2 θ + θ̇ 2 ) + Idy (φ̇ 2 + ψ̇ sin2 θ + 2φ̇ ψ̇ sin θ ) .
2
Eq. 2.40
Considerando pequenos os deslocamentos θ e ψ e que φ̇ = Ω
Td =
1
Md (u̇2 + ẇ2 ) + Idx (ψ̇ 2 + θ̇ 2 ) + Idy (Ω2 + ψ̇ sin2 θ + 2Ωψ̇) .
2
Eq. 2.41
O último termo Idy (Ω2 + 2Ωψ̇θ ) representa a energia e o efeito giroscópico do disco girando
a uma velocidade angular Ω, na ordem da soma, respectivamente. Não há energia potencial
associada ao disco, já que por hipótese, o mesmo foi considerado rígido.
Capítulo 2 Revisão Teórica
15
Formulação por Elementos Finitos
Para descrever o deslocamento do centro de um disco rígido é necessário somente um “nó”
e quatro graus de liberdade, dois de translação u e w e dois para rotação θ e ψ, sendo sua coordenada generalizada q representada por
q = [u, w, θ , ψ]T .
Eq. 2.42
Aplicando a equação de Lagrange (Equação 2.13) na Equação 2.40, obtém-se segundo LALANNE e FERRARIS (2001), as matrizes de contribuição de massa e giroscópica abaixo:




M
0
0 0
0 0 0
0
 d



 0 M



0
0
0
0
0
0
d




[Md ] = 
Eq. 2.43
 e [Gd ] = Ω 


 0 0 0 −Idy 
 0
0
I
0
dx




0
0
0 Idx
0 0 Idy
0
2.2.2
Eixo
Para um elemento de eixo com seção circular S e momento de inércia de área I, constantes
em um segmento finito de comprimento L e massa específica ρ, a energia cinética Te pode
ser obtida por extensão da equação do disco (Equação 2.40) através da integração ao longo do
comprimento L como
ρS
Te =
2
Z L
0
ρI
(u̇ + ẇ )dy +
2
2
2
Z L
2
2
Z L
2
(ψ̇ + θ̇ )dy + ρILΩ + 2ρIΩ
0
ψ̇θ dy.
Eq. 2.44
0
Considerando que o momento de inércia de área I da seção S pode ser obtido da expressão
Z
Ix =
S
Z
2
x dS e Iz =
z2 dS,
Eq. 2.45
S
e que o eixo é axi-simétrico e não esta submetido a carga axial, LALANNE e FERRARIS (2001)
demonstra que a energia potencial Ue integrada sobre a seção S ao longo do comprimento do
eixo L é dada por
EI
Ue =
2
"
Z L 2 2
∂ u
0
∂ y2
∂ 2w
+
∂ y2
2 #
dy.
Eq. 2.46
Formulação por Elementos Finitos
O eixo é dividido em n segmentos de comprimento ∆y de seção constante. Cada nó da
seção elementar tem quatro graus de liberdade. Para pequenos deslocamentos u e w podem-se
Capítulo 2 Revisão Teórica
16
estabelecer relações entre deslocamentos lineares e angulares conforme
θ=
∂w
∂u
e ψ =− .
∂y
∂y
Eq. 2.47
Figura 2.2: Modelo de viga discretizada.
Cada elemento, é em essência, um modelo de um pequeno sólido deformável, no qual
alguns graus de liberdade substituem os infinitos graus de liberdade do sistema contínuo (Figura
2.2). Dentro de cada elemento o deslocamento das coordenadas {u} e {w} é aproximado por
uma combinação linear de funções de forma N(y) e deslocamentos nodais {γ}
{u} = [N1 ]{γu},
Eq. 2.48
{w} = [N2 ]{γw}.
Eq. 2.49
e
As funções de forma são matrizes de uma linha, [N1 ] e [N2 ], cujos elementos variam de acordo
com o modelo de viga utilizado. Estas funções não dependem da variável tempo e por isto, as
velocidades podem ser expressas por
{u̇} = [N]{q̇}.
Eq. 2.50
De modo geral, pode-se escrever a energia cinética T como
1
T = {γ̇}T [M]{γ̇},
2
Eq. 2.51
onde a matriz de massa é obtida da expressão
Z
[M] =
V
ρN T NdV.
Eq. 2.52
Capítulo 2 Revisão Teórica
17
Substituindo adequadamente as funções de forma na Equação 2.44, obtém-se a expressão
Te = Tem + Tes + Teg + ρILΩ2 ,
Eq. 2.53
onde:
1
1
Tem = γ u̇T M1 γ u̇ + γ ẇT M2 γ ẇ,
2
2
1
1
Tes = γ u̇T M3 γ u̇ + γ ẇT M4 γ ẇ,
2
2
T
Teg = Ωγ u̇ M5 γw.
Eq. 2.54
Eq. 2.55
Eq. 2.56
As matrizes M1 e M2 são matrizes consistentes clássicas de massa, M3 e M4 representam a inércia rotatória e M5 o efeito giroscópico, segundo LALANNE e FERRARIS (2001). Finalmente
ao se aplicar a equação de Lagrange (Equação 2.35) na Equação 2.53 obtém-se
d ∂ Te
∂ Te
= ([M] + [Ms ]) γ̈ + [G]γ̇.
−
dt ∂ γ̇
∂γ
As matrizes assim obtidas são


156
0
0
−22L
54
0
0
13L


 0

156
22L
0
0
54
−13L
0




2
2
 0
22L
4L
0
0
13L −3L
0 




2
2

−22L
0
0
4L
−13L
0
0
−3L
ρSL 

.
[M] =


420  54
0
0
−13L 156
0
0
22L 


 0
54
13L
0
0
156 −22L
0 




2
2
 0

−13L
−3L
0
0
−22L
4L
0


2
2
13L
0
0
−3L
22L
0
0
4L

36
0
0
−3L −36

 0
36
3L
0


 0
3L 4L2
0


0
4L2
ρI 
 −3L 0
[Ms ] =
30L 
0
0
3L
 −36

 0
−36 −3L
0


 0
3L −L2
0

−3L 0
0
−L2
0
0
3L
36
0
0
3L
0
0
−3L
3L
Eq. 2.58





−3L −L2
0 


2
0
0
−L 
.

0
0
3L 

36 −3L
0 


2
−3L 4L
0 

0
0
4L2
−36
Eq. 2.57
0
Eq. 2.59
Capítulo 2 Revisão Teórica








ρIΩ 

[G] =
15L 







0
18
−36 −3L
0
0
36
−3L
36
0
0
−3L
−36
0
0
3L
0
0
−4L2 −3L
0
0
0
3L
4L2
0
0
−3L −L2
0
36
3L
0
0
−36
3L
−36
0
0
3L
36
0
0
3L
0
0
L2
−3L
0
0
0
3L
−L2
0
0
−3L
4L2
0


−3L 


2
L 


0 
.

0 

3L 


2
−4L 

0
Eq. 2.60
A energia de deformação do eixo pode ser obtida de modo análogo, a partir da Equação
2.46
1
1
Ue = γuT K1 γu + γwT K2 γw
Eq. 2.61
2
2
onde, segundo LALANNE e FERRARIS (2001), K1 e K2 são as matrizes clássicas de rigidez.
Para levar em conta o cisalhamento da seção da viga durante a flexão, pode-se introduzir um
fator de correção, proveniente da equação diferencial da viga Timoshenko, conforme CLARK
(1972), definindo a relação
12EI
.
Eq. 2.62
Gs Sr L2
A seção reduzida Sr = K 0 S, pode ser obtida da integração do perfil de tensão na área sob cisaa=
lhamento conforme BATHE (1996). Para para uma seção circular
K0 =
6(1 + ν)
w 0.9.
7 + 6ν
Eq. 2.63
Pode-se considerar Sr w S. O módulo de cisalhamento Gs para um material isotrópico é definido
por
E
,
2(1 + ν)
onde ν é o coeficiente de Poisson do material.
Gs =
Eq. 2.64
Aplicando a equação de Lagrange (Equação 2.13) na Equação 2.61, tem-se
∂Ue
= [Kc ]γ.
∂γ
Eq. 2.65
Capítulo 2 Revisão Teórica
19
A matriz obtida é

12
0
0
−6L

 0
12
6L
0


 0
6L (4 + a)L2
0


 −6L 0
0
(4 + a)L2
EI

[Kc ] =
(1 + a)L3 
0
0
6L
 −12

 0
−12
−6L
0


 0
6L (2 − a)L2
0

−6L 0
0
(2 − a)L2
2.2.3
−12
0
0
0
−12
6L
0
−6L (2 − a)L2
6L
0
0
12
0
0
0
12
−6L
0
6L
−6L (4 + a)L2
0
0
−6L






0


2
(2 − a)L 
.

6L



0



0

(4 + a)L2
Eq. 2.66
0
Mancais
Na modelagem dos mancais, somente forças de deslocamento são consideradas neste elementos. Os coeficientes de rigidez e amortecimento são considerados conhecidos e podem ser
determinados, via de regra, com auxílio de aplicativos fornecidos pelos fabricantes dos mancais. Para mancais hidrodinâmicos os parâmetros em geral, são dependentes da carga radial, do
lubrificante, da folga e da rotação. Entretanto no corrente trabalho, são considerados constantes
por hipótese. O modelo leva em conta termos cruzados para amortecimento viscoso e rigidez,
conforme esquematizado na Figura 2.3.
Figura 2.3: Modelo para rigidez e amortecimento do mancal.
O trabalho virtual das forças que atuam sobre o eixo é dado pela relação
δW = −krδ r − cṙδ r,
Eq. 2.67
Capítulo 2 Revisão Teórica
20
sendo
krδ r = kxx uδ u + kxz wδ u + kzx uδ w + kzz wδ w,
Eq. 2.68
cṙδ r = cxx u̇δ u + cxz ẇδ u + czx u̇δ w + czz ẇδ w,
Eq. 2.69
ou na forma matricial

( ) 
( )
(
)
Fu
kxx kxz
u
cxx cxz
u̇


= −
−
.
Fw
w
ẇ
kzx kzz
czx czz
Eq. 2.70
Formulação por Elementos Finitos
Para os mancais, são desprezadas possíveis reações à momentos fletores e rotações, sendo
as matrizes obtidas diretamente da Equação 2.67
Fu = −kxx u − kxz w − cxx u̇ − cxz ẇ
Eq. 2.71
Fw = −kzz w − kzx u − czz ẇ − czx u̇,
Eq. 2.72
e
como é suposta nula a rigidez do mancal para qualquer rotação θ ou ψ. Assim, as forças
decorrentes também são nulas. Ou seja Fθ = 0 e Fψ = 0, obtendo-se a equação


 









F 
k
0 kxz 0 
u 
c
0 cxz 0 
u̇ 





 u 

 xx



 xx








 F 
 0 0 0 0   θ   0 0 0 0   θ̇ 

θ




= −
−
.




 kzx 0 kzz 0  

 czx 0 czz 0  

F
w
ẇ






w










 





 F 

 ψ 


 ψ̇ 
0
0
0
0
0
0
0
0
ψ
2.2.4
Eq. 2.73
Excitação Tipo Desbalanceamento
Em sistemas rotativos o desbalanceamento residual é uma excitação sempre presente em
maior ou menor grau. A força depende da rotação Ω e segundo LALANNE e FERRARIS
(2001), é dada pela relação
F(Ω) = mre Ω2 ,
Eq. 2.74
onde mre representa diretamente o desbalanceamento residual, cuja unidade comumente empregada é o grama-milímetro (1(g.mm) = 10−6 (kg.m)).
Capítulo 2 Revisão Teórica
21
Figura 2.4: Massa residual mu .
A posição D pode ser obtida em termos de coordenadas generalizadas, para uma posição y
qualquer constante, através da transformação de coordenadas




u
+
r
sin
Ωt


2


D=
,
constante




 w + r cos Ωt 
Eq. 2.75
2
e a velocidade através de
Ḋ =



 u̇ + r2 Ω cos Ωt




.
0




 ẇ − r Ω sin Ωt 
2
Eq. 2.76
A energia cinética associada à massa residual mu , conforme LALANNE e FERRARIS (2001)
pode ser aproximada através da relação
Tu ∼
= mu r2 Ω(q̇1 cos Ωt − q̇2 sin Ωt).
Eq. 2.77
Aplicando a equação de Lagrange (Equação 2.35) na expressão para energia cinética (equação
Equação 2.77), para uma massa mu situada no eixo z em t = 0 (Figura 2.4), segundo LALANNE
e FERRARIS (2001), obtém-se
d ∂ Tu ∂ Tu
−
= −(mu r)Ω2
dt ∂ γ̇
∂γ
"
sin Ωt
cos Ωt
#
.
Eq. 2.78
A transformada de Fourier para a excitação periódica (desbalanceamento) a ser considerada
no vetor de força {F(Ω)} para as coordenadas u e w correspondentes à i-ésima coordenada
Capítulo 2 Revisão Teórica
22
generalizada excitada tem a forma

 





...
...
















2




i(m
r)Ω
u




u












 w   −(m r)Ω2 
u
.
⇒
{F(Ω)} ⇒




0
θ
























0
ψ





 






 

...
...
2.2.5
Eq. 2.79
Montagem das Matrizes Globais do Sistema
As matrizes globais de massa [M], rigidez [K] e amortecimento [C], do sistema discretizado em elementos finitos, são obtidas através da montagem de seus elementos individuais,
sobrepostas nos nós em comum, conforme CHILDS (1993).
Figura 2.5: Condições de equilíbrio para 2 graus de liberdade por nó.
Da Figura 2.5, as condições de equilíbrio são
f1 = f11 , f2 = f21 , f3 = f31 + f12 , f4 = f41 + f22 , f5 = f32 , f6 = f42 .
Eq. 2.80
Observando os índices das forças, procedendo as mesmas somas sobre as matrizes individuais
Capítulo 2 Revisão Teórica
de massa dos elementos [M1 ] e [M2 ], obtém-se a matriz global

0
m114
m113
m111 m112

1
1
1
1
 m
0
m24
m23
 21 m22
 1
 m31 m132 m133 + m211
m134 + m212 m213

[M] = 
m144 + m222 m223
 m141 m142 m143 + m221

 0
m233
m232
0
m231

m243
m142
0
0
m241
23
0


0 


2
m14 
.

m224 

m234 

2
m44
Eq. 2.81
O mesmo procedimento pode ser estendido para a matriz de rigidez [K]. Tem-se uma matriz
de interface entre elementos, do número de graus de liberdade do nó, no exemplo acima 2 × 2.
Para quatro graus de liberdade, tem-se uma matriz de superposição 4 × 4 conforme o esquema
da Figura 2.6.
Figura 2.6: Montagem da matriz global, para 6 elementos 8 × 8.
Os termos das matrizes de rigidez e amortecimento dos mancais são adicionadas diretamente aos termos do nó do elemento de eixo onde o mesmo está aplicado.
Capítulo 2 Revisão Teórica
2.3
24
O Material Viscoelástico
Nesta seção é revisado o comportamento do material viscoelástico e o modelo utilizado para
descrever seu comportamento, já que o mesmo é empregado para o projeto de neutralizadores
dinâmicos.
Materiais viscoelásticos, devido sua alta capacidade de dissipar energia, são utilizados com
êxito no controle e isolamento de vibrações e ruído. Este tipo de material tem o comportamento
elástico definido pela deformação imposta e também pressupõe a existência de uma função de
relaxação da tensão devido a sua capacidade de escoamento viscoso.
Suas características dinâmicas, o módulo de cisalhamento (Gr ) e fator de perda (η), variam
com a temperatura e com a freqüência, principalmente.
2.3.1
Relação Constitutiva Generalizada
Uma das principais características de materiais viscoelásticos é que a tensão σ (t) e suas
derivadas dependem da deformação ε(t) e suas derivadas, sendo que, de modo generalizado,
pode-se escrever (2 ) a relação
b0 σ (t) + b1
d κ1 σ (t)
d κn σ (t)
d β1 ε(t)
d βm ε(t)
+
...
+
b
=
a
ε(t)
+
a
+
...
+
a
.
n
o
m
1
dt κ1
dt κn
dt β1
dt βm
Eq. 2.82
Segundo JONES (2001), com este modelo não é possível representar corretamente o comportamento dinâmico do material viscoelástico, numa ampla faixa de freqüência, com poucos termos
da relação acima (Equação 2.82) e utilizando somente valores inteiros para β e κ.
No entanto, operadores integro-diferenciais podem ser definidos para quaisquer β e κ reais
(3 ).
No domínio do tempo, a expressão diferencial generalizada em termos de β ∈ R | 0 < β < 1,
assume uma forma pouco usual para aplicação direta (Riemann-Liouville)
d β σ (t)
1
d
=
β
Γ(1 − β ) dt
dt
Z t
0
σ (τ)
dτ.
(t − τ)β
Eq. 2.83
A função gama que generaliza a função fatorial para números reais, Γ(n + 1) = n! é definida
como
Γ(x) =
Z ∞
t (x−1) e−t dt.
Eq. 2.84
0
O modelo de derivada fracionária encontra sua real aplicação no domínio da freqüência. Apli2 Para
materiais elásticos clássicos, no caso unidimensional, tensão e deformação estão relacionados pela lei de
Hooke: σ , Eε, onde E é o módulo de elasticidade do material (Young).
3 A diferenciação também pode ser definida para valores complexos de β ou κ.
Capítulo 2 Revisão Teórica
25
cando a transformada de Fourier sobre a Equação 2.82 reduzida aos primeiros dois termos,
obtém-se a expressão para o módulo de elasticidade complexo dada pela relação
Ē(Ω) =
a0 + a1 (iΩ)β
,
b0 + b1 (iΩ)κ
Eq. 2.85
assumindo que κ = β , obtém-se o modelo de derivada fracionária de quatro parâmetros
Ē(Ω) =
ao + a1 (iΩ)β
.
b0 + b1 (iΩ)β
Eq. 2.86
Este modelo, como apresentado por PRITZ (1996) e JONES (2001), representa de forma satisfatória o comportamento dinâmico de uma grande gama de materiais viscoelásticos utilizados
em engenharia.
A expressão para o módulo de cisalhamento complexo, obtida de modo totalmente similar
a Equação 2.86, é dada por
Ḡ(Ω) =
a0 + a1 (iΩ)β
.
b0 + b1 (iΩ)β
Eq. 2.87
A representação complexa para o módulo de cisalhamento, a uma determinada freqüência e
temperatura pode ser redefinida como
Ḡ(Ω) = Gr (1 + iη(Ω)),
Eq. 2.88
sendo o “módulo de cisalhamento dinâmico” propriamente dito, dado pela parte real
Gr (Ω) = Re(Ḡ(Ω)),
Eq. 2.89
e o fator de perda do material pela relação entre as partes imaginária e real (4 ).
η(Ω) =
Im(Ḡ(Ω))
.
Re(Ḡ(Ω))
Eq. 2.90
A mesma representação se aplica ao módulo de elasticidade dinâmico Ē(Ω).
2.3.2
Efeito da Temperatura
O módulo de cisalhamento e o fator de perda dos materiais viscoelásticos variam com a
temperatura de operação. Em geral, o módulo de cisalhamento aumenta com a freqüência e
diminui com a temperatura.
Para materiais ditos “plásticos” a temperatura de “amolecimento” Ts é razoavelmente alta,
4η
representa a tangente do argumento do número complexo Ḡ(Ω)(o ângulo de defasagem).
Capítulo 2 Revisão Teórica
26
enquanto que para um elastômero típico é muito mais baixa, em geral, abaixo da temperatura
ambiente. Acima da temperatura de amolecimento Ts , o módulo de cisalhamento Gr cai rapidamente, enquanto o fator de perda η sobe atingindo seu máximo, voltando a cair em seguida.
Este comportamento estabelece a região de transição. A maioria dos materiais viscoelásticos
utilizados em engenharia, tem uma faixa de temperatura de trabalho que compreende três zonas
distintas:
• Região vítrea (zona III), em baixas temperaturas. o material apresenta baixo fator de
amortecimento e o módulo de elasticidade / cisalhamento atinge seu máximo. Nesta
faixa, o material deforma pouco, e tem por isso boa estabilidade estrutural (5 ).
• Região de transição (zona II). O módulo de elasticidade / cisalhamento diminui e o fator
de perda aumenta até atingir o máximo. Neste ponto a variação do módulo é máxima. O
material nesta região é aplicado em neutralizadores, para aproveitar o alto valor do fator
de perda. Não tem boa estabilidade estrutural.
• Região “Rubberlike” (zona I), em altas temperaturas. O módulo de elasticidade / cisalhamento e fator de perda apresentam baixos valores, o material tem estabilidade estrutural
baixa.
Segundo ESPÍNDOLA e SILVA (1992), alguns materiais possuem propriedades à temperatura
de trabalho, correspondentes as das regiões acima, sendo por isso denominados propriamente
como: Tipos I, II e III (Figura 2.7).
Borracha f = 10 (Hz)
3
Gr
eta
2
log(Gr (MPa)),log(eta)
1
Vitreo
Transicao
0
−1
−2
"Rubberlike"
−3
−4
−60
−40
−20
0
20
40
60
80
100
120
140
Temperatura (°C)
Figura 2.7: Variação com temperatura (freqüência constante).
5 Baixa
estabilidade estrutural traduz-se por uma elevada fluência sob tensão constante.
Capítulo 2 Revisão Teórica
27
Uma forma comum de apresentação dos valores das propriedades de materiais viscoelásticos é o nomograma como apresentado por JONES (2001). Nele são apresentadas simultaneamente variações para temperatura e freqüência. A escala encontrada para o eixo das abscissas
é denominada freqüência reduzida log[ f α(T )]. A variação da temperatura é encontrada sob a
forma de linhas iso-temperatura. No eixo das ordenadas y à direita se entra com a freqüência
de trabalho, segue-se horizontalmente até interceptar a linha iso-temperatura correspondente à
temperatura de trabalho requerida. Neste ponto procura-se verticalmente pelas curvas de fator de perda η e pelo módulo de cisalhamento G cujos valores estão na escala logarítmica das
ordenadas à esquerda.
Abaixo, a Figura 2.8 apresenta um nomograma típico de um material viscoelástico utilizado
em engenharia e comercializado pela empresa E-A-R Specialty Composites como “isodamp C1002”.
Isodamp C−1002
Temperatura (C)
40
30
20
10
0
−10
−20
C.Perda(pu) Cisa.(Pa)
10
1.00
0.10
0.01
1E+04
Cisa.
C.Perda
8
1E+03
7
1E+02
6
1E+01
10
10
10
1E+00
5
10 −3
10
Frequencia (Hz)
9
−1
10
1
10
3
10
5
10
7
10
9
10
Frequencia Reduzida
Figura 2.8: Nomograma material Isodamp C-1002.
2.3.3
Efeito da Freqüência
Para materiais elásticos sólidos, o efeito da freqüência é pequeno, nos materiais viscoelásticos porém, são muito evidentes . Qualitativamente o efeito é inverso ao da temperatura, ou
seja, um aumento na freqüência tem um efeito similar a uma diminuição na temperatura.
Por outro lado, diferença significativa ocorre na amplitude das escalas. Enquanto uma
mudança de região devido a temperatura pode ocorrer em no máximo algumas centenas de
graus, uma alteração correspondente em freqüência pode demandar varias ordens de magnitude.
Capítulo 2 Revisão Teórica
28
Uma escala em freqüência pode ir tipicamente de 10−8 até 108 (Hz) (6 ).
2.3.4
Outros Efeitos
Amplitude da deformação dinâmica: Os parâmetros do material viscoelástico sofrem variação em função da deformação. Para grandes deformações cíclicas, devido ao elevado fator de
perda, é provocada uma geração interna de calor e portanto um aumento na temperatura do material. Esta variação é muito difícil de modelar. Em geral, o módulo de cisalhamento se reduz e
o fator de perda aumenta com o aumento da deformação. Se a deformação aumenta ainda mais,
o fator de perda volta a diminuir e o módulo de cisalhamento tende à estabilidade (NASHIF et
al. (1985)).
Pré-carga estática: Importante quando o material trabalha na região “Rubberlike”. O módulo de cisalhamento aumenta, enquanto o fator de perda cai com o aumento da pré-carga
(NASHIF et al. (1985)).
2.3.5
Efeito Combinado de Temperatura e Freqüência
Quando tanto a freqüência quanto a temperatura variam, como geralmente acontece na
prática, é possível combinar os efeitos da variação de ambos em uma só variável, a freqüência
reduzida. Esta combinação é baseada no fato que o valor do módulo da rigidez complexa, em
qualquer freqüência f1 e qualquer temperatura de referência T1 arbitrários, é idêntico a um outro,
em outra freqüência f2 e temperatura T2 de forma a satisfazer a relação
Ḡ( f1 , T1 ) = Ḡ ( f2 α(T2 )) ,
Eq. 2.91
a função α(T2 ) pode ser determinada experimentalmente e plotada num gráfico do tipo log α(T )
versus T1 .
As temperaturas em questão são absolutas (K). Uma das funções α(T ) mais utilizadas,
denominada “shift factor” ou fator de deslocamento é a WLF (Williams-Landel-Ferry), onde
θ1 e θ1 são constantes à determinar, que dependem do material e T0 é uma temperatura de
referência, arbitrariamente escolhida. O fator de deslocamento é definido pela equação
log α(T ) = −θ1
(T − T0 )
.
(θ2 + T − T0 )
Eq. 2.92
Determinando-se adequadamente os parâmetros da Equação 2.86, obtém-se as relações:
6A
escala no tempo pode corresponder a anos.
Capítulo 2 Revisão Teórica
b0 = 1 e
a0
b0
=
E0
b0
29
= E0 , logo
a1
b1
=
E1
b1
= E∞ .
Com o fator de deslocamento α(T ) determinado, pode-se obter as expressões finais para os
módulos de elasticidade e cisalhamento em função da freqüência reduzida Ωα(T ):
Ē(Ω) =
E0 + E∞ b1 (iΩα(T ))β
1 + b1 (iΩα(T ))β
Eq. 2.93
Ḡ(Ω) =
G0 + G∞ b1 (iΩα(T ))β
,
1 + b1 (iΩα(T ))β
Eq. 2.94
onde G0 representa o valor assintótico do módulo para freqüências muito baixas e G∞ para
freqüências muito altas.
A potência de um número imaginário puro pode ser representada por uma função trigonométrica complexa: iβ = cos β π2 + i sin β π2 . Aplicando esta relação na Equação 2.94, obtém-se
G0 + G∞ b1 (Ωα(T ))β cos β π2 + i sin β π2
,
Ḡ(Ω) =
1 + b1 (Ωα(T ))β cos β π2 + i sin β π2
Eq. 2.95
o que possibilita a separação da parte real Gr (Ω) e do fator de perda η(Ω) conforme
G0 + (G0 + G∞ ) b1 (Ωα(T ))β cos β π2 + G∞ b21 (Ωα(T ))2β
Gr (Ω) =
,
1 + 2b1 (Ωα(T ))β cos β π2 + b21 (Ωα(T ))2β
Eq. 2.96
(G∞ − G0 ) b1 (Ωα(T ))β sin β π2
.
G0 + (G0 + G∞ ) b1 (Ωα(T ))β cos β π2 + G∞ b21 (Ωα(T ))2β
Eq. 2.97
η(Ω) =
2.4
Sistemas com Um Grau de Liberdade
Nesta seção são revisados os modelos para sistemas de um grau de liberdade do tipo massamola-amortecimento, abordando os modelos de amortecimento viscoso e viscoelástico e definições de grandezas relacionadas, como introdutório para o controle de vibração.
2.4.1
Modelo Viscoso
Aplicando a equação do movimento (Equação 2.37) para um sistema de um grau de liberdade (Figura 2.9), o sistema de equações diferenciais resume-se à expressão
mẍ(t) + cẋ(t) + kx(t) = f (t).
Eq. 2.98
Capítulo 2 Revisão Teórica
30
Figura 2.9: Sistema massa mola 1gl, modelo viscoso.
Considerando o sistema em regime permanente, pode-se aplicar diretamente a transformada
de Fourier e obter a expressão
(−Ω2 m + iΩc + k)X(Ω) = F(Ω).
Eq. 2.99
A função resposta em freqüência (FRF), receptância, é definida pela relação
X(Ω)
1
=
.
2
F(Ω) −Ω m + iΩc + k
H(Ω) ,
Também podem ser definidas ainda expressões para:
• Mobilidade, γ(Ω) = iΩH(Ω).
• Inertância, A(Ω) = −Ω2 H(Ω).
As funções complementares inversas, também características do sistema são:
• Rigidez dinâmica, K(Ω) ,
F(Ω)
X(Ω)
• Impedância mecânica, Z(Ω) ,
• Massa dinâmica, M(Ω) ,
= −Ω2 m + iΩc + k.
K(Ω)
iΩ .
K(Ω)
.
−Ω2
São definidas ainda as relações auxiliares para:
• Freqüência natural não amortecida, Ωn ,
• Amortecimento crítico, cc , 2mΩn .
• Razão de amortecimento, ζ ,
• Razão de freqüências, ε =
Ω
Ωn .
c
cc .
q
k
m.
Eq. 2.100
Capítulo 2 Revisão Teórica
31
Utilizando as relações auxiliares obtém-se a expressão adimensional para a receptância H(Ω)
Had (Ω) =
H(Ω)
1
=
.
1/k
(1 − ε 2 + 2iζ ε)
Eq. 2.101
Receptância Adimensional |Hadi|
1GL − Viscoso qsi=0.05
20
15
|Hadi| (dB ref=1)
10
5
0
−5
−10
−15
−20
0
50
100
150
200
250
300
350
400
450
500
frequencia (Hz)
Figura 2.10: Função resposta em freqüência adimensional.
A Figura 2.10 apresenta a receptância para o modelo de amortecimento viscoso (Figura 2.9)
de um grau de liberdade com m =1 (kg), ζ = 0.05 e k=1E6 (N/m).
Rigidez Dinamica Adimensional |Kadi|
1GL − Viscoso qsi=0.05
20
15
|Kadi| (dB ref=1)
10
5
0
−5
−10
−15
−20
0
50
100
150
200
250
300
350
400
450
500
frequencia (Hz)
Figura 2.11: Rigidez dinâmica adimensional.
Como pode ser observado na Figura 2.11, que apresenta a rigidez dinâmica, o sistema não
Capítulo 2 Revisão Teórica
32
oferece resistência à vibração quando excitado na sua freqüência natural. Por esta razão, a
resposta característica do sistema apresenta valores elevados naquela freqüência.
2.4.2
Modelo Viscoelástico
O arranjo da Figura 2.12 mostra um sistema de um grau de liberdade, tipo massa-elemento
resilente para o modelo de amortecimento viscoelástico. Fazendo o diagrama de corpo livre e
aplicando a segunda lei de Newton ∑ F = mẍ(t), pode-se obter a expressão
f (t) − LG(Ω)x(t) = mẍ(t).
Eq. 2.102
Figura 2.12: Sistema massa-elemento resilente 1gl.
No domínio do tempo, Ω representa a freqüência angular da excitação f (t) do tipo senoidal.
Levando a Equação 2.102 para o domínio da freqüência através da transformada de Fourier,
obtém-se a relação
H(Ω) =
X(Ω)
1
=
.
2
F(Ω) −Ω m + LG(Ω)
Eq. 2.103
A Equação 2.101 é mais geral no domínio da freqüência, sendo válida para qualquer Ω.
Uma expressão adimensional equivalente à Equação 2.101 é
Had (Ω) =
H(Ω)
1
=
.
1/LGr (Ωn ) −ε 2 + r(Ω)(1 + iη(Ω))
Definindo-se:
• Relação de cisalhamentos, r(Ω) =
• Razão das freqüências, ε =
Gr (Ω)
Gr (Ωn ) .
Ω
Ωn .
A freqüência natural do sistema é obtida através da equação
r
LGr (Ωn )
Ωn ,
.
m
Eq. 2.104
Capítulo 2 Revisão Teórica
33
Para uma massa de 1 (kg) e com os dados do modelo de derivada fracionária do material viscoelástico empregado definidos como sendo o fator de forma do material L = 0.064, a temperatura
de referência T0 = 273 (K), a temperatura de operação T1 = 293 (K), o módulo cisalhamento
inferior G0 = 1.53E6 (N/m2 ), o módulo cisalhamento superior, G∞ = 1.11E8 (N/m2 ), a potência da derivada β1 = 0.0807, o deslocamento de freqüência b1 = 11.34E − 2, o coeficiente
temperatura 1, θ1 = 15.1, e o coeficiente temperatura 2 θ2 = 171.
Receptância Adimensional |Hadi|
1GL − Viscoelastico
20
15
|Hadi| (dB ref=1)
10
5
0
−5
−10
−15
−20
0
50
100
150
200
250
300
350
400
450
500
frequencia (Hz)
Figura 2.13: Função resposta em freqüência adimensional.
A Figura 2.13 apresenta a receptância adimensional para o modelo de amortecimento viscoelástico (Figura 2.12) de um grau de liberdade.
Para determinar a rigidez K destes materiais, pode-se definir o parâmetro L como um fator
de forma geométrico, definido pela relação da área carregada e a espessura ou comprimento de
um elemento viscoelástico sujeito à tração - compressão (a) ou cisalhamento (b) (Figura 2.14).
Para o cisalhamento (b) o fator de forma L é definido por
L=
A
.
h
Eq. 2.105
Para tração - compressão (a), pode-se considerar o módulo de elasticidade aproximadamente
E∼
= 3G, logo
A
A
Eq. 2.106
K(a) = 3 G e K(b) = G.
h
h
Capítulo 2 Revisão Teórica
34
Figura 2.14: Área de carregamento: a) axial, b) cisalhamento.
Para certas aplicações do tipo tração-compressão, onde h é muito pequeno em relação a A,
ou onde a forma do material viscoelástico não permite deformação lateral ou ainda a mesma é
limitada, não é possível utilizar diretamente o módulo de elasticidade E para calcular a rigidez
K. Nestes casos SNOWDON (1959) sugere um módulo aparente definido pela equação
Ea = (1 + βe S2 )E,
Eq. 2.107
onde βe é uma constante numérica que depende da forma e S é a razão de uma área carregada pela área total livre. Para uma seção circular, quadrada ou moderadamente retangular é
aproximadamente igual a dois βe ∼
= 2, segundo NASHIF et al. (1985).
2.5
Controle de Vibração
Nesta seção são apresentados conceitos básicos de controle de vibração, utilizando o modelo de dois graus de liberdade.
Uma forma de controle do nível de vibração consiste em modificar a estrutura do sistema
primário representado pela equação do movimento (Equação 2.37), de tal forma que essa alteração provoque uma diminuição da amplitude do deslocamento resultante {q}.
Dependendo da circunstância, pequenas alterações nos parâmetros podem conseguir grandes reduções nas amplitudes. Por exemplo, para um sistema trabalhando em ressonância (7 ) é
possível alterar a rigidez k, a massa m, o amortecimento c, ou fazer uso de neutralizadores dinâmicos de vibração. Muitas vezes é impossível alterar os parâmetros do sistema ou a excitação é
do tipo banda larga de freqüência de forma que uma modificação estrutural (através da alteração
dos parâmetros físicos não é suficiente para a solução do problema. Nestes casos a aplicação de
neutralizadores dinâmicos viscoelásticos resultam em um controle comprovadamente eficaz.
7 As
freqüências de excitação e natural do sistema coincidem Ω = Ωn .
Capítulo 2 Revisão Teórica
2.5.1
35
Rigidez Dinâmica na Base
Quantifica o deslocamento na base devido a uma excitação f (t) aplicada na base.
Modelo Viscoelástico
Figura 2.15: Modelo viscoelástico no domínio da freqüência.
De forma equivalente ao modelo viscoso, a rigidez dinâmica na base, para o modelo viscoelástico da Figura 2.15, pode ser obtida de através da relação
Kb (Ω) ,
2.5.2
(LG(Ω))(−mΩ2 )
F(Ω)
=
.
Xb (Ω)
−mΩ2 + LG(Ω)
Eq. 2.108
Sistemas com Dois Graus de Liberdade
Analisando a curva da Figura 2.11 se pode verificar, pela resposta em freqüência do sistema,
que quando a freqüência de excitação se aproxima da freqüência natural (em geral rigidez dinâmica pequena), a resposta do sistema é amplificada. Mesmo para excitações f (t) com pequena
amplitude, quando a freqüência de excitação e a freqüência natural se tornam iguais (Ω = Ωn ) ,
a resposta terá níveis elevados se o amortecimento do sistema primário é baixo. Por outro lado,
na base de um sistema de um grau de liberdade, a curva da Figura 2.16 mostra um aumento da
rigidez, justamente na freqüência natural. Pode-se generalizar que todo sistema em qualquer
uma das freqüências naturais, apresenta uma impedância mecânica baixa. Por outro lado, todo
sistema de um grau de liberdade, na freqüência natural, apresenta uma impedância mecânica
elevada na base. Então, acrescentando ao sistema primário, um novo sistema projetado para
Capítulo 2 Revisão Teórica
36
que a freqüência natural seja apropriada, pode-se injetar uma impedância elevada, limitando a
amplitude da resposta do novo sistema, reduzindo assim o nível de vibração.
Rigidez Dinamica na Base |Kb|
1GL − Viscoso qsi=0.05
140
130
|Kb| (dB ref=1)
120
110
100
90
80
70
60
0
50
100
150
200
250
300
350
400
450
500
frequencia (Hz)
Figura 2.16: Gráfico da rigidez dinâmica Kb (Ω) na base.
Este é o princípio de funcionamento do neutralizador dinâmico de vibração.
Princípio de Controle em Banda Estreita
Acoplando um sistema de controle ao sistema primário, agrega-se mais um grau de liberdade e uma nova freqüência natural surge na função resposta em freqüência resultante. Desprezando os amortecimentos ca e c1 (auxiliar e primário), a função de transferência na massa m1 ,
obtida como anteriormente, pode ser definida pela relação (receptância), conforme a equação
H11 (Ω) =
(−Ω2 ma + ka )
X1 (Ω)
.
=
F1 (Ω) [−Ω2 m + (k + ka )](Ω2 ma + ka ) − ka2
Eq. 2.109
Como se pode observar na Equação 2.109, o numerador tende a zero, quando sintonizada
a freqüência do sistema auxiliar, ou seja Ω = Ωa , e a resposta do sistema primário é anulada.
Esta condição, que aparentemente é desejável, apresenta dois problemas importantes:
1. A fadiga do neutralizador, que ocorre em geral logo após sua instalação, acarretando dano
estrutural.
2. A controle é totalmente ineficaz quando a excitação apresenta espectro de banda larga,
devido a presença de novas freqüências naturais (linha tracejada na Figura 2.17).
Capítulo 2 Revisão Teórica
37
Acrescentando-se amortecimento ao neutralizador o controle pode ser projetado em uma ampla
faixa de freqüência.
Para um sistema primário massa-mola sem amortecimento HARTOG (1956) demonstrou
ser possível projetar de forma ótima um neutralizador dinâmico com amortecimento viscoso,
segundo o método dos pontos fixos. As curvas para vários valores de amortecimento ca passam
sempre por dois “pontos fixos” (Figura 2.17).
Segundo o método dos pontos fixos, a razão ótima entre a freqüência natural do neutralizador (Ωa ) e a freqüência natural do sistema primário (Ωn ), conforme a relação
αo =
Ωa
,
Ωn
Eq. 2.110
é definida pela equação
1
,
Eq. 2.111
1+µ
na qual µ representa a relação da massa do neutralizador (auxiliar) e do sistema primário sendo
αo =
definida através da relação
ma
.
Eq. 2.112
m
Geralmente se adota um valor entre 10 a 20 (%). Já a razão de amortecimento ótimo ζo pode
µ=
ser obtida através da equação
1
ζo =
αo
s
3µ
.
8(1 + µ)
Eq. 2.113
Receptância |H| 1GL+neutralizador − Considerando amortecimentos
Variacao amortecimento ca − modelo viscoso
−50
sem neut. qsi=0.0001
com neut. qsia=0.01
−60
com neut. qsia=0.05
com neut. qsia=0.15
−70
com neut. qsia=0.25
|H| (dB ref=1)
−80
Pt. Fixos
−90
−100
−110
−120
−130
−140
0
50
100
150
200
250
300
350
frequencia (Hz)
Figura 2.17: Pontos fixos para vários amortecimentos.
Capítulo 2 Revisão Teórica
38
Controle Ótimo
Em um sistema real, como os módulos de cisalhamento (elasticidade) G(Ω) do material
viscoelástico variam com a freqüência (para uma determinada temperatura) e também dado o
amortecimento presente no sistema primário não existem “pontos fixos”. Neste caso é praticamente impossível estabelecer as mesmas relações analíticas do método para os parâmetros
ótimos do neutralizador.
Para obter parâmetros ótimos em sistemas com materiais viscoelásticos, é possível utilizar
alguma técnica de otimização não linear, tendo como função objetivo a ser minimizada, o valor
máximo do módulo da função resposta em freqüência, numa faixa onde se encontra a amplificação a ser controlada, por exemplo. O vetor projeto, que recebe as variáveis a serem otimizadas,
precisa conter tão somente a freqüência natural do neutralizador Ωa .
A Figura 2.18 mostra o desempenho de dois modelos de neutralizador, com material viscoelástico e viscoso com fator de perda equivalente constante (ηcte = η(Ωa )), projetados com
parâmetros ótimos obtidos através de uma técnica de otimização não linear para um sistema
primário do tipo massa-elemento resilente.
Sistema 2GL − Material Viscoelastico wn = 11.402559(Hz)
mi = 0.15 Beta = 0.25 wa = 10.107146 wv = 10.174224(Hz)
35
sem neutralizador
otim. viscoelastico (wa)
30
otim. viscoso (wv)
25
|Hadi| (dB ref=1)
20
15
10
5
0
−5
−10
0
2
4
6
8
10
12
14
16
18
20
frequencia (Hz)
Figura 2.18: Desempenho dos modelos viscoelástico e viscoso.
Os resultados numéricos demonstram o desempenho superior do material viscoelástico embora ambos confirmem o sucesso do projeto ótimo. O neutralizador de modelo viscoso serve
apenas para comparação pois é de difícil construção prática.
Capítulo 2 Revisão Teórica
2.6
39
Parâmetros Equivalentes Generalizados
Nesta seção é reapresentada uma teoria originalmente apresentada por ESPÍNDOLA e
SILVA (1992), e que permite mostrar modelos dinamicamente equivalentes empregados no projeto ótimo de neutralizadores dinâmicos viscoelásticos.
Um sistema no domínio da freqüência, como o da Figura 2.19, pode ser modelado através
de uma abordagem por parâmetros equivalentes.
Figura 2.19: Modelo com material viscoelástico no domínio da freqüência.
Como apresentado por CRUZ (2004), pode-se pode definir um modelo equivalente viscoso
através das expressões
ce (Ω) , Re(Zb ),
Eq. 2.114
me (Ω) , Re(Mb ).
Eq. 2.115
Como pode ser mostrado, ambos modelos tem a mesma impedância dinâmica na base. Fazendo
o diagrama de corpo livre, aplicando a segunda lei de Newton e a transformada de Fourier (ℑ)
sobre o modelo equivalente (Figura 2.19 b) se obtém
F(Ω) = iΩme (Ω)V (Ω) + ce (Ω)V (Ω),
Eq. 2.116
V (Ω) = ℑ (ẋ(t)) .
Eq. 2.117
F(Ω)
,
V (Ω)
Eq. 2.118
onde
Da definição de impedância dinâmica
Zb (Ω) =
obtém-se
Zb (Ω) = iΩme (Ω) + ce (Ω).
Eq. 2.119
Utilizando as relações auxiliares já estabelecidas, sendo que todos os parâmetros são relativos ao
Capítulo 2 Revisão Teórica
40
sistema auxiliar e sendo Ωa sua freqüência natural, pode-se obter as expressões adimensionais
r(Ω)η(Ω)ε 3
(ε 2 − r(Ω))2 + (r(Ω)η(Ω))2
r(Ω) ε 2 − r(Ω) 1 + η(Ω)2
me (Ω) = −m 2
(ε − r(Ω))2 + (r(Ω)η(Ω))2
ce (Ω) = mΩa
Eq. 2.120
Eq. 2.121
Obtendo-se expressões adequadas para r(Ω) e η(Ω), é possível utilizar este conceito para diversos modelos de neutralizadores como mostrado em BAVASTRI et al. (2007).
Um sistema primário qualquer (Figura 2.20 a) ao qual foram acrescentados dispositivos
de controle de vibração pode ser modelado utilizando parâmetros equivalentes generalizados
(Figura 2.20 b).
Figura 2.20: Equivalência de um sistema primário no domínio da freqüência.
É fácil mostrar que a resposta em freqüência H(Ω) de um sistema de dois graus de liberdade pode ser obtida de forma generalizada a partir do diagrama de corpo livre, aplicando-se a
segunda lei de Newton no domínio da freqüência, conforme a Figura 2.21
H(Ω) =
X(Ω)
1
=
.
F(Ω) −Ω2 (m1 + me (Ω)) + iΩce (Ω) + LG(Ω)
Eq. 2.122
A expressão acima é totalmente similar a de um grau de liberdade e mais ainda, é dinamicamente equivalente à expressão completa para dois graus de liberdade H11 (Ω).
Capítulo 2 Revisão Teórica
41
Figura 2.21: Sistema dois graus de liberdade no domínio da freqüência.
A dinâmica do sistema composto (primário + neutralizador) foi escrita apenas em função da
coordenada do sistema primário. Isto é viabilizado através do emprego do modelo de parâmetros equivalentes generalizados para o sistema auxiliar. Essa técnica é importante pois permite
o projeto de neutralizadores utilizando apenas os parâmetros modais do sistema primário.
Na Tabela 2.1, são apresentados dois modelos de neutralizadores através de seus parâmetros equivalentes generalizados:
Tabela 2.1: Tabela de modelos equivalentes.
Modelo
Viscoso
Viscoelástico
2
2
r(Ω){ε 2 −r(Ω)[1+η(Ω)2 ]}
{ε −[1+(2ζ ε) ]}
me (Ω) −m (ε 2 −1)2 +(2ζ ε)2 −m (ε 2 −r(Ω))2 +(r(Ω)η(Ω))2
ce (Ω)
Ωa
4
ε
mΩa (ε 2 −1)2ζ2 +(2ζ
ε)2
q
k
m
3
r(Ω)η(Ω)ε
mΩa (ε 2 −r(Ω))
2 +(r(Ω)η(Ω))2
q
Lv Gr (Ωa )
m
Ficam definidos para o neutralizador:
• Massa equivalente, me (Ω) (kg).
• Constante de amortecimento viscoso equivalente, ce (Ω) (Ns/m).
• Freqüência natural, Ωa (rad/s).
• Massa do neutralizador, m (kg)
• Constante de mola, k (N/m)
• Constante de amortecimento viscoso, c (Ns/m).
• Fator de aplicação da carga sobre o material viscoelástico, Lv (m).
Capítulo 2 Revisão Teórica
42
• Parte real do módulo de cisalhamento complexo Re Ḡ(Ω) , Gr (Ω) (N/m).
Im(Ḡ(Ω))
• Fator de perda material viscoelástico Re(Ḡ(Ω)) , η(Ω).
• Razão cisalhamentos primário, r(Ω) =
• Razão de freqüências, ε(Ω) =
Gr (Ω)
Gr (Ωa ) .
Ω
Ωa .
• Razão de amortecimento modelo viscoso, ζ =
2.7
c
2mΩa .
Múltiplos Graus de Liberdade, Sistemas Girantes
A equação dinâmica de um sistema girante com n graus de liberdade, desprezada a matriz
de circulação [H] na Equação 2.37, é dada por
[M]{q̈(t)} + ([C] + [G(Ωrpm )]){q̇(t)} + [K]{q(t)} = { f (t)},
Eq. 2.123
sendo M, C e K matrizes simétricas constantes e reais.
Supondo que o sistema girante se encontre em regime permanente a uma dada rotação
(Ω = Ωrpm ), a Equação 2.123 toma a forma
[M]{q̈(t)} + [[C] + [G]]{q̇(t)} + [K]{q(t)} = { f (t)}.
Eq. 2.124
Segundo EWINS (2000), a consideração de um modelo viscoso não proporcional para a
matriz de amortecimento [C] e/ou da matriz giroscópica [G] na Equação 2.124 leva a solução
para um problema de 2n autovalores, o qual pode ser equacionado em um espaço 2n dimensional. Este espaço é denominado “espaço de estado” e é definido por
(
)
{q(t)}
{y(t)} =
.
{q̇(t)}
Eq. 2.125
2n×1
Desta forma a Equação 2.124 fica definida através da equação
h
i
h
i
[[C] + [G]] [M] {ẏ(t)} + [K] [0] {y(t)} = { f (t)}.
Eq. 2.126
Fazendo uso da seguinte analogia
h
i
[M] [0] {ẏ(t)} +
h
i
[0] −[M] {y(t)} = {0},
Eq. 2.127
Capítulo 2 Revisão Teórica
43
é possível obter uma expressão com dimensão 2n × 2n no espaço de estado




(
)
{ f (t)}
[K] [0]
[C1 ] [M]
 {y(t)} =
 {ẏ(t)} + 

,
{0}
[0] −[M]
[M] [0]
Eq. 2.128
na qual [C1 ] = ([C] + [G]) constante para uma dada rotação. Simplificando o sistema de equações pode ser reescrito-escrito na seguinte forma:
[A]{ẏ(t)} + [B]{y(t)} = {N(t)},
Eq. 2.129
na qual:

• [A] = 

• [B] = 
[C1 ] [M]
• {N(t)} =
2n×2n
[0]
[0] −[M]
(
{ f (t)}
{0}
,

[M] [0]
[K]


,

2n×2n
)
.
2n×1
O sistema da Equação 2.129 está apto para ser colocado como um problema de autovalores
generalizado com j = 1, 2n, do tipo
[B]{θ } j = λ j [A]{θ } j ,
Eq. 2.130
[B][θ ] = [λ ][A][θ ],
Eq. 2.131
permitindo a montagem do conjunto
sendo [λ ] a matriz diagonal de autovalores e [θ ] a matriz de autovetores do sistema no espaço de
estado. Supondo que no espaço de configurações se proponha a solução livre {q(t)} = {φ }est
e no espaço de estado {y(t)} = {θ }est , sendo s uma variável complexa igual a −λ e define-se
para uma dada rotação
(
{y(t)} =
{q(t)}
{q̇(t)}
)
(
⇒ {θ } =
{φ }
s{φ }
)
.
Eq. 2.132
Como a matriz [A] não é simétrica devido a presença de [G] em [C1 ], é necessário resolver o
conjunto de problemas adjuntos de autovetores [ψ]:
[B]T [ψ] = [λ ][A]T [ψ],
Eq. 2.133
Capítulo 2 Revisão Teórica
44
onde:
(
{ψ} =
{ϕ}
)
s{ϕ}
2.7.1
Eq. 2.134
Ortonormalização
A ortonormalização dos autovetores à esquerda e à direita, para uma dada rotação, conforme
ESPÍNDOLA (1992), pode ser obtida dividindo cada j-ésima coluna pela raiz quadrada dos
√
elementos correspondentes da diagonal da matriz r a jr = [ψ]T [A][θ ], a j j , obtendo-se as
matrizes [Ψ] e [Θ]. Depois desta operação, ESPÍNDOLA (1992) demonstra que
[Ψ]T [A][Θ] =
r
[Ψ]T [B][Θ] =
1r = [I],
Eq. 2.135
r
Eq. 2.136
λ jr .
O produto da j-ésima linha (coluna transposta) de [Ψ], a matriz [B] e a k-ésima coluna de
[Θ],considerando-se novamente s = −λ , é definido através da relação
(
(
)T 
)
)
(
h
i
{ϕ}
[K] [0]
{φ }
[K]{φ
}
k


= {ϕ}Tj −λ j {ϕ}T
,
s{ϕ}
s{φ
}
λk [M]{φ }k
[0]
−[M]
j
k
Eq. 2.137
{Ψ}Tj [B] {Θ}k = {ϕ}Tj [K]{φ }k − λ j λk {ϕ}T [M]{φ }k = b j δ jk ,
Eq. 2.138
com δ jk sendo definido pelo produto de Kronecker, ou seja δi6= j = 0 e δi= j = 1.
O produto da j-ésima linha (coluna transposta) de [Ψ], a matriz [A] e a k-ésima coluna de
[Θ], é dado pela expressão
(
)
(
)
(
)T 
h
i [C ]{φ } − λ [M]{φ }
{ϕ}
[C1 ] [M]
{φ }
1
k
k
k


= {ϕ}Tj −λ j {ϕ}T
,
s{ϕ}
s{φ
}
[M]{φ }k
[M]
[0]
j
k
Eq. 2.139
{Ψ}Tj [A] {Θ}k = {ϕ}Tj [C1 ]{φ }k − λk {ϕ}Tj [M]{φ }k − λ j {ϕ}Tj [M]{φ }k = a j δ jk ,
Eq. 2.140
{Ψ}Tj [A] {Θ}k = −(λ j + λk ){ϕ}Tj [M]{φ }k + {ϕ}Tj [C1 ]{φ }k = a j δ jk .
Eq. 2.141
As características das matrizes [A] e [B] conduzem à solução de autovalores complexos em pares
conjugados do tipo λ = δ ± iν. Para j 6= k, mas aplicando λk = λ j∗ na Equação 2.141, obtém-se
a j = 0, e se pode escrever
− λ j + λ j∗ {ϕ}Tj [M]{φ }k + {ϕ}Tj [C1 ]{φ }k = {0}.
Eq. 2.142
Capítulo 2 Revisão Teórica
45
Considerando-se as seguintes definições
{ϕ}Tj [M]{φ }k , m j
Eq. 2.143
{ϕ}Tj [C1 ]{φ }k , c1 j .
Eq. 2.144
e
Substituindo a Equação 2.143 e a Equação 2.144 na Equação 2.142 obtém-se a relação
2δ j =
c1 j
.
mj
Eq. 2.145
Para j 6= k, mas aplicando λk = λ j∗ na Equação 2.138 obtém-se b j = 0, logo
{ϕ}Tj [K]{φ }k − λ j λ j∗ {ϕ}T [M]{φ }k = {0},
Eq. 2.146
considerando a seguinte definição
{ϕ}Tj [K]{φ }k , k j ,
Eq. 2.147
e substituindo as definições na Equação 2.146 obtém-se a expressão
δ j2 + ν 2j =
kj
= Ω2j .
mj
Eq. 2.148
Definindo o fator de amortecimento ζ como usualmente
ζj =
c1 j
,
2m j Ω j
Eq. 2.149
e substituindo as definições (Equação 2.145) e (Equação 2.149) na Equação 2.148, aparece
naturalmente a expressão
λ j = ζ j Ω j + iΩ j
q
1 − ζ j2 .
Eq. 2.150
A proposta de solução {θ }est , com s = −λ , leva a resposta temporal “livre” do sistema no
espaço de estado para
2n
y(t) =
∑ C j θ j es jt ,
Eq. 2.151
j=1
sendo C j uma constante arbitrária a ser determinada a partir das condições iniciais do sistema
y(t = 0).
Como mencionado anteriormente, os autovalores λ j = δ j + iν j são complexos, ocorrendo
aos pares conjugados. Ao aplicar a j-ésima condição inicial referente ao j-ésimo autovetor, a
resposta livre do sistema terá a participação apenas de um modo, e a j-ésima resposta temporal
correspondente é composta de duas parcelas. Para obter uma resposta real (possível) no tempo,
Capítulo 2 Revisão Teórica
46
∗ (2.151):
as duas parcelas devem ser: θ j2 = θ j1
∗
(C j1 θ j es j t +C j2 θ j∗ e∗s j t ) = eδ j t (C j1 θ j eiν j t +C j2 θ j∗ e−iν j t ),
Eq. 2.152
mostra-se necessário para que resposta seja real que C j2 = C∗j1 , logo
eδ j t (C j θ j eiν j t +C∗j θ j∗ e−iν j t ) = eδ j t (C j θ j +C∗j θ j∗ ) cos ν j t + i(C j θ j −C∗j θ j∗ ) sin ν j t
Eq. 2.153
e
eδ j t (D j cos ν j t + E j sin ν j t) = e−ζ j Ω j t (D j cos Ω j
q
q
1 − ζ j2t + E j sin Ω j 1 − ζ j2t), Eq. 2.154
com D j = (C j θ j +C∗j θ j∗ ), E j = (C j θ j −C∗j θ j∗ ) números reais.
A senóide (D j cos Ω̄ j t + E j sin Ω̄ j t) é multiplicada pela envoltória exponencial e−ζ j Ω j t , cujo
sinal define a estabilidade do sistema. MEIROVITCH (1990) define a estabilidade de sistemas
lineares da seguinte forma:
1. A solução linearizada é estável se os autovalores s = −λ são não positivos, isto é, partes
reais nulas ou negativas.
2. A solução é assintoticamente estável se todos os autovalores s = −λ tem partes reais
negativas.
3. A solução é instável se houver ao menos um autovalor s = −λ com parte real positiva.
A estabilidade na região acima da rotação crítica, como mostrado por CHILDS (1993), depende do amortecimento introduzido no sistema e mostra que mecanismos de amortecimento
“interno” podem levar à instabilidade, enquanto que os introduzidos através dos mancais tendem a aumentar a região de estabilidade do sistema. O projeto proposto leva amortecimento ao
sistema através de um mancal.
2.7.2
Diagrama de Campbell
Considerando que em um sistema girante seus elementos estão expostos ao efeito giroscópico (8 ), as matrizes que dependem da velocidade (Ωrpm ) tem de ser atualizadas para cada rotação e o problema de autovalores deve resolvido novamente. Criando um gráfico das freqüências
características do sistema em função da rotação (Ω j × Ωrpm ) com Ω j representando as freqüências naturais para uma rotação Ωrpm , obtém-se o diagrama de Campbell.
8 Por
hipótese, é desconsiderada a variação dos parâmetros dos mancais com a rotação.
Capítulo 2 Revisão Teórica
47
Traçando uma linha a 45ř (Ω j = Ωrpm ) ou seja, rpm × rpm, os cruzamentos com as curvas
de freqüência característica indicam a localização de rotações críticas, freqüências características para excitação do tipo desbalanceamento. Normalmente o diagrama apresenta também
linhas rpm × 2rpm e rpm× ∼ 0.5rpm que representam excitações do tipo desalinhamento e
instabilidade do filme de óleo (“oil whirl” em máquinas com mancal de deslizamento).
Para excitações cujas freqüências são independentes da rotação as curvas do diagrama de
Campbell devem ser cortadas por uma reta vertical paralela ao eixo das ordenadas. As freqüências características resultantes destas intercessões são denominadas “freqüências naturais” LALANNE e FERRARIS (2001).
RotorDin − Campbell
LALANE p125
500
450
400
Freqüência (Hz)
350
300
250
200
150
100
50
0
0
5000
10000
15000
20000
25000
30000
Rotação (rpm)
Figura 2.22: Diagrama de Campbell.
A Figura 2.22 apresenta um diagrama de Campbell típico onde se observam as retas 1 ×
rpm, 2 × rpm, 0.5 × rpm, e uma reta vertical de freqüências naturais para Ωrpm = 15000 (rpm).
O Diagrama apresenta duas freqüências naturais. Além da rotação do eixo, ocorre a precessão,
que pode ocorrer em sentido contrário ao da rotação, quando é dito reverso ou “backward” ou no
mesmo sentido da rotação, quando a precessão é dita direta ou “forward”. A separação destas
retas ocorre devido ao efeito giroscópico. Estes fenômenos estão descritos em CHILDS (1993),
ERICH (2004), LALANNE e FERRARIS (2001) e VANCE (1988) entre outros.
2.7.3
Solução simplificada
Segundo ESPÍNDOLA e BAVASTRI (1997), para sistemas rotativos nos quais a excitação
é proveniente somente do desbalanceamento, é possível simplificar a solução do sistema de
Capítulo 2 Revisão Teórica
48
equações para achar diretamente as rotações críticas considerando apenas excitação de desbalanceamento (Equação 2.37). Neste caso Ωrpm = Ω então se pode escrever [G] = Ω[G1 ], já que
[G] neste caso, é função explícita de Ωrpm .
Na equação geral do movimento (Equação 2.37), considerando a excitação { f (t)} somente
devido ao desbalanceamento, tem-se a relação
[M]{q̈(t)} + (Ω[G1 ] + [C]){q̇(t)} + [K]{q(t)} = { f (t)},
Eq. 2.155
com [G1 ] representando a matriz giroscópica (Equação 2.43 e Equação 2.60) externando o produto pela rotação: [G] = Ω[G1 ].
A Equação 2.155 pode ser levada ao domínio da freqüência através da transformada de
Fourier, obtendo-se assim a expressão
−Ω2 [M] + iΩ2 [G1 ] + iΩ[C] + [K] {Q(Ω)} = {F(Ω)}.
Eq. 2.156
Considerando
(
{y(t)} =
{q(t)}
)
{q̇(t)}
Eq. 2.157
e substituindo a variável (Equação 2.157) na (Equação 2.155), obtém-se
h
i
h
i
[Ω[G1 ] + [C]] [M] {ẏ(t)} + [K] [0] {y(t)} = { f (t)}.
Eq. 2.158
Os termos em Ω2 podem ser agrupados convenientemente
−Ω2 [[M] − i[G1 ]] + iΩ[C] + [K] {Q(Ω)} = {F(Ω)}.
Eq. 2.159
Definindo a matriz [M̂] = [[M] − i[G1 ]] se obtém a equação
−Ω2 [M̂] + iΩ[C] + [K] {Q(Ω)} = {F(Ω)}.
Eq. 2.160
Como já demonstrado em 2.7, a equação do movimento pode ser representada no espaço de
estado através da relação




(
)
[C] [M̂]
[K] [0]
{ f (t)}

 {ẏ(t)} + 
 {y(t)} =
,
{0}
[M̂] [0]
[0] −[M̂]
sendo definidas as matrizes:

• [A] = 
[C] [M̂]
[M̂] [0]

.
Eq. 2.161
Capítulo 2 Revisão Teórica

• [B] = 
[K]
[0] −[M̂]
(
• {N(t)} =
[0]
{ f (t)}
{0}
49

.
)
.
O sistema matricial de equações diferenciais completo (Equação 2.161), pode ser reescrito de
forma compacta através da expressão
[A]{ẏ(t)} + [B]{y(t)} = {N(t)}.
Eq. 2.162
Utilizando as definições apresentadas na seção (2.7) o problema de autovalores generalizado
para n graus de liberdade e com j = 1, 2n fica definido através da relação
[B]{θ } j = λ j [A]{θ } j ,
Eq. 2.163
[B][θ ] = [λ ][A][θ ],
Eq. 2.164
[B]T [ψ] = [λ ][A]T [ψ].
Eq. 2.165
sendo o conjunto definido por
e o problema adjunto associado
As matrizes ortonormalizadas [Θ] e [Ψ] e a matriz [λ ] complexas, que definem o sistema
primário rotativo são obtidas utilizando-se o software “RotorDin” (ver seção J.1).
Importante notar que assumindo Ωrpm = Ω foi possível, resolvendo apenas uma única vez
um problema de autovalores, representar o sistema primário através de um conjunto de parâmetros modais a serem utilizados no projeto ótimo de neutralizadores. Isto é importante porque
o projeto ótimo de neutralizadores dinâmicos viscoelásticos que é apresentado no capítulo 4,
está relacionado com a excitação do tipo desbalanceamento. Portanto o modelo “simplificado”
é adequado para esta utilização, pois é preciso e computacionalmente leve.
2.7.4
Solução no Domínio da Freqüência
Partindo-se do sistema de equações diferenciais no espaço de estado (Equação 2.162), aplicando a transformada de Fourier se obtém
[iΩ[A] + [B]] {Y (Ω)} = {N(Ω)}.
Eq. 2.166
Capítulo 2 Revisão Teórica
50
Propondo a seguinte transformação de variáveis no domínio do tempo
{y(t)} = [θ ]{p(t)},
Eq. 2.167
e a levando para o domínio da freqüência através da transformada de Fourier se obtém
{Y (Ω)} = [Θ]{P(Ω)}.
Eq. 2.168
Substituindo a transformação proposta (Equação 2.167) e pré-multiplicando pelos autovetores
[Ψ]T já ortonormalizados, obtém-se no espaço de estado
iΩ[Ψ]T [A][Θ] + [Ψ]T [B][Θ] {P(Ω)} = [Ψ]T {N(Ω)} .
Eq. 2.169
As matrizes abaixo, como mostrado anteriormente (em 2.7.1), são diagonais e se reduzem a:
• [Ψ]T [A][Θ] = [r 1r ] = [I].
• [Ψ]T [B][Θ] =
r
λ jr .
Aplicando as relações acima em (2.169), obtém-se
[D(Ω)]{P(Ω)} = [Ψ]T {N(Ω)} ,
Eq. 2.170
[D(Ω)] = iΩ[I] + r λ jr ,
Eq. 2.171
com
assim, a resposta no espaço modal do espaço de estado resulta
{P(Ω)} = [D(Ω)]−1 [Ψ]T {N(Ω)}
Como pode ser constatado a matriz [D(Ω)] é diagonal, portanto, facilmente invertida. Cada
elemento da diagonal de [D(Ω)]−1 é obtido através da inversão simples do elemento da diagonal
de [D(Ω)].
−1
[D(Ω)]
=
r
1
)r .
(
iΩ + λ j
Eq. 2.172
* Caso a matriz [D(Ω)] seja uma matriz mal condicionada (9 ) sua inversão pode ser obtida
através de algoritmos de pseudo-inversão, como mostrado por BAZáN e BAVASTRI (1996).
avaliação do condicionamento da matriz [A] pode ser o valor da relação cond(A) = kAk A−1 , se
cond(A) 1, a matriz [A] é dita mal condicionada.
9 Uma
Capítulo 2 Revisão Teórica
51
Da transformação de variáveis (Equação 2.167 ) e de estado (Equação 2.157), tem-se
(
)
(
)
{F(Ω)}
{Q(Ω)}
−1
T
= [Θ][D(Ω)] [Ψ]
.
Eq. 2.173
iΩ{Q(Ω)}
{0}
Por definição (Equação 2.100), a função de transferência é a relação entre deslocamento e força
generalizados (saída / entrada), assim
[α(Ω)] , [Θ][D(Ω)]−1 [Ψ]T .
Eq. 2.174
A matriz complexa [α(Ω)] (Equação 2.174), é denominada matriz de receptância no espaço
de estado. Observando a forma da Equação 2.173, deduz-se pelas ordens das matrizes que
[α(Ω)] tem a forma

[α(Ω)] = 
[α11 (Ω)] [α12 (Ω)]
[α21 (Ω)] [α22 (Ω)]

.
Eq. 2.175
Logo, a relação entre deslocamento e excitação, que define a matriz de receptância no espaço
de configuração é
{Q(Ω)} = [α11 (Ω)]{F(Ω)}.
Eq. 2.176
Para se obter uma função resposta em freqüência, utiliza-se a relação
2n
αks (Ω) =
∑
Θk j Ψ s j
.
iΩ + λ j
Eq. 2.177
j=1
Com s representando o ponto de excitação, k o ponto de medição, e n o número de modos.
O sistema que precisa ser resolvido é dado pela Equação 2.173, sendo que o vetor {F(Ω)}
atualizado para cada rotação, é composto pela força (Equação 2.74) aplicada às coordenadas generalizadas correspondentes à localização da seção onde a excitação do tipo desbalanceamento
é aplicada.
A solução da Equação 2.176 leva à solução no domínio da freqüência da resposta {Q(Ω)}
como amplitude da resposta ao desbalanceamento.
Para efeito de comparação, são apresentadas figuras de resposta. A Figura 2.23 apresenta o
gráfico de resposta ao desbalanceamento obtida através do software RotorDin para a geometria
e demais condições do sistema, apresentada em LALANNE e FERRARIS (2001)(p.125). A
Figura 2.24 apresenta a resposta obtida pelo autor citado.
Capítulo 2 Revisão Teórica
52
RotorDin − Resposta ao Desbalanceamento
Amplitude − Amplitude − LALANE p125
−2
10
coord.(z) y= 0.5 (m)
−3
10
−4
Ampl. (m)
10
−5
10
−6
10
−7
10
−8
10
0
5000
10000
15000
20000
25000
30000
Rotacao (rpm)
Figura 2.23: Resposta ao desbalanceamento, RotorDin.
Figura 2.24: Resposta ao desbalanceamento, LALANNE e FERRARIS (2001).
Uma vez montadas as matrizes [A], [B] e [λ ] adequadamente, seja para o modelo geral (Equação 2.123), ou para o “modelo simplificado” (Equação 2.155) e obtida a matriz receptância
[α(Ω)], é possível calcular no domínio da freqüência a resposta {Q(Ω)}, para um dado vetor de
força {F(Ω)}. Por exemplo, uma força de resposta plana corresponde a uma constante unitária
à ser introduzida na posição correspondente à coordenada excitada.
53
3
MODELO DO SISTEMA COMPOSTO
Neste capítulo é apresentada a equação do movimento do sistema composto, rotor mais
neutralizadores dinâmicos viscoelásticos. Uma vez definidas estas equações, é utilizada uma
formulação equivalente à utilizada no trabalho desenvolvido por BAVASTRI (1997)(cap. IV)
para o projeto ótimo de neutralizadores viscoelásticos.
São acrescentados fisicamente ao sistema primário p neutralizadores, fixados sobre um
mancal flexível propriamente projetado para este fim. Os dispositivos de controle atuarão somente em deslocamento, nas coordenadas X e Z, desprezando por hipótese os possíveis efeitos
nos graus de liberdade de rotação.
Nas coordenadas da seção correspondente do sistema à aplicação do neutralizador, agregamse as p massas e amortecimentos equivalentes generalizados, utilizando os modelos da Equação
2.121 e Equação 2.120, com massas me (Ω) e amortecimentos ce (Ω) respectivamente. As matrizes de massa e amortecimento equivalentes resultantes, [Me ] e [Ce ] assim obtidas, contém
elementos somente nas coordenadas generalizadas de deslocamento. Por hipótese, eventuais
momentos são também desprezados. Este tipo de acréscimo somente é possível porque o modelo equivalente não acrescenta grau extra de liberdade, como mostrado na seção (2.6).
Também são apresentados resultados obtidos através de uma versão especialmente modificada do software RotorDin, na qual são introduzidos os parâmetros equivalentes generalizados
diretamente nas matrizes de massa e amortecimento do modelo do eixo, no ponto de fixação
dos neutralizadores, com o objetivo de comparar os resultados com aqueles obtidos no software
de projeto.
3.1
Implementação
Para determinar a resposta do sistema composto, monta-se novamente o sistema de equações utilizando as variáveis de estado como definidas na Equação 2.157. Propondo uma solução
no domínio da freqüência do tipo
{Y (Ω)} = [Θ]{P(Ω)},
Eq. 3.1
Capítulo 3 Modelo do Sistema Composto
54
a equação de movimento no espaço modal do espaço de estado do sistema primário é facilmente
montada e calculada graças a utilização dos parâmetros equivalentes generalizados. Matematicamente, a equação do sistema composto fica definida pela relação
iΩ[Ψ]T [Ā][Θ] + [Ψ]T [B̄][Θ] {P(Ω)} = [Ψ]T {N(Ω)} .
Eq. 3.2
Supondo que o sistema primário “sem neutralizadores” seja definido pela equação
iΩ[Ψ]T [A][Θ] + [Ψ]T [B][Θ] {P(Ω)} = [Ψ]T {N(Ω)} .
Eq. 3.3
Por definição, tem-se que
[Ā] = [A] + [Ã],
Eq. 3.4
[B̄] = [B] + [B̃]
Eq. 3.5
com:

[C] [M̂]
[A] = 

[M̂] [0]

,
Eq. 3.6

[Ceq (Ω)] [Meq (Ω)]
[Ã] , 
[Meq (Ω)]
e

[B] = 

[B̃] , 
[K]
[0]
[0] −[M̂]
[0]
,
[0]
Eq. 3.7

,
[0]
[0] −[Meq (Ω)]
Eq. 3.8

.
Eq. 3.9
Assim, substituindo Equação 3.6 e Equação 3.7 na Equação 3.8 e também a Equação 3.9 e
Equação 3.4 na Equação 3.5, obtém-se a seguinte relação
iΩ[Ψ]T [Ā][Θ] + [Ψ]T [B̄][Θ] = iΩ [I] + [Ψ]T [Ã][Θ] + r λ jr + [Ψ]T [B̃][Θ] . Eq. 3.10
Para um sistema com grande número de graus de liberdade é conveniente truncar as matrizes modais [Ψ] e [Θ] para um número de autovetores 2n̂ que represente o sistema na faixa de
freqüências de interesse.
[Ψ̂] = [Ψ]2n×2n̂
Eq. 3.11
[Θ̂] = [Θ]2n×2n̂
Eq. 3.12
e
Capítulo 3 Modelo do Sistema Composto
55
Os autovetores truncados formam um sub-espaço modal a ser utilizado para determinar uma
solução aproximada do sistema composto. O novo sistema, agora truncado para os 2n̂ primeiros
autovetores, no sub-espaço modal do espaço de estado, fica
[D̂]{P̂(Ω)} = [Ψ̂]T {N(Ω)} ,
Eq. 3.13
com
h
ii
h ˆ T
T
r
ˆ
[D̂]2n̂×2n̂ = iΩ [I] + [Ψ̂] [Ã][Θ̂] +
λ jr + [Ψ̂] [B̃][Θ̂] .
Eq. 3.14
A resposta no espaço modal do espaço de estado é definida através da equação
(
)
{F(Ω)}
T
,
{P̂(Ω)}2n̂×1 = [D̂(Ω)]−1
2n̂×2n̂ [Ψ̂]2n̂×2n
{0}
Eq. 3.15
2n×1
denominado-se Pj (Ω) a j-ésima coordenada principal generalizada.
Partindo da definição de variável de estado (Equação 2.157) e da transformação de variáveis
(Equação 3.1), obtém-se a relação
(
)
(
)
{F(Ω)}
{Q(Ω)}
= [Θ̂][Ď(Ω)][Ψ̂]T 2n×2n
{0}
iΩ{Q(Ω)}
,
Eq. 3.16
2n×1
2n×1
com [Ď(Ω)] = [D̂(Ω)]−1 . Desenvolvendo o produto, já que [Ď(Ω)] não é diagonal obtém-se
αks (Ω) =
2n̂
2n̂
∑
∑
Ďi j (Ω)Θ̂ki Ψ̂s j .
Eq. 3.17
i=1 j=1
Note-se que a resposta não é completa, já que se trata de um sub-espaço modal de estado 2n̂.
3.2
Resposta em Freqüência
É possível obter a resposta do sistema para rotação nula Ω = 0, ou para qualquer outra
rotação constante, obtendo-se os parâmetros modais Ψ, Θ e λ correspondentes, a partir do “caso
geral” (ver seção 2.7) diagrama de Campbell completo. Na Figura 2.22 a partir de uma reta
paralela ao eixo das ordenadas na rotação desejada Ωrpm = cte.
A função resposta em freqüência no espaço de estado é obtida diretamente através da Equação 2.173 ou Equação 2.177 para o sistema primário, e da Equação 3.16 ou Equação 3.17 para
o sistema composto.
Esta implementação permite comparar as respostas numéricas obtidas com as medidas atra-
Capítulo 3 Modelo do Sistema Composto
56
vés do teste de resposta (“bump test”) com martelo instrumentado (ver seção 5.1).
3.3
Inclusão dos Neutralizadores Diretamente no Modelo do
Eixo
Utilizando a versão modificada do software RotorDin é possível entrar com os características dinâmicas do material viscoelástico, a massa dos neutralizadores e suas coordenadas de
aplicação.
Como foi visto na descrição do modelo do sistema girante (seção 2.7), as matrizes globais
de massa rigidez e amortecimento são montadas e a partir destes, com a equação de movimento,
as características dinâmicas do rotor podem ser obtidas (diagrama de Campbell, resposta ao
desbalanceamento entre outras).
Através dos parâmetros equivalentes generalizados, é possível modificar diretamente as
matrizes de massa e amortecimento, introduzindo o modelo de parâmetro equivalente nas coordenadas generalizadas correspondentes à fixação dos neutralizadores.
Através da utilização dos parâmetros equivalentes generalizados não é acrescentado nenhum novo grau de liberdade, então, para cada freqüência Ω ao qual o sistema está sujeito, no
caso em que Ωrpm = Ω, é necessário o cálculo da massa me (Ω) e do amortecimento equivalentes ce (Ω). Estas grandezas devem ser somadas na posição correta às matrizes globais de massa
[M] e amortecimento [C] do sistema girante primário.

0
...
0

 ... . . .
...
...


 0 ... me1 (Ω)
0


[M̃(Ω)] = [M] +  ... ...
...
...

 0 0
...
me2 (Ω)


..
 ... ...
.
...

0 0
...
0
0
0


... 


0 


...  = [M] + [Me (Ω)]

0 


... 

0
Eq. 3.18
Capítulo 3 Modelo do Sistema Composto
57
e

0
...
0

.
 ... . .
...
...


 0 ... ce1 (Ω)
0


[C̃(Ω)] = [C] +  ... ...
...
...

 0 0
...
ce2 (Ω)


..
 ... ...
.
...

0 0
...
0
0
0


... 


0 


...  = [C] + [Ce (Ω)].

0 


... 

0
Eq. 3.19
Todo o sistema deve ser resolvido novamente para cada Ω (freqüência), utilizando porém,
as novas matrizes [M̃(Ω)] e [C̃(Ω)] aplicadas à Equação 2.123.
3.3.1
Diagrama de Campbell
O diagrama de Campbell utilizando as novas matrizes de massa [M̃(Ω)] e amortecimento
[C̃(Ω)] necessita ser obtido de maneira ligeiramente diferente da apresentada anteriormente
pelo fato de se ter dois sistemas cujos parâmetros variam. O sistema girante que depende da
rotação Ω = Ωrpm e as propriedades do material viscoelástico que dependem da freqüência Ω e
da temperatura T (considerada constante neste trabalho).
Para uma rotação qualquer constante Ωrpm , é resolvido em uma faixa de freqüência apropriada ao número rotações características que se deseja analisar, uma seqüência de problemas
de autovalor montados com as matrizes (Equação 3.18) e (Equação 3.19) na equação geral
(Equação 2.123), isto é necessário devido a presença do material viscoelástico cujas características dependem da freqüência Ω. Através deste conjunto de soluções (freqüências naturais),
é montado um diagrama de Campbell dito interno, ou seja, válido somente para rotação atual
Ωrpm .
Os pontos de intercessão das curvas de freqüência natural com uma linha de inclinação unitária (Ω j = Ω) são consideradas freqüências naturais do sistema composto para a dada rotação.
Depois de percorrida a faixa de rotação desejada, obtém-de o diagrama de Campbell propriamente dito, de maneira similar ao apresentado por FERREIRA (2005). Esta metodologia é
extremamente custosa em termos computacionais.
Para se conseguir uma boa precisão na obtenção do Campbell interno é necessário diminuir
o passo de freqüência na região de influência do neutralizador, onde cálculo dos pontos de
intercessão é mais sensível. A região considerada foi de −40% a +20% em torno da freqüência
natural do neutralizador.
Capítulo 3 Modelo do Sistema Composto
58
Campbell interno
Ωj
loop rotacão Ωrpm

loop f reqüência Ω

 atualiza [M̃(Ω)] e [C̃(Ω)]


 calcula [B][θ ] = [λ ][A][θ ]

Campbell interno
Campbell final
Ωj
RotorDin − Campbell Interno
LAVIB 1D N2 rpm:10.00
RotorDin − Campbell
LAVIB 1D N2
500
350
450
300
400
250
350
300
Freqüência (Hz)











Freq. (Hz)

250
200
150
150
100
100
50
50
DiagramaCampbell
200
0
0
100
200
300
400
500
600
0
0
2000
Freq. (Hz)
4000
6000
8000
10000
12000
Rotação (rpm)
Ωrpm
Ωrpm = 6000 (rpm)
No diagrama interno de Campbell, percebe-se a influência dos neutralizadores no sistema
composto, na região próxima da freqüência natural dos mesmos.
Gráficos
Os gráficos das simulações foram obtidos para o modelo experimental de dinâmica de rotores do LAVIB conforme geometria apresentada nas figuras (Figura 4.3 e Figura 4.5). seção
4.5, seção 4.9, utilizando os parâmetros do material neoprene (ver Tabela 4.1) na temperatura
de trabalho de 298 (K) com neutralizadores de 0.3 (kg). Estes parâmetros foram obtidos pelo
grupo PISA/CNPq no LVA-PISA da UFSC conforme LOPES et al. (2004).
O gráfico da Figura 3.1 apresenta o diagrama de Campbell interno do sistema para uma
rotação Ωrpm = 10 (rpm).
RotorDin − Campbell Interno
LAVIB 1D N2 rpm:10.00
500
450
400
350
Freq. (Hz)
300
250
200
150
100
50
0
0
100
200
300
400
500
600
Freq. (Hz)
Figura 3.1: Diagrama de Campbell interno.
O detalhe da freqüência natural (Figura 3.2) mostra a obtenção das freqüências naturais
14000
Capítulo 3 Modelo do Sistema Composto
59
introduzidas pelos neutralizadores que é muito sensível aos parâmetros do neutralizador e do
material viscoelástico.
RotorDin − Campbell Interno
LAVIB 1D N2 rpm:10.00
85
80
75
Freq. (Hz)
70
65
60
55
50
45
40
30
40
50
60
70
80
90
100
Freq. (Hz)
Figura 3.2: Detalhe do diagrama de Campbell interno.
Com o sistema a 12000 (rpm) observa-se no detalhe (Figura 3.3) o efeito giroscópico separando as freqüências para a precessão “backward” e “forward”.
RotorDin − Campbell Interno
LAVIB 1D N2 rpm:12000.00
90
85
80
75
Freq. (Hz)
70
65
60
55
50
45
40
30
40
50
60
70
80
90
100
Freq. (Hz)
Figura 3.3: Detalhe do diagrama de Campbell interno a 12500 (rpm).
Depois de obter as freqüências naturais para cada rotação, obtém-se o diagrama final.
Capítulo 3 Modelo do Sistema Composto
60
RotorDin − Campbell
NDV − LAVIB 1D N2
300
250
Freq. (Hz)
200
150
100
50
0
0
2000
4000
6000
8000
10000
12000
Rotação (rpm)
Figura 3.4: Campbell final com neutralizadores.
A Figura 3.4 apresenta o diagrama de Campbell do sistema composto, onde no eixo das
ordenadas para rotação nula, verificam-se as freqüências naturais decorrentes do efeito da aplicação dos neutralizadores e sobre a reta 1 × rpm, as rotações críticas.
Resposta ao Desbalanceamento
A resposta ao desbalanceamento pode ser obtida de ma-
neira análoga à mostrada na seção 2.2.4, Equação 3.15 ou Equação 3.16, porém como os parâmetros equivalentes dependem da freqüência (neste caso igual a rotação) é necessário atualizar
as matrizes de massa e amortecimento e resolver o problema de autovalores para cada freqüência, atualizando também a matriz giroscópica.
Resposta em Freqüência
Neste caso utiliza-se novamente a formulação para o modelo
geral, fixando a matriz giroscópica para uma rotação Ωrpm qualquer, no caso do sistema em repouso, Ωrpm = 0, mas atualizando as matrizes de massa e amortecimento na faixa de freqüência
desejada devido a presença dos parâmetros equivalentes generalizados.
No exemplo da Figura 3.5, a excitação do sistema foi do tipo resposta plana, foi aplicada a
0.6 (m) da extremidade e a resposta obtida a 0.15 (m).
Capítulo 3 Modelo do Sistema Composto
61
RotorDin − Resposta em Frequencia
Inertancia − Amplitude − LAVIB 1D N2
60
coord.(z) y= 0.15 (m)
40
Iner. dB ref = 1(m/Ns²)
20
0
−20
−40
−60
−80
0
50
100
150
200
250
300
350
400
Freq. (Hz)
Figura 3.5: Resposta com neutralizadores de neoprene e Ωrpm = 0.
Observa-se a correspondência entre o diagrama de Campbell (Figura 3.4, linha vertical na
rotação zero) e a resposta em freqüência .
Para poder fazer uma comparação com o método proposto “simplificado”, a ser apresentado
no capítulo 4 (Figura 4.19), o material viscoelástico foi trocado para borracha butílica, e o sistema mantido nas condições anteriores. O gráfico obtido é mostrado na Figura 3.6 e representa
uma função resposta em freqüência com o sistema em repouso, excitado com delta de Dirac a
0.6 (m) e a resposta obtida a 0.15 (m).
RotorDin − Resposta em Frequencia
Inertancia − Amplitude − LAVIB 1D N2
60
coord.(z) y= 0.15 (m)
40
Iner. dB ref = 1(m/Ns²)
20
0
−20
−40
−60
−80
0
50
100
150
200
250
300
350
400
Freq. (Hz)
Figura 3.6: Resposta com neutralizadores de borracha butílica.
Observação: A temperatura de trabalho do material foi considerada constante a 298 (K).
62
4
PROJETO DO NEUTRALIZADOR PARA DESBALANCEAMENTO
Neste capítulo são definidas:
• A massa dos neutralizadores.
• Uma função objetivo para o projeto ótimo de tais dispositivos de controle.
• Através da formulação proposta é apresentada uma implementação de um código numérico para o projeto ótimo de neutralizadores viscoelásticos aplicados a sistemas rotativos.
• São apresentados os resultados de algumas simulações numéricas realizadas sobre um
modelo da bancada de dinâmica de rotores do LAVIB.
• A partir das freqüências ótimas obtidas é mostrado um projeto de construção de neutralizadores para controle da primeira rotação crítica do rotor da bancada.
Neste projeto são empregados os parâmetros modais obtidos através do software RotorDin para
o modelo “simplificado”, ou seja, considerando apenas excitação do tipo desbalanceamento.
4.1
Massa do Neutralizador
Como mostrado por BAVASTRI (1997) (ver seção A.1.1), por semelhança às expressões
obtidas para o modelo de um grau de liberdade é possível calcular a massa ótima do neutralizador de vibrações utilizando uma consideração modal.
Para um sistema rotativo simples, a massa dos neutralizadores é estimada através de uma
aproximação. Esta aproximação considera o modelo do rotor em repouso (Ωrpm = 0), compatível com o de uma viga simplesmente apoiada, e através de uma equivalência, é estimada a
massa do neutralizador. Neste trabalho, a massa dos neutralizadores projetados para o primeiro
modo de flexão, é de aproximadamente de 5 (%) da massa do sistema primário ou seja, cerca
de 10 (%) da sua massa equivalente, como sugerido por HARTOG (1956).
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
4.2
63
Função Objetivo
Neste trabalho, propõe-se empregar uma técnica de otimização não linear, através de uma
função objetivo f {x}n×1 de Rn → R, sendo {x} o vetor projeto (seção 4.3 e B.1) que é composto
pela freqüência natural dos neutralizadores.
A definição da função objetivo utilizando as coordenadas principais P(Ω) (Equação 3.15),
é feita basicamente pelas seguintes razões:
1. As operações matriciais tem ordens de grandeza reduzidas, 2n̂, importante para sistemas
muito discretizados. O número de operações com ponto flutuante (FLOP) envolvidas
no produto de matrizes é aproximadamente proporcional ao cubo da ordem da matriz
(NISHTALA et al. (2004)).
2. O resultado em termos de coordenadas principais não depende da coordenada de avaliação da resposta, mas da excitação. Teoricamente se pode otimizar o sistema através da
norma 2 da matriz resposta em freqüência [α(Ω)] de ordem elevada ou através de uma
função resposta, onde ainda é preciso definir uma coordenada de avaliação k.
Para uma otimização de Rn → R sobre as coordenadas principais reduzidas P̂(Ω), pode-se
obter um vetor de máximos absolutos para cada par de coordenadas principais correspondente
ao modo que se deseja minimizar, na faixa de freqüências estabelecida, retornando então ao
otimizador a norma 2 deste vetor (1 ).
x
P̂(Ω, x) , xT = [Ωa1 , ..., Ωan ] .
fob j (x) = Ωm1a<Ω<Ω
2
Eq. 4.1
onde Ωa1 , ..., Ωan são as freqüências ótimas dos neutralizadores.
A faixa de freqüências Ω1 < Ω < Ω2 , deve ser estabelecida de tal forma que os modos que
se desejam controlar estejam contemplados. Tais modos podem ser obtidos do programa de
dinâmica de rotores, já que os modos de flexão são inerentes às rotações críticas.
A excitação pode ser constante para uma dada rotação ou quadrática para resposta ao desbalanceamento, devendo ser aplicada distante dos nós dos modos que se deseja controlar (ver
anexo G.1) uma vez que a otimização se dá sobre a resposta.
Uma análise sobre a composição da solução da resposta em freqüência, mostra que as coordenadas principais assim como os autovalores, ocorrem aos pares conjugados, conforme a
1 Norma
2 ou euclidiana de um vetor é definida como:kT k ,
q
2
∑i |Ti |
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
64
expressão
[Θ̂]2n×2n̂ {P̂(Ω)}2n̂×1 = {Q(Ω)}2n×1 ,



 Q1




Qi




 Q 
2n
2n×1
=




(Θ11 P1 + Θ12 P2 ) + ... + (Θ12n̂−1 P2n̂−1 + Θ12n̂ P2n̂ )
Eq. 4.2




(Θi1 P1 + Θi2 P2 ) + ... + (Θi2n̂−1 P2n̂−1 + Θi2n̂ P2n̂ )





 (Θ P + Θ P ) + ... + (Θ
2n1 1
2n2 2
2n2n̂−1 P2n̂−1 + Θ2n2n̂ P2n̂ )
, Eq. 4.3
2n×1
{Qi } = {Σnj=1 Li j }
Eq. 4.4
Li j = Θik Pk + Θi(k+1) P(k+1) , k = 2( j − 1) + 1.
Eq. 4.5
com
Sendo assim, na i-ésima coordenada, para cada j-ésima freqüência a ser controlada, deve-se
considerar a contribuição da soma do par correspondente Li j .
A otimização proposta é aplicada sobre um modelo do sistema primário baseado em parâmetros modais. Assim genericamente, pode-se aplicar a formulação tanto para minimizar uma
resposta ao desbalanceamento quanto para uma rotação constante qualquer. O resultado do
projeto ótimo depende dos parâmetros modais correspondentes à condição do sistema primário.
4.3
Implementação
R 2 ). O código
Toda a implementação foi desenvolvida utilizando-se o software livre Scilab(
R acrescentando funcionalidades
do RotorDin foi portado do originalmente escrito em Matlab,
de exportação dos parâmetros modais para o modelo rotação constante e resposta ao desbalanceamento, inclusão de parâmetros equivalentes generalizados modelo viscoelástico nas matrizes
de massa e amortecimento, entre outras.
O código para projeto ótimo dos neutralizadores dinâmicos é totalmente novo e foi desenvolvido para o modelo no espaço de estado “simplificado” e parâmetros equivalentes generalizados para material viscoelástico. Como entrada são necessários os parâmetros modais, parâmetros do modelo de derivada fracionária do material viscoelástico utilizado no neutralizador
(Equação 2.120 e Equação 2.121), massa dos neutralizadores e configurações como:
• faixa de freqüência,
2 Scilab
Home Page http://www.scilab.org
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
65
• posição dos neutralizadores,
• modos de interesse,
• demais parâmetros de controle.
Como resultado é obtido um gráfico comparativo da resposta ao desbalanceamento (ou em
freqüência) com e sem neutralizadores, e a freqüência ótima dos mesmos.
O algoritmo de otimização não linear escolhido é o BFGS (Broyden - Fletcher - Goldfarb Shanno), da família dos métodos de quase-newton (ver seçao B.1), como descrito em ARORA
(2004). As restrições de desigualdade foram implementadas como funções de barreira para que
as freqüências ótimas se situem em uma faixa pré-determinada.
Figura 4.1: Diagrama de otimização.
A Figura 4.1 apresenta um diagrama simplificado da implementação da metodologia de
projeto proposta.
4.4
Detalhamento do processo de projeto
O primeiro passo é definir todos os parâmetros de projeto. Isto deve ser feito baseado no
resultado de um cálculo da resposta ao desbalanceamento realizado previamente, utilizando o
RotorDin, por exemplo. Sabendo-se a localização das rotações críticas definem-se:
• Os modos a serem controlados.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
66
• O número de neutralizadores utilizados, as massas, e posições dos mesmos.
• O número de autovalores do sub-espaço modal.
• A posição e tipo de excitação, se delta de Dirac ou quadrático (desbalanceamento).
• Os limites das barreiras de freqüência para cada modo desejado.
• As freqüências iniciais para cada neutralizador.
Definidos todos os parâmetros, a simulação numérica segue os seguintes passos:
1. A partir da versão modificada do aplicativo RotorDin, são exportados os parâmetros modais do sistema primário girante na forma das matrizes complexas Ψ, Θ e λ (ver seção
2.7.1) obtidas para o modelo “simplificado”, já que neste trabalho o projeto ótimo dos
neutralizadores considera a excitação de desbalanceamento apenas.
2. As matrizes são então truncadas para o número de autovetores desejado (2n̂), sendo este
parâmetro escolhido pelo usuário tendo em vista a faixa de freqüência de atuação do
controle.
3. Dentro de uma estrutura de repetição que varia a freqüência angular (Ω) dentro da faixa
de interesse, são calculadas ponto a ponto, o módulo de cisalhamento dinâmico Gr (Ω) e
o fator de perda η(Ω) do material viscoelástico.
4. Os parâmetros equivalentes me (Ω) e ce (Ω) tem o valor numérico calculado.
5. As matrizes [Ã] e [B̃] são obtidas.
6. Através da equação (3.14) é obtida a matriz [D].
7. Monta-se também o vetor de excitação {F(Ω)} constante ou quadrático {F(Ω2 )}.
ˆ T truncada e transposta, e
8. Multiplica-se a matriz inversa [D]−1 pela de autovetores [Ψ]
pelo vetor de excitação {F(Ω)}. Obtendo-se o assim o vetor de coordenadas principais
reduzidas {P̂(Ω)}.
9. Para uma freqüência Ω em questão, são verificadas as amplitudes do valor absoluto referente às freqüências que se desejam controlar. Caso os valores atuais sejam superiores
aos anteriores, estes são atualizados.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
67
10. Ao final da análise sobre a faixa de freqüência de controle (estrutura de repetição), na
freqüência superior da faixa definida, estarão armazenados os valores máximos do módulo das coordenadas principais escolhidas para o controle. O valor final da função objetivo é obtido através da aplicação da norma 2 sobre o vetor de máximos.
11. Em termos de otimização, são adicionadas eventuais penalidades pela violação das barreiras guia das freqüências ótimas dos neutralizadores.
12. O loop de otimização termina quando a função objetivo atinge seu menor valor, então o
vetor projeto conterá as freqüências ótimas dos neutralizadores.
Figura 4.2: Diagrama detalhado de otimização.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
68
Após concluída a otimização das freqüências dos neutralizadores (Figura 4.2) o software de
projeto apresenta a resposta ao desbalanceamento com e sem neutralizadores. Eventualmente,
de posse dos parâmetros modais do sistema primário e da freqüência ótima do neutralizador
é possível obter somente o gráfico da resposta (com e sem neutralizadores fixados ao sistema
primário).
De posse das freqüências ótimas e definida a massa do neutralizador, é possível determinar
o parâmetro L utilizando-se a relação da Equação 2.104. Geralmente a espessura h do material
é definida. Através da Equação 2.105 obtém-se a área de aplicação do carregamento A, estando
assim o projeto do neutralizador definido. No caso se se trabalhar com o material vulcanizado,
é necessário definir uma espessura dentro dos limites do processo de fabricação, ficando assim
a espessura h a ser calculada.
4.5
Sistema Primário - Identificação
Os neutralizadores são projetados para o sistema primário montado no LAVIB (Figura 4.3),
destinado ao estudo da dinâmica de rotores, também denominado “rotor rig”. A montagem
é acionada por um motor alimentado através de inversor de freqüência, que permite variar a
velocidade numa ampla faixa, tipicamente de 0-7000 (rpm). O sistema primário é montado de
forma que sua primeira rotação crítica se situe próxima aos 3500 (rpm). Esta configuração foi
obtida simulando a geometria no RotorDin e variando o número de discos, suas posições e a
posição dos mancais.
Figura 4.3: Bancada de teste do LAVIB.
A configuração projetada conforme o esquema da Figura 4.3 e Figura 4.5, é composta por:
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
69
1. Eixo:
• Um eixo cilíndrico de aço com 25 (mm) de diâmetro e 1000 (mm) de comprimento
(9). Para o material do eixo é adotado para o módulo de Young o valor de E=207 (GPa)
e ν =0.3 para o coeficiente de Poisson. Sua massa específica é considerada 7850 (kg/m3 ).
2. Discos:
• Um acoplamento elástico que é considerado como um disco de ferro fundido com
80 (mm) de diâmetro e 30 (mm) de espessura (2), localizado na origem do sistema
de coordenadas. A massa específica é considerada 7850 (kg/m3 ). Este elemento
desacopla o motor elétrico do sistema mecânico.
• Um disco de alumínio de 240 (mm) de diâmetro e 9 (mm) de espessura (7), situado a 475 (mm) da extremidade do lado acionado. A massa específica é adotada
2700 (kg/m3 ).
• Os dois discos, um de proteção e outro para medição são construídos em aço e tem
diâmetro de 80 (mm) e 10 (mm) de espessura (4 e 5), estando localizados a 250 e
325 (mm) respectivamente. Massa específica considerada 7850 (kg/m3 ).
3. Mancais:
• Dois mancais de rolamento com rigidez elevada (ver anexo C.1) kxx = kzz ∼
= 1E9 (N/m)
(3 e 8), localizados a 80 e 730 (mm) apóiam o eixo.
• O mancal é composto por um rolamento auto-compensador de esferas, inserido em
um suporte de ferro fundido nodular quadrilateral, ao qual são fixados os neutralizadores. O suporte por sua vez, é mantido estático por meio de molas helicoidais
ligadas a uma estrutura metálica estática, localizado a 380 (mm) da extremidade.
As propriedades de rigidez e amortecimento das molas foram determinadas experimentalmente (ver anexo H.1) e valem respectivamente kxx = kzz = 55.1E3 (N/m)
e cxx = czz = 18 (Ns/m). A massa do falso mancal de 0.83 (kg) é adicionada pelo
aplicativo RotorDin, na matriz de massa [M] (Equação 2.53) na coordenada correspondente do sistema primário.
Para o amortecimento dos mancais de rolamento foi adotado uma amortecimento viscoso c de
aproximadamente cxx = czz = 5 (Ns/m), o que corresponde a um fator de aproximadamente
0.25%, ζ = 0.0025 do valor equivalente para um sistema de um grau de liberdade c = 2mΩn ζ .
O valor é considerado plausível e foi atribuído com base em medições preliminares. Este valor é
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
70
adotado para todas as simulações. São considerados nulos por hipótese, os termos cruzados das
matrizes de coeficientes de rigidez e amortecimento dos mancais kzx = kxz = 0 e czx = cxz = 0.
A massa total do sistema girante foi calculada em 7.1 (kg)
Como o amortecimento introduzido pelos mancais rígidos de rolamento é muito pequeno, a
amplificação é muito grande e o sistema pode ficar instável como mostrado por CHILDS (1993),
sendo perigoso ultrapassar a rotação crítica sem algum dispositivo de proteção ou controle. No
caso do rotor em questão, foram utilizados batentes de material anti-fricção (latão) para limitar o
deslocamento máximo (Figura 4.4), já que em algumas ocasiões, sem o dispositivo, o eixo ficou
irremediavelmente danificado (deformado plasticamente) após uma tentativa de ultrapassar a
rotação crítica.
Figura 4.4: Detalhe do dispositivo de segurança.
Para o modelo numérico, como mostrado na seção 2.2.2, os elementos de viga de seção
circular possuem dois nós, com quatro graus de liberdade, dois de translação e dois de rotação,
modelo Euler-Bernoulli corrigido para levar em conta o efeito do cisalhamento transversal (LALANNE e FERRARIS (2001)). Foi utilizada uma malha de 25 nós, conforme mostra a Figura
4.5 e que proporciona precisão adequada para os primeiros modos, como discutido em LALANNE e FERRARIS (2001). Os discos são considerados por hipótese, solidamente fixados
sobre o eixo e não deformáveis. Observe-se que a presença das tiras de borracha do acoplamento flexível (2 - Figura 4.3), que acrescenta amortecimento no sistema, foi deliberadamente
desprezada.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
71
RotorDin − Geometria
LAVIB 1D N2
0.3
0.2
z (m)
0.1
0.0
−0.1
−0.2
−0.3
0.0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
y (m)
Figura 4.5: Geometria do sistema primário do LAVIB.
O modelo apresentado na Figura 4.5 foi obtido através do software RotorDin.
4.6
Sistema Auxiliar - Montagem
Quatro neutralizadores de 150 (g) (1 - Figura 4.6) são fixados dois a dois nas direções x e z
sobre o suporte do falso mancal de rolamento (2 - Figura 4.6). Este tipo de mancal permite um
giro limitado dentro do seu suporte e a livre rotação do eixo. A fixação dos neutralizadores é
feita através de parafusos M5. Por hipótese, o modelo matemático dos neutralizadores é efetivo
somente na direção do eixo principal.
Figura 4.6: Detalhe da fixação dos neutralizadores.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
4.7
72
Simulação Numérica Sem Neutralizadores
Com todos os dados geométricos e os parâmetros do materiais do sistema primário, é possível utilizando o RotorDin, obter o diagrama de Campbell, como mostrado na Figura 4.7 e na
Figura 4.8. Os neutralizadores ainda não estão aplicados ao sistema primário.
RotorDin − Campbell
LAVIB 1D N2
350
300
Freqüência (Hz)
250
200
150
100
50
0
0
2000
4000
6000
8000
10000
12000
14000
Rotação (rpm)
Figura 4.7: Diagrama de Campbell do sistema primário.
RotorDin − Campbell
LAVIB 1D N2
67
66
Freq. (Hz)
65
64
63
62
61
60
3400
3500
3600
3700
3800
3900
4000
4100
4200
Rotação (rpm)
Figura 4.8: Diagrama de Campbell do sistema primário (detalhe).
Verifica-se que a primeira rotação crítica ocorre quando o sistema gira a aproximadamente
3800 (rpm), da intercessão das linhas primeira horizontal e 1 × rpm.
A resposta a um desbalanceamento residual de 200 (g.mm), aplicado no disco de alumínio
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
73
(0.475 (m)) calculada pelo aplicativo RotorDin na coordenada do disco de medição (0.325 (m))
se encontra na Figura 4.9:
RotorDin − Resposta ao Desbalanceamento
Deslocamento − Amplitude − LAVIB 1D N2
−30
coord.(z) y= 0.325 (m)
−40
−50
Desl. dB ref = 1(m)
−60
−70
−80
−90
−100
−110
−120
−130
0
50
100
150
200
250
300
350
400
Freq. (Hz)
Figura 4.9: Resposta do sistema primário ao desbalanceamento.
4.8
Exportação das Matrizes
As matrizes complexas Ψ, Θ e λ que representam os parâmetros modais para o modelo simplificado (desbalanceamento) foram obtidas para utilização no programa de cálculo do neutralizador através de uma função de exportação para arquivo texto especialmente implementada no
R (ver anexo G.1). Cada elemento das matrizes é separado
software RotorDin (versão Scilab)
em termo real e imaginário, isto significa que são exportadas seis matrizes.
4.9
Parâmetros do Modelo de Derivada Fracionária do Material Viscoelástico
A borracha butílica empregada no projeto do neutralizador foi gentilmente cedida e caracterizada pelo professor Eduardo M. de Oliveira Lopes, da Universidade Federal do Paraná, UFPR
e apresenta os seguintes parâmetros do modelo de derivada fracionária (Equação 2.94):
• Temperatura referência T0 = 273(K), temperatura de trabalho T = 298(K), coeficiente
temperatura θ1 = 15.1, coeficiente temperatura θ2 = 171.
• Módulo cisalhamento inferior G0 = 1.53E6(N/m2 ), módulo cisalhamento superior G∞ =
1.116E8(N/m2 ), inclinação β = 0.396, deslocamento b1 = 1.34E − 2.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
C.Perda(pu) Cisa.(Pa)
10
8
10
−20
−30
1E+05
Gi=1.53E+06
Gf=1.12E+08
Bt=0.3960
Ds=1.34E−02
t1=15.10
t2=171.00
TR=273.00
Cisa.
C.Perda
1E+04
7
1E+03
1.00 6
10
1E+02
0.10 5
10
1E+01
10
0.01 4
10 −3
10
Frequencia (Hz)
9
Borracha Butilica
LAVIB 1D N2
40 30 20 10 0
−10
Temperatura (C)
74
1E+00
0
10
3
10
6
10
9
10
12
10
Frequencia Reduzida
Figura 4.10: Nomograma do material viscoelástico.
A Figura 4.10 apresenta o nomograma correspondente às características dinâmicas do material viscoelástico empregado no presente trabalho.
4.10
Simulação Numérica Com Neutralizadores
Para a simulação numérica, foram utilizadas as matrizes de parâmetros modais do sistema
primário, já previamente exportadas e os parâmetros do modelo do material viscoelástico já
definidos (4.9). Os neutralizadores foram localizados na seção 11 (Figura 4.5), correspondente
à localização do mancal flutuante (3 - Figura 4.6, 6 - Figura 4.3).
Os seguintes parâmetros foram adotados:
• O modo correspondente à primeira rotação crítica foi selecionado para o controle.
• Foi calculada uma massa de 0.3 (kg) para os neutralizadores (2 × 150 (g)/coord).
• Foram considerados os oito primeiros modos para o sub-espaço modal do espaço de estado (2n̂ = 16).
• A faixa de freqüência de projeto considerada foi de 35 a 73 (Hz), com 200 pontos intermediários.
• A freqüência inicial proposta ao otimizador foi de 58.9 (Hz).
• A barreira de freqüências definida foi de 33.4 a 74.8 (Hz).
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
75
Para o código numérico de cálculo do neutralizador foi solicitada a otimização para o controle
do primeiro modo, adotados os parâmetros supra-citados. A resposta ao desbalanceamento
obtida com e sem neutralizadores é mostrada na Figura 4.11. A localização da resposta coincide
com a posição dos “proximitors” (sensores de deslocamento, 5 - Figura 4.3)
Controle − Resposta ao Desbalanceamento
Deslocamento − LAVIB 1D N2 − butilica 298(K)
−40
sem neutralizadores
wo(1) = 60.2 (Hz)
exc.=13 med.=9
−50
−60
Desl (dB) ref=1(m)
−70
−80
−90
−100
−110
−120
−130
0
50
100
150
200
250
300
350
400
Frequencia (Hz)
Figura 4.11: Resposta ao desbalanceamento com e sem controle.
A freqüência ótima calculada foi de Ωn = 60.2 (Hz).
O gráfico em escala linear abaixo (Figura 4.12) permite avaliar diretamente a redução de
amplitude da resposta ao desbalanceamento para o primeiro modo de flexão.
Controle − Resposta ao Desbalanceamento 200(gmm)
Deslocamento − LAVIB 1D N2 − butilica 298(K)
6
sem neutralizadores
wo(1) = 60.2 (Hz)
exc.=13 med.=9
5
Desl (m) * (1000)
4
3
2
1
0
0
50
100
150
200
250
300
350
400
Frequencia (Hz)
Figura 4.12: Resposta linear com e sem controle.
Existe uma correspondência entre o diagrama de Campbell e a resposta ao desbalancea-
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
76
mento. Para sistemas simétricos quanto ao parâmetros do mancal como o “rotor rig” simulado,
conforme LALANNE e FERRARIS (2001), MUSZYNSKA (2005) e GENTA (2005) não há
amplificação na precessão reversa “backward” para excitação do tipo desbalanceamento. Portanto, para cada coincidência da reta de excitação devido ao desbalanceamento (rpm × rpm) e
as freqüências naturais do sistema para precessão direta “forward”, há uma amplificação correspondente a uma rotação crítica. Isto não significa porém, que o sistema não apresenta rotação
crítica para a precessão reversa, depende do tipo de excitação.
4.11
Construção Física do Neutralizador
Através da Equação 2.104 e do módulo de cisalhamento (elasticidade) do material para
a freqüência ótima, obtém-se o fator de forma L. Com o fator de forma L determinado e a
partir da espessura do material e do número de elementos em paralelo que se supõe conhecidos
ou atribuídos, obtém-se a área necessária do material viscoelástico para se obter a freqüência
natural do neutralizador dinâmico (Equação 2.105 ou Equação 2.106). Optou-se pela construção
do neutralizador com três elementos resilentes trabalhando em cisalhamento.
Com a massa do neutralizador e as propriedades do material conhecidas, obtém-se o volume/massa do neutralizador. A forma do mesmo depende muito das características do sistema.
Neste caso, optou-se pela geometria cilíndrica pela facilidade de usinagem e construção do
mesmo. O material escolhido foi aço inoxidável.
De amostras de borracha butílica de 4 (mm) de espessura previamente caracterizadas e cedidas pelo LAVIB, projetou-se o elemento viscoelástico (Figura 4.13).
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
77
1. núcleo central.
2. borracha butílica (material resilente).
3. massa de sintonização.
Figura 4.13: Geometria do projeto físico do neutralizador de 150 (g).
Através de um ensaio de massa aparente conforme descrito em (D.1),mediu-se uma freqüência de 58 (Hz). O neutralizador montado pode ser visualizado na Figura 4.6, número 1.
Figura 4.14: Modelo sólido do neutralizador.
O resultado do projeto do neutralizador está representado através da Figura 4.14. Observase a configuração dos elementos resilentes (estruturas mais escuras). Os mesmos foram projetados para um carregamento do tipo cisalhamento puro.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
4.11.1
78
Outros Materiais
Para o mesmo sistema girante descrito acima (seção 4.5) foram utilizados outros materiais
viscoelásticos para comparação dos resultados. As condições das simulações são:
• Massa do neutralizador 2 × 0.150 (kg)
• Temperatura referência : 273 (K) (*exceto Isodamp 286 (K)).
• Temperatura de trabalho: 298 (K).
A Tabela 4.1, sumariza os resultados obtidos quanto à freqüência natural dos neutralizadores na
linha Ωo .
Tabela 4.1: Freqüências ótimas para diversos materiais.
Material
Butílica Neoprene Butílica 20% Borr. Natural *Isodamp
G0 (N/m2 )
1.53E+06
3.57E+06
1.80E+06
2.79E+06
6.4E+05
G∞ (N/m2 ) 1.11E+08
1.79E+08
2.51E+08
8.16E+08
1.0E+09
b1
1.34E-02
2.46E-03
7.52E-04
3.55E-04
8.1E-03
β
0.396
0.435
0.479
0.297
0.55
θ1
15.1
6.57
5.85
9.74
72.1
θ2
171
68
92
148
696
Aten. (dB)
34
32
35
27
26
Ωo (rad/s)
378.2
396.5
392.7
400
307.2
Ωo (Hz)
60.2
63.1
62.5
63.5
48.9
Figura
4.11
4.15
4.16
4.17
4.18
Gráficos dos Resultados Obtidos com Outros Materiais
Na Figura 4.15 se observa o desempenho do neutralizador utilizando material neoprene.
Percebe-se devido ao menor amortecimento se comparado à borracha butílica, o aparecimento
de duas freqüências adjacentes e a boa redução na amplitude.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
79
Controle − Resposta ao Desbalanceamento
Deslocamento − LAVIB 1D N2 − neoprene 298(K)
−40
sem neutralizadores
wo(1) = 63.1 (Hz)
exc.=13 med.=9
−50
−60
Desl (dB) ref=1(m)
−70
−80
−90
−100
−110
−120
−130
0
50
100
150
200
250
300
350
400
Frequencia (Hz)
Figura 4.15: Desempenho do neoprene.
Abaixo (Figura 4.16) se observa o excelente desempenho obtido com a borracha butílica
com 20% de EPDM, com bom amortecimento, pequenas amplificações adjacentes e ótima atenuação.
Controle − Resposta ao Desbalanceamento
Deslocamento − LAVIB 1D N2 − butilica 20% 298(K)
−40
sem neutralizadores
wo(1) = 62.5 (Hz)
exc.=13 med.=9
−50
−60
Desl (dB) ref=1(m)
−70
−80
−90
−100
−110
−120
−130
0
50
100
150
200
250
300
350
400
Frequencia (Hz)
Figura 4.16: Desempenho da borracha butílica 20%.
A borracha natural apresenta um resultado satisfatório (Figura 4.17), com boa atenuação,
mas devido ao pequeno amortecimento a atuação do neutralizador se dá numa banda muito
estreita. Uma variação na freqüência natural do sistema primário ou na temperatura de trabalho
podem afetar o desempenho do dispositivo.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
80
Controle − Resposta ao Desbalanceamento
Deslocamento − LAVIB 1D N2 − borracha natural 298(K)
−40
sem neutralizadores
wo(1) = 63.5 (Hz)
exc.=13 med.=9
−50
−60
Desl (dB) ref=1(m)
−70
−80
−90
−100
−110
−120
−130
0
50
100
150
200
250
300
350
400
Frequencia (Hz)
Figura 4.17: desempenho da borracha natural.
O resultado obtido com o material C-1002 para a faixa de freqüência analisada, não é tão
bom em termos de redução da amplitude, mas devido ao alto amortecimento, se mantém praticamente constante mesmo se o sistema primário variar a freqüência natural. Também é menos
sensível à temperatura, podendo trabalhar numa ampla faixa.
Controle − Resposta ao Desbalanceamento
Deslocamento − LAVIB 1D N2 − isodamp C−1002 298(K)
−40
sem neutralizadores
wo(1) = 48.9 (Hz)
exc.=13 med.=9
−50
−60
Desl (dB) ref=1(m)
−70
−80
−90
−100
−110
−120
−130
0
50
100
150
200
250
300
350
400
Frequencia (Hz)
Figura 4.18: desempenho do material EAR isodamp C-1002.
4.11.2
Resposta em Freqüência
Utilizando as matrizes (parâmetros modais) exportadas a partir do modelo geral (não para
desbalanceamento), é possível também calcular a resposta em freqüência para o sistema em
repouso, Ωrpm = 0 ou para uma rotação Ωrpm qualquer. No exemplo abaixo, o sistema foi
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
81
excitado a 0.6 (m), da extremidade e a resposta obtida a 0.15 (m) (ver diagrama da Figura 5.2)
com neutralizadores de borracha butílica de 2 × 0.15 (kg) por direção (ver seção 4.5 e 4.9). O
resultado foi obtido com o software de projeto ótimo proposto. Neste caso o tipo de excitação
utilizado foi o unitário (resposta plana) na direção z.
Controle − Resposta em Frequencia
inertancia − LABVIB 1D N2 − butilica 298(K)
60
sem
wo(1) = 60.2 (Hz)
40
iner (dB) ref=1(m/Ns²)
20
0
−20
−40
−60
−80
0
50
100
150
200
250
300
350
400
freq. (Hz)
Figura 4.19: Resposta em freqüência com neutralizadores e Ωrpm = 0.
Para comparação dos custos computacionais envolvidos nos métodos utilizados, para o
mesmo problema de solução da resposta em freqüência do sistema girante em questão (ver seção
4.5) e cujas matrizes são da ordem de 200 × 200, a implementação dos parâmetros equivalentes
diretamente nas matrizes do modelo do sistema girante gastou aproximadamente 6800 (s)(1.9
horas) enquanto o mesmo problema foi resolvido em apenas 40 (s) utilizando a metodologia
proposta.
4.11.3
Variação dos Parâmetros do Material Viscoelástico com a Temperatura
Os parâmetros do material viscoelástico variam com a temperatura e podem afetar o desempenho do neutralizador. No exemplo abaixo aumentou-se em 20 (◦C) a temperatura de trabalho
do exemplo da Figura 4.19. Nestas condições, a diminuição do amortecimento (curva da Figura
4.10), provocou uma piora na performance e o aparecimento do segundo pico de freqüência
natural introduzido pelo neutralizador. Este tipo de comportamento está diretamente ligado ao
ponto de trabalho e do tipo de material viscoelástico utilizado.
Capítulo 4 Projeto do Neutralizador Para Desbalanceamento
82
Controle − Resposta em Frequencia
Inertancia − LAVIB 1D N2 − butilica 318(K)
60
sem neutralizadores
wo(1) = 60.2 (Hz)
exc.=16 med.=5
40
Iner (dB) ref=1(m/Ns²)
20
0
−20
−40
−60
−80
0
50
100
150
200
250
300
350
400
Frequencia (Hz)
Figura 4.20: Resposta na temperatura de 318 (K).
A Figura 4.20 mostra que a variação de temperatura tira o neutralizador de seu ponto ótimo
de trabalho. Assim se na temperatura de trabalho a redução de amplitude obtida foi de 25 (dB)
após a sintonização, devido à variação de temperatura caiu para 20 (dB), perdendo-se 5 (dB).
Cada material viscoelástico tem sua própria variação característica com a temperatura. Esperase que borrachas e neoprenes sejam menos “sensíveis” às variações de temperatura por estarem
classificados como materiais tipo I, ou seja Gr (Ω) ∼
= cte.
83
5
VALIDAÇÃO EXPERIMENTAL
A validação experimental realizada no LAVIB contou com uma adequada infraestrutura.
Uma bancada de medição com isolamento, um “kit” para dinâmica de rotores, equipamentos
de medição tais como, acelerômetros, martelo instrumentado, “proximitors”, analisador de vibração e placa e software de aquisição de sinais. A bancada foi configurada conforme seção
4.
Foram realizados dois tipos de ensaios distintos: resposta em freqüência com o sistema em
repouso (Ωrpm = 0) e medição de órbita com o sistema passando pela primeira rotação crítica.
5.1
Resposta em Freqüência do Sistema em Repouso
Este ensaio é realizado aplicando um impulso mecânico através de um martelo instrumentado e medindo a resposta com um analisador de espectro. Conhecendo-se a excitação e a
resposta, obtém-se a função resposta em freqüência.
Figura 5.1: Medição da resposta em repouso.
Capítulo 5 Validação Experimental
84
Figura 5.2: Esquema de medição FRF em repouso.
A medição realizada com o sistema primário em repouso (Ωrpm = 0) consistiu na medição
da resposta através de um acelerômetro (2 - Figura 5.1, ver anexo I.3) fixado ao eixo, do impacto
realizado com martelo (1 - Figura 5.1, ver anexo I.2) instrumentado conforme arranjo da Figura
5.1 e Figura 5.2. O analisador de sinais utilizado foi um HP 3560A (ver anexo I.1) ajustado
para a faixa de 0 − 500 (Hz) e 800 linhas.
Foram realizadas duas medições, a primeira somente sobre o sistema primário sem os neutralizadores e a segunda com os neutralizadores fixados em sua posição final, sobre o falso
mancal. Foi aplicado somente um neutralizador de 150 (g) por direção de deslocamento. A
função resposta em freqüência foi obtida para Inertância, a uma temperatura ambiente de aproximadamente 15 (◦C).
Resposta em Frequencia − Repouso
exc.=0.60(m) − resp.=0.16(m)
20
sem neutr.
com 1 neutr/coord.
inertancia (dB) ref=1(m/Ns²)
10
0
−10
−20
−30
−40
0
50
100
150
200
250
300
350
400
450
500
frequencia (Hz)
Figura 5.3: Resultado em repouso (inertância).
O resultado completo da simulação numérica é apresentada na Figura 5.4.
Capítulo 5 Validação Experimental
85
Controle − Resposta em Frequencia
Inertancia − LAVIB 1D N2 − butilica 286(K)
60
sem neutralizadores
wo(1) = 60.2 (Hz)
exc.=16 med.=6
40
Iner (dB) ref.=1(m/Ns²)
20
0
−20
−40
−60
−80
0
50
100
150
200
250
300
350
400
450
500
Frequencia (Hz)
Figura 5.4: Resultado completo em repouso.
Neste caso observa-se que o modelo utilizado não conseguiu simular corretamente a geometria do sistema completo, sistema girante, acoplamento e a base de montagem. Essa diferença
pode ser atribuída em parte ao modelo utilizado, já que o mesmo despreza por hipótese, aspectos físicos da montagem completa, principalmente no tocante ao acoplamento e a base de
montagem.
Os resultados medidos conforme a Figura 5.1 foram comparados com aqueles obtidos utilizando o software que implementa a metodologia de projeto proposta (Equação 3.17). Utilizouse o modelo geral com os parâmetros modais para o sistema em repouso, Ωrpm = 0, com e sem
neutralizador, com ênfase para o primeiro modo de flexão.
Controle − Resposta em Frequencia
Inertancia − LAVIB 1D N2 − butilica 286(K)
15
ens. sem
ens. 150g−z
num. sem
10
num. 150g−z
Iner (dB) ref.=1(m/Ns²)
5
0
−5
−10
−15
−20
−25
40
45
50
55
60
65
70
75
80
85
90
Frequencia (Hz)
Figura 5.5: Comparação de resultados em repouso.
Capítulo 5 Validação Experimental
86
A comparação dos resultados obtidos através da simulação numérica e do ensaio, para o
primeiro modo de flexão, mostra uma boa correlação tanto em termos de freqüência quanto
em amplitude. Em freqüência o erro relativo foi de ~3% e em amplitude 5 (dB) aquém do
esperado. Os erros observados podem ser atribuídos em parte ao modelo utilizado, que despreza
por hipótese, alguns aspectos físicos da montagem.
5.2
Resposta ao Desbalanceamento
A resposta ao desbalanceamento foi obtida registrando-se a órbita (Anexo F.1) obtida através da composição dos registros de deslocamento de dois “proximitors” (Anexo E.1) posicionados a 90 (graus) entre si em um suporte fixo (medição relativa)(5 - Figura 4.3).
Figura 5.6: Detalhe da fixação dos transdutores, em destaque eixo z.
O desempenho dinâmico do neutralizador é avaliado comparando uma órbita de referência
em uma rotação baixa, com a órbita na rotação crítica com e sem neutralizadores. Os gráficos
de órbita foram obtidos utilizando-se um sistema de aquisição de sinais tipo “order tracking”
(Anexo I.3). Para este fim, foi utilizando um software desenvolvido em “Labview” (1 ) e destinado a apresentar a órbita somente da freqüência fundamental (de rotação 1 × Ωrpm ), filtrando
eventuais componentes harmônicas presentes (2 × Ωrpm ). O software em questão foi desenvolvido internamente no LAVIB.
1 National
Instruments.
Capítulo 5 Validação Experimental
87
A órbita de referência foi obtida sem neutralizadores na rotação de1800 (rpm). Nesta rotação a escala foi ajustada para uma unidade, conforme a Figura 5.7.
Figura 5.7: Órbita de referência a 1800 (rpm).
Todo o conjunto foi levado até a rotação crítica, quando o deslocamento excessivo levou o
eixo a tocar na proteção mecânica, produzindo ruído elevado e não permitindo a máxima excursão da órbita. Não fosse a proteção, a amplificação seria tão grande que fatalmente inutilizaria
o eixo.
A órbita medida na rotação crítica de aproximadamente 3600 (rpm) com limitação de amplificação e sem controle, foi de aproximadamente 12.5 unidades como mostrado na Figura
5.8.
Capítulo 5 Validação Experimental
88
Figura 5.8: Órbita sem neutralizadores na rotação crítica a 3600 (rpm).
O sistema foi desligado e foram fixados os neutralizadores. O sistema foi então levado
novamente até a rotação crítica. A órbita máxima medida foi de aproximadamente 2.8 unidades
com mostra a Figura 5.9, demonstrando que a amplificação foi bem controlada e permitindo a
operação do sistema na rotação crítica sem que houvesse vibração excessiva, comprovando a
efetividade dos neutralizadores.
Figura 5.9: Órbita com neutralizadores na rotação crítica a 3600 (rpm).
Obteve-se uma redução de aproximadamente 4.5 vezes, ou -13 (dB), condizente com o
resultado da simulação numérica (Figura 4.11).
Capítulo 5 Validação Experimental
5.3
89
Tabela de Resultados
A Tabela 4.1 mostra a comparação dos resultados numéricos e experimentais para os ensaios com o sistema em repouso (realizado com um neutralizador por direção) e dinâmico juntamente com a redução obtida através da aplicação dos neutralizadores. (* Unidades em (pu),
** Referência 1 (m/Ns2 )).
tipo
condição
5.4
cálculo
sem com
Tabela 5.1: Resultados obtidos.
ensaio
unidade
red.
sem
com
red.
rotação
resposta
15
0
15
13
2.5
10.5
(dB) **
Ωrpm = 0
órbita
-
-
-
12.5 *
2.8 *
13
(dB) **
Ωrpm = Ωcrit
Discussão dos Resultados
Mesmo construído de forma artesanal e numa condição mecânica do material viscoelástico
não muito favorável (elemento esbelto em cisalhamento), a freqüência natural medida do neutralizador (ver anexo D.1) ficou muito próxima do valor ótimo obtido numericamente, o que
demonstra a validade do modelo utilizado e da caracterização do material.
A comparação dos resultados da resposta em freqüência para dois modelos consideravelmente diferentes, um obtido pelo método proposto com matrizes “fixas” do modelo “simplificado” (Figura 4.19) e o outro através da inserção dos parâmetros equivalentes diretamente nas
matrizes (Figura 3.6) mostra total correspondência. Por outro lado no que diz respeito ao tempo
computacional, o modelo proposto consome um tempo muito menor (cerca de 150x).
5.4.1
Teste de Resposta
Apesar do resultado medido estar próximo ao simulado (Figura 5.5), foi detectada a presença de amortecimento não previsto que atenuou a resposta nas freqüências naturais mais altas
e também um desvio no valor das freqüências naturais. Algumas hipóteses para tentar explicar
os desvios encontram-se listadas abaixo:
• Na Figura 4.3 observa-se a presença de um acoplamento de tiras de borracha, responsável
por transmitir o conjugado do motor, o qual não é levado em conta e introduz amortecimento não considerado no modelo.
Capítulo 5 Validação Experimental
90
• As tiras de borracha colocadas sobre as molas de sustentação do falso mancal estão modeladas como amortecimento viscoso, por simplicidade, embora sejam de material viscoelástico (Figura 4.6).
• O posicionamento da medição e excitação no software de projeto do neutralizador é discreto, devido à modelagem através de elementos finitos (Figura 4.5).
• O modelo advindo do software de dinâmica de rotores também diferiu no tocante às
freqüências naturais.
• A variação das propriedades do material viscoelástico com a temperatura, que foi considerada constante.
• Por hipótese foi considerado que os neutralizadores atuam somente na direção principal
de funcionamento, porém existe reação também nas outras direções.
• A rigidez dos suportes dos mancais de rolamento também não foi considerada.
Por outro lado, o resultado para a primeira freqüência natural, de análise, pode ser considerado
preciso.
5.4.2
Órbita
Devido à presença dos limitadores mecânicos não foi possível conhecer a amplitude máxima real na rotação crítica e nem o resultado efetivo da redução devido a presença dos neutralizadores.
Sem os neutralizadores e limitadores essa operação é praticamente impossível e normalmente a tentativa resulta em dano permanente do eixo, além do risco de danificar os elementos
de medição e outras partes do rotor.
Mesmo assim, os resultados obtidos foram “expressivos” principalmente no tocante à operação do sistema. Com a aplicação dos neutralizadores o sistema rotativo foi levado acima da
crítica e operado sobre a mesma sem vibração excessiva.
Os resultados obtidos encorajam novos estudos, para dar continuidade nesta linha de pesquisa.
91
6
CONCLUSÕES
Foi proposta uma metodologia geral de projeto ótimo de neutralizadores dinâmicos viscoelásticos para controle da vibração flexional em sistemas girantes, atuando em uma faixa ampla
de freqüências, a qual pode conter uma ou várias rotações críticas. A modelagem do sistema
composto foi desenvolvida no sub-espaço modal do espaço de estado do sistema primário a
partir dos parâmetros modais deste. O projeto considera que a excitação é do tipo desbalanceamento.
Para modelar o sistema primário foram revisados conceitos tais como, modelos de elementos finitos aplicados a dinâmica de rotores, equação de movimento flexional no espaço de
estado, problema de autovalores e o seu problema adjunto, calculo do diagrama de Campbell,
modelo simplificado para excitação do tipo desbalanceamento e resposta em diferentes pontos
do rotor para diferentes excitações no domínio da freqüência.
Para o projeto ótimo dos neutralizadores dinâmicos viscoelásticos foi revisada a metodologia geral desenvolvida pelo grupo PISA/CNPq aplicada a estruturas não girantes. Assim,
foram introduzidos ao longo deste trabalho conceitos como parâmetros equivalentes generalizados, modelo modal do sistema primário e a resposta do sistema composto (sistema primário
+ neutralizadores), função objetivo e vetor projeto, assim como também uma revisão dos modelos utilizados para materiais viscoelásticos, como o modelo de derivada fracionária de quatro
parâmetros.
Foram implementados os códigos numéricos necessários ao projeto dos dispositivos ótimos.
Para tal fim, foi implementado em SCILAB um novo código, que após receber os parâmetros
modais do sistema primário, obtidos com um código numérico próprio chamado RotorDin, é
capaz de obter, de forma ótima, os parâmetros físicos ótimos do conjunto de neutralizadores.
Foram implementados em SCILAB códigos numéricos para o cálculo do diagrama de
Campbell do sistema composto para simples comparação dos resultados obtidos com os códigos desenvolvidos, que usam no seu calculo os parâmetros modais do sistema primário obtido
através do modelo simplificado. Para implementação destes códigos foi necessário revisar a
teoria de dinâmica de rotores usando mancais compostos com material viscoelástico. Aqui,
deve-se considerar que os materiais viscoelásticos são fortemente dependentes da freqüência (e
Capítulo 6 Conclusões
92
da temperatura), o que leva a ter que resolver um diagrama de Campbell dentro de outro para
obter o diagrama de Campbell final.
Um sistema girante real, composto por um rotor simples, dois mancais rígidos, vários discos e um mancal flexível para posterior fixação dos neutralizadores foi modelado no código
RotorDin desenvolvido no LAVIB. Considerando excitação do tipo desbalanceamento, com
base no modelo simplificado, foram obtidos e exportados os parâmetros modais do sistema a
controlar. A partir do código numérico desenvolvido no presente trabalho, foram projetados de
forma ótima um conjunto de neutralizadores dinâmicos viscoelásticos destinados ao controle
do primeiro modo de flexão do rotor acima citado. Desta forma, foram obtidos os parâmetros
físicos ótimos e, para verificar a eficácia destes dispositivos, foram obtidas as curvas de resposta
do sistema girante com e sem a fixação do dispositivo de controle. Como se pode observar, foi
possível obter numericamente uma redução de 23 (dB).
Para validar a metodologia aqui apresentada, os neutralizadores dinâmicos viscoelásticos
foram construídos e montados sobre o sistema girante proposto no modelo numérico. O sistema composto foi então submetido ao teste de resposta em freqüência e medição de órbita a
diferentes rotações e os resultados foram comparados com o sistema sem neutralizadores para
verificar a eficácia destes dispositivos experimentalmente.
As figuras das órbitas obtidas experimentalmente mostram uma redução da resposta vibratória mínima de 13 (dB). A redução medida seria maior se os limitadores pudessem ter sido
retirados durante a operação sem neutralizadores, o que por outro lado poderia levar a ruptura
do rotor ou provocar grandes deformações plásticas o que poderia danificar os sensores. A
função resposta em freqüência obtida do teste de resposta em freqüência com o rotor parado,
comprovou a validade dos modelos adotados ao apresentar resultados compatíveis com os obtidos numericamente. O erro na diferença entre os resultados obtidos com e sem neutralizadores,
numérica e experimentalmente foi de 4.5 (dB) para mais.
Utilizando o código numérico desenvolvido, foi possível simular o sistema composto também através da inclusão do modelo de parâmetros equivalentes generalizados dos neutralizadores dinâmicos, diretamente nas matrizes do modelo do sistema primário. Desta forma, foi
possível obter o diagrama de Campbell e as respostas em freqüência e ao desbalanceamento
do sistema composto. Os resultados obtidos através desta simulação, se mostraram idênticos
quando comparados aos obtidos através da metodologia proposta.
Verificou-se durante as simulações numéricas, que as propriedades do material viscoelástico
empregado nos protótipos dos neutralizadores, são sensíveis à variação da temperatura. No
presente trabalho, foi realizada uma simulação variando a temperatura para mostrar como a
Capítulo 6 Conclusões
93
mesma afeta o desempenho dos dispositivos. Dependendo do material viscoelástico utilizado,
pequenas variações de temperatura podem provocar perdas na redução de vibração já que o
neutralizador deixa de trabalhar no seu ponto ótimo.
Como ficou demonstrado através da metodologia proposta, é possível projetar de forma
ótima, neutralizadores dinâmicos viscoelásticos para sistemas girantes. A aplicação destes dispositivos amplia a faixa de utilização do sistema. Permite não somente a passagem suave pelas
rotações críticas, mas também a operação nestas rotações. Para bancada de dinâmica de rotores
do LAVIB, os dispositivos vão permitir a realização de testes acima da rotação crítica, como
balanceamento de eixos flexíveis, por exemplo.
94
7
SUGESTÕES PARA TRABALHOS FUTUROS
Por envolver diversas áreas de estudo, controle de vibração, material viscoelástico, otimização e dinâmica de rotores às vezes é necessário desenvolver novas metodologias já que a
literatura existente é fragmentada e escassa para este tipo de aplicação.
As sugestões apresentadas abaixo estão baseadas nas dificuldades encontradas no contexto
do desenvolvimento deste trabalho e nas possibilidades ainda em aberto.
• Desenvolver uma metodologia de cálculo da massa dos neutralizadores no modelo de
espaço de estado.
• Considerar a variação dos parâmetros dos mancais com a rotação.
• Desenvolver um estudo da estabilidade do sistema composto.
• Realizar medições e simulações da amplificação devido ao desbalanceamento e outros
tipos de excitação.
• Avaliar as simplificações e hipóteses adotadas, através de validações experimentais.
• Considerar o efeito dos neutralizadores fora da sua direção principal.
• Rever a metodologia empregada na obtenção do diagrama de Campbell do sistema composto.
95
REFERÊNCIAS BIBLIOGRÁFICAS
ALAUZE, C.; HAGOPIAN, J. Der; GAUDILLER, L. Active balancing of turbomachinery:
Application to large shaft lines. Journal of Vibration and Control, n. 7, p. 249–278, 2001.
ARORA, Jasbir S. Introduction to Optimum Design. San Diego, California, USA: Elsevier
Academic Press, 2004.
BATHE, K. J. Finite Element Procedures. New Jersey, USA.: Prentice Hall, 1996.
BAVASTRI, Carlos Alberto. Redução de Vibrações de Banda Larga em Estruturas Complexas
por neutralizadores Vicoelásticos. Tese (Doutorado) — UFSC - Universidade Federal de Santa
Catarina, Florianópolis - SC, 1997.
BAVASTRI, Carlos Alberto; FILHO, Francisco José Doubrawa; ESPíNDOLA, José João
de; LOPES, Eduardo Márcio de Oliveira; VENANCIO, Herbert W. Modelo geral de
neutralizadores dinâmicos para controle de vibrações e ruído: Parâmetos equivalentes
generalizados. Congresso Ibero Latino Americano sobre Métodos Computacionais em
Engenharia, 2007.
BAVASTRI, Carlos Alberto; POLLI, M. L.; PRESEZNIAK, Flávio A. Optimal passive
control of self-excitation vibrations in cutting process using viscoelastic dynamic neutralizers.
Usinagem 2006, 2006.
BAZáN, Fermín Sinforiano Viloche; BAVASTRI, Carlos Alberto. An optimized pseudo-inverse
algorithm (opia) for multi-input multi-output modal parameter identification. Mechanical
Systems and Signal Processing, v. 10, n. 4, p. 365–380, 1996.
CHILDS, Dara. Turbomachinery Rotordynamics. [S.l.]: John Wiley & Sons, Inc., 1993.
CHÁVEZ, Ramiro Germán Díaz. Dinâmica de um Rotor Horizontal em Apoios Elásticos.
Dissertação (Mestrado) — Pontifícia Universidade Católica do Rio de Janeiro, 2003.
CLARK, S. L. Dynamics of Continuous Elements. [S.l.]: Prentice-Hall, Inc., 1972.
CRUZ, Gilberto Amado Méndez. Projeto Ótimo de Neutralizadores Viscoelásticos Baseado
no Modelo à Derivadas Fracionárias. Tese (Doutorado) — UFSC - Universidade Federal de
Santa Catarina, Florianópolis - SC, 2004.
DUNKERLEY, Stanley. On the whirling and vibration of shafts. Philosophical Transactions of
the Royal Society of London, 1894.
ERICH, Fredric F. Handbook of Rotordynamics. 2nd.. ed. Florida, USA: Krieger Publishing
Company, 2004.
ESPÍNDOLA, José João de. Controle de Vibração,. 1992. Notas de Aula.
Referências Bibliográficas
96
ESPÍNDOLA, José João de; BAVASTRI, Carlos Alberto. An efficient concept of
transmissibility for a general equipment isolation system. ASME Design Engineering Technical
Conferences, DETC9́7/VIB-4120, 1997.
ESPÍNDOLA, José João de; LOPES, Eduardo Márcio de Oliveira; BAVASTRI, Carlos Alberto.
Optimum system of viscoelastic vibration absorbers by fractional calculus. 2nd IFAC Workshop
on Fractional Differentiation and its Applications, Porto, Portugal, v. 01, p. 409–415, 2006.
ESPÍNDOLA, José João de; SILVA, Hilton Penha. Modal reduction of vibrations by dynamic
neutralizers: ageneralized approach. IMAC-10 - 10th International Modal Analysis Conference.
Society for Experimental Mechanics, San Diego., v. 02, p. 1367–1373, 1992.
ESPíNDOLA, José João de; LOPES, Eduardo Márcio de Oliveira; BAVASTRI, Carlos Alberto.
Design of optimum systems of viscoelastic vibration absorbers for a given viscoelastic material
based on the fractional calculus model. Modeling and Control of Autonomous Decision Support
Based Systems, 2005.
EWINS, D. J. MODAL TESTING 2 theory, practice and application. [S.l.]: RESEARCH
STUDIES PRESS LTD., 2000.
FERREIRA, Euda M. da Silva. Modelo de Rotores Dinâmicos com mancal utilizando material
viscoelástico. Dissertação (Mestrado) — Universidade Federal do Paraná, 2005.
FRAHM, H. Device for Damping Vibration of Bodies. October 1909. U.S. Patent No. 989959.
GENTA, Giancarlo. Dynamics of Rotating Systems. Fst. New York, USA: Springer
Science+Business Media, Inc., 2005.
HARTOG, J. P. Mechanical Vibrations. 4th. ed. [S.l.]: Dover Publications INC., 1956. 1985 republication.
ASME/IGTI TURBO EXPO 98. Adaptavive Vibration Control of Industrial Turbomachinery.
JEFFCOTT, H. The lateral vibration of loaded shafts in the neighborhood of a wirling speed the effect of want of balance. Philosophical Magazine, v. 37, n. 6, p. 304–314, 1919.
JONES, David I. G. Viscoelastic Vibration Damping. [S.l.]: John & Wiley Sons LTD., 2001.
KITIS, L. Vibration reduction over a frequency range. Journal od Sound and Vibration, v. 89,
n. 4, p. 559–569, 1983.
KRENEV, Boris G.; REZNIKOV, Leonid M. Dynamic Vibration Absorbers Theory and
Tecnical Appications. West Sussex, England: John Wiley & Sons Ltd, 1993.
LAKES, R. S. Viscoelastic Solids. [S.l.]: CRC Press LLe., 1999.
LALANNE, Michel; FERRARIS, Guy. Rotordynamics prediction in engineering. 2nd. ed.
New York, USA: John Wiley & Sons, Inc., 2001.
LOPES, Eduardo Márcio de Oliveira; BAVASTRI, Carlos Alberto; NETO, J. M. Silva;
ESPíNDOLA, José João de. Caracterização dinâmica integrada de elastômeros por derivadas
generalizadas. III Congresso Nacional de Engenharia Mecânica, 2004.
Referências Bibliográficas
97
LÓPEZ, Ezequiel José. Dinâmica de Rotores: Modelo Matemático - Realização Prática.
Dissertação (Mestrado) — Universidad Nacional del Comahue, 2002.
MARDLE; PASCOE, Sean. An overview of genetic algorithms for the solution of optimisation
problems. CHEER, Computers in Higher Education Economics Review, 1999.
MEIROVITCH, L. Dynamics and Control of Structures. [S.l.]: Wiley Interscience, 1990.
MUSZYNSKA, Agnieszka (Agnes). Rotordynamics. [S.l.]: Taylor & Francis, 2005.
MYKLESTAD, N. O. A new method of calculating natural modes of uncoupled bending
vibration of airplane wings and other types of beams. Journal of Aeronautical Science, v. 2,
n. 11, p. 153–162, 1944.
NASHIF, Ahid D.; JONES, David I. G.; HENDERSON, J. P. Vibration Damping. [S.l.]: John
Wiley & Sons, 1985.
NISHTALA, Rajesh; VUDUC, Richard; DEMMEL, James; YELICK, Katherine. When cache
blocking sparse matrix vector multiply works and why. In: Proceedings of the PARA’04
Workshop on the State-of-the-art in Scientific Computing. Copenhagen, Denmark: [s.n.], 2004.
PRASAD, Dr. F. H. Determination of stiffness of roller bearings – an alternative approach. IEI - MC, v. 84, n. 4, p. 186–192, JAN 2004. Disponível em:
<http://www.ieindia.org/publish/mc/0104/jan04mc9.pdf>.
PRITZ, T. Analysis of four-parameter fractional derivative model of real solid materials.
Journal of Sound and Vibration, v. 195, n. 1, p. 103–115, 1996.
PROHL, Melvin A. A general method for calculating the critical speeds of flexible rotors.
Applied Mechanics, v. 12, n. 3, p. 142–148, 1945.
RANKINE, W. Centrifugal whirling of shafts. The Engineer, April 1869.
ROGERS, L. Operators and fractional derivatives for viscoelastic constitutive equations.
Journal of Rheology, v. 27, p. 351–372, 1983.
SANTOS, Hideraldo Luis Vasconcelos dos. Dinâmica de Rotores: Modelo numérico
considerando mancais viscoelásticos. 2003. Trabalho de Conclusão de Curso. (Graduação em
Engenharia Mecânica) - Universidade Federal do Paraná.
SCHOENBERG, Ronald. Optimization with the quasi-newton method. Aptech Systems, Inc.,
2001.
SILVA, Hilton Penha. Controle Modal de Vibrações por Neutralizadores Dinâmicos: Uma
Abordagem Generalizada. Dissertação (Mestrado) — UFSC - Universidade de Santa Catarina,
1991.
SMITH, David M. The motion of a rotor carried by a flexible shaft in flexible bearings.
Procedings Royal Soc. London, v. 142, p. 92–118, 1933.
SNOWDON, J. C. Steady-state behavior of the dynamic absorber. Journal of Acouustical
Society of America, v. 38, n. 8, p. 1096–1103, August 1959.
STODOLA, A. Steam and Gas Turbines. New York, USA: Mc Graw-Hill, 1927.
Referências Bibliográficas
98
TAMMI, Kari. Active control of radial rotor vibrations. Tese (Doutorado) — Helsinki
University of Technology, Finland, may 2007.
VANCE, John M. Rotordynamics of Turbomachinery. New York, USA: Wiley & Sons, Inc.,
1988.
ZEIDAN, Fouad Y. Application of squeeze film dampers. Turbomachinery International,
September/October 1995.
ZHOU, Shiyu; SHI, Jianjun. Active balancig and vibration control of rotating machinery: A
survey. The Shock and Vibration Digest, p. 360–371, September 2001.
99
ANEXO A -- SISTEMA NÃO GIRANTE
Uma vez determinadas as matrizes de massa [M] e rigidez [K] da Equação 2.37, em termos
de coordenadas generalizadas q, obtém-se
[M]{q̈(t)} + [C]{q̇(t)} + [K]{q(t)} = { f (t)}.
Eq. A.1
Considerando nulas a excitação { f (t)} = {0} e a matriz de amortecimento [C] = [0], uma proposta de solução para o sistema poderia ser: {q(t)} = {φ }eiΩt . Aplicando esta proposta e as
considerações na (Equação A.1 ) e definindo: λ , Ω2 , pode-se chegar ao seguinte conjunto de
problemas de autovalores
[K][φ ] = [λ ][M][φ ].
Eq. A.2
Resolvendo-se este conjunto de problemas, pode-se obter as freqüências naturais através do
elementos da matriz diagonal de autovalores [λ ], e os modos de vibrar das colunas associadas
da matriz de autovetores [φ ]. Obtém-se uma matriz de massa modal ortogonalizada (1 ), através
da operação, como mostrado em ESPÍNDOLA (1992):
r
M jr = [φ ]T [M][φ ].
Eq. A.3
Pode-se ainda, ortonormalizar a matriz de autovetores [φ ], dividindo cada j-ésima coluna pela
√
raiz quadrada dos elementos correspondentes da diagonal da matriz r m jr , m j j , obtendo-se
a matriz [Φ].
A matriz de amortecimento pode ser obtida adotando-se um modelo viscoso, proporcional
às matrizes de massa e rigidez, logo
[C] = [α[M] + β [K]] .
Eq. A.4
Sendo assim demonstra-se (ESPÍNDOLA (1992)), que esta matriz de amortecimento, também
pode ser diagonalizada através das matrizes de autovetores
r
C jr = [ΦT ][C][Φ].
1 por
definição, para uma matriz ortogonal [A]: [A]−1 = [A]T ⇒ [A][A]T = [I]
Eq. A.5
ANEXO A -- Sistema Não Girante
100
Por intermédio da transformação
{q(t)} = [Φ]{p(t)},
Eq. A.6
a Equação A.1 , uma vez pré-multiplicada por [Φ]T , pode ser reescrita como
[Φ]T [M][Φ]{ p̈(t)} + [Φ]T [C][Φ]{ ṗ(t)} + [Φ]T [K][Φ]{p(t)} = [Φ]T { f (t)}.
Eq. A.7
Aplicando-se a transformada de Fourier na Equação A.7 tem-se
−Ω2 [Φ]T [M][Φ] + iΩ[Φ]T [C][Φ] + [Φ]T [K][Φ] {P(Ω)} = [Φ]T {F(Ω)}.
Eq. A.8
Pode-se demonstrar (ESPÍNDOLA (1992)) que as matrizes abaixo são diagonais e se reduzem
a
•[Φ]T [M][Φ] = [I] (matriz identidade [r 1r ]).
2ζ j Ω jr .
•[Φ]T [C][Φ] =
r
•[Φ]T [K][Φ] =
r
λ jr (Ω2j =
kj
m j ).
Aplicando as relações na Equação A.8 , pode-se obter
[D(Ω)]{P(Ω)} = [Φ]T {F(Ω)},
Eq. A.9
[D(Ω)] = −Ω2 [I] + iΩ r 2ζ j Ω jr + r λ jr .
Eq. A.10
com
Da transformação de variáveis da Equação A.6 , tem-se
{Q(Ω)} = [Φ][D(Ω)]−1 [Φ]T {F(Ω)}.
Eq. A.11
Por definição (Equação 2.100), a função de transferência é a relação entre deslocamento e força
generalizadas (saída/entrada)
[H(Ω)] , [Φ][D(Ω)]−1 [Φ]T .
Eq. A.12
Para se obter uma função resposta em freqüência, pode-se calcular
N
αks (Ω) =
∑
j=1
jAks
,
2
−Ω + Ω2j + i2ζ j ΩΩ j
Eq. A.13
ANEXO A -- Sistema Não Girante
101
já que a matriz função resposta em freqüência no espaço modal [D(Ω)]−1 é diagonal. Com s
representando o ponto de excitação, k o de medição e jAks a j-ésima constante modal, jAks =
Φk j Φs j , para N modos.
A.1
Introdução dos Neutralizadores no Modelo
Acrescentando-se p neutralizadores ao sistema primário nas coordenadas generalizadas q,
pode-se agregar as p massas e amortecimentos equivalentes generalizados, como massas me (Ω)
e amortecimentos ce (Ω), somando-as apropriadamente na Equação A.10 .
Figura A.1: Viga e modelo equivalente generalizado.
O novo sistema, pode agora, ser representado pela equação
[M̃]{q̈(t)} + [C̃]{q̇(t)} + [K]{q(t)} = { f (t)},
Eq. A.14
[M̃] = [M] + [Me ] e [C̃] = [C] + [Ce ].
Eq. A.15
sendo
As matrizes de massa e amortecimento equivalentes [Me ] e [Ce ] são diagonais e contém elementos somente nas colunas correspondentes à coordenada generalizada q de fixação do neutralizador.
Aplicando novamente a transformação de variáveis (Equação A.6 ) e pré-multiplicando por
[Φ]T , obtém-se
[D̃(Ω)]{P(Ω)} = [Φ]T {F(Ω)},
Eq. A.16
ANEXO A -- Sistema Não Girante
102
[D̃(Ω)] = −Ω2 [I] + [Φ]T [Me ][Φ] + iΩ r 2ζ j Ω jr + [Φ]T [Ce ][Φ] + r λ jr . Eq. A.17
Um sistema muito discretizado pode apresentar centenas e até milhares de graus de liberdade.
Normalmente se trabalha apenas com o primeiros n̂ autovalores e autovetores, que representam
os modos dentro da faixa de freqüência de interesse. As vezes se trabalha com alguns poucos
modos, dentro da faixa de interesse, e dois modos residuais, inferior e superior, que representam
os modos anteriores e posteriores à faixa de interesse. Então, é possível utilizar somente as
primeiras colunas da matriz de autovetores: [Φ̂] = [Φ]n×n̂ para obtenção das novas matrizes
adicionais de massa:
[M̂] = [Φ̂]T [M̃][Φ̂],
Eq. A.18
[Ĉ] = [Φ̂]T [C̃][Φ̂],
Eq. A.19
e de amortecimento
sendo n̂ × n̂ a dimensão final destas matrizes.
Note-se que essas novas matrizes “não” são diagonais, já que o modelo dos neutralizadores
foi acrescentado “sobre” a matriz de autovetores do sistema primário. Isto implica em que
os modos não são ortogonais em relação à matrizes que não estavam presentes quando foram
geradas.
Considera-se com base em SILVA (1991), no entanto, que estas matrizes são “predominantemente” diagonais. Esta consideração simplificadora levou à definição de uma metodologia de
projeto ótimo de neutralizadores, como é mostrado à seguir, mas na realidade, as matrizes não
são rigorosamente predominantes, para quaisquer conjuntos de neutralizadores.
Obtém-se um novo sistema de equações “truncado” para o número de modos, apropriado
para o sistema em estudo
[D̂(Ω)]{P̂(Ω)} = [Φ̂]T {F(Ω)}.
Eq. A.20
O vetor de coordenadas principais pode ser obtido a partir da seguinte equação
{P̂(Ω)} = [D̂(Ω)]−1 [Φ̂]T {F(Ω)},
h
i
ˆ + [M̂e ] + iΩ r 2ζ ˆj Ω jr + [Ĉe ] + r λˆ jr .
[D̂(Ω)] = −Ω2 [I]
Eq. A.21
Eq. A.22
com
[M̂e ] = [Φ̂]T [Me ][Φ̂] e [Ĉe ] = [Φ̂]T [Ce ][Φ̂] .
Eq. A.23
Na Figura A.2, pode-se observar funções resposta em freqüência para as coordenadas principais
{P(Ω)}, dos cinco primeiros modos de um sistema, onde é possível identificar o acoplamento
modal introduzido pela presença de dois neutralizadores projetados para o segundo e terceiros
ANEXO A -- Sistema Não Girante
103
modos.
Sistema MGL − VE mi=0.2
Viga simples − coord. principais
−60
1º modo
2º modo
3º modo
4º modo
5º modo
−70
−80
P (dB ref=1)
−90
−100
−110
−120
−130
Acoplamento
−140
−150
0
500
1000
1500
2000
2500
3000
Freq. (rad/s)
Figura A.2: Acoplamento das coordenadas principais {P(Ω)}.
A partir da Equação A.20 , é possível determinar a resposta do sistema
{Q(Ω)} = [Φ̂][D̂(Ω)]−1 [Φ̂]T {F(Ω)}.
Eq. A.24
Verificando a Equação A.24 , pode-se obter uma expressão matricial para a função resposta em
freqüência, novamente definindo a função de transferência através da relação entre a resposta e
a excitação
[H(Ω)] , [Φ̂][D̂(Ω)]−1 [Φ̂]T .
Eq. A.25
Para uma resposta na coordenada k com excitação em s, a função pode ser dada por
n̂
αks (Ω) =
n̂
∑ ∑ D̆i j (Ω)Φ̂siΦ̂k j ,
Eq. A.26
j=1 i=1
com [D̆(Ω)] = [D̂(Ω)]−1 . A expressão acima (A.26 ), é análoga a obtida normalmente, como
na expressão da Equação A.13 .
A.1.1
Massa dos Neutralizadores
Como mostrado na Equação A.14 , as matrizes de massa [M̃] e amortecimento [C̃] são
resultado da adição das matrizes do sistema primário e das matrizes de massa e amortecimento
ANEXO A -- Sistema Não Girante
104
equivalentes, função da freqüência Ω, do sistema auxiliar (neutralizadores).


0 0
...
0
0


 0 ...
...
0
0 





 0 0 me1 (Ω)
0
0




[M̃] = [M] +  ... ...
...
...
...  = [M] + [Me ]



 0 0
...
m
(Ω)
0


e2


..
 0 0
.
...
0 


0 0
...
0
0
Eq. A.27
e








[C̃] = [C] + 






0
...
0
...
0

0

0 


0 0 ce1 (Ω)
0
0 


... ...
...
...
...  = [C] + [Ce ].

0 0
...
ce2 (Ω) 0 


..
.
0 
0 0
...

0 0
...
0
0
0
...
0
Eq. A.28
As matrizes ortonormalizadas de massa e amortecimento equivalentes, truncadas para o número
de modos n̂, podem ser obtidas como na Equação A.23

p
p
2
m
Φ
ei
∑ meiΦki1Φki2 ...
 ∑
ki1
 i=1
i=1
 p
p

mei Φ2ki2
...
 ∑ mei Φki2 Φki1
∑

[M̂e ] =  i=1
i=1

...
...
...

 p
p

∑ meiΦkin̂Φki1 ∑ meiΦkin̂Φki2 ...
i=1
p

i=1
p











∑ meiΦki1Φkin̂
∑ meiΦki2Φkin̂
i=1
i=1
...
p
∑ meiΦ2kin̂
i=1
,
Eq. A.29
n̂×n̂
com mei = mei (Ω) e i = 1 a p. A expressão para [Ĉe ] é totalmente idêntica.
Uma comparação entre as equações do sistema com e sem p neutralizadores (A.10 e A.17
), revela que a massa adicional a ser somada é dada pela expressão: [Φ̂]T [Me ][Φ̂].
Dessa observação e supondo que a matriz [M̂e ] estivesse perfeitamente diagonalizada, SILVA
(1991) propôs a seguinte relação para a massa adicional ma em sistemas com múltiplos graus
de liberdade, considerando n̄ modos a serem controlados “modo a modo”
p
µ j m j = ∑ mai Φ̂2ki j .
i=1
Eq. A.30
ANEXO A -- Sistema Não Girante
105
Expandindo o somatório, obtém-se
Φ̂2k11 ma1 + Φ̂2k21 ma2 + ...+
Φ̂2kp1 map =
µ1
Φ̂2k12 ma1 + Φ̂2k22 ma2 + ...+
Φ̂2kp2 map =
µ2
...
...
...
...
...
Φ̂2k1n̄ ma1 + Φ̂2k2n̄ ma2 +
...
Eq. A.31
+Φ̂2kpn̄ map = µn̄ ,
com i = 1 a p. O índice ki representa a i-ésima coordenada generalizada do neutralizador.
Para os n̄ modos de interesse a serem controlados e os p neutralizadores, é possível montar
um sistema de equações de dimensão n̄ × p, teoricamente solucionável
2
Φ̂ki n̄×p {ma } p×1 = {µ}n̄×1 ,
Eq. A.32
2
Φ̂ki {ma } − {µ} = {0} .
Eq. A.33
logo
As massas ótimas podem ser obtidas através de alguma técnica de solução por mínimos
quadrados, por exemplo. Este tipo de solução nem sempre conduz à resultados fisicamente
coerentes. Em alguns casos, onde há simetria no sistema, por exemplo, não existe uma solução
realizável.
Nestes casos, pode-se optar por considerar todas as massas dos neutralizadores iguais. Para
o j-ésimo modo, a massa do neutralizador é dada por
µ jm j
ma j =
.
p
∑ Φ̂2ki j
Eq. A.34
i=1
Observação: m j = 1, se os autovetores estiverem ortonormalizados.
São obtidas n̄ massas mai , uma para cada modo a ser controlado. A massa final é a média
n̄
∑ mai
ma =
i=1
n̄
.
Eq. A.35
A obtenção da massa dos neutralizadores através da matriz de autovetores do sistema primário
traz intrinsecamente uma consideração modal. Se a localização do neutralizador for apropriada, há uma minimização da massa calculada.
A massa média modal foi proposta por BAVASTRI (1997), para controle de n̄ modos simultâneos em banda larga de freqüência.
ANEXO A -- Sistema Não Girante
A.2
106
Otimização Multi-Modal
Em um sistema com múltiplos graus de liberdade, é possível calcular neutralizadores visando atenuar um ou mais modos específicos. Existem duas aproximações possíveis:
1.Modo a modo: Neste caso assume-se que o desacoplamento é total. Os neutralizadores
podem então, ser projetados modo a modo, utilizando a técnica dos pontos fixos, por
exemplo. No entanto, a premissa básica, como mostrado na Figura A.2, não pode ser satisfeita e o desempenho do sistema auxiliar pode ser afetado. Esta técnica é recomendada
para sistemas com modos bem separados e definidos.
2.Banda larga: Projetam-se os neutralizadores simultaneamente, considerando a resposta
sistema para n̄ modos, através da utilização de uma técnica de otimização não linear,
buscando a minimização de uma função objetivo. Normalmente se trabalha com um
vetor projeto contendo as freqüências naturais dos neutralizadores e o amortecimento do
neutralizador, uma vez a massa fixada através da Equação A.30 , por exemplo. Quando
o neutralizador é construído com material viscoelástico, somente as freqüências naturais
são suficientes para a determinação das características físicas dos neutralizadores através
da Equação 2.104.
O resultado obtido depende em grande parte, da escolha correta da posição dos neutralizadores
e da função objetivo.
Para a otimização, as matrizes de massa e amortecimento equivalentes, são montadas e
somadas aos parâmetros modais do sistema primário, previamente estabelecidos, para cada
freqüência na faixa adotada e com o número de modos n̂, escolhido.
Dentro da faixa de freqüências é buscado o máximo para a n̄-ésima coordenada principal
modal escolhida a ser controlada. Ao final, retorna ao otimizador, a “norma 2” do vetor de
máximos calculados.
Abaixo a função resposta em freqüência dos cinco primeiros modos de um sistema, ao
qual foram acrescentados dois neutralizadores projetados para o segundo e terceiros modos
utilizando a função objetivo proposta por BAVASTRI (1997).
ANEXO A -- Sistema Não Girante
107
Sistema MGL − VE mi=0.2
Viga simples − FRF k=3 s=3
−60
sem neutralizadores
com 2 neutralizadores
−70
−80
H (dB ref=1)
−90
−100
−110
−120
−130
−140
−150
0
500
1000
1500
2000
2500
3000
Freq. (rad/s)
Figura A.3: Otimização multi-modal.
Algumas questões devem ser observadas:
•O vetor de força escolhido deve representar o tipo de excitação em questão. Uma excitação localizada em apenas um nó, de valor unitário, pode representar um impulso,
excitando diversos modos. Atenção especial ao ponto escolhido, visando evitar nós dos
modos que se pretende controlar. Pode-se utilizar também uma excitação distribuída, com
impulsos em vários ou todos os nós, representado por exemplo ação do vento sobre cabos.
•A escolha das coordenadas de aplicação dos neutralizadores é fundamental nos aspectos
de minimização de massa e redução da amplitude de vibração minimizada
•As propriedades do material viscoelástico devem estar adequados à faixa de freqüência
que se deseja controlar, caso contrário pode ser virtualmente impossível a realização física
dimensional dos mesmos.
A.3
Funções Objetivo
Uma função objetivo a ser minimizada pela técnica não linear, deve retornar um número
que represente uma quantidade da resposta vibratória do sistema composto dentro da faixa de
freqüência de interesse, para um dado conjunto de neutralizadores. As variáveis que compõem o
vetor projeto, são os parâmetros dos neutralizadores (as freqüências naturais no caso viscoelástico). Como resultado, os parâmetros obtidos devem levar o sistema composto a uma resposta
vibratória minimizada.
ANEXO A -- Sistema Não Girante
108
O sistema primário é representado através de suas matrizes de autovalores [λ ], constantes
de amortecimento [ζ ] (ou fatores de perda [η]) e autovetores normalizados à esquerda [Ψ] e à
direita [Θ]. Estas matrizes podem ser obtidas de distintas maneiras, normalmente através de
algum software de elementos finitos ou alguma técnica de identificação (análise modal), da
estrutura que se deseja modificar.
BAVASTRI (1997) propôs a minimização da função objetivo “norma 2”, do vetor de máximos modais das coordenadas principais reduzidas {P̂(Ω)}, em uma certa faixa de freqüências. A otimização é máxima como pretendido, sobre as coordenadas principais, o que pode
ser comprovado pela análise qualitativa do gráfico (A.2), já que os picos tem a mesma amplitude (HARTOG (1956)). O resultado na função resposta em freqüência pode ser ligeiramente
diferente.
109
ANEXO B -- ALGORITMOS DE OTIMIZAÇÃO
B.1
“Quasi-Newton”
Segundo SCHOENBERG (2001) os extremos de uma função f (x) são caracterizados pela
anulação da sua derivada f 0 (x) = 0. Em funções com múltiplas variáveis pode não haver solução
direta. Newton propôs uma solução iterativa: a partir de uma condição inicial, faz-se uma
aproximação quadrática da função através da série de Taylor para uma máximo local xm , faz-se
a diferenciação e iguala-se a zero, calculando o ponto de máximo local. Faz-se deste ponto, a
nova condição inicial e assim sucessivamente
1
f (x) = f (xm ) + f 0 (xm )(x − xm ) + f 00 (xm )(x − xm )2 .
2
Eq. B.1
Derivando em relação à (x − xm ) e igualando a zero (0)
f 0 (x) = f 0 (xm ) + f 00 (Ωxm )(x − xm ) = 0,
Eq. B.2
x = xm − [ f 00 (xm )]−1 f 0 (xm ),
Eq. B.3
obtém-se
faz-se novamente xm = x, até que seja atingido uma tolerância |xm − x| 6 tol.
Generalizando:
xm+1 = xm − δm ,
Eq. B.4
δm = [ f 00 (xm )]−1 f 0 (xm ) = Hm−1 gm ,
Eq. B.5
com
δm é denominado direção, sendo um vetor que representa um segmento do caminho do
ponto inicial até o ponto solução, com o inverso do “Hessiano” Hm determinando o ângulo e o
gradiente gm seu tamanho.
A função para o cálculo do Hessiano geralmente não está determinada, porém pode ser
estimada:
ANEXO B -- Algoritmos de Otimização
110
Considerando sm = xm+1 − x = δm e ηm = gm+1 − gm , então o próximo Hessiano pode ser
determinado através da relação Hm+1 sm = ηm que representa a razão da variação do gradiente
para a variação dos parâmetros e é denominada condição de “quasi-Newton”. Existem várias
formulações para a solução deste conjunto de equações.
Este tipo de otimização é adequado para sistemas cujo máximo/mínimo encontre-se relativamente próximo ao valor inicial com tendência de atingi-lo, de outra forma, poderá convergir
para um mínimo local não desejado se a função não for mono-tônica.
B.2
Algoritmo Genético
Segundo MARDLE e PASCOE (1999), representa uma classe de técnicas de otimização
evolutiva, utilizada em modelos não lineares complexos com vários mínimos locais. Segue um
conceito de evolução através da utilização de gerações estocásticas de populações de solução
obtidas através de uma função objetivo. Devido à sua natureza probabilística, a solução encontrada pode não representar o mínimo ou máximo global exato, mas em geral está muito
próximo.
Baseia-se no princípio Darwiniano de sobrevivência dos mais aptos. Uma população inicial
é criada, contendo um número pré-definido de indivíduos (soluções), cada qual representado por
um gene contendo as informações das variáveis, e o valor da função objetivo. Os indivíduos
mais aptos, com maiores valores máximos/mínimos da função objetivo, são selecionados para
serem a base da próxima geração, que ainda contará com um esquema de “mutação” randômica
visando garantir que hajam indivíduos capazes de buscar outra solução, caso a maioria deles,
estiver próxima a um máximo/mínimo local.
A técnica consiste em quatro estágios. A evolução que mede a aptidão (valor da função objetivo) de cada indivíduo solução. A seleção, que randomicamente seleciona os indivíduos da
população que servirão de base para a próxima geração. A reprodução (“crossover”), que toma
dois indivíduos e os combina para gerar dois novos indivíduos da próxima população e a mutação que randomicamente modifica os indivíduos através de alterações nos genes. O número de
indivíduos da população deve representar a complexidade do problema. O processo de geração
se repete até que atinja um número suficiente pré-determinado. A solução é representada pelo
elemento mais apto da última geração.
ANEXO B -- Algoritmos de Otimização
111
Figura B.1: Comparação de técnicas de otimização.
Uma combinação muito boa é a utilização do algoritmo genético para uma aproximação
segura ao mínimo/máximo global (Figura B.1), com o resultado sendo aplicado como condição
inicial para o método “quasi-Newton”, que garante a localização exata do ponto de inflexão.
112
ANEXO C -- MANCAIS DE ROLAMENTO
C.1
Rigidez
Segundo ERICH (2004), mancais de rolamento são utilizados em muitas classes de máquinas rotativas onde são requeridos pequeno tamanho e alta capacidade de carga da rotação zero
até a nominal. Devido à sua construção, nos rolamentos sempre há contato de entre as pistas e o
elemento de rolamento, não apresentam amortecimento interno intrínseco e devem ser projetados com folga radial. O deslocamento entre os centros das pistas é constituído por três parcelas
distintas (Figura C.1):
A folga radial δc (constante).
A compressão dos elementos girantes δh (depende da carga).
A ovalização da pista devido à distribuição de carga δo (depende da carga).
A rigidez é obtida da razão da carga aplicada Fa pelo deslocamento total
kr =
Fa
.
(δc + δh + δo )
Eq. C.1
Figura C.1: Deflexão devido compressão dos elementos.
Para um rolamento típico, (ERICH (2004)) apresenta um valor de rigidez de aproximadamente 300E6 (N/m), outro autor (PRASAD (2004)) apresenta valores da ordem de 1E9 (N/m).
113
ANEXO D -- FREQÜÊNCIA NATURAL DO NEUTRALIZADOR
D.1
Medição
Após a construção do neutralizador é necessário medir sua performance. Como este tipo de
sistema é bastante amortecido, a obtenção da função resposta em freqüência através de métodos
convencionais não é simples.
EWINS (2000) sugere a utilização de um “shaker” para excitar o neutralizador através
de um transdutor de força, medindo a aceleração resultante através de um acelerômetro. O
resultado obtido em msN−2 é denominado massa dinâmica.
Como excitação para o “shaker” foi utilizado a função “sweep”, uma senóide com variação
linear lenta de freqüência no tempo.
Figura D.1: Montagem para medição da massa aparente.
Como pode ser observado na figura do exemplo acima (D.1 1 ) a montagem é bastante simples, o pino suporte do neutralizador é fixado ao “shaker”(3) através de um transdutor de força
1A
WEG Indústrias S.A. Motores gentilmente cedeu laboratório e equipamentos para medição.
ANEXO D -- Freqüência Natural do Neutralizador
114
(2). Ainda sobre o pino suporte é fixado o acelerômetro (1) de medição.
2 )” obteve-se a
R
Utilizando as funções de controle e medição do software “Pulse(B&K
maior massa dinâmica na freqüência de aproximadamente 58 (Hz).
Figura D.2: Gráfico massa aparente obtido para um neutralizador de 150 (g).
Os parâmetros do material viscoelástico e do neutralizador podem ser obtidos na seção (4).
Para realização das medições foram utilizados os seguintes equipamentos:
•Four channel pulse data acquisition unit: B&K modelo 3560 C.
•Acelerômetro: B&K modelo 4395 Deltatron.
•Software de análise: Pulse Labshop B&K v9.
•Amplificador de potência B&K modelo Power Amplifier Type 2719 (180 VA).
•Shaker (atuador) B&K modelo PM Vibration Exciter Type 4808 (112 N).
Como resultado da medição, obteve-se a freqüência de 58 (Hz), muito próxima à freqüência
ótima projetada 60..2 (Hz). Observe-se a razão de amortecimento de 33%.
2 B&K
= Brüel & Kjaer.
115
ANEXO E -- TRANSDUTORES DE DESLOCAMENTO
E.1
“Proximitor”
Segundo MUSZYNSKA (2005), os transdutores de deslocamento (“proximitors”) são as
melhores ferramentas para o monitoramento contínuo do comportamento de rotores através da
medição sem contato de vibrações relativas à elementos estáticos. O princípio de funcionamento
é baseado na modificação do campo magnético devido à correntes induzidas nas proximidades
do sensor. A tensão de saída é linearmente proporcional à distância do transdutor e a superfície
do material observado. Tipicamente a sensibilidade se encontra próxima dos 8 ( mV
µm ) cobrindo
uma faixa de freqüências de 0 (zero) a 10 (kHz). Para que se tenha precisão na medição a
superfície observada requer tratamento superficial (polimento).
O instituto americano do petróleo (API - American Pretroleum Institute) recomenda a instalação de “proximitors” na configuração XY (Figura E.1) para monitoração da posição do centro
do eixo.
Figura E.1: Transdutores em configuração ortogonal XY.
116
ANEXO F -- MOVIMENTO ELÍPTICO
F.1
Definição
Segundo GENTA (2005), o ponto P (Figura F.1) se move no plano x − z e sua trajetória
pode ser circular, elíptica ou retilínea em qualquer direção dependendo apenas das condições
iniciais em t = 0. Pode ser interpretada como projeções dos vetores A e B que tem velocidade
angular Ωt. Note-se que a velocidade do ponto P não é constante, mas o período vale T =
2π
Ωt .
O movimento nas direções x e z é representado pelas expressões
x p (t) = A cos Ωt − ψx e z p (t) = B cos Ωt − ψz
Figura F.1: Movimento elíptico (órbita).
Eq. F.1
117
ANEXO G -- PARÂMETROS MODAIS
Para o modelo da Figura 4.5, são apresentados o primeiro modo de vibração, no intuito de
demonstrar a relevância da escolha da posição dos neutralizadores, de acordo com o o modo
que se deseja controlar, e as primeiras 10 linhas das matrizes modais truncadas obtidas para o
R
modelo simplificado. Ambos foram obtidos utilizando a versão 2004 em Scilabdo
do Rotor-
Din.
G.1
Modos de Vibrar
A posição do neutralizador é relevante e deve ser escolhida para ser efetiva aos modos que
se deseja controlar. Perto de “nós” modais neutralizadores não são efetivos. Para o sistema
da Figura 4.5 os neutralizadores foram montados a 380 (mm) da extremidade esquerda do eixo
(cota 0).
RotorDin − Modos de Vibrar
modo=2 f= 64.3(Hz) (FW)
0.15
z (m/1.00E+01)
0.10
0.05
0.00
−0.05
−0.10
−0.15
0.1
0.0
x (m/1.00E+01)
−0.1
0.0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1.0
y (m)
Figura G.1: Forma de vibrar do primeiro modo forward.
A posição escolhida é efetiva para o primeiro modo (Figura G.1).
ANEXO G -- Parâmetros Modais
G.1.1
Matriz [Φ] Complexa e Truncada 10 × 4
.
( 1.97E-03+i3.89E-03) ( 9.44E-04+i4.29E-03) ( 4.10E-03+i1.29E-03) (-3.84E-03-i1.97E-03)
( 3.89E-03-i1.97E-03) ( 4.29E-03-i9.44E-04) (-1.29E-03+i4.10E-03) ( 1.97E-03-i3.84E-03)
(-4.90E-02+i2.48E-02) (-5.40E-02+i1.19E-02) ( 1.62E-02-i5.14E-02) (-2.48E-02+i4.81E-02)
( 2.48E-02+i4.90E-02) ( 1.19E-02+i5.40E-02) ( 5.14E-02+i1.62E-02) (-4.81E-02-i2.48E-02)
( 9.82E-04+i1.94E-03) ( 4.70E-04+i2.14E-03) ( 2.04E-03+i6.42E-04) (-1.91E-03-i9.83E-04)
( 1.94E-03-i9.82E-04) ( 2.14E-03-i4.70E-04) (-6.42E-04+i2.04E-03) ( 9.83E-04-i1.91E-03)
(-4.87E-02+i2.47E-02) (-5.37E-02+i1.18E-02) ( 1.61E-02-i5.13E-02) (-2.47E-02+i4.80E-02)
( 2.47E-02+i4.87E-02) ( 1.18E-02+i5.37E-02) ( 5.13E-02+i1.61E-02) (-4.80E-02-i2.47E-02)
(-1.35E-06-i2.69E-06) (-6.63E-07-i2.96E-06) (-2.88E-06-i8.95E-07) ( 2.69E-06+i1.40E-06)
(-2.69E-06+i1.35E-06) (-2.96E-06+i6.63E-07) ( 8.94E-07-i2.88E-06) (-1.40E-06+i2.69E-06)
G.1.2
Matriz [Ψ] Complexa Truncada 10 × 4
.
( 3.88E-03+i1.97E-03) (-4.22E-03-i9.27E-04) (-1.30E-03-i4.11E-03) (-1.96E-03-i3.82E-03)
(-1.97E-03+i3.88E-03) ( 9.27E-04-i4.22E-03) (-4.11E-03+i1.30E-03) (-3.82E-03+i1.96E-03)
( 2.48E-02-i4.88E-02) (-1.17E-02+i5.31E-02) ( 5.16E-02-i1.63E-02) ( 4.79E-02-i2.46E-02)
( 4.88E-02+i2.48E-02) (-5.31E-02-i1.17E-02) (-1.63E-02-i5.16E-02) (-2.46E-02-i4.79E-02)
( 1.93E-03+i9.79E-04) (-2.10E-03-i4.61E-04) (-6.46E-04-i2.05E-03) (-9.77E-04-i1.90E-03)
(-9.79E-04+i1.93E-03) ( 4.61E-04-i2.10E-03) (-2.05E-03+i6.47E-04) (-1.90E-03+i9.77E-04)
( 2.46E-02-i4.85E-02) (-1.16E-02+i5.28E-02) ( 5.15E-02-i1.62E-02) ( 4.78E-02-i2.45E-02)
( 4.85E-02+i2.46E-02) (-5.28E-02-i1.16E-02) (-1.62E-02-i5.15E-02) (-2.45E-02-i4.78E-02)
(-2.67E-06-i1.37E-06) ( 2.91E-06+i6.28E-07) ( 9.22E-07+i2.89E-06) ( 1.37E-06+i2.69E-06)
( 1.37E-06-i2.67E-06) (-6.28E-07+i2.91E-06) ( 2.89E-06-i9.22E-07) ( 2.69E-06-i1.37E-06)
G.1.3
Matriz [Λ] Complexa Truncada 4 × 4
.
( 1.68E+00+i3.99E+02) (0.00E+00-i0.00E+00) (0.00E+00-i0.00E+00) (0.00E+00-i0.00E+00)
(0.00E+00-i0.00E+00) ( 1.68E+00-i3.99E+02) (0.00E+00-i0.00E+00) (0.00E+00-i0.00E+00)
(0.00E+00-i0.00E+00) (0.00E+00-i0.00E+00) ( 1.70E+00-i4.04E+02) (0.00E+00-i0.00E+00)
(0.00E+00-i0.00E+00) (0.00E+00-i0.00E+00) (0.00E+00-i0.00E+00) ( 1.70E+00+i4.04E+02)
Pode-se observar a esperada ocorrência de pares conjugados.
118
119
ANEXO H -- PARÂMETROS DO MANCAL FLUTUANTE
H.1
Obtenção
Para se obter o valor de amortecimento associado às molas que suportam o mancal flutuante onde são fixados os neutralizadores, foi utilizada a teoria da banda de média potência (ver
EWINS (2000)). Através da curva de resposta em freqüência obtida experimentalmente e mostrada na Figura H.1, é possível calcular o valor do amortecimento considerado constante por
hipótese, c = 18 ( Ns
m ).
Figura H.1: Função resposta medida experimentalmente.
A freqüência natural medida foi de Ω = 41 (Hz), para uma massa de m = 0.83 (kg). A
rigidez associada obtida foi de k = 55.1 (kN/m).
120
ANEXO I -- EQUIPAMENTOS DE MEDIÇÃO
I.1
Analisador de Sinais
Para medição no domínio da freqüência foi utilizado uma analisador de sinais da marca
HP (Hewlett-Packard) modelo 3560A (Figura I.1), portátil, dois canais, FFT1 , fonte de corrente
integrada. As principais características de desempenho são (http://www.metrictest.com):
•Faixa de freqüência: 0.03125 (Hz) a 40 (kHz)
•Precisão em freqüência: ± 8 (Hz)
•Precisão em amplitude: ± 0.2 dB
•Resolução de largura de banda: 0.03125 a 400 (Hz)
•Nível de ruído médio: -100 (dBv)
Figura I.1: Analisador HP 3560A.
1 Fast
Fourier Transform
ANEXO I -- Equipamentos de Medição
I.2
121
Martelo Instrumentado
Para o teste de resposta em freqüência em repouso foi utilizado um martelo instrumentado
R
fabricado pela PCB, modelo 086C04, “Modally TunedImpulse
Hammer w/force sensor” (Fi-
gura I.2). As principais características são:
•Sensibilidade: (±15%) 1.1 (mV/N)
•Faixa de medição: ±4400 (N pico)
•Massa do martelo: 0.16 (kg)
Figura I.2: Martelo PCB 086C04.
I.3
Acelerômetro
Para medição da resposta ao impulso foi utilizado um acelerômetro da marca PCB, modelo
R
352C68 “ceramic shear ICPaccel.
conn” (Figura I.3). Suas principais caracterísitcas são:
•Sensibilidade: (±10%) 10.2 (mV/(m/sš))
•Faixa de medição: ±491 (m/sš pico)
•Resolução de freqüência: (1 a 10000 Hz) 0.0015 (m/sš RMS)
•Faixa de freqüência: (±5%) 0.5 a 10000 (Hz)
•Massa: 2.0 (gramas)
ANEXO I -- Equipamentos de Medição
122
Figura I.3: Acelerômetro PCB 352C68.
I.4
Sistema de Aquisição de Sinal
O sistema de aquisição é baseado em “hardware” e “software” comercializado pela empresa
“National Instruments” (http://www.ni.com). O Hardware é composto por chassis (“frame”)
SCXI 1000 (Figura I.4), módulo de entrada para acelerômetros SCXI1531, módulo de cabeamento SCXI 1180, e uma placa de aquisição de sinais DAQCard 6062E. O Software utilizado é
R
o LabViewna
versão 8.2.
Figura I.4: Acelerômetro PCB 352C68.
As principais características são:
•SCXI1531 (Figura I.5):
ANEXO I -- Equipamentos de Medição
–Oito canais.
–Ganho programado por canal.
–Fonte de corrente 4 (mA), 20 (V ).
–Filtro passa baixas.
Figura I.5: Placa de entrada SCXI 1531.
•DAQCard 6062E (Figura I.6):
–16 entradas analógicas.
∗Resolução 12 bits.
∗Taxa de amostragem 500kS/s(2 ).
∗Faixa de entrada de tensão ±0.05 − ±10 (V ).
–Duas saídas analógicas.
∗Taxa de amostragem 850kS/s.
∗Faixa de tensão ±10 (V )
–8 linhas digitais (TTL/CMOS)
2 500
mil amostragens por segundo.
123
ANEXO I -- Equipamentos de Medição
124
Figura I.6: Placa de aquisição DAQCArd 6062E.
I.5
Proximitor
Transdutor de deslocamento (sensor de proximidade) (Figura I.7), tipo indutivo, marca Vi-
brocontrol, modelo PS 1002 e condicionador de sinal (Figura I.8), modelo GS 5001. Segundo
o fabricante, a sensibilidade para o material Aço SAE 8640 é de 8 (V /mm) e a faixa linear vai
de 0, 5 a 2, 5 (mm) à uma temperatura ambiente de 27 (◦C).
Figura I.7: Sensor PS 1002.
Figura I.8: Condicionador GS 501.
125
ANEXO J -- ROTORDIN
J.1
Introdução:
Em 2002 o aluno Ezequiel José López iniciou o desenvolvimento de um código numérico
R
para o cálculo da dinâmica de rotores em Matab(LÓPEZ
(2002.)). Em 2004 o engenheiro
Hideraldo L. V. dos Santos aperfeiçoou a interface gráfica para permitir uma entrada de dados
confortável e geração de gráficos de saída especializados e incluiu a opção de cálculo com
mancais hidrodinâmicos (SANTOS (2003)). Através do Termo de cooperação 01/2004 - WEG
Indústrias Elétricas e CEFET-PR e a UTFPR foi desenvolvida uma nova versão do código, como
parte do projeto de pesquisa cujos objetivos são desenvolver modelos numéricos utilizando o
Método dos Elementos Finitos, para predizer o comportamento dinâmico de rotores (compostos
por eixo escalonados, discos e diferentes modelos de mancais hidrodinâmicos e rolamentos). O
código numérico foi programado conforme a teoria de sistemas girantes apresentada na seção
2.7, denominado de “RotorDin”.
R
A partir da versão do código de 2004 foi desenvolvida uma nova versão, utilizando Scilab,
R Esta versão foi utilizada neste trabalho para modelauma versão “opensource” do Matlab.
gem da bancada de dinâmica de rotores do LAVIB, com o objetivo de exportar as matrizes de
parâmetros modais em forma de arquivo texto e permitir a inclusão do modelo de parâmetros
equivalentes generalizados dos neutralizadores dinâmicos viscoelásticos diretamente nas matrizes de massa e amortecimento. Foram implementadas as opções de cálculo do diagrama de
Campbell, a resposta ao desbalanceamento e a função resposta em freqüência. Isto permitiu
comparar os resultados numéricos obtidos com a metodologia proposta utilizando o modelo
simplificado.
ANEXO J -- RotorDin
J.2
126
Capacidades
R
As principais capacidades de cálculo e opções disponíveis no RotorDin (versão Matlab,
2004) (1 ) são:
•Entrada de dados gráfica facilitada, no padrão windowsTM , conforme mostrado na Figura
J.1. Opções de cálculo também podem ser acessadas da tela principal.
Figura J.1: Tela principal.
•Visualização da geometria, do modelo de elementos finitos e do posicionamento dos mancais para orientação ao usuário quanto à validade dos dados de entrada, Figura J.2.
1A
versão atual possui opções revistas e adicionais não mencionadas.
ANEXO J -- RotorDin
127
ESQUEMA DO ROTOR
Figura J.2: Geometria e modelo.
•Cálculo do diagrama de Campbell e visualização da estabilidade do sistema, Figura J.3.
A estabilidade é verificada pela cor dos marcadores nas linhas, azul estável, vermelho
instável.
DIAGRAMA DE CAMPBELL − MANCAL VISCOSO
70
FREQUENCIA NATURAL [Hz]
60
50
40
30
20
10
0
0
500
1000
1500
2000
2500
ROTAÇAO [rpm]
3000
3500
4000
Figura J.3: Diagrama de Campbell.
•Cálculo da resposta ao desbalanceamento em módulo e fase, conforme a Figura J.4.
ANEXO J -- RotorDin
128
RESPOSTA AO DESBALANCEAMENTO
−3
10
−4
10
−5
10
| Amp(25) |
−6
10
−7
10
−8
10
−9
10
−10
10
0
500
1000
1500
2000
2500
ROTAÇAO [rpm]
3000
3500
4000
Figura J.4: Resposta ao desbalanceamento.
•Cálculo do modos de vibrar, com orientação quanto à direção de precessão, para uma ou
mais posições e direções, Figura J.5.
MODO 4 − f = 33.6529 Hz − FW
−3
x 10
4
3
2
Z [m]
1
0
−1
−2
−3
−4
5
1.5
−3
x 10
0
1
0.5
X [m]
−5
0
Y [m]
Figura J.5: Modo 1, direção “forward”.
Observação: As figuras apresentadas acima são ilustrativas, utilizadas para demonstrar as
saídas gráficas e podem não ter relação entre si.