Análise Dinâmica de Estruturas com Parâmetros Intervalares

UNIVERSIDADE FEDERAL DE ITAJUBÁ

INSTITUTO DE ENGENHARIA MECÂNICA

PROGRAMA DE PÓS-GRADUAÇÃO EM ENGENHARIA MECÂNICA

TESE DE DOUTORADO

Análise Dinâmica de Estruturas com

Parâmetros Intervalares

Autora: Profa. MSc. Janaína Cunha Vaz Albuquerque

Orientador: Prof. Dr. José Juliano de Lima Junior

Itajubá, 04 de dezembro de 2015




TESE DE DOUTORADO





Curso: Doutorado em Engenharia Mecânica

Área de Concentração: Projeto e Fabricação

Tese submetida ao Programa de Pós-Graduação em Engenharia Mecânica como parte

dos requisitos para obtenção do Título de Doutora em Engenharia Mecânica.

Itajubá, 04 de dezembro de 2015

M.G. – Brasil




TESE DE DOUTORADO





Composição da Banca Examinadora:

Prof. Dr. Fernando José de Oliveira Moreira - EMBRAER Prof. Dr. Antônio Wagner Forti - UNESP Prof. Dr. Nelson Manzanares Filho – IEM/UNIFEI Prof. Dr. Paulo Shigueme Ide - IEM/UNIFEI Prof. Dr. José Juliano de Lima Junior, Presidente - IEM/UNIFEI

Dedicatória

A minha família!

Agradecimentos

Aos meus pais, que apesar de todas as dificuldades enfrentadas, sempre priorizaram a

educação dos filhos.

Ao professor orientador, José Juliano de Lima Junior.

Aos professores Luiz Fernando Barca, Nelson Manzanares Filho, Wlamir Carlos de

Oliveira, Antônio Fernando Moura Santos, pela colaboração e atenção.

Ao professor Genésio José Menon, pela companhia no café.

Ao colega Júlio César Silva de Souza, pela ajuda e parceria.

À Malu, pelo valioso incentivo.

Aos funcionários do Laboratório do IEM da UNIFEI pelo auxílio.

Ao meu marido, pelas palavras de apoio e pelo importante incentivo acadêmico.

A CAPES, pelo apoio financeiro através do programa de bolsas e a FAPEMIG pelo

projeto TEC1670/05.

“Por vezes sentimos que aquilo que fazemos não é senão uma gota de água no mar.

Mas o mar seria menor se lhe faltasse uma gota.”

Madre Teresa de Calcuta

Resumo

VAZ, J. C. (2015), Análise Dinâmica de Estruturas com Parâmetros Intervalares, Itajubá, p.

155, Tese (Doutorado em Projeto e Fabricação), Instituto de Engenharia Mecânica,

Universidade Federal de Itajubá.

O presente trabalho busca avaliar o resultado de problemas de vibração mecânica na

presença de parâmetros de valores intervalares, tais como, massa, comprimento, módulo de

elasticidade, momento de inércia. Um método possibilístico e eficaz computacionalmente é

proposto a fim de se obter a faixa de intervalo da frequência natural de sistemas com

parâmetros intervalares. O problema de autovalor intervalar é formulado pelo Método dos

Elementos Finitos (MEF) cujas matrizes, de rigidez e de massa, são submetidas à perturbação

pela Teoria da Perturbação de Matrizes. Em cada fase da análise, a existência de incertezas na

formulação das matrizes é considerada como perturbação num sistema pseudo-determinístico

capaz de fornecer resultados tecnicamente confiáveis e eficientes. Resultados numéricos são

apresentados usando o programa computacional desenvolvido pela autora em MATLAB®.

Este programa trabalha com estruturas dinâmicas sujeitas a incertezas intervalares. Os

resultados numéricos são validados com a literatura e com o Método de Monte Carlo. Testes

experimentais foram realizados a fim de validar experimentalmente os resultados teóricos.

Palavras-chave

Parâmetros Intervalares, Frequência Natural, Teoria da Perturbação de Matrizes,

Método dos Elementos Finitos, Método de Monte Carlo.

Abstract

VAZ, J. D. C. (2015), Dynamic Analysis of Structures with Interval Parameters, Itajubá,

p. 155, Dr. Thesis – Mechanical Engineering Institute, Federal University of Itajubá.

This study aims to evaluate the result of mechanical vibration problems in presence of

an interval of values for parameters such as mass, length, modulus of elasticity, moment of

inertia. A possibilistic and computationally efficient method is proposed to obtain the limits

of natural frequency of systems with interval parameters. The interval eigenvalue problem is

formulated by the Finite Element Method (FEM) which (stiffness and mass) matrices are

submitted to disturbance by Perturbation Theory of Matrices. At each stage of the analysis,

the existence of uncertainty in matrix formulation is considered as the presence of disorder in

a pseudo-deterministic system capable of providing technically reliable and efficient results.

Numerical results are presented using a tool developed by the author in MATLAB ®. This

program works with dynamic structures with interval uncertainty. The numerical results are

compared with the literature and with a Monte Carlo simulation. Experimental tests were

conducted to validate the results of the proposed method.

Keywords

Interval Parameters, Natural Frequencies, Matrix Perturbation Theory, Finite Element

Method, Monte Carlo Method.

i

SUMÁRIO

SUMÁRIO____________________________________________________________ i

LISTA DE FIGURAS___________________________________________________ iv

LISTA DE TABELAS___________________________________________________ vi

SIMBOLOGIA_________________________________________________________ viii

CAPÍTULO 1__________________________________________________________ 1

INTRODUÇÃO________________________________________________________ 1

1.1 Revisão da Literatura------------------------------------------------------------------------ 5

1.2 Contribuição do Trabalho------------------------------------------------------------------ 9

1.3 Conteúdo-------------------------------------------------------------------------------------- 10

CAPÍTULO 2__________________________________________________________ 12

MÉTODOS NUMÉRICOS_______________________________________________ 12

2.1 Método de Monte Carlo--------------------------------------------------------------------- 13

2.2 Teoria da Perturbação de Matrizes-------------------------------------------------------- 15

2.2.1 Aritmética Intervalar----------------------------------------------------------------- 16

2.2.1.1 Propriedades não comuns entre números reais e intervalares 19

2.2.2 Teoria da Perturbação de Matrizes Aplicada a Problemas de Autovalores

com Parâmetros Intervalares--------------------------------------------------------

23

2.2.2.1 Matriz de perturbação simétrica, positiva semidefinida---------------- 27

2.2.3 Representação Intervalar da Matriz Rigidez Global e Matriz Massa Global- 28

2.2.3.1 Matriz de rigidez intervalar------------------------------------------------ 36

2.2.3.1.1 Matriz de rigidez global intervalar---------------------------- 38

2.2.3.2 Matriz massa intervalar---------------------------------------------------- 39

2.2.3.3 Problema de autovalor intervalar para estruturas dinâmicas---------- 40

2.2.3.3.1 Problema de autovalor------------------------------------------- 40

2.2.3.3.2 Solução do problema de autovalor----------------------------- 41

ii

2.2.3.3.3 Problema de autovalor submetido à matriz de

perturbação --------------------------------------------------------

42

2.2.3.3.4 Frequência natural intervalar----------------------------------- 45

CAPÍTULO 3__________________________________________________________ 48

PROCEDIMENTO DE CÁLCULO DE INCERTEZA DE MEDIÇÃO EM

MEDIÇÕES DIRETAS E INDIRETAS DAS GRANDEZAS FÍSICAS DE UMA

VIGA DE AÇO ________________________________________________________

48

3.1 Incerteza de Medição------------------------------------------------------------------------ 49

3.1.1 Incerteza Padrão----------------------------------------------------------------------- 51

3.1.2 Incerteza Padrão Combinada-------------------------------------------------------- 58

3.1.3 Incerteza Padrão Expandida--------------------------------------------------------- 59

3.2 Incerteza de Medição para o Cálculo das Freqüências Naturais Teóricas------------ 61

3.2.1 Incerteza em Medição Direta-------------------------------------------------------- 61

3.2.2 Incerteza em medição Indireta ----------------------------------------------------- 75

3.3 Incerteza de Medição para o Cálculo das Frequências Naturais Experimentais---- 84

CAPÍTULO 4__________________________________________________________ 89

VALIDAÇÃO_________________________________________________________ 89

4.1 Validação Teórica--------------------------------------------------------------------------- 90

4.1.1 Viga escalonada: módulo de elasticidade intervalar----------------------------- 94

4.1.2 Viga escalonada: seção transversal e momento de inércia intervalares------ 98

4.1.3 Treliça com oito elementos de barra----------------------------------------------- 99

4.1.4 Treliça com três elementos de barra----------------------------------------------- 103

4.2 Validação Experimental-------------------------------------------------------------------- 105

4.2.1 Descrição do ensaio------------------------------------------------------------------ 107

4.2.2 Resultados----------------------------------------------------------------------------- 107

CAPÍTULO 5__________________________________________________________ 110

RESULTADOS________________________________________________________ 110

5.1 Viga------------------------------------------------------------------------------------------- 110

5.2 Treliças---------------------------------------------------------------------------------------- 111

5.2.1 Treliça com três elementos de barras----------------------------------------------- 112

5.2.2 Treliça com dez elementos de barras--------------------------------------------- 114

5.2.3 Treliça com vinte elementos de barras--------------------------------------------- 118

5.3 Pórticos---------------------------------------------------------------------------------------- 121

CAPÍTULO 6__________________________________________________________ 124

CONCLUSÕES E TRABALHOS FUTUROS_____________________________ 124

iii

6.1 Conclusões------------------------------------------------------------------------------------ 124

6.2 Trabalhos Futuros------------------------------------------------------------------------- 125

REFERÊNCIAS BIBLIOGRÁFICAS_____________________________________ 127

ANEXO A ____________________________________________________________ 134

MÉTODO DOS ELEMENTOS FINITOS: FORMULAÇÃO DAS MATRIZES DE

RIGIDEZ E DE MASSA DE ELEMENTOS DE BARRA E DE VIGA__________

134

A.1 Elemento de barra--------------------------------------------------------------------------- 136

A.1.1 Elemento de barra unidimensional------------------------------------------------ 136

A.1.1.1 Elemento de barra inclinado com deslocamentos nodais na

direção do eixo x---------------------------------------------------------- 139

A.2 Elemento de viga---------------------------------------------------------------------------- 142

ANEXO B _____________________________________________________________ 148

PARAMETER SOLUTION VERTEX THEOREM___________________________ 148

ANEXO C _____________________________________________________________ 153

ANEXO C.1 – Paquímetro---------------------------------------------------------- 153

ANEXO C.2 – Balança de Precisão------------------------------------------------- 154

ANEXO C.3 – Acelerômetro Piezoelétrico---------------------------------------- 155

iv

Lista de Figuras

Figura 2.1 – Função densidade de probabilidade acumulada. 15

Figura 2.2 – Representação gráfica das funções g(x) e f(x). 21

Figura 2.3 – Malha formada por subdomínios (ou elementos finitos). 28

Figura 2.4 – Elemento de barra na horizontal 30

Figura 2.5 – Elemento de barra inclinado de um ângulo em relação ao eixo das

abscissas. 31

Figura 2.6 – Elemento de viga na horizontal 32

Figura 2.7 – Elemento de viga inclinado 33

Figura 2.8 – Elemento de viga com dois elementos de comprimentos iguais a L 36

Figura 3.1 – Curva de distribuição gaussiana. 52

Figura 3.2 – Curva de distribuição uniforme (adaptado de GUM, 2008 ). 54

Figura 3.3 – Curva de distribuição triangular (adaptado GUM, 2008). 56

Figura 3.4 – Viga engastada-livre em vibração livre. 61

Figura 4.1 – Viga escalonada em três segmentos. 90

Figura 4.2 – Curvas dos três primeiros autovalores obtidos após 10000 simulações

de Monte Carlo. 93

Figura 4.3 – Treliça com oito elementos de barra. 95

Figura 4.4 – Comparando os limites dos autovalores dentro do intervalo 0 2%

de incerteza obtidos em Qiu et al. (2005) com os autovalores obtidos por

meio de Monte Carlo, usando a matriz de massa apresentada em

Equação (A.21).

102

Figura 4.5 – Treliça com três elementos de barras e módulo de elasticidade

intervalar. 103

Figura 4.6 – Comparação dos parâmetros intervalares e parâmetros aleatórios. 105

v

Figura 4.7 – Régua graduada de aço de seção transversal retangular em vibração

livre. 107

Figura 5.1 – Viga escalonada em n segmentos 111

Figura 5.2 – Curva de inclinação do tempo de execução do algoritmo em função do

aumento do número de elementos da viga da Figura (5.1). 111

Figura 5.3 – Exemplos de treliças. 112

Figura 5.4 – Treliça com três elementos de barras 113

Figura 5.5 – Treliça com 10 elementos de barras. 115

Figura 5.6 – Autovalores resultantes de uma treliça, com parâmetros de entrada

intervalares, variando em função de um fator de incerteza no

intervalo 0 3% .

116

Figura 5.7 – Avaliando o comportamento da superestimação dos resultados com o

aumento do intervalo dos autovalores. 117

Figura 5.8 – Treliça com vinte elementos de barras. 119

Figura 5.9 – Pórtico com 14 elementos de viga e 39 graus de liberdades. 122

Figura A.1 – Malha formada por subdomínios (ou elementos finitos). 134

Figura A.2 – Diagrama de corpo livre de um elemento de barra unidimensional

(fonte: Finit Element Method using Matlab, 1996). 136

Figura A.3 – Pontos nodais de um elemento de barra unidimensional (fonte: Finit

Element Method using Matlab, 1996). 137

Figura A.4 – Pontos nodais de um elemento de barra bidimensional (fonte: Finit

Element Method using Matlab, 1996). 139

Figura A.5 – Elemento de barra bidimensional inclinado (fonte: Finit Element

Method using Matlab, 1996). 140

Figura A.6 – Elemento de viga na horizontal (fonte: Finit Element Method using

Matlab, 1996). 144

Figura A.7 – Funções de interpolação Hermitiana (fonte: Finit Element Method

using Matlab, 1996). 145

vi

Lista de Tabelas

Tabela 3.1 – Grau de liberdade efetivo e o correspondente fator de abrangência. 60

Tabela 3.2 – Valores obtidos após 20 medidas diretas de cada grandeza ( ix ). 63

Tabela 3.3 – Três primeiras frequências naturais coletadas durante os ensaios. 85

Tabela 4.1 – Autovalores da viga da Figura (4.1) com módulo de elasticidade

intervalar

95

Tabela 4.2 – Diferença relativa da Teoria da Perturbação de Matrizes em relação ao

Método de Monte Carlo dos três primeiros autovalores da Tabela (4.1)

96

Tabela 4. 3 – Autovalores da viga da Figura (4.1) com área da seção transversal e

momento de inércia intervalares

98

Tabela 4.4 – Diferença relativa da Teoria da Perturbação de Matrizes em relação ao

Método de Monte Carlo dos três primeiros autovalores da Tabela (4.1)

99

Tabela 4. 5 – Frequências naturais adimensionais da treliça da Figura (4.6) com

módulo de elasticidade intervalar

104

Tabela 4.6 – Diferença relativa das três primeiras frequências naturais adimensionais

apresentadas na Tabela (4.5)

104

Tabela 4.7 – Grandezas física e geométrica da viga de aço ensaiada 108

Tabela 4.8 – Parâmetros intervalares da viga de aço ensaiada 108

Tabela 4.9 – Comparação entre os resultados teóricos e experimentais 109

Tabela 4.10 – Diferença relativa do resultado numérico em relação ao resultado

experimental 109

Tabela 5.1 – Frequências naturais adimensionais de uma treliça simples composta por

três elementos de barras para módulo de elasticidade intervalar e módulo

determinístico.

114

Tabela 5.2 – Diferença relativa dos autovalores para um intervalo em 1% para a Figura

. (5.5).

118

vii

Tabela 5.3 – Intervalo das frequências naturais adimensionais de uma treliça

simples formada por vinte elementos de barras para módulo de elasticidade

intervalar e módulo determinístico.

120

Tabela 5.4 – Intervalo dos autovalores do pórtico da Figura (5.9) 123

viii

Simbologia

Letras Latinas

A área da seção transversal m2

IA matriz intervalar real

[ ]A matriz simétrica do problema de autovalores

a limite inferior do intervalo das distribuições de probabilidade

ia constante

2a largura das curvas de distribuição de probabilidade

b largura da seção transversal da viga m

c limite superior do intervalo das distribuições de probabilidade

ic constante

D matriz espectral ou de autovalores

d vetor deslocamento

d vetor deslocamento nodal

F força N

E módulo de elasticidade longitudinal. GPa

E matriz de perturbação simétrica definida positiva

ix

E(X) valor esperado de uma variável randômica contínua

f (x) função de densidade de probabilidade normalizada

F(X) função densidade de probabilidade acumulada

h altura da seção transversal da viga m

H função de interpolação

i número complexo

k fator de abrangência

K matriz de rigidez N/m

K matriz de rigidez global intervalar N/m

cK matriz de rigidez formada pelos valores médios dos números

intervalares

N/m

RK matriz de rigidez formada pelos valores radiais dos números

intervalares

N/m

eK matriz de rigidez do elemento N/m

e

iK matriz de rigidez determinística do elemento N/m

I momento de inércia de área. m4

I matriz identidade

L comprimento m

uL comprimento útil da viga (régua graduada de inox) usada no estudo

experimental

m

x

tL comprimento total da viga (régua graduada de inox) usada no estudo

experimental

m

[ ]iL matriz de conectividade Boleana do elemento

l limite inferior da matriz intervalar.

M momento fletor

M matriz massa kg

M matriz massa global intervalar kg

cM matriz massa formada pelos valores médios dos números

intervalares

kg

RM matriz massa formada pelos valores radiais dos números intervalares kg

e

iM matriz massa determinística do elemento kg

eM matriz de massa do elemento kg

am massa do acelerômetro. kg

vm massa da viga kg

n número de observações repetidas

p ordem de uma matriz

P carga axial N

0p valor médio da variável aleatórias X

[ ]Q matriz ortogonal ou matriz modal

xi

iQ força N

q carga distribuída N/m

R número aleatório

( )R x quociente de Rayleigh

números reais

( , )i jr x x coeficiente de correlação entre as estimativas ix e jx

s vetor arbitrário

2 ( )is x variância experimental

T matriz de rotação s

t tempo s

u limite superior da matriz intervalar

U energia potencial

TipoAu Incerteza Tipo A.

TipoBu Incerteza Tipo B

( )iu x incerteza padrão

( )cu y incerteza padrão combinada

u deslocamento temporal

u deslocamento temporal

1 2,u u deslocamentos nodais

xii

v deslocamento vertical m

( )V X variância

V força cortante N

( )w x função de teste

X número interval

Ix limite inferior de um número intervalar real

Ix limite superior de um número intervalar real

RX raio de um número intervalar

cX ponto médio de um número intervalar

X raio (ou incerteza) um intervalo

x autovetor

ix observações efetuadas

x coordenada cartesiana

x média aritmética

Xi grandezas de entrada (mensurando)

X variável aleatória

y estimativa da variável resposta Y

y vetor

Y número intervalar

Z número intervalar

xiii

( )z t velocidade m/s

( )z t aceleração m/s2

( )z t deslocamento m

Letras Gregas

1 escalar que multiplica 1x na combinação linear dos autovetores

2 escalar que multiplica 2x na combinação linear dos autovetores

fator de incerteza, que pré multiplica a área da seção transversal,

dentro do intervalo 0 2% .

perturbação das matrizes intervalares

ângulo rad

valor médio

massa específica. kg/m3

deformação específica longitudinal

eff grau de liberdade efetivo

desvio padrão

operador de diferença parcial.

n frequência natural angular rad/s

n frequência natural intervalar angular rad/s

xiv

matriz diagonal dos autovalores

autovalor sujeito a matriz de perturbação

autovalor

fator de amortecimento

frequência natural adimensional

decremento logarítmico

Sobrescritos

T transposta de uma matriz

e índice que representa um elemento no método dos elementos finitos

quantidade determinística.

~ quantidade intervalar

Subscritos

R índice que representa matriz de valores radiais

C índice que representa matriz de valores médios (centrais)

i índice que representa o número do elemento discretizado

max índice do limite superior das frequência natural

mín índice do limite inferior das frequência natural

xv

Abreviaturas

EBE Element by Element

FDA Função Densidade de Probabilidade Acumulada

GDL Grau de Liberdade

MEF Método dos Elementos Finitos.

MEFI Método dos Elementos Finitos Intervalares

PSVT Parameter Solution Vertex Theorem

Siglas

IEM Instituto de Engenharia Mecânica.

1

Capítulo 1

INTRODUÇÃO

O comportamento de estruturas com parâmetros intervalares (ou incertos) têm sido

investigado por diversos autores devido a sua importância nos campos da engenharia

mecânica estrutural. A presença de parâmetros intervalares nos projetos mecânicos pode ser

atribuída às imperfeições físicas, imprecisões do modelo matemático e a complexidade do

sistema. Na prática, busca-se trabalhar com projetos tecnicamente confiáveis, com elevado

grau de segurança e reduzido custo financeiro. Tendo em vista as necessidades do mundo na

atualidade, a avaliação das incertezas nos projetos vem se tornando uma etapa que exige cada

vez mais atenção e cuidados especiais durante a sua elaboração. A fim de se garantir um bom

desempenho ao longo da vida útil do sistema, devem-se levar em consideração as incertezas

presentes nos parâmetros estruturais.

Na prática observa-se que os projetos de engenharia mecânica estrutural sofrem a

influencia de incertezas nos parâmetros em graus variados. Este trabalho está preocupado com

problemas de vibração mecânica estrutural quando sujeita a incertezas em seus parâmetros

físicos, o que produz matrizes de rigidez e de massa de valores variados dentro de um

intervalo. Neste caso, o problema de autovalor é transformado em um problema de

autovalores intervalares.

No campo da engenharia, incerteza é uma palavra importante que deve

sempre ser levada em conta. A análise de incerteza torna-se mais evidente quando o objetivo

dos projetos de engenharia for a segurança estrutural, que é a base para a concepção ou

2

avaliação de qualquer projeto de engenharia. Na prática os materiais de construção não

apresentam propriedades homogêneas e podem apresentar variação das propriedades ao longo

da vida, é o que ocorre com o duralumínio, liga de alumínio usada em estruturas aeronáuticas.

Com o passar do tempo ocorrem mudanças nas propriedades mecânicas do material devido a

perda do tratamento térmico, conhecido como têmpera, ao qual a liga é submetida. Essa perda

ocorre por conta do envelhecimento do material, ocasionado pela temperatura do ambiente.

Esse processo de perda do tratamento térmico do duralumínio produz variação nas

propriedades mecânicas e consequentemente suas grandezas tornam-se variáveis, (Figura 1-

a).

Além das incertezas presentes no material usado na fabricação das estruturas

aeronáuticas, pode-se citar também, as incertezas associadas à rigidez dos mancais que

compõe as linhas de eixos de máquinas rotativas. Neste caso, tem-se que as frequências

naturais das linhas de eixos ficam susceptíveis a incertezas em suas grandezas, em particular a

rigidez dos mancais de apoio. Essa rigidez é um dos parâmetros mais importantes no cálculo

das frequências naturais, podendo haver diferenças significativas nos resultados mesmo na

presença de pequenas variações nos valores da rigidez. A rigidez dos mancais é o resultado

combinado da flexibilidade da película lubrificante, da flexibilidade da carcaça do mancal e

da base de fixação da carcaça do mancal. Todas essas flexibilidades são de difícil

determinação, tendo em vista a complexidade dos fenômenos envolvidos e da geometria das

estruturas. Desta forma, o cálculo de velocidades críticas de operação de linhas de eixo deve

ser executado, tendo-se em mente a incerteza associada à determinação da rigidez dos

mancais e, para finalidades de projeto, seria adequado avaliar o intervalo de possibilidades de

frequências naturais em função das incertezas na determinação da rigidez dos mancais,

(Figura 1- b).

Outra fonte de incerteza em parâmetros mecânicos estruturais é o processo de corrosão

dos materiais de fabricação, pois provoca a deterioração dos materiais como perda uniforme

de massa e consequentemente diminuição da seção transversal, (Figura 1- c).

As incertezas estão presentes também no processo de medição das grandezas físicas. A

exatidão e a precisão de uma grandeza medida estarão sempre limitadas tanto pela

sofisticação do equipamento utilizado como pela habilidade do sujeito que realiza a medida

(Figura 1- d).

3

a) Estrutural aeronáutica.

b) Linhas de eixo de máquinas

rotativas: eixos, rotores, mancais.

c) Corrosão

d) Incertezas em medições

Figura 1 – Exemplos de causas de incertezas.

Vale ressaltar a importância de se garantir a segurança estrutural nos projetos de

engenharia. Uma falha estrutural pode trazer graves problemas para a sociedade, como lesões

humanas, problemas ambientais ou até problemas econômicos. Por esses motivos diversos

estudos têm sido realizados para reduzir os riscos de falhas e garantir maior segurança aos

projetos. A Figura 2 – a, mostra a foto da ponte I-35W localizada nos Estados Unidos, onde o

cálculo correto do projeto poderia ter evitado a quebra da ponte. A Figura2- b, apresenta a

foto de um prédio localizado em Bélem – PA, após o desabamento. O desabamento ocorreu

em janeiro de 2016 e segundo a perícia feita no local uma das possíveis causas do acidente

aponta erro no cálculo estrutural.

4

a) Ponte I-35W após desabamento

b) Prédio desaba em Bélem – PA: erro de

calculo estrutural

Figura 2 – Falhas estruturais.

Os parâmetros físicos utilizados para descrever as estruturas são muitas vezes imprecisos

devido às incertezas físicas e geométricas, ou imprecisões na modelagem do sistema. Estes

parâmetros são, geralmente, identificados como variáveis aleatórias sendo introduzido numa

abordagem estocástica, ou seja, variáveis aleatórias que dependem do tempo. Na literatura,

diferentes métodos são apresentados a fim de solucionar problemas com variáveis aleatórias.

Como exemplo, podem-se citar as técnicas probabilísticas. Entretanto, os métodos

probabilísticos tratam de variáveis aleatórias (estocásticas) cuja densidade de probabilidade

deve ser conhecida. Muitas vezes, têm-se apenas informações empíricas acerca do problema,

sem o conhecimento sobre a densidade de probabilidade, dificultando a avaliação por tais

métodos.

Com base em informações empíricas ou informações coletadas por meio de uma

investigação do problema mecânico, garante-se mais sucesso na avaliação do problema

quando os parâmetros são representados por uma faixa de valores, sem a necessidade de

seguir uma determinada função de distribuição de probabilidade (Pereira, 2002). Neste caso,

busca-se trabalhar com os métodos sob o ponto de vista da aritmética intervalar. Desta

maneira, a incerteza é envolvida pelo tamanho do intervalo e não há necessidade de nenhuma

distribuição de probabilidade para representá-la, simplificando assim o trabalho de operações

e análises.

5

Este trabalho aborda os problemas mecânicos com parâmetros intervalares. Um método,

de solução rápida e eficiente, é apresentado para tratar de problemas dinâmicos intervalares. O

método é validado com os resultados da literatura e com os resultados experimentais. A fim

de obter uma faixa de valores que encontrasse os parâmetros com intervalos estreitos e

confiáveis, utilizaram-se técnicas de probabilidade e estatística apenas para caracterizar a

incerteza das variáveis e obter os limites do intervalo. Vale ressaltar que não há necessidade

de tal operação na prática, pois, com base numa pesquisa de campo é possível obter dados (ou

melhor, os limites das grandezas físicas) a partir de fontes diretas (pessoas) que conhecem,

vivenciam ou tem conhecimento sobre o tema. Em suma, com uma pesquisa empírica é

possível obter o intervalo das variáveis e aplicar a metodologia abordada aqui.

As incertezas nas variáveis de um sistema e em suas relações podem, de forma geral,

serem atribuídas à fontes cognitivas e não-cognitivas:

As fontes não-cognitivas incluem parâmetros do sistema que são tratados como

variáveis aleatórias com distribuição de probabilidade ou dados estatísticos conhecidos.

Nestes casos a resposta do sistema pode ser determinada através da teoria da probabilidade e

de processos aleatórios.

As fontes cognitivas incluem incertezas decorrentes de informações imprecisas ou

vagas, sujeitas a julgamentos ou subjetividade de opiniões. Os axiomas da

probabilidade estatística são limitados para tratar este tipo de incerteza, propiciando

por exemplo a aplicação do Método dos Elementos Finitos (MEF) cujas matrizes são

tratadas com base nas propriedades da matemática intervalar.

Portanto, de acordo com o tipo de incerteza existente em um sistema de engenharia,

diferentes abordagens são utilizadas em sua análise. Assim, duas abordagens são tratadas e

discutidas neste trabalho:

A abordagem probabilística, trata de incertezas numéricas não-cognitivas, na qual as

propriedades são tratadas como variáveis aleatórias. Neste caso é aplicado o MEF

utilizando a simulação de Monte Carlo;

A abordagem possibilística, trata de incertezas numéricas cognitivas, neste trabalho será

representada por um modelo baseado no MEF com aplicação da aritmética intervalar,

6

a fim de que o MEF trabalhe com números intervalares. Neste caso, as matrizes, de

rigidez e de massa, são submetidas a matrizes de perturbação com o objetivo de obter

os limites, inferior e superior, da resposta do sistema.

A técnica de simulação de Monte Carlo baseia-se na geração de um número finito de

amostras de um processo. São gerados vetores de variáveis aleatórias a partir de uma função

de densidade de probabilidade e, a partir deles, realizam-se as simulações do problema. Ao

final do processo, gera-se um resumo estatístico em que são avaliadas as médias e as

variâncias das respostas das simulações. Esta técnica possui como vantagens o fato de ser

simples e absolutamente geral e como desvantagem o número exagerado de simulações

exigido para se obter um resultado com precisão aceitável, o que pode tornar o método

inviável para problemas de grande porte.

A aritmética intervalar, em princípio, é considerada a ferramenta mais simples para tratar

uma incerteza numérica cognitiva, bastando fornecer um intervalo numérico que contém o

valor exato e seguir em frente com as análises, ou seja, não são necessárias informações

estatísticas e nem funções de pertinência das variáveis. As propriedades de subcancelamento e

de subdistributividade da aritmética intervalar fazem com que cuidados especiais sejam

tomados, principalmente na solução de sistemas de equações nos quais os métodos

convencionais não funcionam. Uma diferença notável em aritmética intervalar é observada na

propriedade distributiva x(y + z) = xy + xz que, em geral, não vale para números intervalares,

onde apenas garante-se que: x(y ± z)⊆ xy ± xz, para x, y, z ∈ ℜ e (x ± y)z ⊆ xz ± yz, para x, y,

z ∈ ℜ . Além dos cuidados citados, deve-se ficar atento ao cálculo de expressões na qual a

variável aparece mais de uma vez. Estas expressões podem gerar resultados errados com

superestimação.

As formulações do MEF formam a matriz de rigidez global do sistema por meio de um

procedimento conhecido como assembly, na qual a rigidez dos elementos finitos, que

compartilham um mesmo grau de liberdade, é somada. O MEF foi implementado para

trabalhar com números intervalares onde a matriz de rigidez global e a matriz de massa global

são obtidas usando os princípios da Aritmética Intervalar.

O presente trabalho busca avaliar a solução dos problemas de vibração mecânica na

presença de parâmetros de valores intervalares, tais como, massa, comprimento, módulo de

7

elasticidade, momento de inércia. O resultado teórico é avaliado por meio do MEF associado

à teoria da perturbação de matrizes aplicada a problemas dinâmicos com parâmetros

intervalares (Mehdi, 2005). O método utilizado pelos autores trabalha com problemas de

autovalor intervalar num sistema pseudo determinístico capaz de fornecer resultados

tecnicamente confiáveis de maneira rápida e eficiente. Neste caso, o intervalo das frequências

naturais, de um problema dinâmico, é obtido por meio de um problema de autovalor com

matrizes (de rigidez e/ou de massa) centrais sujeitas a matrizes (de rigidez e/ou de massa) de

perturbação radial. A fim de se obter resultados físicos corretos, são impostas restrições

matemáticas às matrizes de perturbação, tais como: a matriz de massa de perturbação radial

deve ser obtida através de uma combinação linear de matrizes positivas definidas, já a matriz

de rigidez de perturbação radial é obtida por uma combinação de matrizes positivas

semidefinidas.

1.1 REVISÃO DA LITERATURA

Uma questão importante enfrentada na prática da engenharia é como lidar com as

grandezas de valores incertos. As incertezas estão presentes nas grandezas dos sistemas

mecânicos, porém na maioria dos casos trabalha-se com valores médios a partir de grandezas

intervalares. Há casos, como os projetos de veículos automotores, veículos espaciais,

estruturas civis, entre outros, em que se necessita trabalhar com as incertezas, a fim de evitar

que uma aproximação grosseira de um intervalo venha prejudicar a precisão dos resultados.

De acordo com Oberkampf et al. (1999) a incerteza é definida como uma deficiência que

pode ou não ocorrer durante o processo de modelagem de estruturas, devido à falta de

conhecimento do engenheiro acerca do projeto. Elas surgem devido às incompatibilidades e

imprecisões nos parâmetros físicos e geométricos, tais como carga, módulo de elasticidade,

coeficiente de Poisson, comprimento, entre outros. Existem diferentes técnicas para lidar com

incertezas e aleatoriedade (Altus et al., 2005;. Ayyub e Gupta, 1997 ; Biondini et al., 2004;.

Hurtado, 2002; Klir, 1995; Maglaras et al., 1997;. Muhanna e Mullen, 2001; Zimmermann,

2000).

8

Uma maneira simples de representar e executar operações com dados inexatos é utilizar a

aritmética intervalar, também chamada de análise intervalar. Neste contexto, uma incerteza

pode ser representada por um intervalo de números reais, que contém o valor desconhecido

exato do número em questão. Desta maneira, a incerteza é envolvida pelos limites do

intervalo e não há necessidade de nenhuma distribuição de probabilidade para representá-la,

simplificando assim o trabalho de operações e análises. A primeira abordagem da aritmética

intervalar começou com Burkill (1924). Em seguida, um novo estudo sobre o assunto foi

realizado por Young (1931). Anos mais tarde, Sunaga (1958) concretizou o uso da aritmética

intervalar. Mas foi com os estudos da dissertação de doutorado de Moore (1962) que a

aritmética intervalar ganhou uma atenção especial, na qual o autor apresenta um

desenvolvimento moderno da técnica com uma gama de opções de utilização na prática.

Desde então, vários artigos e livros têm abordado este assunto. Congressos, encontros,

softwares e grupos de discussão vêm contribuindo para o aumento do interesse e da

importância deste tema.

A Técnica da Perturbação é essecialmente usada para avaliar a incerteza da resposta

estrutural. Matos (2007) utiliza a técnica na engenharia civil para avaliar a resposta de

sistemas estruturais com parâmetros intervalares, por meio dos valores médios e seus devios

padrão. A metodologia fornece resultados precisos para problemas lineares e variáveis

aleatórias com distribuição de probabilidade normal ou quasi-normal (Eibl e Schmidt-

Hurtienne, 1995.

Diversas formulações têm sido propostas para avaliar eficientemente a incerteza das

variáveis de entrada em diferentes áreas da engenharia, tais como, engenharia civil, mecânica

e aeroespacial (DeLaurentis e Marvris, 2000 ; Du e Sudjianto, 2004;. Haftka et al., 2006;

Huang e Xiaoping, 2006; Rao e Sawyer, 1995; Zou e Mahadevan, 2006). Uma das

metodologias alternativas que tem sido aplicada é a Técnica de Redes Neurais (Ayyub e

Gupta, 1997; Hurtado, 2002; Papadrakakis e Lagaros, 2002).

Muhanna e Mullen (1999, 2001) desenvolveram uma nova formulação do Método dos

Elementos Finitos Intervalar (MEFI) baseada na técnica Element-by-Element (EBE) para o

estudo de incertezas em mecânica estrutural. O método permite o tratamento de incertezas

tanto no carregamento quanto na matriz de rigidez dos elementos. Esta formulação evita as

9

principais fontes de superestimação e gera uma envoltória de soluções com intervalos bastante

estreitos.

Matsumoto e Iwaya (2000) apresentaram duas metodologias para obter as respostas de

uma treliça com incertezas na geometria da estrutura. Na primeira metodologia, o problema é

formulado pelo MEF adaptado a matemática intervalar, introduzido por Moore (1962). Na

segunda metodologia, o MEF é formulado como um problema de otimização, que demonstrou

ser mais eficiente, pois se obteve as respostas que mais se aproximaram da solução verdadeira

obtida por simulação de Monte Carlo.

Conforme mencionado na introdução, o Método de Monte Carlo é uma metodologia

bastante conhecida, de fácil implementação e muito utilizado para validar novas técnicas de

tratamento de incerteza em sistemas mecânicos, porém o seu elevado custo computacional o

torna inviável na prática (Mahadevan e Raghothamachar, 2000; Olsson et al., 2003; Ra-

Rohani e Singh, 2004; Schueller, 2001).

A análise de incertezas na exploração de petróleo é uma área de interesse das companhias

com o objetivo de reduzir o risco exploratório com a consequente otimização dos recursos

financeiros. Tendo em vista este cenário, Pereira (2002) apresentou um estudo com o uso da

aritmética intervalar e da aritmética fuzzy como uma alternativa para avaliar a incerteza em

métodos numéricos para modelagem de bacias. Zadeh (1965) introduziu a teoria dos

conjuntos fuzzy com o objetivo de definir classes de objetos com graduações contínuas de

pertinência, ou associação, compreendidas no intervalo [0,1]. Desta forma, é possível modelar

sistemas complexos que seriam difíceis na teoria dos conjuntos convencionais. Em álgebra

booleana, a noção de valores verdadeiros e falsos está limitada a 1 (um) ou 0 (zero). Em

lógica fuzzy, é possível tratar valores “mais ou menos” verdadeiros ou falsos, definidos por

números reais que variam no intervalo [0,1].

Parâmetros incertos podem causar mudanças significativas nas frequências naturais de

sistemas mecânicos. Tendo em vista a ampla aplicabilidade do assunto na prática, cientistas e

engenheiros vêm estudando problemas de autovalor formado por matrizes intervalares, cujos

elementos são limitados dentro de uma região de possíveis valores. Rohn (1987) estudou

problemas de autovalor intervalar padrão de uma matriz intervalar simétrica e fórmulas

derivadas para autovalores intervalares. Hallot e Bartlett (1987) descobriram que o espectro

10

de autovalores de família de matrizes intervalares depende do espectro de seus limites

definidos. Hudak (1984) investigou maneiras de relacionar a descoberta de Hallot e Bartlett

com os autovalores de uma matriz constante sob certas condições.

Com base nas propriedades de invariância das entradas do vetor característico, Deif

(1991) apresenta um método para resolver problema de autovalor intervalar padrão. Qiu et al.

(1993) estenderam o método de Deif para o problema de autovalor intervalar generalizado.

Porém, a falta de um critério eficiente para julgar as propriedades de invariância dos

autovalores, de acordo com as operações intervalares, apresentam algumas restrições ao

método de Deif que prejudicam a sua aplicabilidade. Com o intuito de superar essa limitação

do método de Deif (1991), Qiu et al. (1995) desenvolveram um método para calcular

autovalor intervalar assumindo semi-definição positiva de erro no par de matrizes

intervalares. Para pequenos erros nas matrizes, Qiu et al. (1996) apresentaram um método de

perturbação intervalar para problemas de autovalor intervalar.

Problemas de autovalor intervalar também foram estudados por Qiu et al. (2005). Os

autores interessados em problemas de vibração estrutural com grandezas intervalares

aplicaram o teorema Parameter Solution Vertex Theorem (PSVT) na solução de

autoproblemas intervalares. Qiu et al. (2005) trabalharam com decomposição positiva de

matrizes de rigidez e de massa.

Adhikari et al. (2009) apresentam dois estudos experimentais que podem ser utilizados

para testar e validar diferentes métodos de quantificação de incerteza estocástica em

diferentes faixas de frequências. Os testes foram rigorosamente controlados e as incertezas

podem ser consideradas conhecidas para fins práticos, permitindo que os estudos na área de

quantificação paramétrica de incerteza em modelos numéricos sejam comparados e validados

com os resultados obtidos experimentalmente pelos autores. O experimento descrito por

Adhikari et al. (2009) utiliza uma viga biengastada com 12 massas colocadas aleatoriamente

ao longo do comprimento da viga. As massas aleatórias correspondem a cerca de 1,6 % da

massa total da viga, o que simula erros aleatórios na matriz de massa. Cem (100) sistemas

dinâmicos nominalmente idênticos foram criados e testados individualmente no laboratório de

Bristol - Laboratory for Advanced Dynamic Engineering.

11

Mehdi (2005) em sua tese de doutorado desenvolveu uma nova metodologia de fácil

aplicação e de elevado desempenho computacional para a análise da resposta dinâmica de um

sistema sujeito a parâmetros intervalares. O método não exige um processo de solução

combinatória (de elevado custo computacional) para o cálculo dos limites da resposta

dinâmica e sugere que este possa ser estabelecido como uma ferramenta padrão na

investigação de problemas gerais de vibração estrutural na presença de grandezas intervalares.

O problema de autovalor intervalar, formulado por meio do MEF, utiliza o conceito de

comportamento monotônico de autovalores cuja matriz de rigidez (simétrica) fica sujeita a

matriz de perturbação (positiva semidefinida), o que leva a um processo computacionalmente

eficiente para determinar os limites das frequências naturais de um sistema mecânico

estrutural. Em seguida, usando os procedimentos para a perturbação de subespaços invariantes

de matrizes, obtêm-se os limites do desvio direcional (inclinação) de cada modo de vibração.

Nesta tese, buscou-se desenvolver uma metodologia capaz de solucionar problemas na

presença de incertezas tanto na matriz global de rigidez quanto na matriz global de massa, já

que na prática, as incertezas presentes em problemas de autovalor aparecem em ambas

matrizes globais.

1.2 CONTRIBUIÇÃO DO TRABALHO

Apesar da importância do tema abordado aqui, este assunto ainda é pouco explorado,

principalmente no Brasil. Com isso, busca-se contribuir com o desenvolvimento de uma

metodologia e um programa computacional para avaliar a resposta de um sistema dinâmico

sujeito a incertezas intervalares. Baseado na teoria de perturbação de subespaços invariantes

de matrizes, Mehdi (2005), obtém-se os limites intervalares de cada frequência natural. Medhi

(2005) apresenta uma metodologia para a solução de problemas sob a influência de incertezas,

apenas, na matriz de rigidez global. Vendo a necessidade de uma metodologia que

apresentasse uma solução para problemas práticos de engenharia, na qual ambas matrizes, de

rigidez de massa, apresentassem incertezas, desenvolveu-se nesta tese uma metodologia capaz

de solucionar problemas com incertezas tanto em parâmetros que compõe a matriz de rigidez

quanto em parâmetros que compõe a matriz de massa global. As matrizes globais de rigidez e

de massa, que compõem o problema de autovalor são perturbadas por matrizes positivas

12

semidefinidas. A metodologia abordada nesse estudo evita um problema, inerente às

operações com números intervalares, conhecido como superestimação de resultados.

Para isto, primeiramente, busca-se colocar os coeficientes intervalares como

multiplicadores da matriz intervalar. A partir daí, a substituição das matrizes intervalares

pelas matrizes de valores médios e pelas matrizes radiais converte o problema, antes não

determinístico, em um problema pseudo-determinístico.

Baseado na metodologia de Mehdi (2005) foi desenvolvido pela autora um programa

computacional, em MATLAB®, formulado pelo MEF e adaptado a matemática intervalar

capaz de fornecer a resposta de problemas dinâmicos sob a influência de incertezas

intervalares.

1.3 CONTEÚDO DO TRABALHO

No capítulo 1, faz-se uma revisão dos trabalhos publicados na área de tratamento de

incertezas de estruturas com parâmetros intervalares e das ferramentas disponíveis na

literatura para a análise estrutural estática e dinâmica capaz de fornecer os limites do intervalo

da resposta de sistema estrutural.

No capítulo 2, são apresentados dois métodos para análise de sistemas dinâmicos com

parâmetros incertos, a saber: Método de Monte Carlo (Fishman, 1996; Liu, 2001), usado para

calibrar e validar novas técnicas de tratamento de incertezas, e a Teoria da Perturbação de

Matrizes aplicada a problemas de autovalores intervalares (Mehdi, 2005), que avalia a média

aritmética e os limites inferior e superior dos parâmetros intervalares. Um código

computacional de elementos finitos adaptado à matemática intervalar, por meio da Teoria da

Perturbação de Matrizes, foi desenvolvido e implementado em MATLAB® para tratar de

sistemas dinâmicos sujeitos a parâmetros intervalares.

No capítulo 3, discute-se sobre o procedimento de avaliação das incertezas presentes

nas medições diretas e indiretas realizadas para o cálculo dos limites do intervalo das

grandezas de uma viga.

13

No capítulo 4, validam-se os resultados obtidos pela autora com a literatura e com dados

experimentais. A resposta dinâmica de vigas e treliças (com parâmetros intervalares)

formuladas pela Teoria da Perturbação de Matrizes é comparada com a resposta obtida por

meio do Parameter Solution Vertex Theorem (Qiu et al. (2005)) e com a resposta obtida

simulação de Monte Carlo. Os resultados numéricos são comparados, também, com os

resultados experimentais coletados a partir de uma viga em balanço em vibração livre.

No capítulo 5, apresentam-se problemas que representam alguns dos sistemas

encontrados na prática. São discutidos e analisados os autovalores e as frequências naturais de

sistemas dinâmicos.

As conclusões da tese são apresentadas no capítulo 6, com comentários sobre alguns dos

resultados numéricos e experimentais, além de sugestões para trabalhos futuros.

ANEXO A: formulação das matrizes de rigidez e de massa de elementos de barra e de

viga.

ANEXO B: no anexo B apresenta-se o teorema desenvolvido por Qiu et al. (2005) -

Parameter Solution Vertex theorem, que é baseado na teoria de otimização.

ANEXO C: informações de catálogo – paquímetro – fabricante: Mitutoyo; balança de

precisão – fabricante: Marte; certificado de calibração – acelerômetro piezoelétrico –

fabricante: PCB.

14

Capítulo 2

MÉTODOS NUMÉRICOS

Neste capítulo são apresentados dois métodos para análise de sistemas mecânicos

com parâmetros incertos, a saber: Método de Monte Carlo (Fishman, 1996; Liu, 2001),

usado para calibrar e validar novas técnicas de tratamento de incertezas, e a Teoria da

Perturbação de Matrizes aplicada a problemas de autovalores intervalares (Mehdi,

2005), que avalia a média aritmética e os limites inferior e superior dos parâmetros

intervalares. Baseado na metodologia de Mehdi foi desenvolvido pela autora um

programa computacional, em MATLAB®, formulado pelo Método dos Elementos

Finitos adaptado a matemática intervalar. Neste estudo, o sistema dinâmico pode estar

sujeito tanto a matriz de rigidez intervalar, quanto a matriz de massa intervalar,

diferentemente de Mehdi (2005), que apresentou a metodologia somente para matrizes

de rigidez intervalar.

15

2.1 MÉTODO DE MONTE CARLO

A técnica de simulação direta de Monte Carlo (Fishman, 1996; Liu, 2001) baseia-

se na geração de um número finito de amostras de um processo. São gerados vetores de

variáveis aleatórias a partir de funções de densidade de probabilidade estipuladas. A

partir deles, realiza-se um grande número de simulações onde são calculadas as

respectivas respostas do processo. Ao final, faz-se um resumo estatístico, onde são

avaliadas as freqüências, médias e variâncias das respostas. A principal vantagem deste

método probabilístico é a sua aplicabilidade geral, e devido a esta característica, tem

sido utilizado para calibrar diversas novas técnicas de tratamento de incertezas. Já a

desvantagem é o número grande de simulações exigido para se obter um resultado com

precisão aceitável, o que pode tornar o método inviável para problemas de grande porte.

Geradores de números aleatórios ou randômicos constituem o cerne do Método

de Monte Carlo. Os computadores não conseguem gerar números aleatórios no estrito

sentido de sua definição, já que estes são, teoricamente, imprevisíveis e irreproduzíveis,

o que deixa de ser verdade se existirem algoritmos para gerá-los. Na verdade, existem

geradores de números pseudo-aleatórios, que são previsíveis e reproduzíveis, mas que

podem gerar séries com uma aparência aleatória (aparentemente, sem padrão),

característica suficiente para realizar simulações como as de Monte Carlo.

Os números aleatórios são gerados no intervalo [0,1] e devem respeitar duas

propriedades básicas:

uniformidade – todos os números gerados devem ter a mesma probabilidade de

ocorrer;

independência – o valor corrente de um número aleatório não tem relação com o

número anterior.

Os diversos compiladores já possuem geradores de números pseudo-aleatórios,

como a função Rnd do Visual Basic, ou a função Rand do FORTRAN. A função RAND

16

do Matlab retorna uma matriz quadrada, que contém os números pseudo-aleatórios,

obtida a partir de uma distribuição uniforme padrão no intervalo [0,1].

Muitos processos físicos, para serem simulados, necessitam da geração de

valores representativos de variáveis que obedecem a uma determinada distribuição

estatística. Isto é feito através da geração de números aleatórios, sendo esta uma etapa

chave na Simulação de Monte Carlo, e particularmente na Simulação Direta de Monte

Carlo.

A distribuição contínua de uma variável aleatória x pode ser descrita por uma

função de distribuição normalizada f(x), tal que a probabilidade de um valor da variável

estar entre x e x + dx é dada por f(x)dx. Se o valor x está num intervalo [a, c], então a

probabilidade total é dada por:

( )

c

a

f x∫ dx = 1 (2.1)

Define-se a função de distribuição acumulada como sendo:

( ) ( )

x

a

F x f x= ∫ dx (2.2)

Com estas definições, torna-se possível gerar um número aleatório R e fazê-lo

igual à ( )F x . O valor representativo para x será obtido a partir de:

( )F x R= (2.3)

Considerando um exemplo trivial, no qual a variável x é uniformemente

distribuída no intervalo [a,c], a distribuição ( )f x será constante e igual a:

1( )f x

c a=

− (2.4)

17

Substituindo a Equação (2.4) na Equação (2.2) obtém-se a Equação (2.5),

conforme Figura (2.1).

1( )

x

a

F xc a

=−

∫ dx x a

c a

−=

− (2.5)

Distribuição uniforme

( )x a

F x Rc a

−= =

−

Figura 2.1 – Função densidade de probabilidade acumulada

Fonte: adaptado de Pereira (2002).

Com base nas Equações (2.5) e (2.3) obtém-se a Equação (2.6).

( )x a c a R= + − (2.6)

A Equação (2.6) relaciona o número aleatório gerado com o valor do parâmetro

em análise, para cada R gerado obtém-se um i

x .

2.2 TEORIA DA PERTURBAÇÃO DE MATRIZES

Mehdi (2005) desenvolveu uma metodologia capaz de avaliar sistemas estruturais

sujeitos a grandezas físicas intervalares como dimensões de geometria, material ou

módulo de elasticidade. Um programa de elementos finitos foi desenvolvido para

aplicar a metodologia e obter os limites da resposta dinâmica de estruturas com

parâmetros intervalares.

18

Em cada etapa da análise, a existência da variação é considerada como a

presença de perturbação no sistema. Inicialmente, um problema de autovalor intervalar

linear é analisado com base no conceito de comportamento monotônicamente crescente

de autovalor para matriz simétrica, sujeita a perturbação definida não negativa. Isso

conduz a um procedimento eficiente computacionalmente para determinar os intervalos

das frequências naturais e dos modos de vibração presentes nos sistemas mecânicos.

A fim de expressar os problemas intervalares em função da perturbação,

introduz-se uma perturbação intervalar, [ 1,1]ε = − , e o intervalo geral é representado

pela soma do valor médio e o valor radial. A transformação de uma problema de

autovalor intervalar em uma perturbação intervalar é facilmente realizada com noções

básicas sobre matemática intervalar (ou aritmética intervalar) que será apresentado a

seguir.

2.2.1 Aritmética intervalar

A maior parte dos problemas científicos envolve incerteza em parâmetros ou

dados inexatos. Uma maneira de representar e executar operações com dados de entrada

inexatos é utilizar a aritmética intervalar, também chamada de análise intervalar. Neste

contexto, uma incerteza em um dado real pode ser representada por um intervalo de

números reais, que contém o valor desconhecido exato do número em questão.

A primeira abordagem de aritmética intervalar começou com Burkill (1924).

Young (1931) realizou um novo estudo a respeito do assunto e mais tarde Sunaga

(1958) concretizou o uso da aritmética intervalar. O desenvolvimento moderno da

aritmética intervalar, com um enfoque global, começou com a tese de doutorado de

Moore (1962). Desde então, vários artigos e livros têm abordado este assunto.

Congressos, encontros, softwares e grupos de discussão vêm contribuindo para o

aumento do interesse e da importância deste tema.

Na análise intervalar, um subconjunto de números reais ℜ da forma X (ver

Equação (2.7)) é chamado de intervalo real fechado ou simplesmente de intervalo. O

19

símbolo (~) representa uma quantidade intervalar. Os limites do intervalo são

representado por: l , limite inferior do intervalo (lower bound) e, u, limite superior do

intervalo (upper bound). O conjunto de todos os intervalos reais fechados é dado por

( )I ℜ .

[ ; ] : , , l u l u l uX x x x x x x x x= = ≤ ≤ ∈ℜ (2.7)

Ponto médio de um intervalo.

2

l uc x x

X+

= (2.8)

Valor radial: representa a incerteza (ou o máximo erro) de um intervalo.

[ ], 1;12

u lR x x

X ε ε −

= = −

(2.9)

Um número intervalar também pode ser representado pelo ponto médio e pelo

raio, ao invés do limite superior e inferior. Esta representação é chamada formulação

centrada e é assim escrita:

C RX X X= + (2.10)

Expressões similares de intervalos, ponto médio e raio, podem ser definidos para

vetores e matrizes intervalares:

2.2.1.1 – Operações básicas

Somente algumas regras algébricas válidas para os números reais permanecem

válidas para números intervalares; algumas regras são respeitadas de forma “fraca”, ou

seja, a igualdade das regras não é verdadeira, valendo apenas as inclusões. Há duas

regras básicas em aritmética intervalar:

20

1. Duas expressões aritméticas equivalentes em aritmética real são equivalentes

em aritmética intervalar, se e somente se as variáveis da expressão aparecem apenas

uma vez de cada lado da igualdade.

2. Se f (x), g(x) são duas expressões aritméticas equivalentes em aritmética real,

então a inclusão f (x)⊆ g(x) é verdadeira, se e somente se todas as variáveis da

expressão aparecem apenas uma vez em f(x) .

A segunda regra quer dizer que o lado esquerdo, f (x), da equação leva ao

intervalo correto e o lado direito, g(x), contém este intervalo. Esta propriedade implica

em formas “fracas” de várias regras aritméticas familiares em aritmética real.

As operações básicas da aritmética intervalar são definidas a seguir para os

seguintes números intervalares: [ ; ] , [ ; ]l u l uX x x Y y y= = .

Adição:

[ , ]+ = + + l l u uX Y x y x y (2.11(a))

Subtração:

[ , ]l u u lX Y x y x y− = − − (2.11.(b))

Multiplicação:

min , , , , max , , , × = l l l u u l u u l l l u u l u u

X Y x y x y x y x y x y x y x y x y (2.11.(c))

Divisão:

( )1 1 1

, , , , 0 ;l u

l u l u

u l u l

X x xX x x y y

Y Y y y y y

= × = × = ∉

(2.11.(d))

21

Observe que a subtração e a divisão de números intervalares não são as

operações inversas da adição e multiplicação, como acontece em operações com

números reais.

2.2.1.2 - Propriedades

2.2.1.2.1 - Propriedades comuns entre números reais e intervalares

A seguir, são apresentadas algumas propriedades comuns entre números reais e

números intervalares.

Comutativa:

X Y Y X× = × (2.12)

X Y Y X+ = + (2.12)

Associativa:

( ) ( ) , [ ; ]l uX Y Z X Y Z Z z z× × = × × = (2.12)

( ) ( )X Y Z X Y Z+ ± = + ± (2.12)

Elemento neutro

0 0X X X+ = + = (2.14)

1 1X X X× = × =

- Outras propriedades comuns:

( )X Y X Y Y X− = + − = − +

( )( )X Y XY− − =

( ) ( )X Y X Y XY− = − = −

22

2.2.1.2.2 - Propriedades não comuns entre números reais e intervalares

Neste item serão apresentadas as propriedades dos números intervalares que os

diferem dos números reais.

• Propriedade da inclusão:

Uma importante propriedade dos números intervalares é a inclusão isotônica

(Neumaier, 1990). Esta propriedade diz que o resultado do cálculo de expressões

intervalares irá sempre incluir o resultado adequado, isto é, obtém-se no máximo uma

superestimação do resultado intervalar apropriado. Matematicamente, a propriedade da

inclusão pode ser escrita como:

& , , , , , ,X Y U V X opU Y opV para op e X Y U e V⊆ ⊆ ⇒ ⊆ ∈ + − × ÷ ∈ ℜ

(2.14)

Para todo conjunto de números reais s, define-se hull do conjunto como o menor

intervalo que engloba todos os números deste conjunto.

[ ]inf , suphull s s s=

(2.15)

na qual:

inf = Infimum, é o maior limite inferior do conjunto s; é definido como uma

quantidade m, de modo que nenhum membro do conjunto seja menor do que m.

sup = Supremum, é o menor limite superior de um conjunto s, é definido como

uma quantidade m, de modo que nenhum membro do conjunto exceda m.

• Propriedade subdistributiva:

23

Uma diferença notável entre os números reais e os intervalares é observada na

propriedade distributiva ( )x y z x y x z× ± = × ± × que, em geral, não vale para números

intervalares, onde apenas garante-se que:

( ) , , , , ,l u l u l uX Y Z X Y X Z para x x y y e z z× ± ⊆ × ± × ∈ ℜ (2.16)

Exemplo da propriedade subdistributiva, seja [ 1;1], [0;1], [ 1;0]X Y Z= − = = − :

( ) ( ) [ 1;1]([0;1] [ 1;0]) [ 1;1]f X X Y Z= × + = − + − = − (2.17)

( ) [ 1;1][0;1] [ 1;1][ 1;0] [ 2;2]g X X Y X Z= × + × = − + − − = − (2.18)

De onde se conclui que: [ 1;1] [ 2;2]− ⊆ −

Figura 2 – Representação gráfica de g(x) e f(x).

A falta da propriedade distributiva é uma fonte frequente de superestimação e

deve-se tomar cuidado em operações de sistemas de equações intervalares. ;

• Propriedade do subcancelamento:

Esta propriedade causa vários problemas na solução de sistemas de equações

intervalares. Esta propriedade diz que a igualdade entre certas expressões não é

satisfeita, podendo apenas afirmar que uma expressão está inclusa na outra. As

principais relações de subcancelamento serão descritas a seguir.

24

( ) ( )

( )

( )

X Y X Z Y Z

X X Z

Y Y Z

− ⊆ + − +

×⊆

×

(2.19)

Exemplo:

Seja [1;2]X = ,

Subtração: 0, 0 ( )− ≠ ∈ − X X X X

[1;2] [1;2] [1;2] [ 1; 2] [ 1;1] 0 ( )X X− = + − − = − ⇒ ∈ − (2.20)

Divisão: 1, 1≠ ∈

X X

X X

[1;2] 1[ ;2] 1

[1;2] 2

X

X

= ⇒ ∈

(2.21)

2.2.1.2 – Cuidados com a aritmética intervalar

Além dos cuidados com as propriedades de subcancelamento e subdistributiva

da aritmética intervalar, deve-se ficar atento ao cálculo de expressões em que aparece

uma mesma variável mais de uma vez. Estas expressões podem gerar resultados errados

com superestimação sendo as principais limitações na aplicação da aritmética intervalar.

Por isso, tem-se a necessidade de desenvolver algoritmos a fim de calcular o intervalo

da solução de um problema com parâmetros intervalares sem que problemas como o da

superestimação prejudique os resultados.

Por exemplo, considerando a função a seguir:

1

( )1

1

f x

x

=

+

para x ≠ 0 (2.22)

Na aritmética convencional, a Equação (2.22) é equivalente à fórmula

simplificada dada em Equação (2.23).

25

( )1

xg x

x=

+ (2.23)

que requer apenas uma divisão ao invés de duas para f (x), mas contém a variável x duas

vezes. Assim, fazendo x = [2,3] e usando aritmética intervalar para as duas expressões,

tem-se: ;

[ ]

[ ]

1 1 1 2 3( 2;3 ) ;

1 1 1 4 3 3 41 1 ; ;

2;3 3 2 3 2

f

= = = = + +

(2.24)

[ ][ ]

[ ]( )

[ ]

[ ][ ]

2;3 2;3 1 1 1( 2;3 ) 2;3 ; ;1

3;4 4 3 22;3 1g

= = = × =

+

Na Equação (2.24) observa-se a propriedade da inclusão ( ) ( )f x g x⊆ . Em g(x)

a variável x aparece duas vezes e, neste caso, g(x) gera como resultado um intervalo

maior que em f(x), confirmando uma superestimação inadequada, maior que o

necessário na expressão onde ocorre repetição de variável.

Figura 2.2 – Representação gráfica das funções g(x) e f(x).

A Figura (2.2) apresenta graficamente os resultados das duas expressões,

observe que g([2;3]) ≠ f ([2;3]) e que a expressão simplificada de g(x), onde a variável x

se repete, superestima de forma significante o valor de f (x).

Para evitar os problemas citados, deve-se realizar um trabalho prévio antes de

chegar na solução do sistema de equações intervalar. Este trabalho consiste no

entendimento total do problema que está sendo proposto e de uma análise criteriosa de

todas as operações intervalares que estão sendo realizadas. O entendimento do problema

26

orienta o caminho mais adequado para trabalhar com as operações intervalares, por

exemplo, a escolha do método numérico. A análise criteriosa das operações intervalares

é realizada para verificar a repetição de variáveis em expressões aritméticas e se o

resultado das expressões realmente representam o intervalo das soluções.

Na metodologia abordada neste estudo, com base nas operações matemáticas e

nas propriedades intervalares citadas, é possível que um problema de autovalor

intervalar não determinístico seja reformulado e convertido em um problema pseudo-

determinístico sujeito a uma matriz de perturbação. Nesta formulação, busca-se

primeiramente colocar os coeficientes intervalares como um multiplicador da matriz,

facilitando as manipulações matemáticas com números intervalares. A partir daí, a

substituição das matrizes intervalares por matrizes de valores médios e por matrizes

radiais converte o problema, antes não determinístico, em um problema pseudo-

determinístico.

Sendo assim, na próxima seção, descreve-se sobre a teoria da matriz de

perturbação aplicada a problemas dinâmicos de autovalores com parâmetros incertos.

2.2.2 Teoria da perturbação de matrizes aplicada a problemas

de autovalores com parâmetros intervalares

O problema linear clássico de autovalor para uma matriz simétrica

([ ] [ ] )TA A= de ordem p é dado por:

[ ]A x xλ= (2.25)

no qual λ são os autovalores e x são os autovetores associados.

Pré-multiplicando a Equação (2.25) por T

x , obtém-se:

[ ] λ=T T

x A x x x (2.26)

27

Rearranjando a Equação (2.26), tem-se:

[ ] λ=T T

x A x x x (2.27)

na qual, a norma Euclidiana de um vetor x poderá ser escrita de acordo com Equação

(2.28).

T

x x x= (2.28)

A Equação (2.27) pode ser transformada em uma razão conhecida como

quociente de Rayleigh, conforme Equação (2.29).

[ ]

( ) =

T

T

x A xR x

x x (2.29)

O quociente de Rayleigh serve para caracterizar o maior ou o menor autovalor

de uma matriz simétrica.

A matriz simétrica é obtida por meio da fatoração espectral, conforme Equação

(2.30) e a mudança de variáveis da Equação (2.29) é obtida através da Equação (2.31).

[ ] [ ] [ ][ ]=T

A Q D Q (2.30)

[ ] T

y Q x= (2.31)

A substituição das Equações (2.30) e (2.31) na Equação (2.29) obtém-se:

[ ] [ ][ ]

[ ] [ ]( )

[ ] [ ][ ]

[ ] [ ]

[ ]

( ) = = =

T T T T T

T T TT T

x Q D Q x x Q D Q x y D yR x

x Q x Q y yx x Q Q (2.32)

28

A Equação (2.32) é o quociente de Rayleigh transformado em uma base

principal com uma matriz ortogonal [ ]Q , matriz de autovetores. A matriz [ ]D

representa a matriz diagonal dos autovalores.

Fazendo as multiplicações matriciais, vem:

[ ]

2 2 2

1 1 2 2

2 2 2

1 2

...( )

...

T

p p

T

p

y y yy D yR x

y y yy y

λ λ λ+ += =

+ + (2.33)

2 2 22 21 1 2

1

12 2 2

1 2

...

( )...

p

p

p

y y y

R xy y y

λ λλ

λ λλ

+ +

= ≥

+ + (2.34)

2 2 21 21 2

2 2 2

1 2

...

( )...

p p

p p

p

p

y y y

R xy y y

λ λλ

λ λλ

+ +

= ≤

+ + (2.35)

Então, o quociente de Rayleigh, para uma matriz simétrica, é limitado entre o

menor, 1λ , e o maior autovalor, p

λ

1 ( )p

R xλ λ≤ ≤ (2.36)

O primeiro autovalor 1λ pode ser obtido por meio de uma minimização irrestrita

do quociente de Rayleigh, ( )R x , Equação (2.37).

[ ]

1

0 0

min ( ) minp p

T

Tx R x R

x x

x A xR x

x xλ

∈ ∈

≠ ≠

= =

(2.37)

A Equação (2.37) pode ser usada para determinar a primeira frequência natural o

primeiro autovalor do sistema. Para isso, selecionamos um vetor de teste x para

representar o primeiro modo de vibração 1x e o substituímos no quociente de

29

Rayleigh. Como o quociente é estacionário, podemos obter com boa exatidão o primeiro

autovalor, ainda que o desvio entre o vetor de teste x e o verdadeiro modo de

vibração 1x seja grande. Obviamente, que o valor da frequência natural será mais

preciso se o vetor de teste x escolhido for próximo ao modo natural verdadeiro 1x .

A fim de obter os autovalores intermediários, restrições adicionais devem ser

impostas no problema de minimização. O segundo autovalor pode ser determinado pela

imposição de uma restrição única, por exemplo, o vetor x deve ser perpendicular a

um vetor arbitrário s . Assim como o modo de vibração x , o vetor de teste s deve

ser escolhido, arbitrariamente, a fim de impor restrições ao quociente de Rayleigh para

que os demais autovalores sejam obtidos. A imposição de perpendicularidade garante

que o produto interno dos vetores seja nulo ( )0T

x s = . Com esta restrição, tem-se

um problema de minimização restrita, cujo limite superior é o segundo menor autovalor

2λ . Assim, adicionando s no problema, obtém-se:

1 20

min ( ) ( )T

x sR x sµ λ

=

= ≤ (2.38)

A Equação (2.38) pode ser provada considerando o vetor x como uma

combinação linear, diferente de zero, do primeiro e do segundo autovetores

normalizados.

1 1 2 2x x xα α= + (2.39)

na qual, x é ortogonal a s , o que impõe uma condição única em 1α e 2α . Para

qualquer combinação linear dos dois primeiros autovetores, tem-se:

[ ]

2 21 1 2 2 1 1 2 2 1 1 2 2

22 2

1 1 2 2 1 1 2 2 1 2

( ) ( )( )

( ) ( )

T

T

x x A x xR x

x x x x

α α α α λ α λ αλ

α α α α α α

+ + += = ≤

+ + + (2.40)

30

Portanto, conforme Equação (2.40), 2λ é obtido com a minimização de ( )R x ,

com 0x ≠ e sujeito a restrição única ( )0T

x s = , na qual o vetor arbitrário s ,

adotado, deve ser capaz de maximizar os mínimos de ( )R x para obter o segundo menor

autovalor, 2λ , conforme Equação (2.41).

[ ]2 max min ( )λ = R x (2.41)

Sujeito à restrição 0T

x s =

Generalizando a Equação (2.41), obtêm-se os demais autovalores aplicando as

restrições adicionais 0, 0T

ix s x= ≠ a ( )R x , na qual 1,..., 1= −i k ,

sendo 2≥k .

max min ( )k R xλ =

(2.42)

(Sujeito as restrições ( )0 , 1,..., 1, 2T

ix s i k k= = − ≥ )

O princípio que rege as Equações (2.41) e (2.42) é o teorema Max-Min (Teschl,

2009) de autovalores para matrizes simétricas.

2.2.2.1 Matriz de perturbação simétrica, positiva semidefinida

Para uma matriz simétrica [ ]A ficará sujeita a uma matriz de perturbação

simétrica, positiva e semidefinida [ ]E que goza da seguinte propriedade.

[ ] ( )0T

x E x ≥ (2.43)

31

Comparando o quociente de Rayleigh da matriz simétrica [ ]A com o quociente

da matriz simétrica perturbada [ ] [ ]A E + , tem-se:

[ ] [ ]

[ ]

T T

T T

x A E x x A x

x x x x

+ ≥ (2.44)

O primeiro autovalor da matriz perturbada pode ser obtido por meio da

minimização irrestrita da Equação (2.44) conforme Equação (2.45).

[ ] [ ]

[ ]

0 0

min min p p

T T

T Tx R x R

x x

x A E x x A x

x x x x∈ ∈

≠ ≠

+ ≥

(2.45)

Portanto,

[ ] [ ]( ) [ ]( )1 1λ λ+ ≥

A E A (2.46)

Os demais autovalores de matrizes perturbadas são dados por:

[ ] [ ]

[ ]

1 10 01,..., 1 1,..., 1

0 0

max min max min T T

T T

T Tx s x s

i k i kx x

x A E x x A x

x x x x= =

= − = −

≠ ≠

+ ≥

(2.47)

Logo,

[ ] [ ]( ) [ ]( )λ λ+ ≥

k kA E A (2.48)

Portanto, os autovalores de uma matriz simétrica sujeita a uma perturbação

positiva semidefinida aumentam monotonicamente a partir dos autovalores da matriz

exata. Similarmente, os autovalores de uma matriz simétrica sujeitos a uma perturbação

negativa e semidefinida decrescem monotonicamente, a partir dos autovalores da matriz

exata.

32

[ ] [ ]( ) [ ]( )λ λ− ≤

k kA E A (2.49)

Este conceito é conhecido como um “comportamento monotônico” dos

autovalores de matrizes simétricas sujeitas à matrizes de perturbação simétrica positiva

(ou negativa) e semidefinida (Bellman, 1960).

2.2.3 Representação intervalar das matrizes de rigidez global e

massa global

Os problemas dinâmicos discutidos neste trabalho são resolvidos pelo Método

dos Elementos Finitos. A matriz de rigidez global [ ]K e a matriz de massa global

[ ]M são obtidas por meio das matrizes de rigidez dos elementos discretizados em

elementos finitos, conforme desenvolvido no Anexo A.

A análise de problemas dinâmicos em vibração livre consiste, basicamente, em

resolver a equação diferencial dada a seguir, Equação (2.50).

[ ] [ ] 0 0( ) ( ) 0 , (0) , (0) (0)M z t K z t z z z z z+ = = = (2.50)

sendo,

( )z t : vetor de deslocamento, m

( )z t : vetor de velocidade, m/s

( )z t : vetor de aceleração, m/s2

( )[ ][ , ]e e

i i i iK l u K = (2.61)

( )[ ][ , ]e e

i i i iM l u M = (2.62)

33

na qual e

iK é a matriz de rigidez intervalar do th

i elemento; e

iM é a matriz de

massa intervalar do thi elemento;

il é o limite inferior e

iu o limite superior do número

intervalar [ , ]l u , que pré-multiplica as matrizes determinísticas de rigidez [ ]e

iK e de

massa [ ]e

iM do elemento.

Trabalhando com o intervalo [ , ]l u como um multiplicador das matrizes dos

elementos é possível realizar as manipulações matemáticas necessárias sem alterar as

características físicas dos elementos, como frequência natural, modos de vibração e as

propriedades de positividade e simetria das matrizes, de rigidez e de massa. Esta forma

paramétrica deve ser utilizada para preservar fortemente os limites do intervalo e evitar

a superestimação do intervalo da resposta do sistema intervalar. A incerteza na rigidez

de cada elemento é assumida como sendo independente. Para uma substrutura com uma

incerteza global intervalar, as Equação (2.61) e (2.62) são usadas para montar a matriz

de rigidez desta subestrutura.

2.2.3.1 Matriz de rigidez intervalar

A matriz de rigidez global [ ]K da estrutura pode ser vista como uma soma

linear das matrizes de rigidez local, [ ]e

iK , de cada elemento, que compõe uma malha

discretizada por elementos finitos.

[ ]1

[ ][ ] [ ]p

e T

i i i

i

K L K L=

=∑ (2.63)

na qual, [ ]i

L é a matriz de conectividade do elemento e [ ]e

iK é a matriz de rigidez do

elemento.

As entradas na matriz de conectividade dependem da posição de um elemento

em particular dentro da malha de elementos finitos. Em geral, a matriz de conectividade

[ ]i

L é formada por elementos iguais a 0 ou 1. Estas matrizes são essenciais e

desempenham um papel importante no desenvolvimento numérico do sistema, já que

através da multiplicação entre matrizes faz a montagem da matriz global (Bathe, 1976).

A seguir é apresentado um exemplo de montagem da matriz global usando as matrizes

34

1[ ]L e 2[ ]L . Adota-se como exemplo um elemento de viga discretizado em dois

elementos de comprimentos iguais a L, conforme Figura 2.8.

Figura 2.8 – Elemento de viga com dois elementos de comprimentos iguais a L

1 1 2 2v vθ θ

2 2

1 3

2 2

4 4

12 6 12 6

6 4 6 2[ ]

12 6 12 6

6 2 6 4

e

L L

L L L LEIK

L LL

L L L L×

−

− = − − −

−

2 2 3 3v vθ θ

2 2

2 3

2 2

4 4

12 6 12 6

6 4 6 2[ ]

12 6 12 6

6 2 6 4

e

L L

L L L LEIK

L LL

L L L L×

−

− = − − −

−

1

6 4

1 0 0 0

0 1 0 0

0 0 1 0[ ]

0 0 0 1

0 0 0 0

0 0 0 0

L

×

=

(2.64)

2

6 4

0 0 0 0

0 0 0 0

1 0 0 0[ ]

0 1 0 0

0 0 1 0

0 0 0 1

L

×

=

(2.65)

[ ]2

6 4 4 4 6 4 1 1 1 2 2 26 61

[ ] [ ] [ ] [ ][ ] [ ] [ ][ ] [ ]e T e T e T

i i i

i

K L K L L K L L K L× × ××

=

= = +∑ (2.66)

35

1 1 2 2 3 3v v vθ θ θ

( ) ( )

( ) ( )

2 2

1

2 2 23 2

22

6 6

12 6 12 6

6 4 6 2

12 6 12 12 6 6 612[ ]

4 4 266 626

12 612 6

6 46 2

L L

L L L L

L L L LEIK

L L LLL L LLL

LL

L LL L×

−

− − − + − + − =

+ −− +

−− − −

1

1

2

2

3

3

v

v

v

θ

θ

θ

(2.67)

2.2.3.1.1 Matriz de rigidez global intervalar

Conforme mencionado, uma matriz de rigidez global determinística pode ser

vista como uma soma linear das matrizes de rigidez intervalar (não determinístico) de

cada elemento, de uma estrutura discretizada por elementos finitos. A montagem da

matriz global intervalar na presença de parâmetros intervalares é realizada de tal modo

que as características físicas dos elementos e as propriedades matemáticas das matrizes

são preservadas.

Do mesmo modo que se obtém a Equação (2.63) para [ ]K , obtém-se também, a

Equação (2.68) para K com parâmetros intervalares.

[ ] [ ]1

pe T

i i i

i

K L K L=

= ∑ (2.68)

[ ] ( )[ ] [ ]1

[ , ]p

e T

i i i i i

i

K L l u K L=

= ∑ (2.69)

( )[ ][ ] [ ]1

[ , ]p

e T

i i i i i

i

K l u L K L=

= ∑ (2.70)

( )1

[ , ]p

e

i i i

i

K l u K=

= ∑ (2.71)

36

na qual e

iK é a matriz de rigidez determinística do elemento posicionada na

coordenada global.com dimensão igual a da matriz global.

2.2.3.2 Matriz de massa intervalar

Assim como descrito para a matriz de rigidez global, a matriz de massa global

determinística da estrutura pode ser vista como uma soma linear das matrizes de massa

global de cada elemento que compõe uma malha de elementos finitos.

[ ]1

[ ][ ] [ ]p

e T

i i i

i

M L M L=

=∑ (2.72)

na qual [ ]e

iM representa a matriz de massa do elemento no sistema de coordenada

global.

A matriz de massa global intervalar é dada por:

1

[ ][ ] [ ]p

e T

i i i

i

M L M L=

= ∑ (2.73)

[ ] ( )[ ] [ ]1

[ , ]p

e T

i i i i i

i

M L l u M L=

= ∑ (2.74)

( )[ ][ ] [ ]1

[ , ]p

e T

i i i i i

i

M l u L M L=

= ∑ (2.75)

( )1

[ , ]p

e

i i i

i

M l u M=

= ∑ (2.76)

na qual e

iM é a matriz de massa determinística do elemento posicionada na

coordenada global com dimensão igual a da matriz global.

37

2.2.3.2.1 Obtenção da matriz de massa de um elemento de

barra inclinado de corpos em vibração

Segundo Hutton (2004) o elemento de barra inclinado sob vibração mecânica

deve ser tratado de maneira diferente à trabalhada em corpos estáticos. Quando a

estrutura está em vibração a matriz de massa do elemento deve ser desenvolvida a partir

de elementos com deslocamentos nodais nas direções axial e transversal, conforme

Figura (2.9).

De acordo com Kassimali (2011) no método dos elementos finitos, as relações

envolvidas na obtenção das matrizes de rigidez e de massa do elemento são baseadas

assumindo variações de deslocamento nos elementos. Tal variação de deslocamento é

referida como função de deslocamento ou função de interpolação. Uma função de

deslocamento descreve a variação da componente do deslocamento ao longo do eixo

centroidal do elemento.

Considere um elemento de barra de uma treliça plana submetido aos

deslocamentos 1 2 3 4, ,u u u eu , conforme Figura (2.9). O elemento analisado se desloca

nas direções x e y, e neste caso, é necessário definir duas funções de deslocamento xu e

yu , conforme mostrado na Figura (2.9). Sendo que os deslocamentos em x e y do nó b

são 1u e 2u , respectivamente, e os deslocamentos no nó e são 3u e 4u . As funções de

deslocamentos xu e

yu são as coordenadas que localizam o deslocamento do centroide

G do elemento.

38

Figura 2.9 – Deslocamento do elemento no sistema de coordenadas locais (fonte: Matrix

Analysis of Structures, Kassimali, 2011)

As funções de deslocamentos podem ser dadas na forma de polinômios,

conforme Equação (2.77).

0

( ) , 0n

i i i

i

u x a x a=

= ≠∑

(2.77)

Figura 2.10 – Deslocamento do elemento na direção x (fonte: Matrix Analysis of

Structures, Kassimali, 2011)

A função xu , ver Figura (2.10), para o elemento de uma treliça plana é dada na

forma de um polinômio linear.

0 1( )u x a a x= +

(2.78)

na qual 0a e 1a são constantes determinadas pela aplicação de duas condições de

contorno:

39

Em 0x = 1xu u=

(2.79)

Em x L= 3xu u=

(2.80)

Tem-se então que:

0 1a u=

(2.81)

3 13 1 1 1

u uu u a L a

L

−= + → =

(2.82)

Substituindo as Equações (2.81) e (2.82) na (2.78), obtém-se:

3 11x

u uu u x

L

− = +

(2.83)

ou

1 31

x

x xu u u

L L

= − +

(2.84)

Figura 2.11 – Deslocamento do elemento na direção y (fonte: Matrix Analysis of

Structures, Kassimali, 2011)

A função de deslocamento yu , para o deslocamento na direção y, ver Figura

(2.11), pode ser determinada de maneira similar, isto é, por meio de um polinômio

linear, conforme Equação (2.85).

2 41

y

x xu u u

L L

= − +

(2.85)

40

Funções de interpolação

As funções de deslocamento (ou funções de interpolação) apresentadas na

Figura (2.12), desenvolvida nas Equações (2.83) e (2.84) podem ser expressos

alternativamente por (2.86) e (2.87), respectivamente.

1 1 3 3xu H u H u= + (2.86)

2 2 4 4yu H u H u= + (2.87)

na qual 1 2 1x

H HL

= = − e 3 4

xH H

L= = , conforme indicado na Figura (2.12)

Na forma matricial fica:

1

1 3 2

32 4

4

( )

( , ) ( ) 0 ( ) 0 ( )

( , ) ( )0 ( ) 0 ( )

( )

x

y

u t

u x t H x H x u t

u x t u tH x H x

u t

=

(2.88)

De maneira simplificada:

[ ] ( , ) ( ) ( )T

u x t H x u t= (2.89)

Figura 2.12 – Funções de interpolação do elemento de uma treliça plana (fonte: Matrix

Analysis of Structures, Kassimali, 2011).

41

Energia Cinética

De acordo com a energia cinética do elemento, a matriz massa de um elemento

inclinado de uma treliça plana pode ser dada por:

[ ][ ]

[ ]

2

0

0

0

1 ( , )( ) ( )

2

1( ) ( ) ( ) ( ) ( )

2

1( ) ( ) ( )

2

1( ) ( )

2

L

L T T

L T

t e

u x tT t m x dx

t

m x u t H x H x u t dx

m x u t u t dx

u t M u t

∂ =

∂

=

=

=

∫

∫

∫

(2.90)

na qual

[ ] [ ] [ ]0

( ) ( ) ( )Le T

M m x H x H x dx= ∫ (2.91)

[ ]0

01

1 1 0 00

00 1 0

0

Le

x

L x x x

L L LM A dx

x x x

L L Lx

L

ρ

−

− − = −

∫ (2.92)

Após multiplicar as matrizes [ ]( )H x e [ ]( )T

H x e integrar a matriz resultante de

0 a L obtém-se a matriz dada pela Equação (2.93).

[ ]

2 0 1 0

0 2 0 1

1 0 2 06

0 1 0 2

e ALM

ρ

=

(2.93)

na qual ρ é a nassa específica, A é a área da seção transversal e L é o comprimento do

elemento de barra. Neste caso, L possui valor fixo, válido somente em casos onde o

elemento de barra não mude de direção.

42

A Figura (2.13) apresenta as funções de interpolação de um elemento de barra

inclinado usado na transformação do sistema de coordenadas local no sistema global.

Figura 2.13 – Funções de interpolação do elemento de uma treliça plana (fonte: Matrix

Analysis of Structures, Kassimali, 2011).

Os elementos de barra inclinados possuem comprimentos diferentes dos

elementos na vertical ou horizontal e, neste caso, devem-se fazer algumas considerações

de maneira que os comprimentos dos elementos fiquem em função de suas coordenadas

nodais 1 2 3 4, ,X Y X e Y , conforme Figura (2.14). Portanto, o comprimento do elemento

de barra da Equação (2.93) após a reformulação é dado na Equação (2.94).

Figura 2.14– Deslocamentos 1 2 3 4, ,v v v e v de um elemento de barra inclinado (adaptado

de Matrix Analysis of Structures, Kassimali, 2011).

43

( ) ( )2 2

4 2 3 1L Y Y X X′ = − + − (2.94)

4 23 1

3 1

, ( ) 0Y Y

arct com X XX X

θ −

= − ≠ −

(2.95)

Para 3 1 0X X− = e 4 2Y Y> , então:

090θ =

Para 3 1 0X X− = e 4 2Y Y< , então:

090θ = −

Agora, a matriz [ ]e

M é dada por Equação (2.96).

[ ]

[ ] ( ) ( )2 2

4 2 3 1

2 0 1 0

0 2 0 1

1 0 2 06

0 1 0 2

2 0 1 0

0 2 0 1

1 0 2 06

0 1 0 2

e

e

ALM

AM Y Y X X

ρ

ρ

′′ =

′ = − + −

(2.96)

De acordo a Figura (2.13) obtêm-se as expressões da Equação (2.97).

1 1 2

2 1 2

3 3 4

4 3 4

cos

cos

cos

cos

F Q Q sen

F Q sen Q

F Q Q sen

F Q sen Q

θ θ

θ θ

θ θ

θ θ

= −

= +

= −

= +

(2.97)

Na forma matricial:

44

1 1

2 2

3 3

4 4

cos 0 0

cos 0 0

0 0 cos

0 0 cos

F Qsen

F Qsen

F Qsen

F Qsen

θ θ

θ θ

θ θ

θ θ

− =

−

(2.98)

ou

[ ] F T Q= (2.99)

na qual [ ]T é a matriz de transformação.

A matriz de massa de um elemento de barra inclinado no sistema global de

coordenadas é dado por:

[ ] [ ]e t e

M T M T ′ = (2.100)

sendo [ ]e

M ′ a matriz de massa de um elemento dada na Equação (2.96).

2 2 2 2

2 2 2 2

2 2 2 2

2 2 2 2

2cos 2 0 cos 0

0 2cos 2 0 cos

6 cos 0 2cos 2 0

0 cos 0 2cos 2

e

sen sen

sen senALM

sen sen

sen sen

θ + θ θ + θ

′ θ + θ θ + θρ = θ + θ θ + θ

θ + θ θ + θ

ou

( ) ( )

2 2 2 2

2 2 2 22 2

4 2 3 1 2 2 2 2

2 2 2 2

2cos 2 0 cos 0

0 2cos 2 0 cos

6 cos 0 2cos 2 0

0 cos 0 2cos 2

e

sen sen

sen senAM Y Y X X

sen sen

sen sen

θ + θ θ + θ

θ + θ θ + θρ = − + − θ + θ θ + θ

θ + θ θ + θ

(2.101)

45

2.2.3.3 Problema de autovalor intervalar para estruturas

dinâmicas

2.2.3.3.1 Problema de autovalor

Em estruturas dinâmicas, a equação de equilíbrio para a vibração livre de um

sistema com múltiplos graus de liberdade é definida como um conjunto de equações

diferenciais lineares ordinárias de segunda ordem com coeficientes constantes,

conforme Equação (2.102).

[ ] [ ] ( ) ( ) 0M z t K z t+ = (2.102)

Condições iniciais:

0

0

(0)

(0)

z z

z z

=

=

(2.103)

Assumindo um movimento harmônico para o deslocamento temporal ( )z t ,

substituindo a Equação (2.104) na (2.77) seguida de manipulações algébricas

pertinentes, obtém-se:

( ) pi tz t x e

ω= (2.104)

[ ] [ ] 2 0pK x M xω− = (2.105)

[ ] [ ] 2

pK x M xω= (2.106)

Para 2

pω λ= ,

[ ] [ ]K x M xλ= (2.107)

46

[ ]( ) [ ] 0K M xλ− = (2.108)

Pré-multiplicando a Equação (2.108) por [ ]1−

M

[ ] [ ]( ) [ ] 1 1

[ ] 0M K M x Mλ− −

− = (2.109)

[ ] [ ]( ) 1[ ] 0M K I xλ

−

− = (2.110)

Em sistemas dinâmicos, o problema de autovalor generalizado é dado pela

Equação (2.110). As frequências naturais são representadas por nω e os modos de

vibração são dados pelo vetor x .

Seja A uma matriz real intervalar e [ ] [ ]

1[ ]

−

=A K M um membro da matriz

intervalar ([ ]A A ∈ ) ou em termos das componentes ( ij ija a∈ ), então, o problema de

autovalor intervalar é dado por:

[ ] [ ]( ) [ ]( )0 ,A I x A Aλ − = ∈ (2.111)

2.2.3.3.2 Solução do problema de autovalor

A solução de um problema de autovalores intervalares é definida como um

conjunto abrangente de valores reais λ , tal que, para qualquer membro da matriz

intervalar, a solução é um membro do conjunto de solução. Portanto, a solução para o

problema de autovalores intervalares pode ser expresso matematicamente como:

[ ] [ ] [ ]( ) , | : 0l uA A A I xλ λ λ λ λ ∈ = ∀ ∈ − =

(2.112)

47

2.2.3.3.3 Problema de autovalor submetido à matriz de

perturbação

O problema de autovalor intervalar generalizado para estruturas dinâmicas pode

ser obtido fazendo a substituição das matrizes de rigidez e massa global intervalar

(Equações (2.71) e (2.76)) no problema de autovalor generalizado, Equação (2.106),

obtém-se então:

2( )pK x M xω =

(2.113)

( ) ( ) 2

1 1

[ , ] ( ) [ , ]p p

e e

i i i p i i i

i i

l u K x l u M xω= =

=

∑ ∑ (2.114)

na qual x representa os modos de vibração intervalar e pω , a frequência natural

intervalar.

O problema de autovalor intervalar, Equação (2.114), pode ser transformado em

um problema pseudo-determinístico, sujeito a uma matriz de perturbação, através da

introdução da matriz de rigidez e massa central e da matriz de rigidez e massa de

perturbação radial, como apresentado nas Equações (2.115) e (2.116).

1

1

2

2

peC i i

i

i

peC i i

i

i

u lK K

u lM M

=

=

+ =

+ =

∑

∑

(2.115)

na qual CM e CK representam as matrizes de valores médios.

[ ]1

1

( ) , 1,12

( )2

peR i i

i i i

i

peR i i

i i

i

u lK K

u lM M

ε ε

ε

=

=

− = = −

− =

∑

∑

(2.116)

na qual RK e RM

representam as matrizes radiais.

48

Deve-se notar que a teoria da perturbação estuda o comportamento de um

sistema sujeito a pequenas alterações em torno das suas variáveis de projeto original.

Neste caso, para expressar os problemas intervalares usando a teoria da perturbação,

adiciona-se [ 1,1]ε = −i

as matrizes radiais, que representa pequenas alterações sofridas

pelos parâmetros de projeto. Conforme mostrado na Equação (2.10), assim como os

números intervalares, as matrizes intervalares podem ser apresentadas em função da

soma das matrizes de valores médios com as matrizes radiais, conforme Equações

(2.117) a (2.128).

Para a matriz de rigidez intervalar tem-se que:

C RK K K = + (2.117)

1

( )2

peC i i

i i

i

u lK K Kε

=

− = + ∑ (2.118)

[ ]1

1,12

peC i i

i

i

u lK K K

=

− = + − ∑ (2.119)

1 1

,2 2

p pe eC i i i i

i i

i i

u l u lK K K K

= =

− − = + − ∑ ∑ (2.120)

De acordo com a aritmética intervalar, Equação (2.11), tem-se:

1 1

,2 2

n ne eC Ci i i i

i i

i i

u l u lK K K K K

= =

− − = − + ∑ ∑ (2.121)

1 2

nel C i i

i

i

u lK K K

=

− = −

∑ (2.122)

1 2

peu C i i

i

i

u lK K K

=

− = +

∑ (2.123)

49

,l uK K K =

(2.124)

na qual lK é a matriz formada pelos valores de limite inferior e

uK é a matriz

formada pelos valores de limite superior dos números intervalares.

Seguindo a mesma linha de raciocínio da matriz de rigidez intervalar, a matriz

de massa intervalar pode também ser representada pela soma da matriz de perturbação

radial com a matriz de valores médios.

C RM M M = + (2.125)

1

( )2

peC i i

i i

i

u lM M Mε

=

− = + ∑ (2.126)

1 1

,2 2

p pe eC Ci i i i

i i

i i

u l u lM M M M M

= =

− − = − + ∑ ∑ (2.127)

, = l u

M M M (2.128)

na qual lM é a matriz formada por valores de limite inferior e uM é a matriz

formada por valores de limite superior dos números intervalares.

Substituindo as Equações (2.117) e (2.125) na Equação (2.113) obtém-se:

( ) ( ) 2[ ] [ ] [ ] [ ]C R C R

pK K x M M xω+ = + (2.129)

Assim, conforme Equação (2.129), a determinação dos limites das frequências

naturais e dos modos de vibração na presença de incertezas é interpretada,

matematicamente, como um problema de autovalor centrado nas matrizes

([ ] [ ]C CK e M ) e sujeito a perturbações radiais representados pelas matrizes

50

([ ] [ ]R RK e M . Esta perturbação é obtida pelo somatório linear das matrizes

determinísticas do elemento ([ ] [ ]e e

i iK e M ) multiplicada por um coeficiente intervalar

iε .

2.2.3.3.4 Frequência natural intervalar

A fim de se obter resultados físicos corretos é necessário que sejam impostas

restrições ao problema de autovalor intervalar. Estas restrições estão presentes,

intrinsecamente, no problema de autopar intervalar. Elas resultam em matrizes de

perturbação radial ([ ] [ ]R RK e M ), que são combinações lineares de matrizes simétricas

positivas dimensionadas por números reais intervalares. Essa característica da matriz de

perturbação radial deve ser considerada no desenvolvimento de qualquer problema de

autovalor intervalar para limitar as frequências naturais.

( )

,λ

=

T

T

x K xK M

x M x (2.130)

( )

[ ] [ ],

[ ] [ ]

T C R

T C R

x K K xK M

x M M xλ

+ = +

(2.131)

( )

,,

,λ

=

T Tl u

T Tl u

x K x x K xK M

x M x x M x (2.132)

Trabalhando com o conceito de caracterização de mínimos e máximos de

autovalores compostos por matrizes simétricas, Equações (2.37) e (2.42), pode-se obter

os autovalores de um sistema dinâmico intervalar, conforme demonstrado nas Equações

(2.133) e (2.135)

51

Para o primeiro autovalor:

( )

1

0

[ ] [ ], min

p

T C R

T C Rx R

x

x K K xK M

x M M xλ

∈

≠

+ = +

(2.133)

( )

1

0

,, min

,p

T Tl u

T Tl ux R

x

x K x x K xK M

x M x x M xλ

∈

≠

=

(2.134)

Para os próximos autovalores:

( )

01,..., 1

0

[ ] [ ], max min

T

i

T C R

k T C Rx s

i kx

x K K xK M

x M M xλ

=

= −

≠

+ = +

(2.135)

( )

01,..., 1

0

,, max min

,T

i

T Tl u

k T Tl ux s

i kx

x K x x K xK M

x M x x M xλ

=

= −

≠

=

(2.136)

Substituindo e expandindo os termos do lado direito das Equações (2.133) e

(2.135), tem-se:

1 1

1 1

[ ] ( ) [ ][ ] [ ] 2 2

[ ] ( ) [ ]2 2

p pT Te ei i i iT C R i i i

i i

p pT C RT Te ei i i i

i i i

i i

u l u lx K x x K x

x K K x

u l u lx M M x x M x x M x

ε

ε

= =

= =

+ − + +

=+ − + +

∑ ∑

∑ ∑

(2.137)

Uma vez que [ ]CK é uma matriz simétrica positiva semidefinida

( [ ] 0T C

x K x ≥ ), então, [ ]CK estará sujeita a uma matriz de perturbação simétrica

52

positiva semidefinida RK , ( )0

T Rx K x ≥

. Já a matriz [ ]CM , simétrica positiva

definida, ( [ ] 0T C

x M x > ), estará sujeita a uma perturbação também simétrica

positiva definida RM , ( )0

T Rx M x >

.

De acordo com as propriedades da matemática intervalar , Equação (2.14), os

limites dos autovalores são obtidos por:

max , maxλ

=

T u

k T l

x K xK M

x M x (2.138)

min , minλ

=

T l

k T u

x K xK M

x M x (2.139)

( )

1 1

1 1

[ ] [ ]2 2

max ,

[ ] [ ]2 2

p pT Te ei i i i

i i

i i

k p pT Te ei i i i

i i

i i

u l u lx K x x K x

K Mu l u l

x M x x M x

λ= =

= =

+ − +

= + − −

∑ ∑

∑ ∑ (2.140)

( )

1 1

1 1

[ ] [ ]2 2

min ,

[ ] [ ]2 2

pnT Te ei i i i

i i

i i

k pnT Te ei i i i

i i

i i

u l u lx K x x K x

K Mu l u l

x M x x M x

λ= =

= =

+ − −

= + − +

∑ ∑

∑ ∑ (2.141)

Portanto, os limites (superior e inferior) das frequências naturais de um

problema de autovalor intervalar não determinístico podem ser obtidas por meio de um

problema de autovalor pseudo-determinístico, facilitando a solução do problema quando

na presença de parâmetros incertos, conforme Equações (2.142) e (2.143).

53

( ) ( ) 2

min

1 1

[ ] ( ) [ ]p p

e e

i i i i

i i

l K x u M xω= =

=∑ ∑ (2.142)

( ) ( ) 2

max

1 1

[ ] ( ) [ ]p p

e e

i i i i

i i

u K x l M xω= =

=∑ ∑ (2.143)

De acordo com as Equações (2.142) e (2.143), conclui-se que os limites

superiores exato das frequências naturais de um problema dinâmico, na presença de

incertezas, são obtidos trabalhando-se com os valores superiores dos elementos que

formam a matriz de rigidez juntamente com os valores inferiores dos elementos que

compõem as matrizes massa. Seguindo a mesma linha de raciocínio, os limites

inferiores das frequências são obtidos por meio de um problema de autovalor formado

por matrizes de rigidez de valores mínimos e matrizes de massa de valores máximos.

54

Capítulo 3

PROCEDIMENTO DE CÁLCULO DE INCERTEZA

DE MEDIÇÕES DIRETAS E INDIRETAS DAS

GRANDEZAS FÍSICAS DE UMA VIGA DE AÇO

Neste capítulo serão apresentadas noções básicas sobre a avaliação das incertezas

de medição das grandezas físicas de uma viga de aço. Com base no intervalo da incerteza

das grandezas calculadas, obtém-se o intervalo da incerteza das frequências naturais

inerentes a viga em estudo, quando em vibração livre.

Conforme mencionado, o objetivo deste trabalho é apresentar uma metodologia

numérica eficaz, com baixo custo computacional, capaz de gerar resultados intervalares a

partir de grandezas de entrada também intervalares. Portanto, é necessário que se obtenha

o intervalo de valores possíveis de serem atribuídas às grandezas físicas da viga de aço.

Neste estudo, compara-se o intervalo de valores que pode ser atribuído as frequências

naturais experimentais com o intervalo das frequências naturais calculadas

numericamente. O Guide to the expression of the Uncertainty in Mesureament (GUM,

2008) apresenta critérios e regras para expressar e combinar diferentes componentes de

incertezas de medição. Todo o procedimento de cálculo para se estimar as incertezas é

baseado no GUM (2008), o que nos permitiu estimar as incertezas presentes nas seguintes

55

grandezas da viga: comprimento total, comprimento útil, largura, altura, volume, área,

massa específica, módulo de elasticidade, momento de inércia de área e também permitiu

avaliar as incertezas presentes nos dados coletados pelos acelerômetros. A análise

probabilística é uma forma de caracterizar a incerteza tanto nos dados coletados pelos

sensores, como nos parâmetros do modelo numérico.

3.1 INCERTEZA DE MEDIÇÃO

A presença de incertezas nos projetos pode ser atribuída às imperfeições físicas,

imprecisões do modelo matemático e a complexidade do sistema e essas incertezas

devem ser quantificadas nas etapas de elaboração dos projetos de engenharia.

Segundo GUM (2008) o conceito de incerteza como um atributo quantificável é

relativamente novo na história da medição, embora o erro e análise do erro tenham uma

longa participação da prática da ciência da medição ou metrologia. É atualmente

reconhecido que, quando todos os componentes conhecidos ou suspeitos do erro tenham

sido avaliados e as correções apropriadas tenham sido aplicadas, há ainda uma incerteza

remanescente acerca da correção do resultado apresentado, isto é, uma dúvida acerca de

quão bem o resultado da medição representa o valor da quantidade sendo medida.

Na prática, existem muitas fontes possíveis de incerteza em uma medição, tais

como as listadas a seguir, (GUM, 2008):

a) definição incompleta do mensurando;

b) realização imperfeita da definição do mensurando;

c) amostragem não representativa - a amostra medida pode não representar o mensurando

definido;

d) conhecimento inadequado dos efeitos das condições ambientais sobre a medição ou

medição imperfeita das condições ambientais;

e) erro de tendência pessoal na leitura de instrumentos analógicos;

f) resolução finita do instrumento ou limiar de mobilidade;

g) valores inexatos dos padrões de medição e materiais de referência;

56

h) valores inexatos de constantes e de outros parâmetros obtidos de fontes externas e

usados no algoritmo de redução de dados;

i) aproximações e suposições incorporadas ao método e procedimento de medição;

j) variações nas observações repetidas do mensurando sob condições aparentemente

idênticas.

A origem da incerteza de medição pode ser atribuída à existência de erros

sistemáticos e aleatórios. Os erros sistemáticos geralmente são de causas identificáveis e

possíveis de serem corrigidos, portanto, em princípio, pode-se eliminar a influência

desses efeitos nas fontes de incertezas. Os erros aleatórios são o resultado de influências

externas e internas, não controladas, que provocam o surgimento de erros não repetitivos.

Os efeitos aleatórios podem ser quantificados por análise estatística.

A medição direta é aquela cuja indicação resulta naturalmente da aplicação do

sistema de medição sobre o mensurando. Há apenas uma grandeza de entrada envolvida.

Exemplos de medições diretas efetuadas neste trabalho: medição de comprimento,

largura, altura e massa. Já a medição indireta envolve a determinação do valor associado

ao mensurando a partir da combinação de duas ou mais grandezas por meio de expressões

matemáticas. Exemplo de medições indiretas realizadas neste trabalho: medição de

volume, área, massa específica, módulo de elasticidade e momento de inércia de área.

Na medição direta, os efeitos aleatórios de cada uma das fontes de incertezas

(comprimento, largura, altura, massa da viga e massa do acelerômetro) foram estimados

baseados em parâmetros estatísticos e não estatísticos. O efeito aleatório é estimado por

meio da incerteza padrão de cada fonte de incerteza, indicado pela faixa de dispersão em

torno do valor central equivalente a um desvio padrão. A avaliação da incerteza padrão

pode ser classificada em tipo A e tipo B. Ambos são baseados em distribuições de

probabilidade (normal, retangular, entre outras.). A avaliação tipo A é realizada por uma

análise estatística de uma série de observações ou medidas, sob as mesmas condições de

operação . Neste estudo, para uma avaliação do tipo A, efetuou-se 20 medições diretas

das seguintes grandezas: comprimento, largura, altura e massa. Já a avaliação do tipo B é

realizada por procedimentos não estatísticos, a partir de um julgamento científico baseado

em informações relevantes sobre os instrumentos usados durante o ensaio e sobre o

processo de medição. Neste trabalho foram utilizados os seguintes equipamentos de

57

medição: paquímetro, balança, miniacelerômetro e analisador de sinais. Para efetuar uma

avaliação do tipo B, coletaram-se informações de catálogos e certificados de calibração.

Segundo GUM (2008) uma completa avaliação da incerteza de medição pressupõe

a estimativa da incerteza padrão, da incerteza padrão combinada e da incerteza

expandida. A incerteza padrão u(xi) de um dado efeito aleatório corresponde à estimativa

equivalente a um desvio padrão da ação deste efeito sobre a indicação. A incerteza

combinada uc(xi) de um processo de medição é estimada considerando a ação simultânea

de todas as fontes de incerteza e ainda corresponde a um desvio padrão da distribuição

resultante. A incerteza expandida U(xi) associada a um processo de medição é estimada a

partir da incerteza combinada multiplicada por um fator de abrangência (k) apropriado e,

isto, reflete a faixa de dúvidas presente na medição para nível de confiança de 95%. Na

engenharia é comum trabalhar com níveis de confiança de 95% e para isso é necessário

que uc seja multiplicado por k.

Nos próximos tópicos serão apresentadas alguns procedimentos de avaliação de

incertezas de medição em medições diretas efetuadas no cálculo das grandezas:

comprimento, largura, altura e massa e, em medições indiretas realizadas no cálculo do

módulo de elasticidade, do momento de inércia de área, da massa específica e da área da

seção transversal da viga. Serão discutidos também, algumas da distribuições de

probabilidades mais usuais, tais como, normal, uniforme, trapezoidal e triangular.

3.1.1 Incerteza padrão

A medição de uma grandeza envolve um conjunto de operações para estimar o seu

valor verdadeiro. Esta grandeza pode ser medida diretamente (como o comprimento, a

altura, a massa, o tempo) ou indiretamente (como o volume, a velocidade, a massa

específica etc.). Conforme mencionado, em uma medição indireta supõem-se uma relação

matemática que permita determinar o valor de uma grandeza desconhecida a partir dos

valores de N outras grandezas conhecidas ( 1 2, ,...,N

X X X ), denominadas grandezas de

entrada, através de uma relação funcional f, como mostra a Equação (3.1).

58

1 2( , ,..., )N

Y f X X X= (3.1)

As grandezas de entrada, por sua vez, podem também ser consideradas

mensurandos que dependem de outras grandezas. Seus valores e respectivas incertezas

podem ser obtidos a partir de uma única observação ou de repetidas observações, de

dados fornecidos pelos fabricantes dos instrumentos, da experiência do observador, da

literatura, de medições realizadas anteriormente, de padrões de calibração, de materiais

de referência ou de certificados de calibração. Assim, a incerteza associada ao resultado

da medição deve levar em consideração as incertezas individuais que afetam o processo

de medição. A função f deve ser interpretada como aquela função que contém todas as

grandezas, incluindo todas as correções e fatores de correções que possam contribuir com

um componente significativo da incerteza para o resultado de medição (Santos, 2011).

Os componentes da incerteza podem ser analisados por diferentes modos de

avaliação como a do tipo A e do tipo B. Os dois tipos de avaliação são baseados em

distribuições de probabilidade e os componentes de incerteza resultantes de qualquer tipo

são quantificados por variâncias ou desvios padrão.

Avaliação do tipo A: é realizada estatisticamente a partir de medições repetidas

de uma dada grandeza, assumindo uma distribuição normal de probabilidade ou outra

qualquer.

Distribuição Normal ou Gaussiana: existem muitas funções de distribuição de

probabilidade, porém a distribuição normal representa um bom modelo para uma série de

dados experimentais em diversas áreas do conhecimento humano, conforme mostra a

curva da Figura (3.1).

Figura 3.1 – Curva de distribuição gaussiana.

59

A curva de distribuição normal tem o aspecto de um sino, sendo simétrica em

relação à média aritmética ( µ ). Pode-se afirmar com um nível de confiança de 68,27%

que a média está no intervalo µ σ± . Este nível aumenta para 95,45% no intervalo

2µ σ± e para 99,74% no intervalo 3µ σ± .

A média aritmética é uma medida de tendência central de um conjunto de dados

isolados ou amostral. Sabe-se que a média aritmética, Equação (3.2), é a melhor

estimativa disponível do valor esperado de uma grandeza X que varia aleatoriamente e

para a qual n observações independentes foram obtidas sob as mesmas condições de

medição.

1

ni

i

xx

n=

=∑ (3.2)

na qual x é a média aritmética, 1 2, ,...,n

x x x são as observações efetuadas, e n é o

número de observações repetidas.

A variabilidade dos valores das n observações realizadas da grandeza X ocorre

devido aos efeitos aleatórios, gerando um alargamento da curva de distribuição de

probabilidade e aumentando a imprecisão do método analítico. Assim, a variância

experimental, 2 ( )i

s x , das observações é dada por:

2 2

1

1( ) ( )

1

n

i i

i

s x x xn =

= −−∑ (3.3a)

A raiz quadrada da variância experimental 2 ( )i

s x é o desvio padrão experimental

dado por ( )i

s x .

A incerteza padrão ( )i

u x , de sua estimativa i

x , equivalente à variância da média

2 ( )i

s x . A raiz quadrada da variância da média é chamada de desvio padrão experimental

da média dada por ( )i

s x .

22 2

1

( ) 1( ) ( ) ( )

( 1)

ni

i i i

i

s xu x s x x x

n n n =

= = = −−∑ (3.3b)

60

Avaliação do tipo B: nesta avaliação a incerteza é obtida a partir de um

julgamento científico baseado em todas as informações relevantes disponíveis sobre o

instrumento e o processo de medição, tais como, a experiência do profissional ou o

conhecimento geral sobre os materiais e os instrumentos, além das especificações do

manual do fabricante.

Se a estimativa i

x é tomada de um manual do fabricante ou outra fonte técnica a

sua incerteza, cotada deste manual, é estabelecida como um múltiplo de um desvio

padrão, ou seja, a incerteza padrão ( )i

u x é simplesmente o valor cotado dividido pelo

multiplicador e a variância estimada 2 ( )i

u x é a raiz deste quociente (GUM, 2008). Sendo

assim, se um mensurando possui uma incerteza (cotada do manual) ao nível de confiança

de 68,26%, então a incerteza padrão é igual ao desvio padrão, ( )i

u x σ= .

Se a incerteza declarada de i

x é um parâmetro ao qual está associado um dado

nível de confiança de 90, 95 ou 99%, o cálculo da incerteza padrão só será efetuado se a

distribuição de probabilidade caracterizada pela estimativa do mensurando e sua incerteza

for conhecida. Neste caso, a incerteza padrão é a incerteza citada dividida pelo fator de

abrangência apropriado para a distribuição adotada. Tal fator é encontrado nas tabelas de

distribuição de probabilidades (Draper e Smitth, 1996). Se não for especificado o tipo de

distribuição pode-se assumir uma distribuição normal.

Quando o intervalo é bem definido, mas a distribuição de erros não é bem

conhecida, deve-se admitir uma distribuição diferente da normal, tal como a retangular,

trapezoidal ou triangular (GUM, 2008). A seguir são apresentadas as distribuições de

probabilidades mais utilizadas nas avaliações do tipo B.

Distribuição retangular: é utilizada quando não há nenhum conhecimento sobre

os possíveis valores do mensurando dentro do intervalo.

Em alguns caso é possível estimar apenas os limites, superior e inferior, para i

X e

estabelecer que a probabilidade para que o valor de i

X pertença ao intervalo [ ai , as ] é

100% e a probabilidade para que o valor i

X esteja fora desse intervalo é zero. Se não

61

houver conhecimento específico de possíveis valores de i

X dentro do intervalo, pode-se

assumir que i

X esteja em qualquer outro ponto dele. Neste caso, assume-se uma

distribuição uniforme ou retangular, Figura (3.2), e consequentemente o seu grau de

liberdade é infinito (Link, 1997).

Figura 3.2 – Curva de distribuição uniforme (adaptado de GUM, 2008 ).

Neste caso, i

x (estimativa ou valor esperado de i

X ) é o ponto médio do intervalo

e pode ser expressado pela Equação (3.4)

( )

2

i si

a ax

+= (3.4)

A variância associada é dada por:

22 ( )( )

12

s ii

a au x

−= (3.5)

Se a diferença entre os limites, −s i

a a , for designada por 2a , tem-se que a

Equação (3.5) torna-se:

22 ( )

3=

i

au x (3.6a)

E a incerteza padrão correspondente é dada por:

62

( )3

=i

au x (3.6b)

Dentre as variáveis deste tipo de distribuição estão à linearidade, a resolução e a

sensibilidade dos instrumentos de medição, encontrados em manuais do fabricante.

Distribuição Triangular: em muitos casos, é mais realista esperar que valores

perto dos limites sejam menos prováveis do que os que estejam perto do ponto médio. É,

então, razoável substituir a distribuição retangular simétrica, por uma distribuição

trapezoidal simétrica, tendo lados inclinados iguais (um trapezóide isósceles), uma base

de largura 2s i

a a a− = e um topo de largura 2 aβ , com 0 1β≤ ≤ (GUM, 2008).

Na medida em que 1β → esta distribuição trapezoidal se aproxima da

distribuição retangular, vista anteriormente, enquanto que, para 0β = , torna-se uma

distribuição triangular.

Neste caso, a variância associada é dada por:

2 22

22

(1 ): ( )

6

0: ( )6

β

β

+=

= =

i

i

aDistribuição trapezoidal u x

aDistribuiçãotriangular com u x

(3.7a)

E a incerteza padrão correspondente é dada por:

2(1 ): ( )

6

0: ( )6

β

β

+=

= =

i

i

aDistribuição trapezoidal u x

aDistribuiçãotriangular com u x

(3.7b)

Numa distribuição de probabilidade triangular a probabilidade de que um valor de

iX esteja dentro do intervalo [ ai , as] é igual a 1, para todos os pontos, e a probabilidade

63

de que i

X esteja fora deste intervalo é essencialmente zero. A Figura (3.3) mostra a curva

da função densidade de probabilidade triangular.

Figura 3.3 – Curva de distribuição triangular (adaptado GUM, 2008).

Em suma, devem ser coletadas informações que permitam estimar a incerteza

associada a cada fonte de erro. Recomenda-se apresentar o valor associado aos limites de

variação da fonte de incertezas em sua unidade natural e identificar o tipo de distribuição

de probabilidade envolvida (normal, retangular, triangular ou outra). Em função do tipo

de distribuição será definido o divisor utilizado para converter o valor conhecido na

incerteza padrão. Para distribuições normais este valor geralmente é unitário no caso da

avaliação de incerteza tipo A, ou coincide com o fator de abrangência utilizado na fonte

de informação quando a avaliação tipo B é considerada. Conforme apresentado na

descrição acima, os divisores das distribuições de probabilidade retangular e triangular

são 3 e 6 .

3.1.2 Incerteza padrão combinada

A incerteza padrão combinada pode ser calculada a partir das incertezas padrões

individuais, para cada uma das variáveis que interferem no processo de medição, através

de uma lei conhecida como “Lei de Propagação de Incertezas”. A incerteza, assim

determinada, é definida como incerteza padrão combinada e é designada por c

u . Neste

64

caso, a incerteza padrão combinada, ( )c

u y , é dada pela raiz quadrada positiva, conforme

mostra Equação (3.8).

2

2

1 1 1

( ) ( ) 2 ( ) ( , )N N N

N

c i j i j

i i j ii i j

f f fu y u x u x r x x

x x x= = = +

∂ ∂ ∂= +

∂ ∂ ∂ ∑ ∑ ∑ (3.8)

na qual y é a estimativa da variável resposta Y , i

x é a estimativa da variável i

X ,

2 ( )i

u x é variância associada a i

x , para todo i variando de 1 a N , N é o número de

variáveis que afetam a resposta Y , 2 ( )ju x é a incerteza associada à fonte de erro

representada pela estimativa jx e ( , )i jr x x é o coeficiente de correlação entre as

estimativas ix e jx .

A Lei de Propagação de Incerteza, porém, somente pode ser aplicada quando o

modelo matemático que relaciona a variável de resposta da medição com as variáveis que

afetam o seu comportamento for conhecido. A Equação (3.8), referenciada como a Lei de

Propagação de Incerteza, é baseada numa aproximação da série de Taylor de primeira

ordem de 1 2( , ,..., )NY f x x x= . As derivadas parciais em função de cada variável ix com

1,...,i N= que aparecem na Equação (3.8) são denominadas coeficientes de

sensibilidade. A grandeza desses coeficientes descreve a contribuição de cada fonte de

incerteza no valor final da incerteza de medição. Por sua vez, o segundo termo da referida

equação expressa a correlação existente entre duas fontes de incertezas ix , jx com i j≠ .

O coeficiente r , ver Equação (3.9), fornece uma medida do grau de correlação

entre as variáveis ix e jx .

2 2

1

2 2

1 1

( ) ( )

( , )

( ) ( )

M

ik i jk j

ki j

M M

ik i jk j

k k

x x x x

r x x

x x x x

=

= =

− −

=

− −

∑

∑ ∑ (3.9)

na qual M representa o número de valores atribuídos às variáveis i

x e j

x e i

x e j

x

representam, respectivamente, as médias aritmética dos M valores atribuídos à i

x e j

x .

65

O coeficiente de correlação varia de -1 a 1. Quando esse valor se aproxima dos

extremos significa que as variáveis i

x e j

x são altamente correlacionadas. Por outro lado,

se o coeficiente é zero significa que não há correlação entre as variáveis, assim o segundo

termo da Equação (3.8) desaparece conforme a Equação (3.10)

2

2 2

1

( ) ( )N

c i

i i

fu y u x

x=

∂=

∂ ∑ (3.10)

Portanto, se as estimativas i

x e j

x são independentes entre si, o coeficiente de

correlação é igual a zero diminuindo-se, assim, o número de cálculos necessários para

determinar a incerteza padrão combinada.

3.1.3 Incerteza padrão expandida

A incerteza padrão combinada calculada através da Lei de Propagação de

Incertezas apresenta uma probabilidade de abrangência de 68,27% sendo muito pequena

para a maioria das aplicações em engenharia. Assim sendo, o Comitê Internacional de

Pesos e Medidas propõe descrever a incerteza de medição através de intervalos que

representam os valores esperados para os erros de medição, com uma probabilidade

maior, geralmente de 95,45%. Esta incerteza recebe o nome incerteza expandida (U ) e

pode ser estimada pela Equação (3.11)

cU k u= (3.11)

na qual c

u é a incerteza padrão combinada e 0k > é o fator de abrangência.

O fator de abrangência k deve levar em conta, além do nível de confiança

desejado, o número de graus de liberdade efetivos associados ao caso para o intervalo

y U− a y U+ . Geralmente, k assume valores entre 2 e 3, mas, para aplicações especiais

k pode assumir valores fora dessa faixa

66

O grau de liberdade efetivo eff

ν é obtido através da fórmula de Welch-

Satterthwaite, conforme mostra a Equação (3.12).

4

4

1

ceff N

i

i i

u

uν

ν=

=

∑ (3.12)

na qual i

ν é o grau de liberdade associado a cada fonte de incerteza (ou variável de

entrada); c

u é a incerteza combinada; i

u e a incerteza padronizada associada à i-ésima

fonte de incerteza e N é o número total de fontes de incertezas analisadas.

Com a Equação (3.12) é possível calcular o número de graus de liberdade efetivo,

permitindo que seja obtido o valor de “ k ” para nível de confiança de 95% através da

Tabela (3.1).

Tabela 3.1 – Grau de liberdade efetivo e o correspondente fator de abrangência.

effν 1 2 3 4 5 6 7 8 10 12 14 16

95k 13,97 4,53 3,31 2,87 2,65 2,52 2,43 2,37 2,28 2,23 2,20 2,17

effν 18 20 25 30 35 40 45 50 60 80 100 ∞

95k 2,15 2,13 2,11 2,09 2,07 2,06 2,06 2,05 2,04 2,03 2,02 2,00

Para valores fracionários de eff

ν , interpolação linear pode ser usada se eff

ν > 3.

Alternativamente o valor de 95k corresponde ao valor de eff

ν imediatamente inferior na

tabela pode ser adotado.

Procedimento de avaliação da incerteza de medição para um caso geral:

1. Determinar o modelo matemático que relaciona a grandeza de entrada com a

saída, 1 2( , ,..., )=n

y f x x x ;

67

2. Identificar todas as correções a serem feitas ao resultado de medição;

3. Listar componentes sistemáticos da incerteza associada a correções e tratar efeitos

sistemáticos não corrigidos com parcelas de incerteza;

4. Atribuir valores de incertezas e distribuição de probabilidades com base em

conhecimentos experimentais práticos ou teóricos;

5. Calcular a Incerteza Padronizada, i

u , para cada componente de incerteza;

6. Calcular a Incerteza Combinada, ( )c

u y ;

7. Calcular a Incerteza Expandida, U .

3.2 INCERTEZA DE MEDIÇÃO PARA O CÁLCULO DAS

FREQUÊNCIAS NATURAIS TEÓRICAS

3.2.1 Incerteza em medição direta

Uma régua graduada de aço inox de seção transversal retangular e uniforme

representa uma viga engastada-livre (Figura (3.4)) colocada em vibração livre na direção

do eixo x, a fim de analisar a influência dos variáveis (intervalares) de entrada na resposta

do sistema. Mais detalhes da viga de aço encontram-se no Capítulo 4, item 4.2.1

68

a) Ensaio da viga engastada

b) Seção transversal retangular da viga ensaiada

Figura 3.4 – Viga engastada-livre em vibração livre.

Após 20 medições diretas da aceleração na direção x da viga, calculou-se a

incerteza padrão de cada uma das seguintes grandezas geométricas: base ( b ) e altura ( h )

da seção transversal, comprimento total ( tL ), comprimento útil ( u

L ), massa da viga ( vm )

e massa do acelerômetro ( am ).

69

Incertezas dos instrumentos de medição: paquímetro (fabricante: Mitutoyo; 150

mm; resolução 0,02 mm); paquímetro (fabricante: Mitutoyo; 300mm; resolução 0,05

mm); balança (fabricante: Marte; resolução 41,5.10− g ), ver Anexo B.

Tabela 3.2 – Valores obtidos após 20 medidas diretas de cada grandeza (i

x ).

(mm)b (mm)h (mm)uL (mm)tL (g)vm (g)am

26,20 0,98 250,50 305,80 54,6182 0,4341

26,20 0,90 250,30 305,85 54,6185 0,4342

26,16 0,96 250,55 305,85 54,6193 0,4341

26,12 1,00 250,00 305,50 54,6199 0,4342

26,10 0,90 250,30 306,10 54,6196 0,4341

25,98 0,94 250,60 305,85 54,6195 0,4341

26,00 0,94 250,50 306,10 54,6194 0,4342

26,00 0,94 250,20 306,20 54,6194 0,4341

25,98 0,90 250,25 306,00 54,6197 0,4341

25,98 0,98 250,50 306,10 54,6198 0,4341

26,00 0,94 250,40 306,10 54,6197 0,4342

26,00 0,96 250,45 306,10 54,6196 0,4341

26,22 0,96 250,50 306,10 54,6194 0,4342

26,24 1,00 250,55 305,85 54,6196 0,4342

26,20 1,00 250,50 305,85 54,6191 0,4342

26,10 0,96 250,30 305,85 54,6190 0,4342

26,00 0,94 250,55 306,20 54,6191 0,4341

25,96 0,96 250,30 305,80 54,6192 0,4342

25,98 0,96 250,55 305,85 54,6192 0,4342

25,96 0,96 250,50 305,80 54,6191 0,4341

x 26,07 0,95 250,42 305,94 54,6191 0,4341

2 ( )i

s x 0,1015 0,0305 0,1540 0,17791 4,19617 410−× 5,130 510−

×

70

Comprimento total ( tL )

Fontes de incerteza:

Incerteza Tipo A (TipoA

u ): 20 medições aleatórias

Incerteza Tipo B (TipoB

u ): 0,08 mm± (exatidão do paquímetro (300 mm) obtido do

catálogo do fabricante), ver Anexo B.

Estimativas dos efeitos aleatórios:

a) Incerteza padrão Tipo A: avaliada estatisticamente a partir de 20 medições repetidas do

comprimento total. A distribuição de probabilidade que rege a Incerteza Tipo A é a

normal, logo, a incerteza padrão Tipo A é obtida dividindo ( )i

u x por um valor unitário.

Média aritmética:

1==

∑n

i

i

x

xn

(3.13)

20n =

ix = valores apresentados na quarta coluna da Tabela (3.1)

305,94 mmx =

Desvio padrão experimental da variável i

x :

2 2

1

1( ) ( )

1

n

i i

i

s x x xn =

= −−∑ (3.14)

Desvio padrão da média ou (incerteza padrão):

2 ( ) ( )( ) ( ) i i

i i

s x s xu x s x

n n= = = (3.15)

71

0,17791( ) 0,03978 mm

20i

u x = =

( ) 0,03978 mmTipoA i

u u x= =

b) Incerteza padrão Tipo B (assume-se aqui distribuição retangular, logo a incerteza

padrão é obtida dividindo metade da exatidão do instrumento por 3 ).

Exatidão do paquímetro (300 mm) = 0,08 mm

Incerteza do paquímetro devido à exatidão = Inst

I

0,080,04 mm

2 2inst

exatidãoI = = = (3.16)

0,04

3 3

0,023 mm

InstTipoB

TipoB

Iu

u

= =

=

(3.17)

c) Incerteza padrão combinada

2 2

2 20,03978 0,023

0,046 mm

c TipoA TipoB

cLt

cLt

u u u

u

u

= +

= +

=

(3.18)

d) Incerteza padrão expandida

95%Lt cLtU k u=

Número de graus de liberdade efetivo é cálculado por meio da equação de Welch-

Satterwaite:

4

4

1

( )

( )

ceff N

i

i i

u y

u yν

ν=

=

∑ (3.19)

72

Para incerteza Tipo A: 1i

nν = −

Para incerteza Tipo B: i

ν = ∞

4

4 4

0,04633,97

0,03978 0,023

19

Ltν = =

+∞

(3.20)

De posse de Lt

ν e com as informações da Tabela (3.1), obtém-se:

95% 2,07k = (3.21)

Portanto:

2,07 0,046

0,095 mm

Lt

Lt

U

U

= ×

= (3.22)

Comprimento útil ( uL )






catálogo do fabricante)


a) Incerteza padrão Tipo A:

Média aritmética:

20n = medições

ix = valores apresentados na terceira coluna da Tabela (3.2)

250, 42mmx =

73

( ) 0,1540( ) 0,03443 mm

20

ii

s xu x

n= = = (3.23)

( ) 0,03443 mmTipoA i

u u x= =




0,080,04 mm

2 2inst


0,04

3 3

0,023 mm

InstTipoB

TipoB

Iu

u

= =

=

(3.25)


2 2

2 20,03443 0,023

0,042 mm

c TipoA TipoB

cLu

cLu

u u u

u

u

= +

= +

=

(3.26)


95%Lu cLuU k u=

4

4 4

0,04239,94

0,03443 0,023

19

Luν = =

+∞

(3.27)

De posse de Lu


95% 2,06k = (3.28)

Portanto:

74

2,06 0,042

0,086 mm

Lu

Lu

U

U

= ×

= (3.29)

Largura da seção transversal ( b )






catálogo do fabricante).



Média aritmética:

20n = medições

ix = valores apresentados na primeira coluna da Tabela (3.2)

26,07 mmx =

( ) 0,1015( ) 0,02269 mm

20

ii

s xu x

n= = = (3.30)

( ) 0,02269 mmTipoA i

u u x= =




0,030,015 mm

2 2inst


75

0,015

3 3

0,00866mm

InstTipoB

TipoB

Iu

u

= =

=

(3.32)


2 2

2 20,02269 0,0087

0,024 mm

c TipoA TipoB

cb

cb

u u u

u

u

= +

= +

=

(3.33)


95%b cbU k u=

4

4 4

0,02424,94

0,02269 0,0087

19

bν = =

+∞

(3.34)

De posse de b


95% 2,11k = (3.35)

Portanto:

2,11 0,024

0,051 mm

b

b

U

U

= ×

= (3.36)

Altura da seção transversal ( h )




76



catálogo do fabricante).



Média aritmética:

20n = medições

ix = valores apresentados na segunda coluna da Tabela (3.2)

0,9540 mmx =

( ) 0,0305( ) 0,006821 mm

20

ii

s xu x

n= = = (3.37)

( ) 0,006821 mmTipoA i

u u x= =




0,030,015 mm

2 2inst


0,015

3 3

0,00866mm

InstTipoB

TipoB

Iu

u

= =

=

(3.39)


2 2

2 20,006821 0,0087

0,011 mm

c TipoA TipoB

ch

ch

u u u

u

u

= +

= +

=

(3.40)

77


95%h chU k u=

42

4 4

0,0111,29 10

0,006821 0,0087

19

hν = = ×

+∞

(3.41)

De posse de h


95% 2,014k = (3.42)

Portanto:

2,014 0,011

0,022 mm

h

h

U

U

= ×

= (3.43)

Massa da viga


Incerteza padrão Tipo A (TipoA


Incerteza padrão Tipo B (TipoB

u ): 0,30 mg± (miligrama) obtido do certificado de

calibração da balança – Fabricante: Marte Balanças e Aparelhos de Precisão Ltda, ver

Anexo B.



Média aritmética:

20n = medições

ix = valores apresentados na quinta coluna da Tabela (3.2)

78

54,6193 gx =

45( ) 4,19617.10

( ) 9,38293 10 g20

ii

s xu x

n

−

−= = = × (3.44)

5( ) 9,38293 10 gTipoA iu u x−

= = ×


padrão é obtida dividindo inst

I por 3 ).

Incerteza do instrumento obtida do certificado de calibração (Balança – Marte) =

0,30 mg±

0,30,15 mg

2inst

I = = (3.45)

4

5

1,5 10

3 3

8,66 10 g

InstTipoB

TipoB

Iu

u

−

−

×= =

= ×

(3.46)


2 2

5 2 5 2

4

(9,38293 10 ) (8,66 10 )

1,3 10

c TipoA TipoB

cmv

cmv

u u u

u

u g

− −

−

= +

= × + ×

= ×

(3.47)


95%mv cmvU k u=

4 4

5 4 5 4

(1,3 10 )65,16

(9,38293 10 ) (8,66 10 )

19

mvν

−

− −

×= =

× ×+

∞

(3.48)

De posse de mv


79

95% 2,04k = (3.49)

Portanto:

( )4

4

2,04 1,3 10

2,6 10 g

mv

mv

U

U

−

−

= × ×

= ×

(3.50)

Massa do acelerômetro





u ): 0,30 mg± (obtido do certificado de calibração da

balança – Fabricante: Marte Balanças e Aparelhos de Precisão Ltda).



Média aritmética:

20n = medições

ix = valores apresentados na quinta coluna da Tabela (3.2)

0,4342 gx =

55( ) 5,130.10

( ) 1,147(1) 10 g20

ii

s xu x

n

−

−= = = × (3.51)

5( ) 1,147(1) 10 gTipoA iu u x−

= = ×



I por 3 ).

Incerteza do instrumento obtida do certificado de calibração (Balança – Marte) =

0,30 mg±

80

0,30,15 mg

2inst

I = = (3.52)

4

5

1,5 10

3 3

8,66 10 g

InstTipoB

TipoB

Iu

u

−

−

×= =

= ×

(3.53)


2 2

5 2 5 2

5

(1,147 10 ) (8,66 10 )

8,7 10 g

c TipoA TipoB

cma

cma

u u u

u

u

− −

−

= +

= × + ×

= ×

(3.54)


95%ma cmaU k u=

5 44

5 4 5 4

(8,7 10 )6,4 10

(1,147 10 ) (8,66 10 )

19

maν

−

−

− −

×= = ×

× ×+

∞

(3.55)

De posse de eff


95% 2,00k = (3.56)

Portanto:

4

4

2,00 (8,7 10 )

1,7 10 g

ma

ma

U

U

−

−

= × ×

= × (3.57)

81

3.2.2 Incerteza em medição indireta

A medição indireta envolve a determinação do valor associado ao mensurando a

partir da combinação de duas ou mais grandezas (de entrada) por meio de expressões

matemáticas. Quando essas grandezas de entrada se comportam de forma desvinculada,

do ponto de vista estatístico, estas variáveis (grandezas) são ditas independentes ou não

correlacionadas, e seu coeficiente de correlação (ver Equação (3.9)) é zero. Assume-se

aqui que as grandezas combinadas não possuem dependência estatística já que foram

obtidas com diferentes instrumentos de medição e não foram medidas sob as condições

de operação.

Grandezas medidas indiretamente: módulo de elasticidade ( E ), momento de inércia

de área ( I ), área da seção transversal ( A ) e massa específica ( ρ ).

Módulo de elasticidade

Por medição indireta, determinou-se o módulo de elasticidade (E) através do ensaio

dinâmico. Partindo da expressão usada no cálculo da frequência natural, tem-se:

2

n

k

mω = , 21

d nsω ω= − (3.58)

2(2 )n

ef

kf

mπ = (3.59)

na qual a rigidez da viga é dada por 33 /u

k EI L= , sendo uL o comprimento útil da viga e

3 /12=I bh , o momento de inércia de área; a massa efetiva total que age na extremidade

de uma viga engastada-livre é dada por ( )ef a eq

m m m= + , sendo, am , a massa do

acelerômetro, (33 /140)eq v

m m= a massa equivalente da viga , v

m , a massa da viga e d

ω ,

a frequência natural amortecida. Portanto, rearranjando a Equação (3.59) se obtém a

expressão matemática, Equação (3.60), para o cálculo do módulo de elasticidade.

82

A partir da Equação (3.60) calculou-se a incerteza padrão expandida do módulo de

elasticidade da viga de aço inox.

( )2 2

32 33

13 140 2

n

a v u

fE m m L

I

π

π

∆ = + +

(3.60)

na qual a razão / 2π∆ representa o fator de amortecimento (ς ), o decremento

logarítmico, dado por 2 1[ ( )] /π∆ = − nf f f , sendo 11,00nf Hz= a primeira frequência

natural e, 2 110,76 11,26f Hz e f Hz= = , as frequências correspondentes a largura de

banda associados aos pontos de meia potência.

a) Incerteza padrão:

Grandezas de entrada: largura da seção transversal ( b ); altura da seção transversal

( h ); massa do acelerômetro ( am ); massa da viga ( vm ); comprimento útil ( uL ).

Grandeza de saída: Módulo de Elasticidade (E).

Expressão matemática que correlacionam às grandezas: Equação (3.60)

Incerteza padrão:

largura da seção transversal: cb

u

altura da seção transversal: ch

u

massa do acelerômetro: acm

u

massa da viga: vcm

u

comprimento útil: ucL

u

frequências naturais, nf , 1 2f e f : incerteza padrão é desprezível

b) Incerteza padrão combinada: grandezas estatisticamente independentes.

A Equação (3.61) é usada para estimar a incerteza padrão combinada de grandezas

estatisticamente independentes, neste caso, as grandezas foram medidas em condições

diferentes de trabalho, por diferentes instrumentos de medição.

83

2

2 2

1

( ) ( )N

c i

i i

fu y u x

x=

∂=

∂ ∑ (3.61)

2 2 22 2

a v ucE cb ch cm cm cL

a v u

E E E E Eu u u u u u

b h m m L

∂ ∂ ∂ ∂ ∂ = + + + +

∂ ∂ ∂ ∂ ∂ (3.62)

( )

( )

2 2

3

2 3

2 2

4 2 3

2 2 4 3

2 331

140 23

12

2 11 33 0,14274.341 10 5,4619 10 1 0, 25042

140 2(2,607 10 ) (9,5 10 )3

12

n

a v u

fEm m L

b b h

E

b

E

π

π

π

π

− −

− −

∂ ∆ = − + + ∂

∂ = − × + × + ∂ × ×

∂

∂

126,857 10b

= − ×

( )

( )

2 2

3

4

2 2

4 2 3

2 4 4

2 333 1

140 23

12

2 11 33 0,14273 4,341 10 5,4619 10 1 0,25042

140 2(2,607 10 )(9,5 10 )3

12

n

a v u

fEm m L

h bh

E

h

E

π

π

π

π

− −

− −

∂ ∆ = − + + ∂

∂ = − × + × + ∂ × ×

∂

∂

145,644 10h

= − ×

( )

( )

2 2

3

3

2 2

3

2 4 3

13

21

23

12

2 11 0,14271 0, 25042

2(2,607 10 )(9,5 10 )3

12

1,343 10

n

u

a

a

a

fEL

m bh

E

m

E

m

π

π

π

π− −

∂ ∆ = + ∂

∂ = + ∂ × ×

∂= ×

∂

84

( )

( )

2 2

3

3

2 2

3

2 4 3

12

2331

140 23

12

2 1133 0,14271 0, 25042

140 2(2,607 10 )(9,5 10 )3

12

3,165 10

n

u

v

v

v

fEL

m bh

E

m

E

m

π

π

π

π− −

∂ ∆ = + ∂

∂ = + ∂ × ×

∂= ×

∂

( )

( )

2 2

2

3

2 2

4 2 2

2 4 3

2 333 1

140 23

12

2 11 33 0,14273 4,341 10 5,4619 10 1 0, 25042

140 2(2,607 10 )(9,5 10 )3

12

n

a v u

u

u

fEm m L

L bh

E

L

E

π

π

π

π

− −

− −

∂ ∆ = + + ∂

∂ = × + × + ∂ × ×

∂

∂

122,142 10uL

= ×

( )

( )

( )

( )

22

12 5 16

22

14 5 19

22

13 8 12

22

12 7

( 6,857 10 )(2, 4 10 ) 2,708 10

( 5,644 10 )(1,1 10 ) 3,854 10

(1,343 10 )(8,7 10 ) 1,365 10

(3,165 10 )(1,3 10 ) 1,693 10

a

v

cb

ch

cm

a

cm

v

Eu

b

Eu

h

Eu

m

Eu

m

−

−

−

−

∂ = − × × = ×

∂

∂ = − × × = ×

∂

∂= × × = ×

∂

∂= × × = ×

∂

( )

11

22

12 5 15(2,142 10 )(4, 2 10 ) 8,093 10

ucL

u

Eu

L

− ∂

= × × = × ∂

6,21 GPacEu = (3.63)

85

c) Incerteza padrão expandida:

Número de graus de liberdade efetivo:

4

4 4 44 4

a v u

a v u

cEE

cm cm cLcb cha v u

b h m m L

u

E E EE Eu u uu u

m m Lb h

ν

ν ν ν ν ν

= = ∞

∂ ∂ ∂∂ ∂ ∂ ∂ ∂∂ ∂

+ + + +

(3.64)

sendo que o grau de liberdade individual das grandezas após 20 medições é

20 1 19a v ub h m m Lν ν ν ν ν= = = = = − = .

( )

( ) ( ) ( ) ( ) ( )

49

4 4 4 4 412 5 14 5 13 8 12 7 12 5

6, 21 1019,02

( 6,857 10 )(2, 4 10 ) ( 5,644 10 )(1,1 10 ) (1,343 10 )(8,7 10 ) (3,165 10 )(1,3 10 ) (2,142 10 )(4, 2 10 )

19 19 19 19 19

Eν− − − − −

×= =

− × × − × × × × × × × ×+ + + +

De acordo com a Tabela (3.1) para 19,01Eν = temos 95% 2,14k = .

95%E cEU k u=

( )92,14 6,21 10EU = × ×

13,29 GPaEU = (3.65)

Momento de inércia ( I )

3

12

bhI =

(3.66)

a) Incerteza padrão

Incerteza padrão das grandezas de entrada:

Incerteza padrão de b = 5

2,4 10cbu−

= ×

Incerteza padrão de h = 5

1,1 10chu−

= ×

86

b)Incerteza padrão combinada: grandezas estatisticamente independentes.

A Equação (3.67) é usada para estimar a incerteza padrão combinada de grandezas

estatisticamente independentes como é o caso do momento de inércia de área, na qual

tanto a base quanto a altura, da seção transversal da viga, foram medidas por um mesmo

instrumento (paquímetro – 300mm) nas mesmas condições de trabalho. Neste caso, a

incerteza combinada para grandezas estatisticamente dependente é calculada através da

Equação (3.67).

2

2 2

1

( ) ( )N

c i

i i

fu y u x

x=

∂=

∂ ∑ (3.67)

2 2

cI cb ch

I Iu u u

b h

∂ ∂ = +

∂ ∂

3

12

I h

b

∂=

∂= 7,14

1110

−×

2

312

I bh

h

∂=

∂= 5,88

910

−×

( )( )( ) ( )( )( )2 2

11 5 9 57,14 10 2,4 10 5,88 10 1,1 10cIu

− − − −= × × + × ×

14 4 6,64 10 mcIu

−= × (3.68)



4

4 4

cII

cb ch

b h

u

I Iu u

b h

ν

ν ν

=∂ ∂

∂ ∂

+

(3.69)

87

( )

( )( )( )( )

( )( )( )( )

414

4 411 5 9 5

6,64 10

7,14 10 2, 4 10 5,88 10 1,1 10

20 1 20 1

Iν

−

− − − −

×=

× × × ×

+− −

Iν = 21,1

Substituindo as variáveis calculadas em seções anteriores sobre medições diretas na

Equação (3.19), obtém-se Iν = 21,1 e o corresponde 95% 2,1256k = (ver Tabela (3.1)).

95%I cIU k u=

( )142,1256 6,64 10IU −= × ×

13 41, 41 10 mIU−

= × (3.70)

Área da seção transversal ( A )

A bh= (3.71)



Incerteza padrão de b = cbu

Incerteza padrão de h = chu


2 2

cA cb ch

A Au u u

b h

∂ ∂ = +

∂ ∂ (3.72)

49,5 10A

hb

−∂= = ×

∂

22,607 10A

bh

−∂= = ×

∂

88

( ) ( )( ) ( ) ( )( )2 2

4 5 2 59,5 10 2, 4 10 2,607 10 1,1 10cAu− − − −

= × × + × ×

7 2 2,88 10 mcAu−

= × (3.73)



4

4 4

cAA

cb ch

b h

u

u uν

ν ν

=

+

(3.74)

4

4 4

cAA

cb ch

b h

u

A Au u

b h

ν

ν ν

=∂ ∂

∂ ∂ +

( )

( )( )( )( )

( ) ( )( )( )

47

4 44 5 2 5

2,88 10

9,5 10 2, 4 10 2,607 10 1,1 10

20 1 20 1

Aν

−

− − − −

×=

× × × ×

+− −

19,3Aν =


Equação (3.74), obtém-se 19,3Aν = e o corresponde 95% 2,14k = (ver Tabela (3.1)).

95%A cAU k u=

72,14 2,88.10AU−

= ×

7 26,16 10 mAU−

= × (3.75)

Massa específica ( ρ )

. .

v

t

m

b h Lρ = (3.76)

89



Incerteza padrão de vm = vmu

Incerteza padrão de tL = tLu

Incerteza padrão de b = cbu

Incerteza padrão de h = chu


2 2 2 2

v tc cm cL cb ch

v

u u u u um Lt b h

ρ

ρ ρ ρ ρ ∂ ∂ ∂ ∂ = + + +

∂ ∂ ∂ ∂ (3.77)

511,32 10

. .v tm b h L

ρ∂= = ×

∂

4

22,36 10

. .

v

t t

m

L b h L

ρ∂= − = ×

∂

5

22,77 10

. .

v

t

m

b b h L

ρ∂= − = ×

∂

6

27,60 10

. .

v

t

m

h b h L

ρ∂= − = ×

∂

( )( )( ) ( )( )( ) ( )( )( ) ( )( )( )2 2 2 2

5 2 4 1 5 5 6 51,32 10 5,462 10 2,36 10 3,059 10 2,77 10 2,4 10 7,60 10 1,1 10c

uρ

− − − −= × × + × × + × × + × ×

4

3

kg1,020 10

mcuρ

= × (3.78)



90

4

44 4 4

v t

c

cm cL cb chv

mv Lt b h

u

u u u um Lt b h

ρ

ρν

ρ ρ ρ ρ

ν ν ν ν

=

∂ ∂ ∂ ∂ ∂ ∂ ∂ ∂

+ + +

(3.79)

37,97ρ

ν =


Equação (3.19), obtém-se 37,97A

ν = e o corresponde 95% 2,064k = .

95%ρ ρ=

cU k u

( )42,064 1,020 10Uρ

= × ×

4

3

kg2,105 10

mU

ρ= × (3.80)

3.3 INCERTEZA DE MEDIÇÃO PARA O CÁLCULO DAS

FREQUÊNCIAS NATURAIS EXPERIMENTAIS

Frequências naturais experimentais





u ) obtida do certificado de calibração do miniacelerômetro

– modelo: 352C22; serial 94793; Fabricante: PCB Piezoeletronics. Ver Anexo B.

Incerteza padrão expandida no intervalo (10 - 99) Hz: 1,5 %± .

Incerteza padrão expandida no intervalo (100 - 1999) Hz: 1,0 %± .

91

Tabela 3.3 – Três primeiras freqüências naturais coletadas durante os ensaios.

Ensaio 1f (Hz) 2f (Hz) 3f (Hz)

1 11,0 70,5 195,0

2 11,0 70,5 195,0

3 11,0 70,5 195,5

4 11,0 70,5 195,5

5 11,0 70,5 195,5

6 11,0 70,5 195,5

7 11,0 70,5 195,5

8 11,0 70,5 195,0

9 11,0 70,5 195,5

10 11,0 70,5 195,5

11 11,0 70,5 195,5

12 11,0 70,5 195,5

13 11,0 70,5 195,5

14 11,0 70,5 195,5

15 11,0 70,5 195,5

16 11,0 70,5 195,5

17 11,0 70,5 195,5

18 11,0 70,5 195,5

19 11,0 70,5 195,5

20 11,0 70,5 195,5

na qual , 1,2,3i

f i = , representam as frequências naturais.


20n = medições.

ix = valores apresentados na Tabela (3.3)

Média aritmética : primeira frequência natural 1ω :

11,0 Hzx =

92

( ) 0 HzTipoA i

u u x= = (3.81)

Média aritmética : segunda frequência natural 2ω :

70,5 Hzx =

( ) 0 HzTipoA i

u u x= = (3.82)

Média aritmética : terceira frequência natural 3ω :

195,425 Hzx =

( ) 0,18317 HzTipoA i

u u x= = (3.83)



I por 3 ).

b1) Incerteza combinada do instrumento: Primeira Freqüência Natural

(0,015 11) Hz 0,165 Hz± × = ±

0,1650,0825 Hz

2inst

I = =

0,0825

3 3

0,0476 Hz

InstTipoB

TipoB

Iu

u

= =

=

(3.84)

Incerteza padrão combinada:

2 2

2 2

1

1

(0) (0,0476)

0,0476 Hz

c TipoA TipoB

cf

cf

u u u

u

u

= +

= +

=

(3.85)

Incerteza padrão expandida

1 95% 1f cfU k u=

Para 1fν = ∞ tem-se: 95% 2,00k =

1

1

2,00 (0,0476)

0,09526 Hz

f

f

U

U

= ×

= (3.86)

93

b2) Incerteza combinada do instrumento: Segunda Freqüência Natural

(0,010 70,5) Hz 0,705 Hz± × = ±

0,7050,3525

2inst

I Hz= =

0,3525

3 3

0,2035 Hz

InstTipoB

TipoB

Iu

u

= =

=

(3.87)

Incerteza padrão combinada

2 2

2 2

2

2

(0) (0, 2035)

0, 2035 Hz

c TipoA TipoB

cf

cf

u u u

u

u

= +

= +

=

(3.88)


2 95% 2f cfU k u=

Para 2fν = ∞ tem-se: 95% 2,00k =

2

2

2,00 (0, 2035)

0,4070 Hz

f

f

U

U

= ×

= (3.89)

b3) Incerteza combinada do instrumento: Terceira Freqüência Natural

(0,010 195,42) Hz 1,9542Hz± × = ±

0,195420,97710

2inst

I Hz= =

0,97710

3 3

0,56413 Hz

InstTipoB

TipoB

Iu

u

= =

=

(3.90)

Incerteza padrão combinada

2 2

2 2

3

3

(0,18317) (0,56413)

0,56563 Hz

c TipoA TipoB

cf

cf

u u u

u

u

= +

= +

=

(3.91)


94

3 95% 3f cfU k u=

Para 3fν = ∞ tem-se: 95% 2,00k =

3

3

2,00 (0,56563)

1,1313 Hz

f

f

U

U

= ×

= (3.92)

95

Capítulo 4

VALIDAÇÃO

Este capítulo tem como objetivo validar os resultados numéricos, obtidos pela autora, com a

literatura e experimentalmente. A resposta dinâmica de vigas e treliças (com parâmetros

intervalares) formuladas pela Teoria da Perturbação de Matrizes é comparada com a resposta obtida

por meio do Parameter Solution Vertex Theorem (Qiu et al., 2005) e validada, também, através da

simulação de Monte Carlo, concluindo a validação teórica do método apresentado neste trabalho.

Baseado em um dos métodos numéricos mais utilizados em análise estrutural, um algoritmo de

elementos finitos foi desenvolvido e aprimorado a fim de trabalhar com a Teoria da Perturbação de

Matrizes aplicada a problemas de autovalores intervalares. É apresentada uma validação

experimental, na qual o intervalo dos parâmetros intervalares de uma viga em balanço de aço (que

segue uma distribuição uniforme de probabilidade) é obtido por meio de métodos probabilísticos.

As soluções destes problemas, com os mesmos parâmetros de entrada dos autores da literatura, são

apresentadas a seguir em forma de tabelas e gráficos.

96

4.1 VALIDAÇÃO TEÓRICA

Os resultados teóricos (frequências naturais) de vigas e treliças gerados no MATLAB®, via

Teoria de Perturbação de Matrizes aplicado a problemas dinâmicos intervalares, são comparados

aos resultados da literatura obtidos via Parameter Solution Vertex, conforme apresentado no Anexo

B, e aos resultados obtidos por simulação de Monte Carlo. O MATLAB é um software poderoso

para a solução de métodos numéricos na engenharia, possui inúmeros recursos e funções com o

objetivo de facilitar a vida do engenheiro, fornecendo recursos práticos para a solução dos

problemas na mecânica estrutural. A função Eig (A,B) do MATLAB, usada para cálculos de

autovalores e autovetores, é aplicada pela autora na solução de sistemas lineares. A função Eig

(A,B) para matrizes simétricas positivas definidas calcula os autovalores usando a decomposição de

Cholesky da matriz B. Se A e B não são simétricas, então, a função calcula os autovalores por meio

do algoritmo QZ ou QR, na qual QZ é uma matriz ortogonal e QR uma matriz triangular superior.

4.1.1 Viga escalonada: módulo de elasticidade intervalar

Inicia-se o capítulo com a validação dos autovalores intervalares de uma viga escalonada em

três segmentos, conforme apresentado em Qiu et al (2005), ver Figura (4.1). No primeiro caso, cada

um dos segmentos da viga apresenta-se com o módulo de elasticidade intervalar e os demais

parâmetros como, área da seção transversal, momento de inércia de área, massa específica e os

comprimentos de cada segmento são apresentados com valores fixos. Neste caso, os parâmetros

intervalares são dados dentro dos seguintes intervalos: 1 [199,7 , 200,3] GPaE = ,

2 [199,8, 200, 2] GPaE = e 3 [199,9, 200,1] GPaE = .

Figura 4.1 – Viga escalonada em três segmentos.

As grandezas físicas adotadas para a viga da Figura (4.1) são:

Massa específica: 37.800 kg / m iρ = , 1, 2,3i =

Comprimento de cada segmento: 0, 4 mi

L =

Área da seção transversal: 2 2 2 2 2 2

1 2 31,44 10 m , 1,0 10 m , 0,64 10 mA A A− − −

= × = × = ×

97

Momento de inércia de área: 4 4 4 4 4 4

1 2 30, 2 10 , 0,1 10 , 0,05 10I m I m I m− − −

= × = × = ×

A Tabela (4.1) apresenta os autovalores da viga da Figura (4.1) obtidos pela autora em

comparação com os resultados obtidos por Qiu et al. (2005) e por simulação de Monte Carlo. Os

valores obtidos pela autora concordam com os obtidos por Qiu et al. (2005) com uma precisão de

até seis casas decimais. Observa-se também boa concordância dos resultados da autora em

comparação com os valores obtidos após 10 mil simulações de Monte Carlo.

Tabela 4.1 – Autovalores da viga da Figura (4.1) com módulo de elasticidade intervalar.

AUTORA Teoria da Perturbação de

Matrizes

Método de Qiu et al. (2005) Parameter Vertex Solution

Algorithm Método de Monte Carlo

Inferior Superior

Inferior Superior

Inferior Superior

1λ 3,645921510× 3,655699

510× 3,645921510× 3,655699

510× 3,646009510× 3,655650

510×

2λ 7,327840 610× 7,343228 610× 7,327840 610× 7,343228 610× 7,328810 610× 7,343021 610×

3λ 4,817373 710× 4,826370 710× 4,817373 710× 4,826370 710× 4,818002 710× 4,825716 710×

4λ 2,577919 810× 2,583228 810× 2,577919 810× 2,583228 810× 2,578083 810× 2,582989 810×

5λ 8,910574 810× 8,928784 810× 8,910574 810× 8,928784 810× 8,911294 810× 8,927422 810×

6λ 2,738051 910× 2,741501 910× 2,738051 910× 2,741501 910× 2,738226 910× 2,741339 910×

A Tabela (4.2) apresenta as diferenças relativas dos três primeiros autovalores (apresentados

na Tabela (4.1)) obtidos pelo método da autora em relação ao método de Monte Carlo.

Tabela 4.2 – Diferença relativa da Teoria da Perturbação de Matrizes em relação ao Método de

Monte Carlo dos três primeiros autovalores da Tabela (4.1).

DIFERENÇA RELATIVA 100

Monte Carlo Autoraq q

q Monte Carloq

rλ

λ λ

λ

−= × (%) ; 1,2,3q =

(da Teoria da Perturbação de Matrizes em relação ao Método de Monte Carlo)

Inferior

Superior

1rλ 0,0024 -0,0013

2rλ 0,0132 -0,00282

3rλ 0,0131 -0,0136

Conforme mencionado no Capítulo 2, a técnica de simulação de Monte Carlo baseia-se na

geração de um número finito de amostras a partir de funções de densidade de probabilidade. De

98

acordo com a literatura, apesar do elevado custo computacional inerente à técnica de Monte Carlo,

diversos autores têm usado a técnica para validar novos métodos de tratamento de incertezas. Neste

estudo, adotou-se a função de densidade de probabilidade para uma distribuição uniforme.

Realizou-se 10 mil simulações permitindo-se obter as curvas, conforme Figura (4.2). Os limites

inferior e superior de cada um dos autovalores (adicionado ao eixo das ordenadas) estão indicados

na Figura (4.2).

0 1000 2000 3000 4000 5000 6000 7000 8000 9000 100003.646

3.647

3.648

3.649

3.65

3.651

3.652

3.653

3.654

3.655

3.656x 10

5

X: 3384

Y: 3.656e+005

Simulação de Monte Carlo

número de simulações

λ1

X: 6825

Y: 3.646e+005

0 1000 2000 3000 4000 5000 6000 7000 8000 9000 100007.326

7.328

7.33

7.332

7.334

7.336

7.338

7.34

7.342

7.344x 10

6

X: 3828

Y: 7.343e+006



λ2

X: 9880

Y: 7.328e+006

0 1000 2000 3000 4000 5000 6000 7000 8000 9000 100004.817

4.818

4.819

4.82

4.821

4.822

4.823

4.824

4.825

4.826

4.827x 10

7

X: 2907

Y: 4.826e+007



λ3

X: 482

Y: 4.818e+007

Figura 4.2 – Curvas dos três primeiros autovalores obtidos após 10 mil simulações de Monte Carlo.

99

4.1.2 Viga escalonada: seção transversal e momento de inércia

intervalares

No segundo caso, os autovalores da Figura (4.1) serão avaliados tendo como parâmetros

intervalares a área da seção transversal e o momento de inércia, cujos valores são os seguintes:

Área da seção transversal:

2 2 2

1

2 2 2

2

2 2 2

3

1, 426 10 ; 1, 454 10 m

0,99 10 ; 1,01 10 m

0,634 10 ; 0,646 10 m

A

A

A

− −

− −

− −

= × ×

= × ×

= × × (4.1)

Momento de inércia de área:

4 4 4

1

4 4 4

2

4 4 4

3

0,1998 10 ; 0, 2002 10 m

0,0999 10 ; 0,1001 10 m

0,04995 10 ; 0,05005 10 m

I

I

I

− −

− −

− −

= × ×

= × ×

= × × (4.2)

As demais grandezas são adotadas com valores fixos como o módulo de elasticidade

( 200 GPai

E = , 1,2,3i = ), a massa específica e os comprimentos dos segmentos da viga. Os

últimos parâmetros citados possuem valores iguais aos do primeiro caso.

De acordo com a Tabela (4.3), os autovalores encontrados pela autora são os mesmos

obtidos por Qiu et al. (2005). E estes valores concordam com boa precisão com os obtidos por

simulação de Monte Carlo.

Tabela 4.3 – Autovalores da viga da Figura (4.1) com área da seção transversal e momento de

inércia intervalares.

AUTORA

Teoria da Perturbação de

Matrizes

Método de Qiu et al. (2005)

Parameter Vertex Solution

Algorithm

Método de Monte Carlo

Inferior Superior

Inferior Superior

Inferior Superior

1λ 3,612794510× 3,689557

510× 3,612794510× 3,689557

510× 3,616628510× 3,686738

510×

2λ 7,257447 610× 7,415161 610× 7,257447 610× 7,415161 610× 7,290289 610× 7,422973 610×

3λ 4,770923 710× 4,873816 710× 4,770923 710× 4,873816 710× 4,809789 710× 4,895430 710×

4λ 2,553292 810× 2,608389 810× 2,553292 810× 2,608389 810× 2,575256 810× 2,619870 810×

5λ 8,824674 810× 9,016557 810× 8,824674 810× 9,016557 810× 8,887011 810× 9,046509 810×

6λ 2,711208 910× 2,768894 910× 2,711208 910× 2,768894 910× 2,717279 910× 2,768069 910×

100

A Tabela (4.4) apresenta a diferença relativa dos três primeiros autovalores (apresentados na

Tabela (4.3)) obtidos pelo método da autora em relação ao método de Monte Carlo.

Tabela 4.4 – Diferença relativa da Teoria da Perturbação de Matrizes em relação ao Método de

Monte Carlo dos três primeiros autovalores da Tabela (4.3).



q Monte Carloq

rλ

λ λ

λ

−= × (%) ; 1,2,3q =

(da Teoria da Perturbação de Matrizes em relação ao Método de Monte Carlo)

Inferior

Superior

1rλ 0,1060 -0,0765

2rλ 0,4505 0,1052

3rλ 0,8081 0,4415

4.1.3 Treliça com oito elementos de barra

O terceiro caso em estudo é uma treliça com 8 (oito) elementos de barras e um total de 8

(oito) graus de liberdade, ver Figura (4.3). As áreas das seções transversais dos elementos de barra

de números 1,2,3,4 e 6 são adotados com intervalos que variam com a mudança gradativa de um

fator de incerteza ( β ), conforme o seguinte intervalo:

,I c c c c

iA A A A Aβ β = − + (4.3)

na qual o valor médio da área da seção transversal para 1,2,3,4,6i = é 4 22,0 10 mc

A−

= × e o fator

de incerteza varia dentro do intervalo 0 2%β≤ ≤ .

A área da seção transversal dos demais elementos de barra é: 4 2

5 7 8 1,0 10A A A m−

= = = × . O

módulo de elasticidade adotado é 200 GPa e massa específica é de 7 800 kg/m3 .

101

Figura 4.3 – Treliça com oito elementos de barra.

A Figura (4.4) mostra o acréscimo (no limite superior) e o decréscimo (no limite inferior)

dos autovalores à medida que o fator de incerteza assume valores crescentes dentro do intervalo

0 2%β≤ ≤ . Os gráficos apresentados na primeira coluna da Figura (4.4) são obtidos pelo Método

da Perturbação de Matrizes (Autora) e por Parameter Solution Vertex Theorem (Qiu et al., 2005) e

os gráficos da segunda coluna são obtidos por simulação de Monte Carlo. A primeira coluna mostra

que os autovalores encontrados pela autora são iguais aos encontrados por Qiu et al (2005) e por

isso elas se sobrepõem. Como já era esperado, o intervalo dos autovalores aumentam

monotonicamente com o aumento de β no intervalo de 0 a 2%.

Os limites inferior e superior obtidos após 10 mil simulações de Monte Carlo apresentam

concordância com os apresentados na primeira coluna da Figura (4.4). A grande vantagem do

método da autora em comparação com o de Monte Carlo, além da praticidade de aplicação na

prática, é a velocidade com que os resultados são gerados. Rapidez de simulação. Os resultados

apresentados pela autora foram obtidos após dois segundos de simulação, enquanto que os

resultados calculados por Monte Carlo foram obtidos após cinco minutos de simulação (trabalhando

em um computador com processador Intel (R) Core i3 - 1,40 GHz).

102

Autora e Método de Qiu et al.(2005)


0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 27.65

7.7

7.75

7.8

7.85

7.9

7.95

8

8.05

8.1

8.15x 10

5 Autovalores

β %

λ1

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 27.65

7.7

7.75

7.8

7.85

7.9

7.95

8

8.05

8.1

8.15x 10

5 Autovalores obtidos por simulação de Monte Carlo

β %

λ1

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 26.35

6.4

6.45

6.5

6.55

6.6

6.65

6.7

6.75x 10

6 Autovalores

β %

λ2

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 26.35

6.4

6.45

6.5

6.55

6.6

6.65

6.7

6.75x 10


β %

λ2

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 28.6

8.7

8.8

8.9

9

9.1

9.2x 10

6 Autovalores

β %

λ3

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 28.6

8.7

8.8

8.9

9

9.1

9.2x 10


β %

λ3

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 22.4

2.42

2.44

2.46

2.48

2.5

2.52

2.54x 10

7 Autovalores

β %

λ4

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 22.4

2.42

2.44

2.46

2.48

2.5

2.52

2.54x 10


β %

λ4

Figura 4.4 – Comparando os limites dos autovalores dentro do intervalo 0 2%β≤ ≤ de incerteza

obtidos pela autora e Qiu et al. (2005) com os autovalores obtidos por meio de Monte Carlo, usando

a matriz de massa apresentada na Equação (A.21).

103

4.1.4 Treliça com três elementos de barra

O estudo seguinte se refere a uma treliça com três elementos de barras inclinadas e três graus

de liberdade, conforme mostrado em Mehdi (2005), ver Figura (4.5). A fim de reduzir o custo

computacional da análise do sistema, Mehdi (2005) trabalhou com matrizes de massa concentrada,

que é obtida aproximando-se o elemento da treliça a um sistema do tipo massa-mola. A formulação

da matriz de massa usada por Mehdi (2005) reduz o custo computacional, porém, como

consequência, tem-se a redução nos valores das frequências naturais quando comparadas a resposta

de um sistema com matriz de massa consistente.

O problema apresentado por Mehdi (2005) possui incertezas nos módulos de elasticidade e

os valores intervalares adotados para cada uma das barras são:

[ ] [ ] [ ]1 2 30,9; 1,1 , 0,7; 1,3 , 0,8; 1,2E E E E E E= × = × = ×

(4.4)

Figura 4.5 – Treliça com três elementos de barras e módulo de elasticidade intervalar.

De acordo com o Teorema de Buckingham, os quatro parâmetros seguintes podem ser

reduzidos a um único parâmetro adimensional ( Ω ): frequência natural, f (1/s), comprimento da

barra, L (m), módulo de elasticidade, E (Pa) e massa específica, ρ (kg/m3). A relação entre esses

parâmetros de forma a torná-los adimensionais é dado pela seguinte expressão:

/

f L

E ρΩ = (4.5)

104

A validação da solução do problema da Figura (4.5) com a literatura (Mehdi, (2005)) é

apresentado na Tabela (4.5). A tabela seguinte apresenta as três primeiras frequências naturais

adimensionais, mínimas e máximas ( Ω ).

Tabela 4.5 – Frequências naturais adimensionais da treliça da Figura (4.5) com módulo de

elasticidade intervalar.

AUTORA Mehdi (2005)

Método de Monte Carlo

Inferior Superior

Inferior Superior

Inferior Superior

1Ω 0,5661 0,6964 0,5661 0,6964 0,5338 0,7271

2Ω 0,8910 1,0936 0,8910 1,0936 0,8368 1,1399

3Ω 1,2188 1,4897 1,2188 1,4897 1,1361 1,5475

A Tabela (4.6) apresenta a diferença relativa das três primeiras frequências naturais

adimensionais (apresentadas na Tabela (4.5)) obtidas pelo método da autora em relação ao método

de Monte Carlo.

Tabela 4.6 – Diferença relativa das três primeiras frequências naturais adimensionais apresentadas

na Tabela (4.5).



q Monte Carloq

rΩ

Ω − Ω= ×

Ω (%) ; 1,2,3q =

Inferior

Superior

1rΩ

- 6,0510 4,2223

2rΩ

- 6,4771 4,0618

3rΩ

- 7,2793 3,7351

4.2 VALIDAÇÃO EXPERIMENTAL

O objetivo aqui é comparar os resultados numéricos (na presença de incertezas padrão do

Tipo A e do Tipo B) com os resultados experimentais obtidos no Laboratório de Vibrações

Mecânicas da UNIFEI, durante a análise dinâmica de uma viga de aço.

Em Matos (2007) os parâmetros aleatórios, que seguem uma distribuição normal de

probabilidade (com aproximadamente 95% de confiabilidade), são convertidos em parâmetros

105

intervalares assumindo uma distribuição uniforme de 100% de confiabilidade ao intervalo da

distribuição normal do mesmo parâmetro, ver Figura (4.6).

O intervalo dos parâmetros da viga ensaiada foi obtido por meio de métodos probabilísticos

num processo de avaliação das incertezas de medição do ensaio. Neste caso, o limite inferior e o

limite superior do intervalo de cada parâmetro seguem uma distribuição uniforme de 100% de nível

de confiança. Os limites inferiores e o superior dos parâmetros são obtidos subtraindo e adicionando

o desvio padrão (δ ) de cada variável aleatória (X) de seus valores médios ( x ). Neste caso, os

parâmetros intervalares são limitados da seguinte maneira: X=[ ],x xδ δ− + , conforme

exemplificado na Figura (4.6).

Figura 4.6 - Comparação dos parâmetros intervalares e parâmetros aleatórios.

Conforme apresentado no Capítulo 2, os resultados numéricos foram obtidos a partir da

teoria da perturbação de matrizes aplicados a um problema de autovalor intervalar. A quantificação

da incerteza teórica requer o cálculo de incertezas padrão e incertezas expandidas de alguns dos

parâmetros utilizados nas Equações (2.142) e (2.143) do Capítulo 2, reapresentadas a seguir.

( ) ( ) 2

min

1 1

[ ] ( ) [ ]p p

e e

i i i i

i i

l K x u M xω= =

=∑ ∑ (4.6)

( ) ( ) 2

max

1 1

[ ] ( ) [ ]n n

e e

i i i i

i i

u K x l M xω= =

=∑ ∑ (4.7)

106

Nas Equações (4.8) e (4.9) é apresentado um exemplo das matrizes intervalares para o

primeiro (1º) elemento de uma viga.

( )( ) ( )

( )

1 1

2 2

1 1 1 1

1 1

1 1 1 1 2 21 1 1 1 1 11 13

1 1 1

12 6 12 6 0 0

6 4 6 2 0 0

12 6 12 6 0 0, ,

, [ ] 6 6 0 02 4,

00 00 0

0 00 00 0

ee

L L

L L L L

L Ll u E l u I

K l u K L LL Ll u L

−

− − − −

= = −

…

…

…

…

(4.8)

( )( ) ( )

1 1

2 2

1 1 1 1

1 13

1 1 1 1 1 2 21 1 1 1 1 11 1

156 22 54 13 0 0

22 4 13 3 0 0

54 13 156 22 0 0, ,

, [ ] 13 22 0 03 4420

00 00 0

0 00 00 0

ee

L L

L L L L

L Ll u A l u L

M l u M L LL Lρ

−

− −

= = − −−

…

…

…

…

(4.9)

na qual os limites, inferior e superior, ( ),i il u , das incertezas de grandezas como: momento de

inércia de área, I , área da seção transversal, A , e módulo de elasticidade, E , foram obtidos

calculando a incerteza expandida desses parâmetros. Calculou-se também a incerteza padrão, após

20 medições de cada uma das seguintes grandezas geométricas: base ( b ) e altura ( h ) da seção

transversal, comprimento total ( tL ) e comprimento útil ( u

L ), que descarta o trecho usado no

engastamento da viga. Por medição indireta, determinou-se o módulo de elasticidade (E) através de

um ensaio dinâmico, como mostrado no Capítulo 3.

4.2.1 Descrição do ensaio

Uma régua graduada de aço de seção transversal retangular e uniforme representa uma viga

em balanço em vibração livre, conforme Figura (4.7).

107

a)Detalhe da montagem

b) Detalhe do analisador de sinais

Figura 4.7 -. Régua graduada de aço de seção transversal retangular em vibração livre.

A viga em balanço é colocada em oscilação mediante a aplicação de um impulso com o

intuito de obter as frequências fundamentais em Hz. O sinal captado pelo miniacelerômetro

(Modelo: 352C22 SN 94793; Serial: 94793; Sensibilidade: 9,45 mV/g; Massa: 0,5g, fixado à

estrutura vibrante, é transmitido ao analisador de sinais (Modelo: SR780; Network Signal Analyzer;

SRS – Stanford Research Systems). O analisador recebe os sinais enviados pelo acelorômetro para

processamento. Um analisador comumente usado é denominado analisador de Transformada

Rápida de Fourier (FFT). Este tipo de analisador recebe sinais analógicos de tensão do

acelerômetro (representados por deslocamento, velocidade, aceleração, deformação ou força) para

determinar as frequências naturais e obtenção dos gráficos: resposta no tempo e resposta de

frequência.

108

As grandezas físicas, da viga de aço em estudo, segue uma distribuição normal de

probabilidade com nível de confiança de 95% e são apresentadas na Tabela (4.7).

Tabela 4.7 – Grandezas física e geométrica da viga de aço ensaiada.

Grandezas Valor médio Incerteza padrão expandida (95% de confiança)

(g)a

m 0,43420 0,00017±

(g)vm 54,61932 ± 0,00026

(mm)u

L 250,415 ± 0,086

(mm)t

L 305,942 ± 0,093

(mm)b 26,069 ± 0,051

(mm)h 0,954 ± 0,022

3(kg / m )ρ 7.010,90 ± 167,60

(GPa)E 176,51 ± 12,33

4(m )I 121,89 10−

× ±120,13 10−

×

2(m )A 52,490 10−

× ±50,058 10−

×

A Tabela (4.8) apresenta as grandezas físicas intervalares, com uma distribuição uniforme de

confiança de 100%, correspondentes às grandezas da Tabela (4.7).

Tabela 4.8 – Parâmetros intervalares da viga de aço ensaiada.

Grandezas Limite inferior Limite superior

(g)a

m 0,43397 0,43432

(g)vm 54,61905 54,61957

(mm)u

L 250,330 250,500

(mm)t

L 305,847 306,038

(mm)b 26,018 26,120

(mm)h 0,932 0,976

3(kg / m )ρ 6843,31 7178,50

(GPa)E 164,18 188,84

4(m )I 121,76 10−

× 122,02 10−

×

2(m )A 52,425 10−

× 52,549 10−

×

109

A comparação dos resultados numéricos com os experimentais é apresentada na Tabela

(4.9). Os resultados experimentais foram obtidos conforme descrito no capítulo 3, item 3.3.

Calculou-se a incerteza expandida das três primeiras frequências naturais, influenciadas pelas

incertezas do tipo A (incertezas aleatórias) e do tipo B (incerteza do instrumento de medição).

Observa-se concordância entre os resultados, já que alguns dos pontos, presentes nos intervalos, são

coincidentes. Segundo Matos (2007), quando pelo um dos valores do intervalo numérico coincidir

com um valor do intervalo experimental, pode-se afirmar que a estrutura encontra-se em bom

estado físico. Em contrapartida, quando os valores não concordam, é possível que haja danos na

estrutura ou o modelo numérico não é adequado ao sistema em estudo.

Tabela 4.9 – Comparação entre os resultados numérico e experimentais.

A Tabela (4.10) apresenta a diferença relativa das três primeiras frequências naturais obtidas

pelo método numérico em relação às frequências obtidas experimentalmente.

Tabela 4.10 – Diferença relativa do resultado numérico em relação ao resultado experimental


Numerico Experimentalq q

fq Numericoq

f fr

f

−= × (%) ; 1,2,3q =

(do resultado numérico em relação ao resultado experimental)

Inferior Superior

1fr 1,482 9,162

2fr -1,035 7,360

3fr -2,573 310−

× 8,314

RESULTADOS

Frequência natural (Hz) NUMÉRICOS EXPERIMENTAIS

Limite inferior Limite superior Limite inferior Limite superior

1f 11,069 12,214 10,905 11,095

2f 69,372 76,544 70,09 70,91

3f 194,285 214,373 194,29 196,55

110

Capítulo 5

RESULTADOS

Neste capítulo, apresentam-se estruturas dinâmicas comuns na prática com o

objetivo de simular casos ainda não encontrados nas referências bibliográficas do assunto

abordado. Serão discutidos e analisados os autovalores e as frequências naturais de

sistemas dinâmicos, com resultados numéricos gerados pela autora.

5.1 Viga

O objetivo aqui é avaliar o tempo de execução do programa desenvolvido em

Matlab em função do aumento do número de elementos uma viga escalonada engastada,

conforme Figura (5.1).

Para este estudo, adotou-se parâmetros adimensionais e unitários, tais como: massa

específica, comprimento do elemento, área da seção transversal e momento de inércia de

área. Para o módulo de elasticidade, optou-se por um valor intervalar e adimensional,

[ ]0,99 , 1,01 , 1,2,3,..., .iE i n= =

111

Figura 5.1 – Viga escalonada em n segmentos.

Devido a inclinação da curva apresentada na Figura (5.2) observa-se que o tempo

de execução do algoritmo desenvolvido nesta tese, aumenta seguindo uma função

logarítmica, com o aumento do número de elementos e do GDL da viga.

50 100 150 200 250 3000

0.5

1

1.5

2

2.5

3

No. Elementos

tem

po (

s)

Figura 5.2 – Curva de inclinação do tempo de execução do algoritmo em função do

aumento do número de elementos da viga da Figura (5.1).

5.2 Treliças

Uma treliça é uma estrutura reticulada composta por elementos unidos nas suas

extremidades, de modo a formar uma estrutura rígida. As treliças são formadas pelo

cruzamento de elementos retos (de madeira, aço, alumínio, etc) cujas extremidades são

ligadas em pontos conhecidos como nós. As treliças são dimensionadas de modo que as

únicas forças atuantes nos elementos sejam de tração ou compressão. As treliças podem

112

ser planas ou tridimensionais. Pontes, suportes para telhados, guindastes e outras

estruturas similares são exemplos comuns de treliças, conforme Figura (5.3).

Treliças Planas: o elemento básico de uma treliça plana é o triângulo. Três barras

unidas por pinos em suas extremidades constitui uma estrutura rígida. Em especial, as

treliças planas se situam em um único plano e geralmente são usadas para sustentar

telhados e pontes.

Treliças tridimensionais: a idealização de uma treliça tridimensional consiste em

barras rígidas, conectadas por rótulas em suas extremidades. A unidade básica estável é

composta por seis barras conectadas em suas extremidades de modo a forma arestas de

um tetraedro.

a) Suportes para telhados b) Guindaste

Figura 5.3 – Exemplos de treliças.

5.2.1 Treliça com três elementos de barras

Uma treliça, Figura (5.4), formada por 3 elementos de barras e 2 graus de liberdade

efetivos é apresentada com imprecisão em seus módulos de elasticidades. As

frequências naturais adimensionais são dadas pela Equação (5.1) e os limites do intervalo

das frequências são obtidos pelas Equações (5.3) e (5.4)

113

/i

i

L

E

ω

ρΩ = (5.1)

Figura 5.4 – Treliça com três elementos de barras.

Conforme já mencionado, o problema de autovalor generalizado para sistemas

dinâmicos é dado por:

[ ] [ ]( ) ( ) 0K M xλ− = (5.2)

na qual 2

nλ ω= , λ representa os autovalores, n

ω , as frequências naturais, x , os modos

de vibração, [ ]K a matriz de rigidez global e [ ]M a matriz de massa global.

O problema de autovalor pseudo-determinístico formulado para a treliça da Figura

(5.4), na presença de incertezas no módulo de elasticidade, somente, em dois dos

elementos de barras, como 1 1 1[ , ] ([0,85 ; 1,1])l uE E E E= =

(primeiro elemento de barra)

e 2 2 2[ , ] ([0, 75 ; 1, 4 ])l uE E E E= =

(segundo elemento) já o terceiro elemento não possui

incertezas, 3 3 3[ , ] ([1 ; 1])l uE E E E= =

.

114

1 3 1 3 3 3

2

1 3 1 3 3 max 4

53 3 1 3

3 3 1 0 0 0

3 3 3 3 3 ( ) 0 1 0 04

0 0 1 03 4

u u u u u

u u u u u

u u u u

E E E E E xA

E E E E E AL xL

xE E E E

ω ρ

+ − − − + − = − +

(5.3)

1 3 1 3 3 3

2

1 3 1 3 3 min 4

53 3 1 3

3 3 1 0 0 0

3 3 3 3 3 ( ) 0 1 0 04

0 0 1 03 4

l l l l l

l l l l l

l l l l

E E E E E xA

E E E E E AL xL

xE E E E

ω ρ

+ − − − + − = − +

(5.4)

A Tabela (5.1) apresenta as frequências naturais adimensionais para a treliça da

Figura (5.4) em uma situação com parâmetros intervalares e uma outra, em que todos os

parâmetros são determinísticos. Lembrando que as frequências naturais foram

adimensionalisadas por meio da Equação (5.1), neste caso, devem-se adotar valores

unitários para os parâmetros 1L A Eρ= = = = .

Tabela 5.1 – Frequências naturais adimensionais de uma treliça simples composta por três

elementos de barras para módulo de elasticidade intervalar e módulo determinístico.

Frequência Natural Adimensional

Módulo intervalar

1 1 1[ , ] ([0,85 ; 1,1])l uE E E E= =

2 2 2[ , ] ([0, 75 ; 1, 4 ])l uE E E E= =

3 3 3[ , ] ([1 ; 1])l uE E E E= =

Módulo determinístico

1

2

3

0,975

1,075

E E

E E

E E

=

=

=

Limite Inferior Limite Superior

1Ω 0,7723 0,9578 0,8752

2Ω 1,3429 1,6830 1,5195

5.2.2 Treliça com 10 elementos de barras

115

Com uma treliça de dez elementos de barras, dois graus de liberdade por nó e um

total de nove GDL’s efetivos, deseja-se avaliar o grau de superestimação com o aumento

o aumento gradativo do intervalo dos parâmetros de área. A estrutura em análise é

apresentada na Figura (5.5).

Figura 5.5 – Treliça com 10 elementos de barras.

Os autovalores apresentados são adimensionais. As grandezas intervalares são

dadas na forma de números fracionários, como é o caso da área da seção transversal dos

elementos 1,2,3,4 e 6, únicas grandezas intervalares, cujos intervalos variam com a

mudança gradativa de um fator de incerteza ( β ), conforme Equação (5.5.)

,I c c c c

iA A A A Aβ β = − +

(5.5)

na qual o valor médio da área da seção transversal para 1,2,3,4,6i = é 1,0cA = e o fator

de incerteza varia dentro do intervalo 0 3%β≤ ≤ .

Por se tratar de uma análise adimensional, os demais parâmetros adotados são

unitários, tais como, módulo de elasticidade, massa específica, comprimento, momento

de inércia.

A avaliação do grau de superestimação é feito comparando os resultados obtidos

pelo método da autora (Método da Perturbação de Matrizes) com os resultados obtidos

pelo Método de Monte Carlo. Observa-se nas Figuras ((5.6a), (5.6b) e (5.6c)) que os

autovalores resultantes crescem à medida que aumenta o intervalo do parâmetro de

entrada do problema. O intervalo resultante dos autovalores obtidos por Monte Carlo são

mais estreitos do que os obtidos pela autora, indicando uma superestimação aceitável,

quando ocorre um aumento da grandeza intervalar.

116

Figura (5.6a) – Autovalores obtidos pela

autora.

Figura (5.6b) - Autovalores obtidos por

Monte Carlo.

0 0.5 1 1.5 2 2.5 30.0375

0.038

0.0385

0.039

0.0395

0.04

0.0405

0.041

0.0415

Beta %

λ1

Figura (5.6c) – Comparando os autovalores apresentados na Fig. (5.6a) com os da Fig.

(5.6b).

Figura 5.6 – Autovalores resultantes de uma treliça, com parâmetros de entrada

intervalares, variando em função de um fator de incerteza no intervalo 0 3%β≤ ≤ .

A Figura (5.7) mostra o comportamento da superestimação dos resultados com o

aumento do intervalo dos autovalores, à medida que o fator de incerteza da área assume

valores crescentes dentro do intervalo 0 3%β≤ ≤ . Observe nas Figuras (5.7a) e (5.7b)

que o aumento da superestimação dos resultados segue uma função linear à medida que

aumenta o intervalo das grandezas de entrada da treliça. A Figura (5.7a) mostra

desenvolvimento da superestimação trabalhando com os limites inferiores das grandezas

117

intervalares e a Figura (5.7b) mostra o comportamento da função linear na presença dos

limites superiores.

0 0.5 1 1.5 2 2.5 30

0.005

0.01

0.015

Beta %

ri =

(la

mbda

1im

c -

lam

bda

1ia

uto

ra)

/ la

mbda

1im

c

Figura (5.7a) – Superestimação obtida

trabalhando com os limites inferiores das

grandezas intervalares.

0 0.5 1 1.5 2 2.5 30

0.005

0.01

0.015

Beta %

rs =

(la

mbda

1s

mc -

lam

bda

1s

au

tora

) /

lam

bda

1sm

c

Figura (5.7b) – Superestimação

trabalhando com os limites superiores das

grandezas intervalares.

Figura 5.7 – Avaliando o comportamento da superestimação dos resultados com o

aumento do intervalo dos autovalores.

No gráfico da Figura (5.7), o eixo da vertical indicado, ,i sr , é a diferença relativa

dos autovalores obtidos pela autora com os obtidos por Monte Carlo, conforme Equação

(5.6) .

, ,

, ,

i s i s

MC autorai s i s

MC

rλ λ

λ

−= (5.6)

na qual i,s

autoraλ representam os autovalores gerados pela Teoria da Perturbação de Matrizes

e i,s

MCλ são os autovalores gerados por Monte Carlo.

118

Exemplo: Os gráficos apresentados nas Figuras (5.6) e (5.7) foram obtidos a

partir dos autovalores intervalares da primeira frequência natural. O exemplo dado na

Tabela (5.2) é obtido a partir dos autovalores gerados para 1%β =

Tabela 5.2 – Diferença relativa dos autovalores para um intervalo em 1%β = para a

Figura (5.5).

Limite inferior Limite superior

Autovalores de

1 1 1,i s

autora autora autoraλ λ λ =

obtidos pela autora

usando a Teoria da

Perturbação de Matrizes

i

1 0,03837autoraλ =

s

1 0,03908autoraλ =

Autovalores de

1 1 1,i s

MC MC MCλ λ λ =

obtidos pelo Método de

Monte Carlo

1

i

MCλ = 0,03857 1 0,03886s

MCλ =

1 1

1

1

MC autora

MC

rλ λ

λ

−=

3

1 5,185 10ir −= ×

3

1 7, 463 10sr −= ×

5.2.3 Treliça com vinte elementos de barras

O segundo caso em estudo envolve uma treliça formada por 20 elementos de

barras, 10 nós e 17 graus de liberdade efetivos, ver Figura (5.8). Cada um dos elementos

de barra possui grandezas físicas e geométricas tais como, massa específica, área da seção

transversal, comprimento de barra e módulo de elasticidade. Os elementos apresentados

na horizontal e na vertical possuem comprimentos iguais, como os elementos de

números: 2, 3, 4, 7, 8, 9, 12, 13, 14, 17, 18 e 20.

119

Neste tópico, tem-se o interesse em apresentar as frequências naturais

adimensionais de uma treliça, por isso não cabe aqui adotar valores das grandezas citadas,

apresentam-se, somente, os intervalos fracionais das grandezas intervalares. Para as

demais grandezas, devem-se adotar valores unitários 1L A Eρ= = = = .

Figura 5.8 – Treliça com vinte elementos de barras.

Na tabela (5.3) apresentam-se as frequências naturais adimensionais da treliça da

Figura (5.8) considerando duas situações diferentes. Na primeira situação, adotam-se

incertezas para os módulos de elasticidades e os módulos intervalares adotados são:

[ , ] ([0,95 ; 1,1])l u

i i iE E E E= =

nos seguintes elementos de barras, 1, 2, 7, 12, 14,15=i .

Na segunda situação, trabalha-se com a treliça sem parâmetros intervalares, neste caso

todos os parâmetros são determinísticos e o valor adotado para as barras

1, 2, 7, 12, 14 15i e= são 1, 025iE E= . Para os demais elementos são adotados valores

unitários, assim como considerado para a treliça da Figura (5.8).

120

Tabela 5.3 – Intervalo das frequências naturais adimensionais de uma treliça

simples formada por vinte elementos de barras para módulo de elasticidade intervalar e

módulo determinístico.

Frequência Natural Adimensional

MÓDULO INTERVALAR

[ , ] ([0,95 ; 1,1])l u

i i iE E E E= =

para 1, 2, 7, 12, 14,15=i

MÓDULO

DETERMINÍSTICO

1, 025iE E=

para 1, 2, 7, 12, 14,15=i

Limite Inferior

Limite Superior

1Ω 0,089294096416283 0,092911617188943 0,091179317224561

2Ω 0,190522742717331 0,200660661479653 0,195722703294325

3Ω 0,304779378690098 0,313520088341272 0,309365182581101

4Ω 0,586813928615123 0,604928832388361 0,596246044799551

5Ω 0,614781014581671 0,636345075436610 0,626020512338296

6Ω 0,864044999138466 0,874298089397304 0,869447484342737

7Ω 0,983916702250502 1,003088393175482 0,994455099348437

8Ω 1,016705183476271 1,026911670162670 1,021371310736960

9Ω 1,076825163527033 1,080481458562693 1,078778274392623

10Ω 1,201611495824762 1,230857996915977 1,217916264335618

11Ω 1,243547006095110 1,266923115964340 1,254015880515910

12Ω 1,306897451183824 1,342594550568073 1,325280798330094

13Ω 1,396403441668350 1,438844021477956 1,420445383378869

14Ω 1,457024873750229 1,478820604627931 1,466219682109757

15Ω 1,470516990347710 1,482181692790799 1,476098228582219

16Ω 1,693041208123021 1,733324422044810 1,712302682624909

17Ω 2,017667566100262 2,059802193324962 2,038205602594668

121

5.3 pórtico

Pórticos são estruturas compostas por elementos lineares (normalmente vigas e

colunas), conectados em suas extremidades de forma a não permitir rotações relativas

(conexões rígidas). Pórticos são capazes de resistir esforços normais, cortantes e,

principalmente, esforços de flexão. Nas edificações, normalmente são utilizados em um

padrão com repetições, resultando em estruturas hiperestáticas.

O caso em estudo nesta seção refere-se a um pórtico com 14 elementos, 15 nós e 39

graus de liberdade efetivos, ver Figura (5.9). Os elementos do pórtico de números 1, 2, 6,

12 e 13 são adotados com incertezas na área da seção transversal e no momento de

inércia, os demais parâmetros são considerados com valores fixos, como a massa

específica 37800 kg/mρ = e o módulo de elasticidade 200 GPaE = . Todos os

elementos possuem comprimentos iguais a 0,2 m.

1º e 12º elemento de viga:

Momento de inércia: 4 44, 0, 0999 10 , 0,[ ] [ 1001 ] m10l uI I I

− −×= ×=

Área da seção transversal: 2 22, 0,99 1[ 0 , 1, 01] 10[ ] ml uA A A

− −× ×= =

2º e 13º elemento de viga:

Momento de inércia: 4 44, 0,1998 10 , 0,[ ] [ 2002 ] m10l uI I I

− −×= ×=

Área da seção transversal: 2 22, 1, 426 10 , 1[ ] [ , 454 10 ] ml uA A A

− −× ×= =

6º elemento de viga:

Momento de inércia: 4 44, 0,04995 10 , 0,0[ 5005 10] [ ] ml uI I I

− −= × ×=

Área da seção transversal: 2 22, 0, 634 10 , 0[ ] [ , 646 10 ] ml uA A A

− −× ×= =

Para os demais elementos: 3, 4, 5, 7, 8, 9, 10, 11 e 14, adotaram-se os seguintes

parâmetros de valores fixos:

Momento de inércia: 4 40,1 1 m0I −×=

Área da seção transversal: 2 2101 mA −×=

122

Massa específica: 37800 kg/mρ =

Módulo de elasticidade: 200 GPaE = .

Figura 5.9 – Pórtico com 14 elementos de viga e 39 graus de liberdades.

Assim como feito nos casos anteriores, serão mostradas frequências naturais para

o pórtico da Figura (5.9) em uma situação na qual alguns dos parâmetros são intervalares

e, em outra situação, na qual todos os parâmetros são de valores fixos. Os valores fixos

foram obtidos por meio da média entre os limites dos parâmetros intervalares. Portanto os

parâmetros determinísticos adotados são:

Para o 1º e o 12º elemento de viga: 4 40,01 m10I−

×= e 2 21 10 mA −×=

.

Para o 2º e o 13º elemento de viga: 4 40, 2 1 m0I−

×= e 2 21, 44 m10A −×=

.

Para o 6º elemento de viga: 4 40, 05 m10I−

×= e 2 20, 64 m10A −×=

.

Para os demais: 3, 4, 5, 7, 8, 9, 10, 11 e 14, os parâmetros adotados foram:

4 40,1 1 m0I −×= , 2 2101 mA −

×= , 37800 kg/mρ = e 200 GPaE =

Na Tabela (5.3) são apresentados os autovalores para o pórtico com parâmetros

intervalares e parâmetros determinísticos.

123

Tabela 5.4 – Intervalo dos autovalores do pórtico da Figura (5.9).

Autovalores

Parâmetros intervalares Parâmetros determinísticos

Limite Inferior Limite Superior

1λ 5,9953668744 510× 6,0137826935 5

10× 6,0045662936 510×

2λ 3,47510348019 610× 3,49523222982 6

10× 3,48516314426 610×

3λ 6,52479786731 610× 6,57035253351 6

10× 6,54747348169 610×

4λ 1,162146631440 710× 1,168281940056 7

10× 1,165209542076 710×

5λ 4,544173616292 710× 4,593217098825 7

10× 4,568792111650 710×

6λ 6,634489763163 710× 6,709407442299 7

10× 6,671856354174 710×

7λ 8,950821154738 710× 9,024151720400 7

10× 8,987416519663 710×

8λ 1,3206006855737 810× 1,3369356349779 8

10× 1,3287551774718 810×

9λ 2,1353188654430 810× 2,1572270914742 8

10× 2,1462483089827 810×

10λ 2,5192569919679 810× 2,5455604858173 8

10× 2,5323516386897 810×

11λ 2,8968728951151 810× 2,9405708216015 8

10× 2,9185910395755 810×

12λ 3,7498304982751 810× 3,7800804912494 8

10× 3,7649137250068 810×

13λ 4,6042133013289 810× 4,6676538425856 8

10× 4,6358258831914 810×

14λ 5,4960180414309 810× 5,5613219566965 8

10× 5,5284712013936 810×

15λ 7,3108465622893 810× 7,3754708830533 8

10× 7,3430992111134 810×

16λ 9,0738638440035 810× 9,1569639401521 8

10× 9,1153435731282 810×

17λ 1,10147237064439 910× 1,11347748038401 9

10× 1,10747965961050 910×

18λ 1,23722073989728 910× 1,25036373970076 9

10× 1,24375236103370 910×

19λ 1,68705672505858 910× 1,70602654046057 9

10× 1,69652838096522 910×

20λ 1,87999570173437 910× 1,89623047878967 9

10× 1,88810767697514 910×

21λ 2,35821044409556 910× 2,37873493434071 9

10× 2,36845753616958 910×

22λ 3,03381937642706 910× 3,06819802648441 9

10× 3,05113972085246 910×

23λ 3,22251612516316 910× 3,25944452451049 9

10× 3,24084007576566 910×

24λ 3,73926431274539 910× 3,78778111550643 9

10× 3,76352193430481 910×

25λ 4,20211104046095 910× 4,26150597423892 9

10× 4,23234786357071 910×

26λ 4,25630059483613 910× 4,31618058211313 9

10× 4,28556983839114 910×

27λ 5,10086811426856 910× 5,16897664243100 9

10× 5,13451129791767 910×

28λ 6,20177347213656 910× 6,27957979811479 9

10× 6,24016395522510 910×

29λ 6,44404396987963 910× 6,48339725184557 9

10× 6,46369702802998 910×

30λ 6,90506093550566 910× 6,95302049322776 9

10× 6,92855810435880 910×

31λ 7,61520284922480 910× 7,68565194468737 9

10× 7,65020858679939 910×

32λ 9,14530811665328 910× 9,24847772302015 9

10× 9,19629484656660 910×

33λ 2,21716325969128 1010× 1,228471235319194 1010× 1,225078277149106 1010×

34λ 1,539181431226008 1010× 1,549099926167509 1010× 1,544108304049372 1010×

35λ 2,272488653930983 1010× 2,283519484706917 1010× 2,278003262758140 1010×

36λ 2,858065763586156 1010× 2,885843311641310 1010× 2,871878764551931 1010×

37λ 3,307753932142422 1010× 3,334788706501518 1010× 3,321242424301232 1010×

38λ 3,933228048294351 1010× 3,962897443669482 1010× 3,947978458384904 1010×

39λ 4,317287636191113 1010× 4,357233362051358 1010× 4,336822160393520 1010×

124

Capítulo 6

CONCLUSÕES E TRABALHOS FUTUROS

6.1 CONCLUSÕES

Neste trabalho apresentou uma metodologia para a avaliação da influencia dos

parâmetros intervalares na resposta dinâmica de um sistema, utilizando o método dos

elementos finitos. Desenvolveu-se um programa computacional, aplicando a

metodologia apresentada, capaz de gerar resultados confiáveis em um curto período de

tempo e com grande possibilidade de aplicação na prática.

Diversos trabalhos científicos apresentam os autovalores calculados com matrizes

de massa desenvolvidas a partir de elementos com deslocamentos unidirecionais.

Porém, segundo Hutton (2004), somente as estruturas estáticas podem ser estudadas

com tais matrizes de massa. Logo, já que o nosso estudo é dedicado a estruturas

dinâmicas, calcularam-se os autovalores (frequências naturais) a partir de matrizes de

massa desenvolvidas para elementos com deslocamentos nas direções axial e

transversal, conforme Equação (2.101) apresentada no Capítulo 2. O estudo dos

autovalores a partir das matrizes globais corretas é de extrema importância para a

obtenção de bons resultados.

A validação teórica foi efetivada comparando a Teoria da Perturbação de

Matrizes com o Parameter Solution Vertex Theorem (Qiu et al.,2005) e com os

125

resultados obtidos por meio das simulações de Monte Carlo. A diferença relativa dos

resultados teóricos em relação aos resultados de Monte Carlo variam de 0,0013% à

7,2%, comprovando que a metodologia desenvolvido não apresentou problemas de

superestimação dos resultados.

No capítulo 5, avaliou-se o tempo de execução do programa desenvolvido em

Matlab, em função do aumento do número de elementos de uma viga. Observou-se que

a inclinação da curva da Figura (5.2) mostra um crescimento do tempo de execução

seguindo uma função logarítmica, com o aumento do número de elementos e dos

GDL’s da viga.

Também no capítulo 5, avaliou-se o grau de superestimação comparando os

resultados obtidos pelo método da autora (Método da Perturbação de Matrizes) com os

resultados obtidos pelo Método de Monte Carlo. As Figuras ((5.6a), (5.6b) e (5.6c))

mostram que os autovalores resultantes crescem à medida que aumenta o intervalo do

parâmetro de entrada do problema. O intervalo resultante dos autovalores obtidos por

Monte Carlo são mais estreitos do que os obtidos pela autora, indicando uma

superestimação aceitável, quando ocorre um aumento da grandeza intervalar. A Figura

(5.7) mostra que o aumento da superestimação dos resultados segue uma função linear

com o aumento do intervalo dos autovalores, à medida que cresce o fator de incerteza,

das grandezas de entrada da treliça, dentro do intervalo 0 3%β≤ ≤ .

Neste trabalho encontram-se também resultados (autovalores e frequências

naturais) com aplicações diretas na prática da engenharia mecânica estrutural (como

treliças e pórticos), a fim de fazer a nossa contribuição à literatura com abordagem aos

parâmetros intervalares.

6.2 TRABALHOS FUTUROS

Para trabalhos futuros sugere-se que seja adaptado o programa computacional

desenvolvido pela autora a sistemas dinâmicos sujeitos a forças externas

intervalares.

126

Outra proposta de trabalho é desenvolver metodologias para sistemas sob

vibração amortecida intervalar.

Sugere-se também, que sejam realizados ensaios experimentais com outras

condições de contorno com o intuito de validar os modelos teóricos.

Por fim, sugere-se que sejam desenvolvidos trabalhos com parâmetros incertos

ocasionados por imprecisões na modelagem do sistema.

Investigar aplicações em materiais compósitos.

Investigar quantificação de incerteza utilizando matemática intervalar.

Analisar estruturas dinâmicas considerando os efeitos do cisalhamento.

Por fim, sugere-se que sejam desenvolvidos trabalhos com parâmetros incertos

ocasionados por imprecisões na modelagem do sistema.

127

REFERÊNCIAS BIBLIOGRÁFICAS

ADHIKARI, S.; FRISWELL, M. I.; LONKA, K.; SARKAR, A. (2009),

“Experimental Case Studies for Uncertainty Quantification in Structural

Dynamics”; Journal of Science Direct, Elsevier, Engenharia Mecânica

Probabilistica, Vol. 24, Issue 4, p. 473 – 492.

AGARWAL, H.; RENAUD, J. E.; PRESTON, E. L.; PADMANABHAN, D. (2004),

Uncertainty quantification using evidence theory in multidisciplinary design

optimization, Reliability Engineering and System Safety, vol. 85(1-3), p. 281-294.

ALTUS, E. , TOTRY, E. , GIVLI, S. (2005), “Optimized functional perturbation

method and morphology based effective properties of randomly heterogeneous

beams”, International Journal of Solids and Structures, vol. 42 (8), p. 2345 –

2359.

AYYUB, B. M.; GUPTA, M. M. (1997), "Uncertainty Analysis in Engineering and

Sciences: Fuzzy Logic, Statistics, and Neural Network approach”, Kluwer

Academic Publishers, Boston.

BATHE, K. J.; WILSON, E. L. (1976), Numerical Methods in Finit Element Analysis,

Practice-Hall, Inc., Englewood Cliffs, New Jersey, USA.

BELLMAN, R. (1960), Introduction to Matrix Analysis, McGrawHill, New York.

128

BIONDINI, F., BONTEMPI, F., MALERBA, P. G. (2004), “Fuzzy reliability

analysis of concrete structures”, Computational Structures, Vol. 82(13-14): 1033-

1052.

BURKILL, J. C. (1924), “Functions of intervals”, Proceedings of the London

Mathematical Society, v. 22, pp. 375-446.

DEIF, A. (1991), Advanced Matrix Theory for Scientists and Engineers, Gordon and

Breach Science Publishers, London, 2nd

.

DELAURENTIS, D.; MARVRIS, D. (2000), “Uncertainty modelling and

management in multidisciplinary analysis and synthesis”. Proceedings of the 38th

AIAA Aerospace Sciences Meeting, paper No. AIAA 2000-0422, Reno, NV,

January, p. 10-13.

DRAPER, N. R. ; SMITH, H. (1996); Applied regression analysis, New York, USA,

Wiley, p. 407.

DU, X.; SUDJIANTO, A. (2004), “A saddle point approximation method for

uncertainty analysis”, Proceedings of the DETC’04 – ASME, International Design

Engineering Technical Conferences and the Computers and Information in

Engineering Conference, Salt Lake City, Utah, USA.

EIBL, J. , SCHMIDT-HURTIENNE, B. (1995), “General outline of a new safety

format”, CEB, Bulletin d’Information, vol. 229, p. 33-48.

FISHMAN, G. S. (1996), Monte Carlo: Concepts, Algorithms, and Applications,

Springer Series in Operations Research, Springer Verlag, 2 ed., ISBN:

038794527X.

GUM (2008), Guide to the expression of uncertainty in measurement, Geneva, p. 131.

HAFTKA, R. T.; ROSCA, R. I.; NIKOLAIDIS, E. (2006), “An approach for testing

methods for modelling uncertainty”, Journal of Mechanical Design, Vol. 128(5),

p. 1038-1049.

129

HALLOT, C. V.; BARTLETT, A. C. (1987), “ On the eigenvalues of interval

matrices”, Tech. Rep. Department of Electric Engineering And Computer

Engineering., University of Massachusetts, Amherst, MA, USA, Also,

Proceedings of the 26th IEEE, Conference on Decision and Control, Los Angeles

, p. 794–799.

HUANG, B.; XIAOPING, D. (2006), “Uncertainty analysis by dimension reduction

integration and saddle point approximations”, Journal of Mechanic Design, Vol.

128(1), p. 26-33.

HUDAK, D. (1984), “Eigenwertproblem Fur Intervall-Matrizen”, Journal of Applied

Mathematics and Mechanics, vol. 64, 503–505.

HURTADO, J. E. (2002), “Analysis of one-dimensional stochastic finite elements

using neural networks, Probabilistic Engineering Mechanics, vol. 17(1): 35-44”.

HUTTON, D. V. (2004), Fundamentals of Finite Element Analysis, The McGraw-Hill

Companies, p. 500.

KASSIMALI, ASLAM (2011), Matrix Analysis of Structures, Cengage Learning

publisher, EUA, 2º ed., p. 656.

KIUREGHIAN, A. D.; DITLEVSEN, O. (2009), “Aleatory or epistemic? Does it

matter?”, Structural Safety, vol. 31, p. 105–112.

KLIR, G. J. (1995), “Analysis of one-dimensional stochastic finite elements using

neural networks”, Probabilistic Engineering Mechanics, vol. 17 (1), p. 35-44”.

KULPA, Z.; POWNUK, A.; SKALNA I. (1998), “Analysis of linear mechanical

structures with uncertainties by means of interval methods”, Computer Assisted

Mechanics and Engineering Sciences, vol. 5, p.443-477.

LIMBOURG, P. ; APONTE, D. E. (2005), An optimization algorithm for imprecise

multi-objective problem functions, Proceedings of IEEE Congress on

Evolutionary Computation, IEEE Press, p. 459 466

130

LINK, W. (1997), Metrologia mecânica: Expressão da incerteza de medição, Instituto

de Pesquisas Tecnológicas, INMETRO, Rio de Janeiro, p.174.

LIU, J. S. (2001), Monte Carlo Strategies in Scientific Computing, Springer Series in

Statistics, Springer Verlag, ISBN: 0387952306.

MAGLARAS, G. , NIKOLAIDIS, E. , HAFTKA, R. T. , CUDNEY, H. H. (1997)

“Analytical experimental comparison of probabilistic methods and fuzzy set based

methods for designing under uncertainty”, Structural Multidisciplinary

Optimization, vol. 13 (2-3), p. 69-80.

MAHADEVAN, S. ; RAGHOTHAMACHAR, P. (2000), “Adaptive simulation for

system reliability analysis of large structures”, Computational Structures, Vol.

77(6), p. 725-734.

MATOS, J. C. , (2007), “Tratamento de Incertezas em Modelos Numéricos na

Engenharia Civil”, Tese de Doutorado, Engenharia Civil, Faculdade de

Engenharia Universidade do Porto - FEUP, Porto, Portugal, p. 253.

MATSUMOTO, M.; IWAYA, E. (2000), “Interval finite analysis to structural

system”, IEEE International Workshop on Robot and Human Interactive

Communication, p. 132-136.

MEHDI, M. (2005), “Dynamic Analisys of Strutuctures With Interval Uncertainty,

Tese de Doutorado”, Departamento de Engenharia Civil, Case Western Reserve

University, USA, p. 96 .

MOORE, R. E. (1962), “Interval Arithmetic and Automatic Error Analysis in Digital

Computing”, PhD dissertation, Stanford University, USA.

MUHANNA, R. L.; MULLEN, R. L. (1999), “Bounds of structural response for all

possible loading”, Journal of Structural Engineering, ASCE, vol. 125 (1), p. 98-

106.

131

MUHANNA, R. L.; MULLEN, R. L. (2001), “Uncertainty in mechanics problems

interval based approach”, Journal of Engineering Mechanics, 127(6):557-566.

NEUMAIER, A. (1990), Interval Methods for Systems of Equations, Cambridge Univ.

Press, New York.

OBERKAMPF, W. L.; DELAND, S.; RUTHERFORD, B. M.; DIEGERT, K. V.;

ALVIN, K. F. (1999), A new methodology for the estimation of total uncertainty

in computational simulation, 40th AIAA / ASME / ASCE / AHS / ASC structures,

Structural Dynamics and Materials Conference. AIAA, vol. 99 (1612), p. 3061 –

3083.

OLIVEIRA, W. C. (2006), “Método dos elementos finitos aplicado na mecânica

estrutural”, Universidade Federal de Itajubá - UNIFEI, Itajubá, Brasil, p. 494.

OLSSON, A. ; SANDBERG, G.; DAHLBLOM, O. (2003), “On Latin Hypercube

Sampling for structural reliability analysis”, Structural Safety, vol. 25 (1), p. 47-

68.

PAPADRAKAKIS, M. ; LAGAROS, N. D. (2002), “Reliability-based structural

optimization using neural networks and Monte Carlo simulation”, Computational

Methods Applied to Mechanical Engineering, Vol. 191, p. 3491–3507.

PEREIRA, S. C. A. (2002), “Tratamento de Incertezas em Modelagens de Bacias”,

Tese de Doutorado, Engenharia Civil, Universidade Federal do Rio de Janeiro -

UFRJ , Rio de Janeiro, p. 312.

QIU, Z.; CHEN, S. H.; ELISHAKOFF, I. (1995), “Natural frequencies of structures

with uncertain-but-non-random parameters”, Journal of Optimization Theory and

Applications , vol. 86, p. 669–683.

QIU, Z.; CHEN, S. H.; ELISHAKOFF, I. (1996), “Bounds of eigenvalues for

structures with an interval description of uncertain-but-non-random parameters”,

Chaos, Solitons & Fractals, vol. 7, p. 425–434.

132

QIU, Z.; CHEN, S. H.; NA, J. X. (1993), “The Rayleigh Quotient method for

computing eigenvalue bounds of vibrational systems with interval parameters”,

Acta Mechanics Solid Sinica, vol.6, p. 309–318.

QIU, Z.; WANG, X.; FRISWELL, M. I. (2005), “Eigenvalue bounds of structures

with uncertain-but-bounded parameters”, Journal of Sound and Vibration , vol.

282, p. 297–312.

RAIS-ROHANI, M.; SINGH, M. N. (2004), “Comparison of global and local

response surface techniques in reliability-based optimization of composite

structures”, Structures and Multidisciplinary Optimization, Vol. 26(5): 333-345.

RAO, S. S., CHEN, L (1998), “Numerical solution of fuzzy linear equations in

engineering analysis”, International Journal of Numerical Methods in

Engineering, vol. 43, p. 391-408.

RAO, S. S.; SAWYER, P. (1995), “Fuzzy finite element approach for analysis of

imprecisely defined systems”, AIAA Journal, vol. 33(2), p. 2364-2370.

ROHN, J. (1987), “Eigenvalues of a symmetric interval matrix”, Freiburger Intervall-

Berichte, vol. 87 (10), p. 67–72.

SANTOS, A. F. M. (2011), “Caracterização de um sistema medição de vibrações de

baixo custo para aplicações gerais” Dissertação de Mestrado, Engenharia

Mecânica, Universidade Federal de Uberlândia - UFU , Uberlândia, Minas Gerais,

p. 152.

SCHUELLER, G. I. (2001), “Computational stochastic mechanics - recent advances”,

Computational Structures, Vol. 79(22-25): 2225-2234.

SEGERLIND, L. J. (1984), “Applied finite elemento analysis”, John Willey & Sons,

USA, p. 411.

SUNAGA, T. (1958), “Theory of an interval algebra and its application to numerical

analysis”, RAAG Memoirs, Gaukutsu Bunkwen Fukeyu-kai, p. 29-46, Tokyo.

133

TESCHL, G. (2009), “Mathematical Method in Quantum Mechanics”, American

Mathematical Society, Rhode Island, USA , vol. 99, p. 117-119.

YOUNG, R. C. (1931), “The algebra of many-valued quantities”. Math. Ann., vol.104,

p. 260-290.

ZADEH, L. A. (1965), “Fuzzy Sets”, Information and Control, Vol. 8: 338-353.

ZIMMERMANN, H. J. (2000), “An application-oriented view of modeling

uncertainty. European”, Journal of Operational Research, vol. 122, p. 190-198.

ZOU, T.; MAHADEVAN, S. (2006), “Versatile formulation for multiobjective

reliability-based design optimization”, Journal of Mechanic Design, Vol. 128(6),

p. 1217-1226.

134

ANEXO A

MÉTODO DOS ELEMENTOS FINITOS:

FORMULAÇÃO DAS MATRIZES DE RIGIDEZ E DE

MASSA DE ELEMENTOS DE BARRA E DE VIGA

Os problemas dinâmicos discutidos neste trabalho são resolvidos pelo Método dos

Elementos Finitos. A base do MEF consiste em dividir o domínio do problema em

subdomínios (Bathe, 1976). Os subdomínios são chamados de elementos finitos e são

conectados aos elementos vizinhos por pontos denominados “pontos nodais” ou simplesmente

“nós”, Figura (A.1). Chama-se de malha o domínio do problema discretizado por elementos

finitos. Uma malha pode ser formada por várias maneiras diferentes, dependendo do tipo e

quantidade de elementos utilizados.

Figura A.1 – Malha formada por subdomínios (ou elementos finitos).

O método dos elementos finitos pode ser dividido em seis passos básicos (Oliveira,

2006).

135

1. Discretizar a região. Localizar e numerar os pontos nodais e especificar suas

coordenadas.

2. Especificar a equação para a aproximação das variáveis físicas. Esta equação é escrita

em termos das variáveis físicas do elemento. Dependendo do grau das funções de

interpolação, define-se o número de graus de liberdade do elemento.

3. Desenvolver o sistema de equações. Usando o método de Galerkin avalia-se a integral

dos resíduos ponderados. Na formulação da energia potencial, a energia do sistema é

escrita em função dos deslocamentos nodais e depois é minimizado com relação a

esses deslocamentos (Segerlind, 1976).

4. Fazer o acoplamento das equações a nível global. Inclui, neste passo, a entrada das

condições de contorno.

5. Resolver o sistema de equações.

6. Calcular as quantidades de interesse relacionadas aos elementos.

Neste trabalho, o sistema de equações desenvolvido é um problema de autovalor

generalizado para estruturas dinâmicas, [ ]( ) [ ] 0K M xλ− = , na qual a matriz de rigidez

global [ ]K e a matriz de massa [ ]M são obtidas por meio das matrizes dos elementos

discretizados. Quando a equação diferencial dos elementos é conhecida, aplica-se o método

dos resíduos ponderados, via formulação de Galerkin, para deduzir a matriz de rigidez dos

elementos. No entanto, caso a equação diferencial do elemento seja desconhecida, aplica-se o

Princípio da Energia Potencial Mínima (Oliveira, 2006). As formulações dos métodos dos

elementos finitos geram as matrizes globais do sistema por um processo conhecido por

assembly (ou soma dos elementos superpostos) onde as contribuições das rigidezes dos

elementos finitos, que compartilham um mesmo grau de liberdade, são somadas. O

desenvolvimento das matrizes de rigidez e de massa de um elemento barra e de viga é

apresentado nos tópicos a seguir.

136

A.1 ELEMENTO DE BARRA

A.1.1 Elemento de barra unidimensional

A equação do movimento de um elemento de barra sob carga axial, Equação (A.2), é

obtida aplicando a segunda Lei de Newton no diagrama de corpo livro mostrado na Figura

(A.2).

Figura A.2 – Diagrama de corpo livre de um elemento de barra unidimensional (fonte: Finit

Element Method using Matlab, 1996).

2

2

( , ) ( , )( ) ( , ) ( , )

u x t P x tA x dx P x t dx P x t

t x

∂ ∂ ρ = + −

∂ ∂ (A.2)

na qual ( , )u x t é o deslocamento axial ao longo do eixo x da barra, t é o eixo temporal, ρ é a

massa específica, ( , )P x t é a carga axial e A(x) é a área da seção transversal da barra.

Da lei de Hooke, tem-se:

( , )( , ) ( ) ( )

u x tP x t E x A x

tx

∂=

∂ (A.3)

na qual E é o modulo de elasticidade e ( , )u x t

tx

∂

∂ é a deformação específica longitudinal do

elemento unidimensional.

Substituindo A.3 em A.2, tem-se:

2

2

( , ) ( , )( ) ( ) ( )

u x t u x tA x dx A x E x

t x x

∂ ∂ ∂ ρ =

∂ ∂ ∂ (A.4)

137

Considerando a área da seção transversal e o módulo de elasticidade constantes ao

longo do comprimento da barra, tem-se:

2 2

2 2

( , ) ( , )( ) ( ) ( )

u x t u x tA x A x E x

t x

∂ ∂ρ =

∂ ∂ (A.5)

Aplicando o Método de Galerkin dos resíduos ponderados obtém-se a Equação (A.6)

Substituindo a Equação (A.4) e a Equação (A.3) na Equação (A.2), obtém-se:

2 2

2 20

( , ) ( , )( ) ( ) ( ) ( ) 0

L u x t u x tw x A x A x E x dx

t x

∂ ∂ρ − =

∂ ∂ ∫

(A.6)

Discretizando a barra, tem-se:

[ ] ( , ) ( , )u x t u x t H d= =

(A.7)

com:

[ ]L x x

HL L

− =

(funções de interpolação)

1

2

ud

u

=

(deslocamento dos nós 1 e 2, ver Figura (A.3))

Figura A.3 – Pontos nodais de um elemento de barra unidimensional (fonte: Finit Element

Method using Matlab, 1996).

Tomando a Equação (A.6), vem:

2 2

2 20 0

( , ) ( , )( ) ( ) ( ) ( ) ( ) 0

L Lu x t u x tA x w x dx A x E x w x dx

t x

∂ ∂ρ − =

∂ ∂∫ ∫

(A.8)

138

Resolvendo a segunda integral, da Equação (A.8), por partes, vem:

2

020 0

( , ) ( , ) ( ) ( , )( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )

L LLu x t u x t w x u x t

A x E x w x dx A x E x dx A x E x w xx x x x

∂ ∂ ∂ ∂− = −

∂ ∂ ∂ ∂∫ ∫

(A.9)

Aplicando a condição do contorno natural.

( , )( ) ( ) ( , ) 0

u x tA x E x P x t

x

∂− =

∂

(A.10)

na equação (A.9) e substituindo o resultado na Equação (A.8), vem:

2

020 0

( , ) ( , ) ( )( ) ( ) ( ) ( ) ( ) ( ) 0

L LLu x t u x t w x

A x w x dx A x E x dx w x P xt x x

∂ ∂ ∂ρ + − =

∂ ∂ ∂∫ ∫

(A.11)

Derivando a Equação (A.7) em relação às coordenadas cartesianase também em

relação ao tempo, obtém-se:

[ ] [ ] ( , ) ( , )u x t u x t d

H d B dx x dx

∂ ∂≅ = =

∂ ∂

(A.12)

[ ] [ ] ( )w x d

H d B dx dx

∂≅ =

∂

[ ] ( , ) ( , )u x t u x t

H dt t

∂ ∂≅ =

∂ ∂

[ ] 2 2

2 2

( , ) ( , )u x t u x tH d

t t

∂ ∂≅ =

∂ ∂

com,

[ ]1 1

BL L

= −

(A.13)

Aplicando a Equação (A.12) na (A.11) , vem:

[ ] [ ] 2

20 0

( , )( ) ( ) ( )

TL L Tu x tA x w x dx A x H H d dx

t

∂ρ = ρ

∂∫ ∫

(A.14)

139

[ ] [ ] 0 0

( , ) ( )( ) ( ) ( ) ( )

L L T Tu x t w xA x E x dx A x E x B B d dx

x x

∂ ∂=

∂ ∂∫ ∫

(A.15)

[ ]( ) ( , ) ( , )w x P x t H P x t= (A.16)

A Equação (A.17) é chamada de matriz de massa do elemento e Equação (A.18) é a

matriz de rigidez do elemento.

[ ] [ ]0

2 1( )

1 26

TLe AL

M A x H H dx ρ

= ρ =

∫ (A.17)

[ ] [ ]0

1 1( ) ( )

1 1

TLe AE

K E x A x B B dxL

− = = −

∫ (A.18)

A.1.1.1 Elemento de barra inclinado com deslocamentos nodais na

direção do eixo x.

O vetor de deslocamento de um elemento de barra bidimensional, conforme Figura

(A.4), é dada por:

1 1 2 2

Ted u v u v= (A.19)

Figura A.4 – Pontos nodais de um elemento de barra bidimensional (fonte: Finit Element

Method using Matlab, 1996).

A matriz de rigidez de uma elemento de barra bidimensional é:

140

1 0 1 0

0 0 0 0

1 0 1 0

0 0 0 0

e AEK

L

− = −

(A.20)

A relação entre o sistema de coordenadas local e o sistema global de um elemento de

barra inclinado, Figura (A.5), chamada de transformação de coordenadas, é apresentada na

Equação (A.21).

1 1

1 1

2 2

2 2

0 0

0 0

0 0

0 0

u uc s

v vs c

u uc s

v vs c

− = −

(A.21)

na qual cosc e s sen= β = β .

Reescrevendo a Equação (A.21), obtém-se:

[ ] e ed T d= (A.22)

Figura A.5 – Elemento de barra bidimensional inclinado (fonte: Finit Element Method using

Matlab, 1996).

Matriz de rigidez do elemento no sistema global de coordenadas é obtido por meio da

energia de deformação, conforme expressado na Equação (A.23).

141

1

2

Te e e

U d K d = (A.23)

Substituindo a Equação (A.22) na Equação (A.23) obtém-se:

[ ] [ ] 1

2

T Te e eU d T K T d = (A.24)

1

2

Te e e

U d K d = (A.25)

Na qual,

[ ] [ ]Te e

K T K T = (A.26)

Substituindo as Equações (A.20) e (A.21) na Equação (A.26) chega-se a matriz de

rigidez do elemento inclinado, conforme expressada na Equação (A.27).

2 2

2 2

2 2

2 2

e

c cs c cs

cs s cs sAEK

l c cs c cs

cs s cs s

− −

− − = − − − −

(A.27)

O vetor dos graus de liberdades nodais é dado por:

1 1 2 2

ed u v u v= (A.28)

A matriz de massa de um elemento inclinado é obtida usando a expressão da energia

cinética, assim como aplicada na matriz de rigidez do elemento. Dessa forma, segue a

Equação (A.29).

[ ] [ ]Te e

M T M T = (A.29)

142

Deve-se deixar claro que a Equação (A.30) é a matriz de massa de um elemento de

barra inclinado desenvolvida para corpos estáticos, diferente da matriz de massa de um

elemento de barra inclinado desenvolvido para corpos em vibração, conforme demonstrado no

Capítulo 2, Seção 2.2.3.2.1.

2 2

2 2

2 2

2 2

2 2

2 2

6 2 2

2 2

e

c cs c cs

cs s cs sALM

c cs c cs

cs s cs s

ρ

=

(A.30)

A.2 ELEMENTO DE VIGA

De acordo com o modelo de viga de Euler-Bernoulli, considera-se as seguintes

hipóteses:

• O comprimento da viga é muito maior que as dimensões da seção transversal, 10L

e≥ ,

sendo L o comprimento da viga e e, a dimensão da seção transversal, como largura ou

diâmetro.

• A viga é constituída de um material homegeneo, elástico linear;

• Planos perpendiculares a linha neutra, permanecem planos e perpendiculares após a

deformação;

• Os efeitos de cisalhamento e da inercia de rotação são desprezados.

A equação de viga de Euler-Bernoulli para a flexão de viga é dada por:

2 2 2

2 2 2

( , ) ( , )( ) ( ) ( ) ( , )

v x t v x tA x E x I x q x t

t x xρ

∂ ∂ ∂+ =

∂ ∂ ∂ (A.30)

na qual ( , )v x t é o deslocamento transversal da viga, ρ é a massa específica, EI é a rigidez da

viga, ( , )q x t é a carga aplicada, t e x são as coordenadas cartesianas, no tempo e no espaço,

respectivamente.

143

Na formulação do método dos elementos finitos, aplicou-se um dos métodos dos

resíduos ponderados, método de Galerkin, na equação da viga, Equação (A.30), obtendo-se:

2 2 2

2 2 20

( , ) ( , )( , ) ( , ) 0

L v x t v x tA EI q x t w x t dx

t x xρ ∂ ∂ ∂

+ − = ∂ ∂ ∂ ∫

(A.31)

na qual L é o comprimento da viga e w(x) é uma função de teste.

A formulação fraca da Equação (A.31) é obtida a partir da sua integração por partes.

2 2 2

2 2 20 0 00

( , ) ( , ) ( ) ( )( ) ( , ) ( ) ( , ) ( ) ( , ) 0

LL L Lv x t v x t w x w x

A w x dx EI dx q x t w x dx V x t w x M x tt x x x

ρ∂ ∂ ∂ ∂

+ − + − = ∂ ∂ ∂ ∂

∫ ∫ ∫

(A.32)

Sendo, ( , )V x t , o esforço cortante:

3

3

( , )( , )

v x tV x t EI

x

∂=

∂

(A.33)

E, ( , )M x t , o momento fletor:

2

2

( , )( , )

v x tM x t EI

x

∂=

∂

(A.34)

Considere o elemento de viga apresentado na Figura (A.5). O elemento de viga possui

dois nós e cada nó possui dois graus de liberdade, sendo que um dos graus de liberdade

refere-se à inclinação do nó e o outro, refere-se a deflexão. Sendo que a inclinação é dada pela

derivada primeira da deflexão em termos de x, dv dxθ = . Devido às quatro variáveis nodais

presentes no elemento de viga da Figura (A.5), assume-se uma função cúbica para ( )v x .

2 3

0 1 2 3( )v x c c x c x c x= + + + (A.35)

Considerando o modelo de viga de Euler-Bernoulli, a inclinação é dada por:

144

2

1 2 3( ) 2 3x c c x c xθ = + + (A.36)

Figura A.6 – Elemento de viga na horizontal (fonte: Finit Element Method using Matlab,

1996).

Avaliação da deflexão e da inclinação dos nós:

Em x=0,

0 1(0)v c v= =

1 1(0) cθ θ= =

(A.37)

Em x=L,

2 3

0 1 2 3 2( )v L c c L c L c L v= + + + =

2

1 2 3 2( ) 2 3L c c L c Lθ θ= + + =

Após algumas manipulações matemáticas chega-se a:

1 1 2 1 3 2 4 2( ) ( ) ( ) ( ) ( )v x H x v H x H x v H xθ θ= + + +

Na qual,

2 3

1 2 3

2 3

2 2

2 3

3 2 3

2 3

4 2

3 2( ) 1

2( )

3 2( )

( )

x xH x

L L

x xH x x

L L

x xH x

L L

x xH x

L L

= − +

= − +

= −

= − +

(A.38)

145

As funções mostradas na Equação (A.28) são chamadas de funções de interpolação

Hermitiana, a Figura (A.6) apresenta a forma gráfica das funções.

Figura A.7 – Funções de interpolação Hermitiana (fonte: Finit Element Method using Matlab,

1996).

A aplicação das funções de interpolação Hermitiana e do método de Galerkin no

segundo termo da Equação (A.32) resulta em:

[ ] [ ]0

l TeK B EI B dx = ∫ (A.39)

Na qual,

[ ][ ]

2

2

d HB

dx= (A.40)

E o correspondente vetor de graus de liberdade nodal é dado por:

1 1 2 2

Ted v vθ θ= (A.41)

Para esta formulação, a matriz de rigidez de um elemento de viga é:

2 2

3

2 2

12 6 12 6

6 4 6 2

12 6 12 6

6 2 6 4

e

L L

L L L LEIK

L LL

L L L L

−

− = − − −

−

(A.42)

146

Para uma análise dinâmica de vigas, a força inercial necessita ser incluída. Neste caso,

a deflexão transversal é uma função de x e t, do espaço e do tempo, respectivamente. A

deflexão é dada por:

1 1 2 1 3 2 4 2( , ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )v x t H x v t H x t H x v t H x tθ θ= + + + (A.43)

A Equação (A.43) indica que as funções de interpolação são usadas para interpolar a

deflexão no domínio do tempo. Fazendo as substituições pertinentes no primeiro termo da

Equação (A.32) obtém-se:

[ ] [ ] 0

L T eA H H dx dρ∫

(A.44)

na qual,

[ ] 1 2 3 4

TH H H H H= (A.45)

Neste caso, a matriz de massa consistente do elemento de viga é dada por:

[ ] [ ]0

2 2

2 2

156 22 54 13

22 4 13 3

54 13 156 22420

13 3 22 4

L Te

e

M A H H dx

L L

L L L LALM

L L

L L L L

ρ

ρ

=

−

− = − − − −

∫ (A.46)

A matriz de rigidez do elemento de pórticos inclui os efeitos de força axial, da força

cisalhante e do momento fletor. Para obter a matriz de rigidez basta somar a matriz de rigidez

do elemento de barra inclinado com a matriz de rigidez do elemento de viga inclinado.

Somando-se os efeitos da inércia provocada pelos esforços transversais, momentos fletores e

forças axiais, a matriz de massa total do elemento de pórtico pode ser obtida através da soma

da matriz de massa do elemento de barra inclinado pela matriz de massa do elemento de viga

inclinado (Oliveira, 2007).

A matriz de transformação é apresentada na Equação (A.47).

147

[ ]

cos 0 0 0 0

cos 0 0 0 0

0 0 1 0 0 0

0 0 0 cos 0

0 0 0 cos 0

0 0 0 0 0 1

sen

sen

Tsen

sen

θ θ

θ θ

θ θ

θ θ

−

= −

(A.47)

Matriz de rigidez do elemento de viga inclinado (elemento de um pórtico plano):

2 2 2 2

2 2 2 2

2 2 2

2 2

12 12 6 12 12 6cos cos cos cos

12 6 12 12cos cos cos

[ ]

i i i i i i

i i i

ei

i

I I I I I IA sen A sen sen A sen A sen sen

L L L L L L

I I I IAsen A sen Asen

L L L L

EK

L

θ θ θ θ θ θ θ θ θ θ

θ θ θ θ θ θ

+ − − − + − − −

+ − − − +

=

2

2

2 2

2 2

2 2

2

6cos cos

6 64 cos 2

12 12 6cos cos

12 6cos cos

4

i i

i i

i i i

i i

I

L

I II sen I

L L

I I IA sen A sen sen

L L L

I IAsen

L L

simétrica I

θ θ

θ θ

θ θ θ θ θ

θ θ θ

−

+ −

+ −

(A.48)

Matriz de massa total do elemento de viga inclinado (elemento de um pórtico plano):

2 2 2 2

2 2 2 2

2 2

2 2

4(39s 35cos ) 16 cos 22 (54s 70cos ) 16 cos 13

4(35s 39cos ) 22 cos 16 cos (70 54cos ) 13 cos

4 13 13 cos 3[ ]

420 4(39s 35cos ) 16

i i

i i

i i i ie ii

en sen L sen en sen L sen

en L sen sen L

L L sen L LALM

en sen

θ θ θ θ θ θ θ θ θ θ

θ θ θ θ θ θ θ θ

θ θρ

θ θ θ

+ − − +

+ − + −

− −=

+ −

2 2

2

cos 22

4(35s 39cos ) 22 cos

4

i

i

i

L sen

en L

simétrica L

θ θ

θ θ θ

+ −

(A.49)

148

ANEXO B

PARAMETER SOLUTION VERTEX THEOREM

No anexo B apresenta-se o teorema desenvolvido por Qiu et al. (2005) -

Parameter Solution Vertex theorem, que é baseado na teoria de otimização.

Para se entender o teorema de Qiu et al (2005) e obter a solução de um problema

de autovalor, é necessário conhecer as técnicas utilizadas para se obter a solução de um

sistema linear, como a Regra de Cramer. Neste caso, se a matriz de coeficientes do

sistema linear é não singular, então o sistema linear tem solução única, dada por

Equação (B.1).

[ ] A x b= ou 1

p

ji j i

j

a x b=

=∑ , 1, 2,...,i p= (B.1)

[ ] 1

x A b−

=

Sendo que a inversa da matriz dos coeficientes é dada por:

[ ][ ]( )

[ ]( )1 1

detA adj A

A

−

= (B.2)

149

Sendo [ ]( )adj A a matriz adjunta de [ ]A e [ ]( )det A , o determinante de [ ]A . A

matriz adjunta [ ]( )adj A é igual à transposta da matriz dos cofatores [ ]C , portanto,

[ ]( ) [ ]T

adj A C= .

Se [ ]C é a matriz dos cofatores de [ ]A que é formada por elementos ij

a , então,

cada elemento ij

c da matriz dos cofatores é obtida por:

( ) [ ]1i j

ijc A

+

= − (B.3)

sendo que [ ]A é o determinante de[ ]A menos a linha i e a coluna j de [ ]A .

Exemplo para se obter o primeiro elemento da matriz dos cofatores [ ]C

(B.4)

Logo, 11c fica assim:

( )

22 21 1

11

2

1

p

p pp

a a

c

a a

+

= −

Para se obter os demais elementos da matriz dos cofatores, deve-se realizar o

mesmo procedimento da Equação (B.4).

[ ]( )

1 11 21 1 1

2 12 22 2 2

1 2

1

det

p

p

p p p pp p

x c c c b

x c c c b

A

x c c c b

=

(B.5)

150

Assim, a entrada na i-ésima linha de x é:

[ ]( )

1 1 2 2

det

i i p pi

i

b c b c b cx

A

+ + +==

(B.6)

na qual b1, b2, b2,..., bp são entradas de b . Os cofatores nessa expressão vêm da i-

ésima coluna de [ ]A e, portanto, permanecem inalterados se trocarmos a i-ésima coluna

de [ ]A por b (pois a i-ésima coluna é suprimida quando calculamos os cofatores).

Como essa troca fornece a matriz ( )iA , o numerador da Equação (B.6) pode ser

interpretado como a expansão em cofatores ao longo da j-esima coluna dej

A . Assim,

A Equação (B.6) é equivalente a:

[ ] 1

1

det( )

p

i j ji

j

x b cA =

= ∑ (B.7)

ou,

( )[ ]( )

( )det

det

i

i

Ax

A

= (B.8)

na qual,

( )1 2

( ) ( )

1 1

1

det ( 1)j p

pi t i

j ji i i ji pi

j i

A b c a a b a=

= = − ∑ ∑ … … (B.9)

e

[ ]( )1 2

( )

1 1

1

det ( 1)j p

pt i

j ji i i ji pi

j i

A b c a a a a=

= = −∑ ∑ … … (B.10)

na qual ( )t i é o número de inversões na permutação das colunas da matriz que se deseja

obter o determinante.

151

Se as matrizes de um sistema linear são intervalares, têm-se que:

[ ]l u

A A A ≤ ≤ e l ub b b≤ ≤

(B.11)

[ ]

11 12 1 11 12 1 11 12 1

21 22 2 21 22 2 21 22 2

1 2 1 2 1 2

, ,

l l l u u up p p

l l l u u up p pl u

l l l u u up p pp p p pp p p pp

a a a a a a a a a

a a a a a a a a aA A A

a a a a a a a a a

= = =

… … …

… … …

ou,

[ ] ,l uA A A A ∈ =

(B.12)

E o vetor b está contido em um intervalo de vetores de limites inferior lb e

superior ub , conforme Equação (B.13).

1 2, ,...,

T

nb b b b=

, 1 2, ,...,T

l l l l

nb b b b=, 1 2, ,...,

Tu u u l

nb b b b=

ou

,l ub b b b ∈ =

(B.13)

Logo, a solução intervalar pode ser escrita como:

l ux x x≤ ≤

(B.14)

na qual x é a solução exata da solução do sistema linear da Equação (B.1) e o vetores

dos limites do intervalo, inferior e superior, são representados, respectivamente, por

lx e u

x .

,l ux x x x ∈ =

(B.15)

Logo, cada elemento do vetor intervalar é dado por:

152

,l u

i i i ix x x x ∈ =

(B.16)

Um sistema linear formado por matrizes intervalares pode ser representado por

um problema de valores extremos, conforme Equação (B.17).

[ ] [ ]

( )[ ]( )

( )

ext,

detextremum , 1, 2,...,

det

i

iA A b b

Ax i p

A ∈ ∈

= =

(B.17)

Os valores extremos (máximo e mínimo) da Equação (B.17), de um vetor

intervalar, pode ser obtido por meio de um problema de otimização global, conforme

Equações (B.18) e (B.19).

[ ] [ ]

[ ] 1

,min , 1,2,...,l

iA A b b

x A b i p−

∈ ∈

= =

(B.18)

e

[ ] [ ]

[ ] 1

,

max , 1,2,...,u

iA A b b

x A b i p−

∈ ∈

= =

(B.19)

na qual [ ]A é uma matriz positiva definida.

De acordo com a Regra de Cramer, as Equações (B.18) e (B.19) podem ser

expressas por:

[ ] [ ]

( )[ ]( )

( )

,

detmin , 1, 2,...,

det

i

l

A A b b

Ax i p

A ∈ ∈

= =

(B.20)

e

[ ] [ ]

( )[ ]( )

( )

,

detmax , 1, 2,...,

det

i

u

A A b b

Ax i p

A ∈ ∈

= =

(B.21)

153

ANEXO C

O anexo C apresenta as características do paquímetro, da balança de precisão e

do acelerômetro piezoelétrico utilizados na obtenção das incertezas de medição dos

dados experimentais, conforme Capítulo 3.

ANEXO C.1 PAQUÍMETRO

Informações de catálogo: PAQUÍMETRO UNIVERSAL

Fabricante: Mitutoyo Processo antidesgaste.

Fabricados em aço temperado de alta resistência, com guias revestidas com

camada de titânio.

Paquímetro - capacidade: 150mm

Código: 530-312B-10

Graduação: 0,02 mm

Exatidão: ±0,03mm

Paquímetro - capacidade: 300 mm

Código: 530-115

Graduação: 0,05 mm

Exatidão: ±0,08mm

154

ANEXO C.2 BALANÇA DE PRECISÃO

Certificado de calibração:

Fabricante: Marte Balanças e Aparelhos de Precisão Ltda Objeto de calibração: material - aço inoxidável Massa nominal: 200 g Erro: 0,33 mg

Incerteza adotada: 0,3mg±

155

ANEXO C.3 ACELERÔMETRO PIEZOELÉTRICO

Certificado de calibração:

Modelo: 352C22 Número de Série: 94793

Descrição: ICP® Accelerometer

Fabricante: PCB Sensibilidade: 9,45 mV/g

PAQUÍMETROS UNIVERSAIS

ANALÓGICOS ! "# $ % & ' ( $ ) * + , - % * - % * . - $ / , ' 0 % - $ % % 0 & 1 0 ) 0 % 2 ' 3 4 $ ) ' 5 ) * *) $ % / 0 % 1 $ ) 0 / - 0 6 0 7 8 * 9: ; < = > ? : @ A @ B = < @ < C D E @ < F @ G H ? I J @ K = < H ?

L # $ % & ' ( $ ) * + , - % * - % * . - $ / , ' 0 % - $ % % 0 & 1 0 ) 0 % 2 ' 3 4 $ ) ' 5 ) * *) $ % / 0 % 1 $ ) 0 / - 0 6 0 7 8 * 9: ; < = > ? : @ A @ B = < @ < C D E @ < F @ G H ? I J @ K = < H ?M N N O O P Q : ; < = > ? : @ A @ B = < @ < C D E @ < F @ G H ? I J @ K = < H ?R S T U V W XN Y N M O O M O OM N N O O N Y N M O O M O O Z N Y N M O OM N N O O [ Q \ 0 1 - ' 1 * 9

] S U X ^ _ V W `

MA

DE

IN

BR

AZ

ILM

AD

E I

N B

RA

ZIL

a - 4 formas de acesso à medição.

Exterior Interior Em ressaltos Em profundidade

Calibration Certificate ^ Per ISO 1 6063-21

Model Number:

Serial Number:

Description:

Manufacturer:

352C22

94793

ICP® Accelerometer Method: Back-to-Back Comparison (AT401-3)

PCB

Sensitivity @ 100.0 Hz

Discharge Time Constant

9.45

(0.963 mV/m/s^)

3.2 seconds

Calibration Data

mV/g Output Bias

Transverse Sensitivity

Resonant Frequency

9.9 VDC

1.7 %

92.0 kHz

Sensitivity Plot

dB

3,0-1

2.0

1.0

0.0

-1.0

-2.0

Temperature: 71 °F (22 °C) Relative Humidity: 4 1 %

- - i - - -

Hz 10,0 100,0

Data Points

1000.0 10000,0

Frequency (Hz) Dev. (%) Frequency (Hz) Dev. (%) Frequency (Hz) Dev. (%)

10.0 -0.2 300.0 0.1 7000.0 1.6

15.0 -0.3 500.0 0.2 10000.0 2.8

30.0 -0.2 1000.0 0.2

50.0 -0.2 3000.0 0.5

REF. FREQ. 0.0 5000.0 1.1

Mounting Siuface: Tungsten Adapter Fastener: Cyanoaci)'late Adiiesive Fixture Orientation: Vertical Acceleration Level (mis)': 10.0 g (98.1 m/s')'

'The acceleration level may be limited by shaker displacement at low frequencies. If the listed level cannot be obtained, the calibration system uses the following formula to set the vibration amplitude; Acceleration Level (g) = 0.010 ,\ (freq)'. 'The gravitational constant used for calculations by the calibration system is; 1 g = 9.80665 m/s'.

Condition of Unit As Found: n/a

As Left: New Unit, In Tolerance

Notes Calibration is NIST Traceable thru Project 822/274086 and PTB Traceable thru Project 1060.

This certificate shall not be reproduced, except in f u l l , without written approval f rom PCB Piezotronics, Inc.

Calibration is performed in compliance with ISO 9001, ISO 10012-1, ANSI/NCSL Z540-1-1994 and ISO 17025.

See Manufacturer's Specification Sheet for a detailed listing of performance specifications.

Measurement uncertainty (95% confidence level with coverage factor of 2) for frequency ranges tested during calibration

1.

2.

3.

4.

5.

are as follows: 5-9 Hz; +/- 2.0%, 10-99 Hz; +/- 1.5%, 100-1999 Hz; +/- 1.0%, 2-10 kHz; +/- 2.5%.

Technician: Susan Lyon Date: 06/12/08

IACCREDlTEPl

CALIBRATION C E R T #1862.02

PAGE 1 of I

®PCB PIEZOTRONICS' VIBRATION DIVISION

Headquarters: 3425 Walden Avenue, Depew, NY 14043 Calibration Performed at: 10869 Highway 903, Halifax, NC 27839

TEL: 888-684-0013 • FAX: 716-685-3886 • www.pcb.com cal48-3296137930.99

III IIIII illllll III IIIII llli IIII IIIII IIIII 11 IIII IIIII

Documents

Análise Dinâmica de Estruturas com Parâmetros Intervalares