Dynamic analysis in the domain of frequency by the method of the harmonic balance

==Resumo==

A análise de estruturas no domínio da frequência vem ganhando campo nas últimas décadas, primeiramente, em virtude do desenvolvimento computacional e estratégias para avaliar a Transformada Discreta de Fourier com mais eficiência. As vantagens em se avaliar a resposta no domínio da frequência é que não há necessidade de conhecer as condições iniciais, e pode se avaliar a resposta para sistemas com propriedades dependentes da frequência com maior precisão, por exemplo, interação solo-estrutura. Este trabalho tem como objetivo avaliar a resposta dinâmica de um sistema estrutural discreto no domínio da frequência mediante as curvas de ressonância por meio do Método do Balanço Harmônico HBM. A resposta é avaliada na fase permanente para ações e/ou cargas harmônicas para sistemas dotados de amortecimento viscoso tipo proporcional e não proporcional. As propriedades do sistema, como massa e rigidez, são tomadas como constantes. Para os modelos desenvolvidos e analisados neste estudo são apresentados os respectivos espectros de resposta em frequência e histórico de deslocamento no tempo por meio do HBM. As respostas do sistema obtidas por intermédio do HBM são comparadas com outros dois métodos também baseados no domínio da frequência (ImFT e Pseudo-Forças), permitindo, então, avaliar as amplitudes e frequências críticas do sistema.''

'''Palavras-chave''': Análise dinâmica, domínio da frequência, espectro de resposta

==Abstract==

The analysis of structures in the frequency domain has been gaining ground in the last decades, first, because of the computational development and strategies to evaluate the Discrete Fourier Transform more efficiently. The advantages of evaluating the frequency response are that there is no need to know the initial conditions, and the response to systems with frequency-dependent properties can be evaluated more accurately, eg, soil-structure interaction. This work aims to evaluate the dynamic response of a discrete structural system in the frequency domain through the resonance curves using the Harmonic Balance Method HBM. The response is evaluated in the permanent phase for actions and / or harmonic loads for systems with proportional and non-proportional viscous damping. The properties of the system, such as mass and rigidity, are taken as constants. For the models developed and analyzed in this study are presented the respective spectra of frequency response and history of displacement in time through the HBM. The system responses obtained through the HBM are compared with two other methods also based on the frequency domain (ImFT and PseudoForces), allowing to evaluate the amplitudes and critical frequencies of the system.

'''Keywords''': Dynamic analysis, response structural, domain frequency, analysis harmonics of structures

==1. Introdução==

A análise de estruturas no domínio da frequência vem ganhando campo nas últimas décadas, primeiramente, em virtude do desenvolvimento computacional. A resposta no domínio da frequência permite obter um espectro de resposta associado, neste caso, às propriedades do sistema, como massa, rigidez, amortecimento e excitação dinâmica. Em alguns casos a resposta em frequência, que é caracterizada pela curva de ressonância, pode ser obtida indiretamente por meio da análise no tempo Galvão [7] e Silva [13], o que leva a um grande esforço computacional para estruturas com vários graus de liberdade.

Estendendo a formulação para o domínio da frequência, é possível avaliar quais são as frequências críticas e suas respectivas amplitudes, porém se sabe que esta amplitude está associada à resposta na fase permanente, dependendo da excitação dinâmica e do modelo de amortecimento adotado. Na Figura 1 é apresentado um sistema dinâmico discreto qualquer e suas respectivas propriedades, também a resposta no domínio frequência e no domínio do tempo. 

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image1.png|528px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding:10px;"| '''Figura 1'''. Exemplo: Histórico de resposta (tempo) e espectro de resposta (Freq.), SDOF
|}


Entre as soluções possíveis ou alternativas para obter a resposta em frequência, temos os métodos baseados na superposição modal e os métodos baseados em linearização equivalente, ambos abordados de uma forma implícita; uma abordagem diferenciada, porém não nova, que consiste em fazer uso do Método do Balanço Harmônico ou Balanço Harmônico Linearizado para obter a resposta no domínio da frequência.

A transformação da equação de equilíbrio da análise dinâmica por meio de elementos finitos para o domínio da frequência, é feita pelo método do balanço harmônico clássico, seguindo a formulação apresentada por Pasquetti [9].

Após a transformação e mudança de variável temporal, a solução da equação de equilíbrio implica duas soluções possíveis: via Newton-Raphson e/ou incremento de frequência e linearização Bernárdez [1] e Childs [3].

==2. Resposta em frequência==

A resposta em frequência está associada a sistemas lineares invariantes no tempo sob excitação harmônica ou periódica, considerando sua resposta no regime permanente. A importância da resposta em frequência consiste no fato de que sinais periódicos podem ser decompostos ou representados por uma combinação de harmônicos.

Os métodos de análise baseados na resposta em frequência são, talvez, os mais utilizados em ambientes industriais, pois permitem avaliar a resposta de um sistema na presença de incertezas em razão da utilização e/ou operação.

A avalição do comportamento dinâmico de estruturas no domínio da frequência dependerá basicamente das propriedades ortogonais dos modos. O conceito de modos normais em sistemas estruturais dinâmicos é associado às frequências naturais (ou fundamentais) de um sistema linear, no qual a existência destes modos possibilita a introdução de coordenadas normais que proporcionam o desacoplamento do sistema. Nesse contexto, podemos citar os trabalhos desenvolvidos por Claret et al. [4] e Ribeiro [10].

Tendo em vista que algumas propriedades são não ortogonalizáveis, tal como o amortecimento não proporcional dos dependentes da frequência, faz-se necessário, então, um método no qual sua resposta seja avaliada no domínio da frequência de forma implícita, ou seja, a resposta e dada por meio de operações matriciais que envolvam a linearização do sistema, obtendo a resposta em coordenadas físicas, como o método das Pseudo-Forças de Ferreira [6] e Soares [14], e o método do balanço harmônico, linearizado por Muñoz et al. [8].

Logo, a resposta do sistema no domínio da frequência é definida em relação aos termos associados à amplitude <math display="inline">U</math> e à fase <math display="inline">\phi</math>. Observa-se que nos métodos de redução da ordem da equação do equilíbrio dinâmico e da linearização do sistema, por mais que se dê a mudança da variável temporal, a resposta do sistema é dada no domínio do tempo. Para a proposta de obter a resposta tanto no domínio da frequência quanto no tempo para sistemas dotados de amortecimento proporcional e/ou não proporcional, faz-se uso do Método do Balanço Harmônico por meio de duas soluções possíveis: via Newton-Raphson e/ou solução incremental de frequência <math display="inline">\omega</math>  e solução linear em termos de amplitudes modais <math display="inline">\left( A,B\right)</math>.

===2.1 Método do balanço harmônico===

Umas das técnicas geralmente utilizadas para obter soluções aproximadas para equações diferenciais é o Método do Balanço Harmônico. Este método permite gerar soluções periódicas aproximadas a partir de uma série de Fourier truncada Pasquetti [9].

O método do balanço harmônico consiste na aproximação dos deslocamentos na fase permanente por meio da combinação de harmônicos, tal como dado pela expressão:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>u\left(t\right)=\sum_{n=0}^N\left[A_nsin\left(n\omega t\right)+\right. </math><math>\left. B_ncos\left(n\omega t\right)\right]</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (1)
|}

onde <math display="inline">\, {A}_{n}</math> e <math display="inline">{B}_{n}</math> são os vetores das componentes das amplitudes modais no enésimo harmônico e <math display="inline">N</math> é o número de harmônicos usados na aproximação, tendo a equação de equilíbrio dinâmica dada pela equação (2)

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>M\mbox{ }\ddot{u}\left(t\right)+C\mbox{ }\dot{u}\left(t\right)+</math><math>Ku\left(t\right)=p</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (2)
|}

Submetendo a estrutura a uma excitação harmônica com uma frequência de excitação <math display="inline">\omega</math>  que multiplica o vetor de forças de referência <math display="inline">p</math>, o equilíbrio dinâmico é dado por:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>M\mbox{ }\ddot{u}\left(t\right)+C\mbox{ }\dot{u}\left(t\right)+</math><math>Ku\left(t\right)=Fsin\left(\omega \mbox{ }t\right)\mbox{ }p</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (3)
|}

onde <math display="inline">F</math> e <math display="inline">\omega</math>  são a amplitude e a frequência da excitação harmônica.

Com base nos trabalhos publicados por Muñoz et al. [8] e Pasquetti [9], sendo o caso de uma análise linear a equação (3), se faz necessário apenas o segundo termo do somatório.

O vetor de deslocamento, portanto, é aproximado por:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>u\left(t\right)=A_1sin\left(\omega \mbox{ }t\right)+</math><math>B_1cos\left(\omega \mbox{ }t\right)</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (4)
|}

Substituindo a equação (4) e suas respectivas derivadas na equação (3) obtém-se:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>M\left(-A_1sin{\left(wt\right)}^2-B_1cos{\left(wt\right)}^2\right)+</math><math>C\left(A_1cos\left(wt\right)w-B_1sin\left(wt\right)w\right)+</math><math>K\left(A_1sin\left(wt\right)+B_1cos\left(wt\right)\right)-</math><math>Fsin\left(wt\right)p=0</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (5)
|}

Simplificando a equação (5) tem-se:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\left(-A_1Mw^2-B_1Cw+A_1K-Fp\right)sin\left(wt\right)+</math><math>\left(-B_1Mw^2-A_1Cw+B_1K\right)cos\left(wt\right)=</math><math>0</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (6)
|}

Considerando <math display="inline">Fp=\, {F}_{0}</math>, a equação (6) pode ser reescrita como um sistema de equações que iguala os coeficientes das funções harmônicas, como:

Os componentes trigonométricos de cada harmônico podem ser reduzidos a:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\left[\begin{array}{cc}
\omega \mbox{ }C & -{\omega }^2M+K\\
-{\omega }^2M+K & -\omega \mbox{ }C
\end{array}\right]\left\{\begin{array}{c}
A_1\\
B_1
\end{array}\right\}=\left\{\begin{array}{c}
F_0\\
0
\end{array}\right\}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (7)
|}

A equação (6) pode ser escrita de uma forma mais compacta, tal como:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overline{K}\left(\omega \right)\cdot \overline{A}=</math><math>\overline{F}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (8)
|}

onde:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overline{K}\left(\omega \right)=\left[\begin{array}{cc}
\omega \mbox{ }C & -{\omega }^2M+K\\
-{\omega }^2M+K & -\omega \mbox{ }C
\end{array}\right]\mbox{ },\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\overline{A}=</math><math>\left\{\begin{array}{c}
A_1\\
B_1
\end{array}\right\},\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\overline{F}=</math><math>\left\{\begin{array}{c}
F_0\\
0
\end{array}\right\}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (9)
|}

sendo a frequência de excitação <math display="inline">\omega</math>  e o vetor de amplitudes modais <math display="inline">\overline{A}</math> as variáveis a serem determinadas.

Muñoz et al. [8] salientam que, apesar de uma análise linear, o sistema de equações é não linear, portanto é necessário o emprego de um método de solução ou de controle para solução de sistemas de equações não lineares.

====2.1.1 Solução via Newton-Raphson====

Como observado, temos a solução da equação de equilíbrio dinâmico reduzida à solução de um sistema de equações lineares; este, porém, apresenta um número de incógnitas superior às equações de equilíbrio fornecidas, recorrendo ao método de Galerkin para ortogonalizar o resíduo no intervalo de <math display="inline">\left[ 0,\frac{2\pi }{\omega }\right]</math>  e obter a solução das equações algébricas.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\begin{array}{c}
\int_0^{\frac{2\pi }{\omega }}R\left(t\right)sinwt\mbox{ }dt=-F-CB_1\omega +KA_1-MA_1{\omega }^2=0\\
\\
\int_0^{\frac{2\pi }{\omega }}R\left(t\right)coswt\mbox{ }dt=-MB_1{\omega }^2+CA_1\omega +KB_1=0
\end{array}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (10)
|}

O sistema resultante <math display="inline">G\left( A,B,\omega \right) =0</math>, apresenta duas equações e duas variáveis dependentes <math display="inline">A</math> e <math display="inline">B</math> com uma variável independente <math display="inline">\omega</math>.

A solução recai em um dado valor de <math display="inline">\omega</math>. Podemos, então, aplicar Newton-Raphson para um par de valores iniciais <math display="inline">{\left( {A}^{0},{B}^{0}\right) }^{T}</math>, convergindo para uma solução aproximada da resposta permanente do sistema,

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\left(A^{k+1},B^{k+1}\right)}^T={\left(A^k,B^k\right)}^T-</math><math>G_{A,B}^{-1}\left(A^k,B^k\right)G\left(A^k,B^k\right)</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (11)
|}

onde a matriz Jacobiana <math display="inline">{G}_{A,B}</math> será:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>G_{A,B}\left(A,B\right)=\left(\begin{array}{cc}
\frac{\partial G_1}{\partial A}\left(A,B\right) & \frac{\partial G_1}{\partial B}\left(A,B\right)\\
\frac{\partial G_2}{\partial A}\left(A,B\right) & \frac{\partial G_2}{\partial B}\left(A,B\right)
\end{array}\right)</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (12)
|}

A amplitude da resposta harmônica é derivada dos coeficientes <math display="inline">\left( A,B\right)</math> 

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\vert u\vert =\sqrt{A^2+B^2}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (13)
|}

enquanto o ângulo fase entre a resposta e a excitação é uma função de <math display="inline">\mathrm{{tan}^{-1}}\,\left( \frac{B}{A}\right)</math>  para um domínio <math display="inline">\left[ 0,\pi \right]</math>

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\begin{array}{c}
A>0\rightarrow \phi =-{tan}^1\left(\frac{B}{A}\right)\\
\\
A<0\rightarrow \phi ={tan}^1\left(\frac{B}{A}\right)-\pi 
\end{array}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (14)
|}

O sinal de menos do ângulo de fase em todo o intervalo indicará um atraso da resposta em relação à excitação.

====2.1.2 Solução incremental e linear====

Como observa-se na equação (9), e em sua forma reduzida na equação (8), uma das soluções possíveis consiste basicamente em realizar incrementos de frequências <math display="inline">\omega</math>. Neste caso a variável independente, onde o vetor de amplitudes modais <math display="inline">\overline{A}</math> é determinado pela inversão da matriz <math display="inline">\overline{K}\left( \omega \right)</math>, sendo, então, <math display="inline">\omega</math>  já estabelecido dentro do ciclo interativo, logo.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>\overline{A}=\overline{K}{\left(\omega \right)}^{-1}\overline{F}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (15)
|}

A curva de ressonância é obtida a partir do vetor de amplitudes modais <math display="inline">\overline{A}</math>, sendo este avaliado a cada incremento de frequência <math display="inline">\omega</math>. O vetor de amplitudes modais é avaliado pela sua normal.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\vert u\vert }_{m\acute{a}x}=\sqrt{A_1{}^2+B_1{}^2}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (16)
|}

É a solução para cada incremento de <math display="inline">\omega</math>, e o mesmo também assume posições distintas dentro do vetor de amplitudes modais <math display="inline">\overline{A}</math>, por exemplo. Este depende das dimensões das matrizes envolvidas, logo:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\vert u\vert }_{\left(i,j\right)max}=\sqrt{A_{i,j}{}^2+B_{2+i,j}{}^2}</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (17)
|}

sendo <math display="inline">i=1\ldots n</math>, para os <math display="inline">n</math> respectivos graus de liberdade, e <math display="inline">j=</math><math>0\ldots {\omega }_{n}</math> avaliado para cada incremento de frequência.

Podemos ainda associar a resposta com o ângulo fase ou simplesmente fase da resposta do espectro, a partir de uma frequência discreta <math display="inline">\omega</math>, como mostrado na equação (14) e/ou na forma matricial.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\phi }_{\left(i,j\right)}={tan}^{-1}\left(\frac{A_{i,j}}{B_{2+i,j}}\right)</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (18)
|}

A fase no método do balanço harmônico consiste em uma parte real <math display="inline">{A}_{i,j}</math> e a parte imaginária de <math display="inline">{B}_{2+i,j}</math>. É possível avaliar então a relação entre amplitude <math display="inline">\left| u\right|</math>  e a fase <math display="inline">\phi</math>  em relação à frequência <math display="inline">\omega</math>, como mostra a equação (18).

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>{\phi }_{\left(i,j\right)}={tan}^{-1}\left(\frac{Re\left[A_{i,j}\right]}{Im\left[B_{2+i,j}\right]}\right)</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (19)
|}

Tem-se, então, que a mudança na amplitude e na fase não é fixa, mas, sim, dependente das características dinâmicas do sistema e da frequência de entrada <math display="inline">\omega</math>.

Para o caso de análise linear, a solução incremental e linear, baseada na inversão da matriz <math display="inline">{\overline{K}}_{\omega }</math> reduz significativamente o tempo de processamento e também dispensa a necessidade de avaliar os valores iniciais <math display="inline">\left( {A}^{0},{B}^{0}\right)</math>. As duas soluções apresentadas foram implementadas gerando, respectivamente, a mesma resposta na frequência. Verifica-se que, caso não se realize o incremento de frequência <math display="inline">\omega</math>, ambas as soluções são indeterminadas.

==3. Matriz de amortecimento==

Em sistemas em que o amortecimento é considerado distribuído pela estrutura de forma semelhante à massa e à rigidez, ele é dito proporcional ou de Rayleigh. Sistemas com amortecimento proporcional utilizam técnicas de transformação modal para desacoplar as equações de equilíbrio dinâmicas. Assim, a resposta em coordenadas físicas pode ser obtida pela superposição da resposta em coordenadas modais, tanto no domínio do tempo quanto no da frequência Clough et al. [5].

Já em sistemas com características mais complexas de amortecimento, a matriz de amortecimento não é mais proporcional à massa e rigidez do sistema. Diz-se, então, que o amortecimento é não proporcional.

O amortecimento não proporcional é representado por uma matriz que não é ortogonal à matriz modal. Uma das técnicas utilizadas para desacoplar as equações de equilíbrio dinâmica, como dado por Rubinstein et al. [11], consiste em utilizar modos complexos, porém o problema de autovalor e autovetor para este caso é duas vezes maior que o convencional, e os modos complexos dificultam a interpretação física da resposta do sistema.

Conforme Santolin [12], pode-se resolver os sistemas dotados de amortecimento não proporcional por meio do método das Pseudo-Forças. Neste método, os termos da matriz de amortecimento responsáveis pelo acoplamento das equações são transferidos para o segundo membro da equação e tratados como Pseudo-Forças. Assim, o primeiro termo permanece desacoplado, como ocorre no sistema de amortecimento proporcional.

==4. Resposta na fase permanente==

A resposta associada ao sistema na fase permanente nos remete ao modelo de excitação dinâmica a ser adotado ou combinado. Em virtude do uso do método do balanço harmônico, para cada incremento de frequência <math display="inline">\omega</math>  interessa obter a respectiva amplitude correspondente à excitação para alcançar a solução do sistema em razão das amplitudes modais.

Uma forma de corresponder à excitação dinâmica dada no domínio do tempo, é por intermédio da periodização do carregamento. O método a ser empregado aqui é o mesmo utilizado por Camargo [2] e Ribeiro [10] na formulação da ImFT, que consiste em fazer uso da transformada discreta de Fourier e periodizar o vetor de cargas de referência. Os parâmetros, associados à excitação no domínio da frequência, necessitam, por exemplo, do tempo estendido <math display="inline">{T}_{0}</math> e o número de termos <math display="inline">N</math>, onde se deseja avaliar a resposta em frequência. O vetor de incremento de frequências <math display="inline">\omega</math>  será, então, associado a uma respectiva amplitude dada pela <math display="inline">FFT\left( p\right)</math>, quando o vetor de amplitudes modais <math display="inline">\overline{A}</math> apresentam a solução do sistema para o equilíbrio dinâmico. A excitação adotada neste trabalho é uma excitação harmônica, dada por uma combinação linear de seno e cosseno com a mesma frequência, tal como <math display="inline">p\left( t\right) =</math><math>b\mathrm{sin}\,\left( \omega t\right)</math>.

==5.  Exemplos==

A seguir são avaliadas as respostas de sistemas estruturais discretos, representados por suas respectivas matrizes de massa, rigidez e amortecimento, estando o sistema excitado por um carregamento e/ou ação harmônica. São apresentados três modelos para análise, variando entre si a matriz de amortecimento e, então, comparando-se a resposta no tempo com os métodos baseados no domínio da frequência supracitados.

Os espectros de respostas foram obtidos pelo método do balanço harmônico, assim como o histórico de deslocamento para ambos os modelos. O critério de comparação do histórico de deslocamento com os demais métodos depende, respectivamente, do tipo de amortecimento, e optou-se para sistemas dotados de amortecimento proporcional comparar a resposta com a ImFT Ribeiro [10], e para sistemas dotados de amortecimento não proporcional com a Pseudo-Forças Santolin [12].

===5.1 Modelo 1===

O sistema estrutural discreto apresentado na Figura 2 é tratado como Shear Building, sendo analisado por Clough et al. [5] e Santolin [12], e será submetido a um carregamento harmônico a fim de obter e verificar quais são as frequências e as amplitudes críticas.

O sistema apresenta os seguintes coeficientes para a matriz de massa <math display="inline">M</math> e rigidez <math display="inline">K</math>:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>M=\left[\begin{array}{ccc}
1,0 & 0 & 0\\
0 & 1,5 & 0\\
0 & 0 & 2,0
\end{array}\right]\mbox{ }ton,\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }\mbox{ }K=</math><math>\left[\begin{array}{ccc}
600,0 & -600,0 & 0\\
-600,0 & 1800,0 & -1200,0\\
0 & -1200,0 & 3000,0
\end{array}\right]kN/m</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (20)
|}

Na matriz de amortecimento <math display="inline">C</math> será assumido um amortecimento de <math display="inline">5\%</math> para os modos 1 e 3, resultando, desta forma, em uma matriz de amortecimento proporcional em coordenadas físicas

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math display="inline">C=\left[\begin{array}{ccc}
2,09 & -0,99 & 0\\
-0,99 & 4,63 & -1,98\\
0 & -1,98 & 7,16
\end{array}\right]kN\cdot s/m</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (21)
|}

No vetor de forças de referência adota-se um carregamento harmônico, dado

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math display="inline">p\left(t\right)=100sin\left(\omega t\right)</math>
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (22)
|}

com incidência sobre o grau de liberdade <math display="inline">{u}_{1}</math>, como mostra a Figura 2.

Serão desenvolvidas neste Modelo 1 duas análises: a primeira com <math display="inline">p\left( t\right)</math>  assumindo <math display="inline">\omega =</math><math>{\omega }_{crit}</math>, e uma segunda análise com base no trabalho publicado por Santolin [12], assumindo <math display="inline">p\left( t\right) =</math><math>100\mathrm{sin}\,\left( \omega t\right)</math>  com <math display="inline">\omega =</math><math>5,0</math>, atuando de modo pertinente a excitar a estrutura.

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
<!--|style="padding:10px;"|[[Image:draft_Zaltron_287705376-image24.png|498px]]-->
|style="padding:10px;"|[[File:Zaltron_Neves_2017a_8229_Fig2.png]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 2'''. Modelo estrutural discreto 1
|}


Já os parâmetros associados à discretização do tempo e frequência de amostragem são: intervalo de tempo <math display="inline">\Delta t=</math><math>0,01</math> segundos, tempo de análise de <math display="inline">3,0</math> segundos. Para obter o espectro de resposta, os parâmetros foram definidos <math display="inline">N=</math><math>1024</math> e a frequência de amostragem <math display="inline">{F}_{s}=100</math>. As respectivas frequências naturais foram obtidas a partir da solução do problema de autovalor e autovetor, pois será utilizada para comparar com as obtidas pelo espectro de resposta, onde <math display="inline">{\omega }_{1}=</math><math>14,52\, rad/s</math>, <math display="inline">{\omega }_{2}=31,05\, rad/s</math> e <math display="inline">{\omega }_{3}=</math><math>46,09\, rad/s</math>.

Na Figura 3 temos a resposta no domínio da frequência e/ou espectro de resposta via HBM. Como é possível verificar, as amplitudes críticas coincidem respectivamente com as frequências naturais do sistema.

Já para a Figura 4 temos o histórico de deslocamento avaliado pela ImFT, também se comparou a resposta no tempo por algoritmo de integração direta via método de Newmark. O propósito desta avaliação e garantir que as soluções desenvolvidas venham ao convergir para a mesma resposta, sendo os resultados satisfatórios.

Verifica-se que as amplitudes obtidas no espectro de resposta, refere-se a amplitude da resposta na fase permanente do sistema, isto quando assume-se <math display="inline">\omega =</math><math>{\omega }_{n}</math>, caso contrário a resposta deverá ser avaliada a partir das amplitudes modais.

Logo, a partir da Figura 5 temos então a apresentação das respostas no domínio da frequência e no tempo para os modelos adaptados de Clough et al. [5] e Santolin [12].

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image25.png|400px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 3'''. Espectro de resposta para o grau de liberdade <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math> 
|}


{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"| [[Image:draft_Zaltron_287705376-image26.png|500px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 4'''. Histórico de resposta (tempo) para <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math> via ImFT com <math display="inline">\mathit{\boldsymbol{\omega }}=</math><math>{\mathit{\boldsymbol{\omega }}}_{\mathit{\boldsymbol{crit}}}</math> 
|}


Observa-se nas Figuras 3 e 4 que as amplitudes são idênticas quando assume-se que <math display="inline">\omega =</math><math>{\omega }_{n}</math>, tendo, assim, as respectivas amplitudes críticas <math display="inline">\left| {u}_{crit}\right|</math>, o que mostra a convergências de valores para ambos os métodos.

A resposta a um carregamento harmônico monta o vetor de forças de referência assumindo <math display="inline">p\left( t\right) =</math><math>100\mathrm{sin}\,\left( \omega t\right)</math>, tal como apresentado por Santolin [12] via ImFT, atuando de modo pertinente a excitar a estrutura, ou seja, agindo durante todo o período definido pela análise, logo:

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image27.png|400px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 5'''. Espectro de resposta via HBM para <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math>, <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{2}}}</math> e <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{3}}}</math> 
|}


A Figura 6 mostra o histórico de resposta para o grau de liberdade <math display="inline">{u}_{1}</math>, comparando-o com o ImFT e HBM para um sistema dotado de amortecimento do tipo proporcional.

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image28.png|500px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 6'''. Histórico de resposta (tempo) para <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math> via ImFTe HBM
|}

===5.2 Modelo 2===

O modelo 2 consiste em um sistema com as mesmas propriedades atribuídas ao modelo 1, porém sendo acrescentado um amortecedor discreto ao primeiro grau de liberdade <math display="inline">{u}_{1}</math>, onde, <math display="inline">{c}_{1}=</math><math>20,0\, kN\cdot s/m</math>, como observa-se na Figura 7. Desta forma, o problema não pode ser mais avaliado utilizando técnicas de superposição modal. A solução agora fica basicamente restringida ao uso do método do balanço harmônico e da Pseudo-Forças.

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image29.png|600px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 7'''. Modelo estrutural discreto 2
|}


As matrizes de massa e rigidez são as mesmas apresentadas no modelo 1; logo a matriz de amortecimento, levando em consideração a inclusão de um amortecedor discreto, é dado por:

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>C=\left[\begin{array}{ccc}
22,09 & -0,99 & 0\\
-0,99 & 4,63 & -1,98\\
0 & -1,98 & 7,16
\end{array}\right]kN\cdot s/m</math> .
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (23)
|}


A inclusão do amortecedor discreto faz com que se some o valor de <math display="inline">{c}_{1}</math> ao termo <math display="inline">{C}_{1,1}</math> da matriz de amortecimento.

Na Figura 8, é possível observar uma mudança no espectro de resposta devido a inclusão deste amortecedor junto ao ponto de excitação, fazendo com que as amplitudes sejam menores que as apresentadas no modelo 1. Também é possível verificar que as faixas de frequências que podem ser assumidas críticas para o sistema estão à esquerda de <math display="inline">{\omega }_{1}</math> ou seja a primeira frequência natural. Já para os demais graus de liberdade esta condição não é aplicável, pois é facilmente verificada junto à <math display="inline">{\omega }_{2}</math> e <math display="inline">{\omega }_{3}</math> as amplitudes apresentam picos, assumindo uma amplitude crítica.

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image31.png|400px]] 
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 8'''. Espectro de resposta via HBM para <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math>, <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{2}}}</math> e <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{3}}}</math>
|}


Na Figura 8 temos a resposta no domínio da frequência e/ou espectro de resposta via HBM para um sistema dotado de amortecimento do tipo não proporcional.

A Figura 9 mostra a resposta no tempo para o grau de liberdade <math display="inline">{u}_{1}</math>, comparando-o com a Pseudo-Forças e HBM para um sistema dotado de amortecimento do tipo não proporcional.

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"| [[Image:draft_Zaltron_287705376-image32.png|500px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 9'''. Histórico de resposta (tempo) para <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math> via Pseudo-Forças e HBM
|}

===5.3 Modelo 3===

No modelo 3 o grau de liberdade <math display="inline">{u}_{1}</math> segue como no modelo 1, porém são adicionados amortecedores discretos <math display="inline">{c}_{2}</math> e <math display="inline">{c}_{3}</math>, como observa-se na Figura 10, onde se pretende avaliar a influência que o amortecedor pode trazer em relação à amplitude máxima avaliada pelo espectro de resposta.

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image33.png|600px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 10'''. Modelo estrutural discreto 3
|}


Os amortecedores discretos <math display="inline">{c}_{2}=10\, kN\cdot s/m</math> e <math display="inline">{c}_{3}=</math><math>20\, kN\cdot s/m</math>, sendo a matriz de amortecimento do sistema, apresentada junto a expressão 24.

{| class="formulaSCP" style="width: 100%; text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;" 
|-
| <math>C=\left[\begin{array}{ccc}
3,1 & 3,1 & 0\\
-3,1 & 18,9 & -5,8\\
0 & -5,8 & 36,4
\end{array}\right]kN\cdot s/m</math> .
|}
| style="width: 5px;text-align: right;white-space: nowrap;" | (24)
|}

Na Figura 11 temos o espectro de resposta via HBM para o modelo 3.

{| style="text-align: center; border: 1px solid #BBB; margin: 1em auto; width: auto;max-width: auto;"
|-
|style="padding:10px;"|  [[Image:draft_Zaltron_287705376-image34.png|400px]]
|- style="text-align: center; font-size: 75%;"
| colspan="1" style="padding-bottom:10px;"| '''Figura 11'''. Espectro de resposta via HBM para <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math>, <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{2}}}</math> e <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{3}}}</math>
|}


Ao observar os espectros de respostas gerados via HBM e comparar as Figuras 8 e 11, verifica-se que os amortecedores discretos a serem inseridos no sistema se tornam mais eficientes sobre o grau de liberdade ou posição onde incide a excitação dinâmica. Tal conclusão é verificada em razão das amplitudes dos deslocamentos geradas pelo espectro de resposta.

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
[[Image:draft_Zaltron_287705376-image35.png|600px]] </div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;font-size:75%;">
'''Figura 12'''. Histórico de resposta (tempo) para <math display="inline">{\mathit{\boldsymbol{u}}}_{\mathit{\boldsymbol{1}}}</math> via Pseudo-Forças e HBM</div>


A Figura 12 mostra a resposta para o grau de liberdade <math display="inline">{u}_{1}</math>, comparando-o com a Pseudo-Forças e HBM para um sistema dotado de amortecimento do tipo não proporcional.

==6.  Conclusões==

Ao se proceder uma análise harmônica, o HBM apresenta resultados satisfatórios ao se avaliar a resposta em uma frequência discreta <math display="inline">\omega</math>. Com isso, podemos obter facilmente a amplitude dos deslocamentos <math display="inline">\left| u\right|</math>  associada às propriedades do sistema.

O HBM permitiu, dentro da mesma solução, obter a resposta tanto no domínio do tempo quanto no domínio da frequência, tanto para sistemas dotados de amortecimento proporcional ou não proporcional. O histórico de deslocamento dado via HBM e comparado com os demais métodos (ImFT e Pseudo-Forças), apresentou resultados satisfatórios tanto em amplitude quanto em fase.

Quanto à inclusão dos dissipadores de energia, neste caso amortecimento discreto, pode-se avaliar que o mesmo tem influência significativa na resposta do sistema em relação à amplitude crítica. Ao avaliar a resposta no domínio do tempo, porém, temos que o melhor modelo de amortecimento a ser adotado é aquele distribuído de forma homogênea no sistema.

==Agradecimentos==

Os autores agradecem à CAPES, PROPEC e UFOP.

<span id='_GoBack'></span>

==Referências==
<div class="auto" style="width: auto; margin-left: auto; margin-right: auto;font-size: 85%;">
[1] Bernádez, S. A. ''Nonlinear normal modes in mechanical systems'': Concept and computational with numerical continuation. 2016. Thesys (Doctor) – ETSII/UPM, Madrid, Espanha, 2016.

[2] Camargo, R. S. Técnicas eficientes de correção de erros na análise dinâmica de estruturas no domínio da frequência. 2008. Dissertação (Mestrado) – PPGEC/UFES, Vitória, ES, Brasil, 2008.

[3] Childs, D. W; Harihara, P. ''Solving problems in dynamics and vibrations using MATLAB''. 1. ed. Texas: Texas A & M University College Station, 2012.

[4] Caret, A. M.; Venâncio-Filho, F. A modal superposition Pseudo-Forças method for dynamics analysis of structural system with non-proportional damping. ''Earthquake Engineering and Structural Dynamics'', n. 20, p. 303-315, 1991.

[5] Clough, R. W.; Penzien, J. ''Dynamics of structures''. 2. ed. NewYork: McGraw-Hill, 1993.

[6] Ferreira, W. G. Análise dinâmica não linear no domínio da frequência em sistemas estruturais com amortecimento não proporcional. 1998. Tese (Doutorado) – COPPE/UFRJ, Rio de Janeiro, RJ, Brasil, 1998.

[7] Galvão, S. A. ''Instabilidade estática e dinâmica de pórticos planos com ligações semi-rígidas''. 2004. Tese (Doutorado) – PUC-Rio, Rio de Janeiro, RJ, Brasil, 2004.

[8] Muñoz, L. P. et al. ''Análise dinâmica não linear de estruturas reticulares submetidas a ações sísmicas no domínio da frequência''. Conference: Ibero-Latin American Congress on Computational Methods in Engineering – CILAMCE, 35., v. 1, Fortaleza, 2014.

[9] Pasquetti, E. ''Métodos aproximados de solução de sistemas dinâmicos não-lineares''. 2008. Tese (Doutorado) – PUC-Rio, Rio de Janeiro, RJ, Brasil, 2008.

[10] Ribeiro, M. V. ''Formulação matricial da análise dinâmica de estruturas no domínio da frequência.'' 1998. Dissertação (Mestrado) – PROPEC/UFOP, Ouro Preto, MG, Brasil, 1998.

[11] Rubinstein, C. W.; Hurty, F. M. 1964. ''Dynamics of structures''. New York: McGraw-Hill, 1964.

[12] Santolin, A. Formulação da análise dinâmica de sistemas de múltiplos graus de liberdade no domínio da frequência em coordenadas físicas. 2006. Dissertação (Mestrado) – PROPEC/UFOP, Ouro Preto, MG, Brasil, 2006.

[13] Silva, A.R.D.; Sistema computacional para análise avançada estática e dinâmica de estruturas metálicas. Ouro Preto, 2009. Tese de (Doutorado ) – PROPEC/UFOP, Ouro Preto, MG, Brasil, 2009.

[14] Soares, JR. D. Análise dinâmica no domínio do tempo e da frequência de modelos não lineares discretizados por elementos finitos. 2002. Dissertação (Mestrado) – COPPE/UFRJ, Rio de Janeiro, RJ, Brasil, 2002.
</div>