Esta é uma pré-visualização de arquivo. Entre para ver o arquivo original
UNIVERSIDADE FEDERAL DE UBERLÂNDIA
FACULDADE DE ENGENHARIA MECÂNICA
Programa de Pós-Graduação em Engenharia Mecânica
INCT de Estruturas Inteligentes em Engenharia
Laboratório de Mecânica de Estruturas Prof. José Eduardo Tannús Reis
Prof. Domingos Alves Rade
2011
MÉTODO DOS ELEMENTOS FINITOS
APLICADOS À ENGENHARIA
MECÂNICA
Introdução ao Método dos Elementos Finitos D.A. Rade
1
CAPÍTULO 1
INTRODUÇÃO AO MÉTODO DOS ELEMENTOS FINITOS
1.1 - Fundamentos do método dos elementos finitos
O método dos elementos finitos (MEF) é uma técnica de análise numérica
destinada à obtenção de soluções aproximadas de problemas regidos por equações
diferenciais. Embora o método tenha sido originalmente desenvolvido para a análise
estática de sistemas estruturais, ele tem sido utilizado no estudo de uma grande
variedade de problemas de Engenharia, nos domínios da Mecânica dos Sólidos,
Mecânica dos Fluidos, Transmissão de Calor e Eletromagnetismo. Devido à sua
eficiência e flexibilidade, além de sua adequação à implementação em computadores
digitais, o MEF tem hoje uma grande difusão tanto no meio acadêmico como no
industrial, estando disponível em grande número de “pacotes” comerciais existentes
no mercado (ANSYS®, NASTRAN®, ABAQUS®, SYSTUS®, COMSOL®, etc.).
Contudo, deve ser lembrado que a utilização eficaz destes programas e a correta
interpretação dos resultados requerem o amplo conhecimento, por parte do
Engenheiro, dos fundamentos do MEF.
A principal motivação para o uso do MEF reside no fato que, devido à
complexidade dos problemas práticos de Engenharia, soluções analíticas em forma
fechada tornam-se inviáveis ou mesmo impossíveis. Assim, devemos recorrer a
técnicas capazes de fornecer soluções numéricas aproximadas. A título de exemplo,
consideremos os problemas de determinação da capacidade de carga de uma placa
contendo enrijecedores e entalhes de formas complexas, ou de determinação da
concentração de poluentes sob condições atmosféricas não uniformes, ou ainda de
caracterização do perfil de velocidades em torno de um aerofólio. Para cada um
destes problemas podemos obter, sem grande esforço, as equações governantes e as
condições de contorno, utilizando princípios elementares da Física. Contudo,
nenhuma solução analítica simples poderá ser obtida quando o problema exibir
geometria e/ou condições de contorno complicadas, o que quase sempre ocorre em
situações práticas. Para contornar esta dificuldade, uma estratégia possível é a
simplificação do problema (em termos de sua geometria e/ou condições de contorno)
de modo a viabilizar a construção de um modelo matemático cuja resolução analítica
seja possível. Contudo, em grande número de casos (talvez na maioria das vezes),
este procedimento tem como conseqüência graves imprecisões nas previsões do
modelo. Uma segunda alternativa consiste em preservar a complexidade do modelo
e empregar técnicas aproximadas de resolução. Esta segunda estratégia, na qual
está inserido o MEF, tem sido cada vez mais viabilizada pela crescente capacidade
de processamento dos computadores digitais.
Em todo problema formulado em domínios contínuos, as incógnitas do
problema, denominadas variáveis de campo (que podem ser grandezas escalares,
como temperaturas ou vetoriais, como deslocamentos) podem assumir valores
independentes em cada ponto do domínio. Conseqüentemente, o problema tem
número infinito de incógnitas, sendo caracterizado como um problema
infinito-dimensional. Este tipo de problema é geralmente modelado por equações
diferenciais parciais, cuja solução analítica é dada por funções que fornecem os
Introdução ao Método dos Elementos Finitos D.A. Rade
2
valores das variáveis de campo em função das coordenadas espaciais para todos os
pontos do domínio.
O MEF é essencialmente um processo de discretização, que visa transformar
um problema infinito-dimensional em um problema finito-dimensional, com número
finito de incógnitas. O método consiste em dividir o domínio sobre o qual o problema
é estudado em várias regiões interconectadas, denominadas elementos. Cada
elemento dispõe de um certo número de pontos (interiores e/ou limítrofes),
denominados nós ou pontos nodais. O conjunto de elementos utilizados na
discretização é denominado malha. Um exemplo é apresentado na Figura 1.1, que
mostra a seção transversal de uma palheta de turbina de geometria complexa,
discretizada em elementos de forma triangular, tendo, em cada vértice, um nó.
Neste exemplo, o problema em questão poderia ser a determinação da distribuição
de temperaturas sobre a seção da palheta, conhecidos o fluxo de calor e as condições
de contorno.
Uma vez definidos os elementos e seus respectivos nós, no interior de cada
elemento são admitidas soluções aproximadas para as variáveis de campo,
expressas como funções arbitrárias dos valores que as incógnitas assumem nos nós
(valores nodais). Estas funções são denominadas funções de interpolação ou funções
de forma. São também impostas condições garantindo a continuidade da solução no
nós compartilhados por vários elementos. As incógnitas do problema, denominadas
graus de liberdade (g.d.l.), passam a ser os valores das variáveis de campo nos
pontos nodais, sendo o número destas incógnitas (agora finito), denominado número
de graus de liberdade do modelo. Dependendo da natureza do problema, após a
discretização, o modelo matemático regente resulta representado por um número
finito de equações diferenciais ordinárias ou de equações algébricas, cuja resolução
numérica conduz aos valores das incógnitas nodais. Uma vez determinadas estas
incógnitas, os valores das variáveis de campo no interior dos elementos podem ser
avaliados empregando as funções de interpolação.
Conforme será visto mais adiante, a precisão da solução obtida depende
essencialmente do número de elementos e do tipo de funções de forma empregadas
na discretização. Sendo satisfeitas algumas condições, admite-se que a solução do
problema discretizado convirja para a solução exata do problema contínuo à medida
que se aumenta o número de incógnitas nodais.
Figura 1.1 – Ilustração da malha de um modelo de elementos finitos
Introdução ao Método dos Elementos Finitos D.A. Rade
3
Em comparação com outras técnicas numéricas, as principais vantagens o
método dos elementos finitos são as seguintes:
• elementos de diferentes formas e tamanhos podem ser associados para
discretizar domínios de geometria complexa.
• a divisão do contínuo em regiões facilita a modelagem de problemas
envolvendo domínios não homogêneos, onde as propriedades físicas variam em
função das coordenadas espaciais.
• o método pode ser todo formulado matricialmente, facilitando sua
implementação computacional.
A implementação do MEF pode sempre ser efetuada em etapas sucessivas, de
forma estruturada. As principais etapas são as seguintes:
1ª) Discretização do domínio. O primeiro passo é a divisão do domínio em
elementos. O tipo e número de elementos a serem utilizados devem ser escolhidos de
modo a representar adequadamente a geometria do problema e caracterizar
convenientemente as variações da solução ao longo do domínio. Alguns tipos de
elementos freqüentemente empregados para a discretização de domínios
unidimensionais, bidimensionais e tridimensionais são ilustrados na Figura 1.2.
Neste aspecto, deve-se observar que problemas unidimensionais são aqueles
definidos em domínios representados por apenas uma coordenada espacial (linhas),
ao passo que problemas bidimensionais e tridimensionais são aqueles definidos em
domínios representados por duas coordenadas espaciais (superfícies) e três
coordenadas espaciais (volumes), respectivamente. Os elementos axissimétricos,
mostrados na Figura 1.2, são elementos utilizados para a discretização de
problemas tridimensionais caracterizados pela existência de simetria geométrica e
de carregamento em relação a um dado eixo. Neste caso, o problema tridimensional
pode ser formulado como um problema bidimensional.
(a)
elementos unidimensionais
(b)
elementos bidimensionais
Introdução ao Método dos Elementos Finitos D.A. Rade
4
(c)
elementos tridimensionais
(d)
elemento axissimétrico
Figura 1.2 – Ilustração de diferentes tipos de elementos
2ª) Escolha das funções de interpolação. Nesta etapa são escolhidas as funções
de interpolação que representam as variáveis de campo no interior de cada
elemento. Freqüentemente, mas nem sempre, funções polinomiais são escolhidas
como funções de interpolação, devido à facilidade que oferecem para derivação e
integração. Os graus dos polinômios utilizados estão relacionados ao número de
incógnitas nodais de cada elemento, devendo também atender a certos requisitos de
continuidade das variáveis de campo a serem satisfeitos nos nós e nas fronteiras
entre elementos imediatamente vizinhos.
3ª) Construção das matrizes elementares. Uma vez escolhidos o tipo e número
de elementos e as funções de interpolação, devemos estabelecer as relações
matriciais expressando o comportamento (relações de causa-efeito), em termos de
propriedades físicas e geométricas, para cada elemento, individualmente. Em outras
palavras, procede-se à formulação em nível elementar. Para tanto, podem ser
utilizados os seguintes processos:
• processo direto, que é baseado no método da rigidez da análise estrutural, através
do qual são obtidas as relações matriciais entre forças e deslocamentos nodais, a
partir das relações de equilíbrio de forças e compatibilidade de deslocamentos.
Procedimento similar pode ser utilizado na modelagem de problemas
unidimensionais de transmissão de calor. Embora o processo direto só seja
conveniente no tratamento de problemas mais simples, ele tem a grande vantagem
de ser de fácil entendimento, permitindo clara interpretação física do significado das
relações matriciais obtidas, interpretação esta que pode ser estendida a problemas
mais complexos.
• processo variacional, que é baseado no Cálculo Variacional e envolve a busca dos
pontos críticos – geralmente pontos de mínimo – de um funcional associado ao
problema estudado. De acordo com este processo, as relações matriciais em nível
Introdução ao Método dos Elementos Finitos D.A. Rade
5
elementar resultam da imposição da condição de estacionaridade do funcional
associado ao problema. Em problemas de Mecânica dos Sólidos, por exemplo, este
funcional pode representar a energia de deformação ou a energia potencial
complementar. O processo variacional, embora mais complicado que o processo
direto sob o ponto de vista teórico, permite estender o MEF a problemas mais
complexos. Contudo, sua aplicabilidade é limitada aos problemas regidos por
princípios variacionais, que estabelecem a existência de funcionais.
• processo dos resíduos ponderados. Este é um procedimento ainda mais versátil
que os dois procedimentos anteriores, sendo baseado integralmente em operações
matemáticas. O processo opera diretamente sobre as equações diferenciais que
governam o problema e prescinde da existência de um funcional ou de um princípio
variacional.
4ª) Montagem das matrizes elementares para obtenção das matrizes
globais. Para caracterizar o comportamento do sistema completo, resultante da
associação dos vários elementos, devemos agrupar as matrizes de cada um dos
elementos de uma forma adequada. Em outras palavras, devemos combinar as
equações matriciais expressando o comportamento dos elementos individuais para
formar as equações matriciais que descrevem o comportamento do sistema em todo
o domínio. Este processo é conhecido como montagem das matrizes globais. No
processo de montagem, impõe-se a condição que em cada nó onde vários elementos
estão interconectados, os valores das variáveis de campo são os mesmos para cada
elemento compartilhando aquele nó.
No final deste processo, as equações matriciais globais devem ser modificadas para
satisfazer as condições de contorno do problema. A ordem das matrizes globais
coincide com o número total de incógnitas nodais. Este número é chamado número
de graus de liberdade do modelo.
5ª) Imposição dos carregamentos externos e das condições de contorno. As
equações matriciais globais devem ser modificadas para satisfazer as condições de
contorno do problema, que expressam o fato que alguns valores das incógnitas
nodais são prescritos. Assim, por exemplo, em problemas de transferência de calor,
os valores da temperatura em alguns pontos do contorno podem ser previamente
conhecidos. Da mesma forma, deve-se alterar as equações globais para leva em
conta que, em alguns nós, cargas externas conhecidas (forças, fluxos de calor, etc.)
são aplicadas. Ao final deste processo, o número total de incógnitas nodais
remanescentes define o chamado número de graus de liberdade do modelo.
6ª) Resolução do sistema de equações. Ao final do processo de montagem das
matrizes globais, o modelo matemático do problema estará representado por um
conjunto de equações, que podem ser lineares ou não lineares, algébricas ou
diferenciais, dependendo da natureza do problema enfocado. Estas equações devem
então ser resolvidas numericamente para a determinação dos valores das variáveis
de campo nos pontos nodais. Neste processo de resolução, procedimentos numéricos
apropriados, implementados sob a forma de rotinas computacionais, devem ser
utilizados.
Introdução ao Método dos Elementos Finitos D.A. Rade
6
7ª) Realização de cálculos complementares. Em várias situações, cálculos
complementares devem ser realizados para a determinação de grandezas
dependentes das variáveis de campo, determinadas na etapa precedente. Assim, por
exemplo, nos problemas de Mecânica dos Sólidos, uma vez determinados os
deslocamentos, cálculos adicionais são necessários para a determinação das
deformações (utilizando as relações deformação-deslocamento) e das tensões
(utilizando as relações tensão-deformação).
1.2 - Domínios de aplicação do MEF
Os problemas passíveis de tratamento pelo MEF podem ser divididos em três
categorias:
• problemas de equilíbrio. Esta é a classe de problemas cuja solução é independente
do tempo. São exemplos os problemas da Mecânica dos Sólidos envolvendo a
determinação de tensões e deformações em elementos estruturais submetidos a
carregamentos estáticos e os problemas da Mecânica dos Fluidos tratando da
determinação de distribuições de pressão, velocidade em regime permanente e os
problemas de Transferência de Calor em regime permanente. Para este tipo de
problema, o processo de discretização através do MEF conduz a um modelo
matemático representado por um conjunto de equações algébricas, que podem ser
lineares ou não lineares.
• problemas de autovalor. Nesta classe de problemas, o modelo matemático obtido é
representado por um conjunto de equações lineares homogêneas, caracterizado pela
dependência em relação a um parâmetro, cuja resolução conduz a um conjunto de
autovalores e autovetores. São exemplos os problemas que tratam da determinação
de freqüências naturais
e modos de vibração de meios sólidos e fluidos, além de
cargas de flambagem de elementos estruturais. No primeiro caso os autovalores
correspondem às freqüências naturais e os autovetores associam-se aos modos
naturais de vibração; no segundo, os autovalores correspondem às cargas de
flambagem e os autovetores dizem respeito aos campos de deslocamentos
correspondentes.
• problemas de propagação. Os problemas de propagação são aqueles em que se
busca caracterizar a evolução das variáveis de campo em função do tempo. É o caso
típico de fenômenos que se desenvolvem em regime transitório. Os seguintes
exemplos podem ser mencionados: determinação do movimento de sistemas
estruturais submetidos a cargas de impacto e determinação de distribuições de
temperatura geradas por fluxos de calor variáveis.
Alguns exemplos de aplicações práticas do MEF são ilustrados a seguir.
Introdução ao Método dos Elementos Finitos D.A. Rade
7
Figura 1.3 - Modelo de EF para análise aeroelástica de um Airbus A-320 (100.000 g.d.l.)
(extraído de J.F. Imbert, "Analyse des Structures par Eléments Finis)
Introdução ao Método dos Elementos Finitos D.A. Rade
8
(b)
Figura 1.4 - Análise térmica por EF de uma tubulação em forma de estrela
(de Desai C.S. e Abel J.F., "Introduction to the Finite Element Method")
Introdução ao Método dos Elementos Finitos D.A. Rade
9
(c)
Figura 1.5 - Análise por EF do escoamento de um fluido ideal em torno de um cilindro
posicionado entre duas placas paralelas
(de Desai C.S. e Abel J.F., "Introduction to the Finite Element Method")
Introdução ao Método dos Elementos Finitos D.A. Rade
10
Figura 1.6 - Análise por EF do resfriamento de um lingote metálico
(fonte: http://www.comsol.com (Comsol Group))
Figura 1.7 - Análise dinâmico de freios automotivos pelo MEF
(fonte: http://www.simulia.com (Dassault Systems))
Introdução ao Método dos Elementos Finitos D.A. Rade
11
Figura 1.8 – Simulação numérica do processo de soldagem pelo MEF
(fonte: http://www.simulia.com (Dassault Systems))
Figura 1.9 – Simulação do comportamento da bexiga humana pelo MEF
(fonte: http://www.comsol.com (Comsol Group))
Introdução ao Método dos Elementos Finitos D.A. Rade
12
Figura 1.10 – Caracterização do fluxo magnético em um gerador pelo MEF
(fonte: http://www.comsol.com (Comsol Group))
Figura 1.11 – Simulação de impacto em palhetas de uma turbina pelo MEF
(fonte: http://www.ansys.com (Ansys, Inc.))
Introdução ao Método dos Elementos Finitos D.A. Rade
13
Figura 1.12 – Simulação de “crash-test” pelo MEF
(fonte: http://www.ansys.com (Ansys, Inc.))
1.3 - Sobre as limitações do método dos elementos finitos
Embora o MEF seja reconhecidamente uma ferramenta útil e eficiente para a
análise de diversos tipos de problemas de Engenharia, é importante estar atento às
limitações do método, lembrando que ele fornece modelos matemáticos aproximados
para representar o comportamento de sistemas físicos. Conforme será evidenciado
no desenvolvimento do curso, em geral, a elaboração de um modelo de EF envolve a
admissão de uma série de simplificações que terão efeito direto sobre a precisão das
previsões do modelo. Deve também ser levado em conta que, na maioria das vezes,
sob considerações de custo e limitações de recursos computacionais, busca-se
estabelecer um compromisso entre a complexidade do modelo e a precisão
considerada satisfatória.
Algumas fontes de incerteza inerentes à modelagem por EF são:
• a não consideração de certos tipos de efeitos físicos, tais como não
linearidades, histerese, amortecimento, etc.
• erros de discretização, devidos à impossibilidade de se obter uma perfeita
representação de domínios de geometria complexa utilizando os tipos de elementos
disponíveis.
• conhecimento impreciso dos valores de alguns parâmetros físicos e/ou
geométricos que são utilizados na elaboração do modelo (ex.: módulo de elasticidade,
densidade, condutividade térmica, viscosidade, etc.)
• dificuldade de modelar efeitos localizados, tais como junções parafusadas e
rebitadas.
• erros oriundos do processo de resolução numérica.
Introdução ao Método dos Elementos Finitos D.A. Rade
14
Visando a validação de modelos de elementos finitos, um procedimento
recomendável é a confrontação de suas previsões com dados provenientes de outros
tipos de análises, em particular, de análises experimentais. Várias técnicas foram
desenvolvidas, principalmente no âmbito da Dinâmica Estrutural, objetivando a
correção sistemática de modelos de elementos finitos a partir da confrontação com
dados experimentais (consultar, por exemplo, o livro de Friswell e Mottershead).
Outro fator a ser considerado é que a elaboração de modelos de EF de
problemas complexos é, na maioria das vezes, um processo interativo, fazendo apelo
ao conhecimento do Engenheiro acerca do próprio método e do problema em estudo.
1.4 - Bibliografia
• Bathe, K. J., Finite Element Procedures in Engineering Analysis, Prentice-Hall,
1982.
• Cook, R.D., Concepts and Applications of Finite Element Analysis, John Wiley &
Sons, 1974.
• Cook, R., Finite Element Modeling for Analysis, John Wiley & Sons, 1995.
• Desai, C.S., Abel, J. F., Introduction to the Finite Element Method, Van Nostrand
Reinhold Company, 1972.
• Imbert, J. F., Analyse des Structures par Eléments Finis, 3ième édition, Cépauès
Editions, 1991.
• Friswell, M.I., Mottershead, J.E., Finite Element Model Updating in Structural
Dynamics, Kluwer Academic Publishers, 1995.
• Moaveni, S., Finite Element Analysis. Theory and Application with ANSYS,
Prentice-Hall, 1999.
•Zienkiewicz, O.C., The Finite Element Method, 3rd edition, McGraw-Hill Book Co.,
1977.
•Zienkiewicz, O.C., Morgan, K, Finite Elements and Approximation, John Wiley &
Sons, 1983
•Huebner, K.H., Thornton, E.A., The Finite Element Method for Engineers, John
Wiley & Sons, 1982
D.A.Rade Formulação do MEF pelo Processo Direto
15
CAPÍTULO 2
FORMULAÇÃO DO MÉTODO DOS ELEMENTOS FINITOS
PELO PROCESSO DIRETO
No Capítulo 1 havíamos mencionado que as relações matriciais traduzindo o
comportamento de cada subdomínio (elemento) podiam ser obtidas empregando um
dos três métodos: Método Direto, Método Variacional e Método dos Resíduos
Ponderados.
Neste capítulo vamos ilustrar, com o auxílio de exemplos simples, a utilização
do Processo Direto, que tem a vantagem de proporcionar facilitar a interpretação
física para as diversas etapas da elaboração de um modelo de elementos finitos.
Embora sua utilização seja viável apenas no tratamento de problemas
unidimensionais simples, os conceitos derivados deste método podem ser estendidos,
com vantagem, a problemas mais complexos.
2.1 – Análise estática de sistemas compostos por associações de molas
lineares.
Um dos problemas mais elementares que podem ser examinados sob o
enfoque do MEF é a análise estática do sistema constituído por um conjunto de
molas lineares, com constantes de rigidez ik , i=1, 2, ..., como aquele exemplificado
na Figura 2.1. Nesta figura, if e iu designam, respectivamente, as forças externas
aplicadas e os deslocamentos dos nós. O problema consiste em determinar os
deslocamentos nodais, dados os valores das forças aplicadas.
Figura 2.1
1k 2k 3k
nó 1 nó 2 nó 3 nó 4
2f 3f 4f 1f
1u 2u 3u 4u
D.A.Rade Formulação do MEF pelo Processo Direto
16
2.1.1 – Obtenção das equações de equilibro em nível elementar
Cada mola é identificada com um elemento finito. Para obter a matriz de
rigidez para um elemento genérico, este elemento é considerado isoladamente,
conforme mostrado na Figura 2.2, onde os índices E e D designam as grandezas
associadas aos nós das extremidades esquerda e direita, respectivamente.
Figura 2.2
Admitindo que todo o conjunto e, portando, cada elemento, encontra-se em
equilíbrio, podemos escrever:
EiDi ff −=
Sabemos ainda que para as molas lineares, o alongamento é proporcional à
força aplicada em suas extremidades, sendo a constante de proporcionalidade o
coeficiente de rigidez. Assim, escrevemos:
( )EiDiiEi uukf −−=
( )EiDiiDi uukf −=
As duas equações acima podem ser postas na seguinte forma matricial:
⎪⎭
⎪⎬⎫⎪⎩
⎪⎨⎧=⎪⎭
⎪⎬⎫⎪⎩
⎪⎨⎧⎥⎦
⎤⎢⎣
⎡
−
−
D
i
E
i
D
i
E
i
ii
ii
f
f
u
u
kk
kk (2.1)
A equação acima pode ainda ser escrita na forma compacta:
E
iu Diu
ik
E
if
D
if
D.A.Rade Formulação do MEF pelo Processo Direto
17
( )[ ] ( ){ } ( ){ }eieiei FUK = i=1,2,... (2.2)
onde ( ){ } [ ] ( ){ } [ ]TDiEieiTDiEiei ffFeuuU == são, respectivamente, os vetores de
deslocamentos nodais e forças nodais em nível elementar e ( )[ ]eiK é denominada
matriz de rigidez elementar. Sobre esta matriz, podemos observar:
• o elemento genérico ( )( )mneiK representa a força de reação aplicada no nó m
quando provocamos um deslocamento unitário no nó n, mantendo bloqueado o nó m.
• é uma matriz simétrica ( )[ ] ( )[ ]Teiei KK = . Levando em conta a interpretação
dada acima, a simetria da matriz de rigidez traduz o Princípio de Reciprocidade de
Maxwell-Betti, aplicável a sistemas mecânicos lineares. Isto significa que a força de
reação que surge no nó m quando provocamos um deslocamento unitário no nó n,
mantendo bloqueado o nó m é idêntica à força de reação que surge no nó n quando
provocamos um deslocamento unitário no nó m, mantendo bloqueado o nó n.
• é uma matriz singular (não inversível). Isto porque, como não foram
introduzidas as condições de contorno em nível elementar, a relação (2.2) deve
contemplar a existência de uma configuração de equilíbrio sem deformação da mola ( )1+= ii uu e sem o aparecimento de forças nos nós ( )01 == +ii ff .
• é uma matriz semi-definida positiva { } ( )[ ] { } { } { }⎟⎠⎞⎜⎝⎛ ≠∀≥ 00 x,xKx eiT
2.1.2 – Obtenção das equações de equilíbrio em nível global
Após a obtenção das equações que descrevem o comportamento de cada
elemento, isoladamente, devemos considerar o fato que os elementos estão, na
realidade, interconectados nos pontos nodais. Fisicamente, a interconexão significa
que deve haver, nos nós compartilhados por mais de uma mola, equilíbrio das forças
e compatibilidade de deslocamentos. Considerando dois elementos vizinhos, i e i+1,
ilustrados na Figura 2.3, estas condições são expressas segundo:
11 ++ =+ iEiDi fff (equilíbrio do nó i+1) (2.3)
11 ++ == iEiDi uuu (compatibilidade de deslocamentos no nó i+1) (2.4)
Para imposição destas condições, vamos primeiramente desenvolver a
equação matricial (2.1) para os elementos i e i+1:
D.A.Rade Formulação do MEF pelo Processo Direto
18
E
i
D
ii
E
ii fukuk =− (2.5)
D
i
D
ii
E
ii fukuk =+− (2.6)
E
i
D
ii
E
ii fukuk 11111 +++++ =− (2.7)
D
i
D
ii
E
ii fukuk 11111 +++++ =+− (2.8)
Somando as equações (2.7) e (2.8) e introduzindo nas três equações
resultantes as equações (2.3) e (2.4), obtemos as relações.
E
iii
E
ii fukuk =− +1 (2.9)
( ) iDiiiiiEii fukukkuk =−++− ++++ 1111 (2.10)
D
i
D
iiii fukuk 1111 ++++ =+− (2.11)
Retornando à notação matricial, as equações (2.9) a (2.11) são dispostas sob a
forma:
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
=
⎪⎭
⎪⎬
⎫
⎪⎩
⎪⎨
⎧
⎥⎥
⎥
⎦
⎤
⎢⎢
⎢
⎣
⎡
−
−+−
−
+
+
+
+
++
++
D
i
i
E
i
D
i
i
E
i
ii
iiii
ii
f
f
f
u
u
u
kk
kkkk
kk
1
1
1
1
11
11
0
0
(2.12)
Figura 2.3
E
iu Diu
ik
E
if
E
iu 1+ Diu 1+
1+ik
E
if 1+ Dif 1+
D
if
nó i nó i+1 nó i+1 nó i+2
D.A.Rade Formulação do MEF pelo Processo Direto
19
Considerando o exemplo ilustrado na Figura 2.1, aplicando sucessivamente o
mesmo procedimento aos dois pares de molas que compartilham os nós 2 e 3,
obtemos a seguinte relação matricial:
⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
=
⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
−
−+−
−+−
−
4
3
2
1
4
3
2
1
33
3322
2211
11
00
0
0
00
f
f
f
f
u
u
u
u
kk
kkkk
kkkk
kk
(2.13)
A equação acima pode ainda ser escrita na forma compacta:
( ) ( ){ } ( ){ }g g gK U F⎡ ⎤ =⎣ ⎦ (2.14)
onde ( ){ } [ ] ( ){ } [ ]1 2 3 4 1 2 3 4eT Tg gU u u u u F f f f f= = são, respectivamente, os
vetores de deslocamentos e forças em nível global e ( )gK⎡ ⎤⎣ ⎦ é denominada matriz de
rigidez global. Sobre esta matriz, podemos fazer as mesmas observações
anteriormente apresentadas para as matrizes de rigidez elementares:
• o elemento genérico ( )g
mn
K⎡ ⎤⎣ ⎦ representa a força de reação aplicada no nó m
quando provocamos um deslocamento unitário no nó n, mantendo bloqueado o nó m.
• é uma matriz simétrica ( ) ( )
Tg gK K⎡ ⎤ ⎡ ⎤=⎣ ⎦ ⎣ ⎦ .
• é uma matriz singular.
• é uma matriz semi-definida positiva { } ( ) { } { } { }( )0, 0T gx K x x⎡ ⎤ ≥ ∀ ≠⎣ ⎦
O processo de obtenção da matriz de rigidez global a partir das matrizes
elementares e denominado montagem da matriz global. Um procedimento mais
geral de montagem, bem adaptado à implementação computacional, pode ser
formulado através da introdução de matrizes de transformação, compreendendo as
seguintes etapas:
1ª) Para cada elemento, estabelecemos relações matriciais expressando as
transformações dos deslocamentos nodais em nível elementar com os deslocamentos
nodais em nível global.. Assim, no exemplo considerado:
D.A.Rade Formulação do MEF pelo Processo Direto
20
[ ]( )
⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
=⎪⎭
⎪⎬⎫⎪⎩
⎪⎨⎧ ×
4
3
2
1
42
u
u
u
u
T
u
u
iD
i
E
i i=1,2,3
2ª) Pré-multiplicamos as matrizes de rigidez elementares por [ ]TiT e a pós-
multiplicamos por [ ]iT , para obter as matrizes elementares expandidas, com a
mesma dimensão da matriz de rigidez global:
( )[ ]( ) [ ]( ) ( )[ ]( ) [ ]( )42222444 ×××× = ieiTiei TKTK i=1,2,3 (2.15)
3ª) Adicionamos as matrizes elementares expandidas para obter finalmente a
matriz de rigidez global. No exemplo:
( ) ( )3
1
g e
i
i
K K
=
⎡ ⎤ ⎡ ⎤=⎣ ⎦ ⎣ ⎦∑ (2.16)
Detalhamos, a seguir, o procedimento de construção da matriz de rigidez
global para o exemplo considerado:
• elemento nº 1:
⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
⎥⎦
⎤⎢⎣
⎡=⎪⎭
⎪⎬⎫⎪⎩
⎪⎨⎧
4
3
2
1
1
1
0010
0001
u
u
u
u
u
u
D
E
( )[ ] [ ] ( )[ ] [ ]
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
−
−
=⎥⎦
⎤⎢⎣
⎡⎥⎦
⎤⎢⎣
⎡
−
−
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
==
0000
0000
00
00
0010
0001
00
00
10
01
11
11
11
11
1111
kk
kk
kk
kk
TKTK eTe
• elemento nº 2:
⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
⎥⎦
⎤⎢⎣
⎡=⎪⎭
⎪⎬⎫⎪⎩
⎪⎨⎧
4
3
2
1
2
2
0100
0010
u
u
u
u
u
u
D
E
D.A.Rade Formulação do MEF pelo Processo Direto
21
( )[ ] [ ] ( )[ ] [ ]
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
−
−=⎥⎦
⎤⎢⎣
⎡⎥⎦
⎤⎢⎣
⎡
−
−
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
==
0000
00
00
0000
0100
0010
00
10
01
00
22
22
22
22
2222 kk
kk
kk
kk
TKTK eTe
• elemento nº 3:
⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
⎥⎦
⎤⎢⎣
⎡=⎪⎭
⎪⎬⎫⎪⎩
⎪⎨⎧
4
3
2
1
3
3
1000
0100
u
u
u
u
u
u
D
E
( )[ ] [ ] ( )[ ] [ ]
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
−
−=⎥⎦
⎤⎢⎣
⎡⎥⎦
⎤⎢⎣
⎡
−
−
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
==
33
3333
33
3333
00
00
0000
0000
1000
0100
10
01
00
00
kk
kkkk
kk
TKTK eTe
( ) ( ) ( ) ( )
1 2 3
g e e eK K K K⎡ ⎤ ⎡ ⎤ ⎡ ⎤ ⎡ ⎤= + +⎣ ⎦ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
−
−+−
−+−
−
=
33
3322
2211
11
00
0
0
00
kk
kkkk
kkkk
kk
(2.17)
2.1.3 – Imposição das condições de contorno
As equações (2.13) devem ser modificadas para levar em conta que os valores
dos deslocamentos nos nós extremos (1 e 4) podem ser prescritos. Admitamos que as
condições de contorno sejam as seguintes:
11 uu = (2.18.a)
44 uu = (2.18.b)
onde 1u e 4u são os valores conhecidos.
Desenvolvendo (2.13) considerando as relações (2.18), obtemos o seguinte
conjunto de equações:
D.A.Rade Formulação do MEF pelo Processo Direto
22
( )
( )
( )
( )
⎪⎪
⎪⎪
⎩
⎪⎪
⎪⎪
⎨
⎧
=+−
+=++−
−=−+
=−
⇒
⎪⎪
⎪⎪
⎩
⎪⎪
⎪⎪
⎨
⎧
=+−
=−++−
=−++
=−
44333
43333222
11232221
12111
44333
34333222
23222111
12111
fukuk
ukfukkuk
ukfukukk
fukuk
fukuk
fukukkuk
fukukkuk
fukuk
(2.19)
Determinamos os deslocamentos nodais incógnitos 2u e 3u resolvendo
simultaneamente a segunda e a terceira equações do sistema (2.19):
⎭⎬
⎫
⎩⎨
⎧
−
−
⎥⎦
⎤⎢⎣
⎡
+−
−+=
⎭⎬
⎫
⎩⎨
⎧⇒
⎭⎬
⎫
⎩⎨
⎧
−
−=
⎭⎬
⎫
⎩⎨
⎧⎥⎦
⎤⎢⎣
⎡
+−
−+ −
433
112
1
322
221
3
2
433
112
3
2
322
221
ukf
ukf
kkk
kkk
u
u
ukf
ukf
u
u
kkk
kkk (2.20)
Os valores das forças de reação nos nós 1 e 4 são finalmente determinados
através da primeira e quarta equações do sistema (2.19):
⎪⎩
⎪⎨
⎧
+−=
−=
43334
21111
ukukf
ukukf
(2.21)
É importante observar que, após a imposição das condições de contorno, a
matriz de rigidez modificada, deve ser inversível para que a operação indicada em
(2.20) possa ser efetuada. Este será sempre o caso quando as condições de contorno
impostas forem suficientes para impedir que o sistema mecânico seja
cinematicamente variável ou, em outras palavras, que possa se movimentar sem
que haja deformação de pelo menos uma das molas.
Um procedimento sistemático para imposição das condições de contorno e
cálculo das forças de reação, bem adaptado à implementação computacional,
consiste em particionar os graus de liberdade em dois conjuntos: graus de liberdade
livres e graus de liberdade impostos. No nosso exemplo, estes dois conjuntos seriam:
• graus de liberdade livres: ( ){ } [ ]2 3 TgU u u=A ; ( ){ } [ ]2 3 TgF f f=A
• graus de liberdade impostos: ( ){ } [ ]1 4 TgU u u=A ; ( ){ } [ ]1 4 TgF f f=A
Após reordenação das equações, o sistema (2.13) pode ser expresso sob a
forma:
( ) ( )
( ) ( )
( ){ }
( ){ }
( ){ }
( ){ }
g gg g
i
g gg g
i ii ii
U FK K
U FK K
⎡ ⎤ ⎧ ⎫ ⎧ ⎫⎡ ⎤ ⎡ ⎤⎣ ⎦ ⎣ ⎦ ⎪ ⎪ ⎪ ⎪⎢ ⎥ =⎨ ⎬ ⎨ ⎬⎢ ⎥⎡ ⎤ ⎡ ⎤ ⎪ ⎪ ⎪ ⎪⎢ ⎥ ⎩ ⎭ ⎩ ⎭⎣ ⎦ ⎣ ⎦⎣ ⎦
A AAA A
A
D.A.Rade Formulação do MEF pelo Processo Direto
23
Do sistema acima tiramos duas equações matriciais:
( ) ( ){ } ( ) ( ){ } ( ){ }g g g g gi iK U K U F⎡ ⎤ ⎡ ⎤+ =⎣ ⎦ ⎣ ⎦AA A A A (2.22)
( ) ( ){ } ( ) ( ){ } ( ){ }g g g g gi ii i iK U K U F⎡ ⎤ ⎡ ⎤+ =⎣ ⎦ ⎣ ⎦A A (2.23)
De (2.22) obtemos os valores dos deslocamentos nodais correspondentes aos
graus de liberdade livres e, em seguida, a partir de (2.23) calculamos as forças de
reação correspondentes aos graus de liberdade impostos:
( ){ } ( ) ( ){ } ( ) ( ){ }( )1g g g g gi iU K F K U−⎡ ⎤ ⎡ ⎤= −⎣ ⎦ ⎣ ⎦A AA A A (2.24)
( ){ } ( ) ( ){ } ( ) ( ){ }g g g g gi i ii iF K U K U⎡ ⎤ ⎡ ⎤= +⎣ ⎦ ⎣ ⎦A A (2.25)
2.2 – Exemplo resolvido utilizando o MATLAB®
Nesta seção apresentamos a resolução do problema mostrado na Figura 2.1,
particularizado para a seguinte situação:
41 20 10k = × N/m 42 30 10k = × N/m 43 50 10k = × N/m
10002 =f N 1 3 0f f= = N
041 == uu
A Figura 2.4 mostra o código implementado em ambiente MATLAB® . A
solução obtida para o problema é a seguinte:
32 3,64 10 mu −= × 33 2,73 10 mu −= ×
N37271 ,f −= N72724 ,f −=
D.A.Rade Formulação do MEF pelo Processo Direto
24
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Programa para análise estática por EF de sistemas constituídos por molas lineares
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Domingos Alves Rade
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ENTRADA DE DADOS
% valores dos coeficientes de rigidez das molas
k_values = 1e4*[ 20 30 10 ]; % valores em N/m
% número de nós
nb_node = 4;
% matriz de conectividade
mat_conect=[1 2 ; 2 3 ; 3 4 ]
% condições de contorno nos gdl impostos
cond_cont=[1 0 ; 4 0];
% forças externas aplicadas nos gld livres
forcas_aplic= [ 2 1000 ; 3 0]; % valores em Newtons
nb_ele=length(k_values);
% CONSTRUÇÃO DAS MAT. ELEMENTARES E MONTAGEM DA MATRIZ DE
% RIGIDEZ GLOBAL
K_global=zeros(nb_node);
for ii=1:nb_ele
K_elementar=k_values(ii)*[1 -1 ; -1 1];
mat_ident=eye(nb_node);
mat_transf=[mat_ident(mat_conect(ii,1),:); mat_ident(mat_conect(ii,2),:)];
K_global=K_global+mat_transf'*K_elementar*mat_transf;
end
%IMPOSIÇÃO DAS CONDIÇÕES DE CONTORNO PELO MÉTODO DO
%PARTICIONAMENTO DA MATRIZ DE RIGIDEZ
% identificação dos gdl livres e gdl impostos
gdl_livres=forcas_aplic(:,1);
gdl_impostos = cond_cont(:,1);
% construção das submatrizes de rigidez
K_ll=K_global(gdl_livres,gdl_livres);
K_li=K_global(gdl_livres,gdl_impostos);
K_ii=K_global(gdl_impostos,gdl_impostos);
% construção dos vetores de forças nos gdl livres e de deslocamentos nos gdl
% impostos
f_liv=forcas_aplic(:,2);
d_imp=cond_cont(:,2);
% CÁLCULO DOS DESLOCAMENTOS NOS GDL LIVRES E FORÇAS DE REAÇÃO % NOS
GDL IMPOSTOS
d_liv=inv(K_ll)*(f_liv-K_li*d_imp)
f_imp=K_li*d_liv+K_ii*d_imp
Figura 2.4
D.A.Rade Formulação do MEF pelo Processo Direto
25
2.3 – Análise estática de treliças planas.
A idéia de modelar estruturas como uma série de elementos finitos surgiu
como uma extensão dos métodos matriciais utilizados para análise de treliças. Estas
estruturas consistem de barras interconectadas em certos pontos considerados como
rótulas perfeitas, de modo que apenas forças axiais podem ser transmitidas às
barras, não sendo considerados efeitos de flexão e cisalhamento. Um exemplo de
treliça é ilustrado na Figura 2.5.
Figura 2.5
Na modelagem, cada barra da treliça é tratada com um elemento finito, e os pontos
de conexão das barras são associados aos nós do modelo de EF.
Para determinar as relações forças-deslocamentos para a treliça completa,
devemos, primeiramente, determinar estas relações para um elemento genérico,
considerado isoladamente, levando em conta sua orientação, e proceder, em seguida,
à montagem das equações matriciais que levam em conta a compatibilidade de
deslocamentos e o equilíbrio de forças nos nós.
Considerando a Figura 2.5, que mostra um elemento genérico i, orientado
segundo um ângulo iα , cujas propriedades relevantes são: comprimento iL , área de
seção transversal iA e módulo de elasticidade iE . Vamos utilizar os conceitos
elementares da Resistência dos Materiais para obter as relações forças-
deslocamentos em nível elementar.
Para cada elemento, definimos dois sistemas de referência ilustrados na
Figura 2.6:
- o sistema de referência global X-Y, comum a todos os elementos.
- um sistema de referência local x-y, cuja orientação é determinada pela
orientação de cada barra.
D.A.Rade Formulação do MEF pelo Processo Direto
26
Figura 2.6
Deve ser notado, na Figura 2.6, que cada um dos nós apresenta duas
componentes independentes de deslocamento nas direções X e Y, havendo, portanto,
dois graus de liberdade por nó. Além disso, as forças aplicadas nos nós são
orientadas na direção da barra.
Da condição de equilíbrio do elemento, temos:
( ) ( )E Di ix xF F= − (2.26)
A relação força-deslocamento, derivada da lei de Hooke para o elemento
considerado permite escrever:
( ) ( ) ( )Ei iD E xi i ix x i i
F L
L u u
E A
Δ = − = − (2.27)
( ) ( ) ( )Di iD E xi i ix x i i
F L
L u u
E A
Δ = − = (2.28)
As duas últimas equações acima podem ser dispostas na seguinte forma
matricial.
Y
X
y
x
iL
i
α
E
D
ii A,E
( ) ( ),D Di ix xu F
( ) ( ),E Ei ix xu F
( ) ( )XEiXEi F,u
( ) ( )YEiYEi F,u
( ) ( ),D Di iY Yu F
( ) ( ),D Di iX Xu F
D.A.Rade Formulação do MEF pelo Processo Direto
27
( )
( )
( )
( )
i i i i E E
i ii i x x
D Di i i i i ix x
i i
E A E A
u FL L
E A E A u F
L L
⎡ ⎤− ⎧ ⎫ ⎧ ⎫⎢ ⎥ ⎪ ⎪ ⎪ ⎪⎢ ⎥ =⎨ ⎬ ⎨ ⎬⎢ ⎥ ⎪ ⎪ ⎪ ⎪−⎢ ⎥ ⎩ ⎭ ⎩ ⎭⎣ ⎦
(2.29)
Devemos observar que, na equação acima os deslocamentos e forças nodais
estão expressas em termos de suas componentes nas direções dos eixos locais. Como
cada um dos elementos possui, em geral, orientação diferente da dos demais, antes
da montagem das equações matriciais globais é necessário expressar as relações
forças-deslocamentos em termos das componentes no sistema de referência global.
Isto é feito introduzindo a seguinte transformação para os deslocamentos, que pode
ser facilmente extraída da Figura 2.4:
( ) ( ) ( )
( ) ( ) ( )
( )
( )
( )
( )
( )
( )
cos sen cos sen 0 0
0 0 cos sencos sen
E
i X
E E E E E
i i i i i i ix X Y x Yi i
D D D D Di ii i i i i i ix X Y x X
D
i Y
u
u u u u u
u u u u u
u
α α α α
α αα α
⎧ ⎫⎪ ⎪⎪ ⎪⎧ ⎫= + ⎪ ⎪⎡ ⎤⎪ ⎪ ⎪ ⎪⇒ =⎨ ⎬ ⎨ ⎬⎢ ⎥⎣ ⎦= + ⎪ ⎪ ⎪ ⎪⎩ ⎭ ⎪ ⎪⎪ ⎪⎪ ⎪⎩ ⎭
ou, em forma compacta:
( )
( )
( )
( )
( )
( )
i
E
i X
E E
i ix Y
D D
i ix X
D
i Y
u
u u
T
u u
u
α
⎧ ⎫⎪ ⎪⎪ ⎪⎧ ⎫ ⎪ ⎪⎪ ⎪ ⎪ ⎪⎡ ⎤=⎨ ⎬ ⎨ ⎬⎣ ⎦⎪ ⎪ ⎪ ⎪⎩ ⎭ ⎪ ⎪⎪ ⎪⎪ ⎪⎩ ⎭
(2.30)
onde:
[ ] ⎥⎦⎤⎢⎣⎡= iiii sencossencosT i ααααα 00 00 (2.31)
Introduzindo (2.30) em (2.29) e pré-multiplicando a equação resultante por [ ]T
i
Tα , obtemos:
D.A.Rade Formulação do MEF pelo Processo Direto
28
[ ]T
i
Tα
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
−
−
i
ii
i
ii
i
ii
i
ii
L
AE
L
AE
L
AE
L
AE
[ ]
i
Tα
( )( )( )( ) ⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
Y
D
i
X
D
i
Y
E
i
X
E
i
u
u
u
u
= [ ]T
i
Tα
( )
( )
E
i x
D
i x
F
F
⎧ ⎫⎪ ⎪⎨ ⎬⎪ ⎪⎩ ⎭
(2.32)
Com base na Figura 2.4, podemos desenvolver o membro do lado direito de
(2.32) sob a forma:
[ ]T
i
Tα
( )
( )
E
i x
D
i x
F
F
⎧ ⎫⎪ ⎪⎨ ⎬⎪ ⎪⎩ ⎭
=
⎥⎥
⎥⎥
⎦
⎤
⎢⎢
⎢⎢
⎣
⎡
i
i
i
i
cos
cos
sen
cos
α
α
α
α
0
0
0
0 ( )
( )
E
i x
D
i x
F
F
⎧ ⎫⎪ ⎪⎨ ⎬⎪ ⎪⎩ ⎭
=
( )( )( )( ) ⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
=
⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
Y
D
i
X
D
i
Y
E
i
X
E
i
i
D
i
i
D
i
i
E
i
i
E
i
F
F
F
F
senF
cosF
senF
cosF
α
α
α
α
(2.33)
Assim, a equação (2.32) assume a seguinte forma, representando a relação
forças-deslocamentos expressa no sistema de referência global:
⎥⎥
⎥⎥
⎥
⎦
⎤
⎢⎢
⎢⎢
⎢
⎣
⎡
−
−−
−−
−−
iiiiii
iiiiii
iiiii
iiiiii
i
ii
sensencossensencos
sencoscossencoscos
sensencossensencos
sencoscossencoscos
L
AE
αααααα
αααααα
αααααα
αααααα
22
22
22
22 ( )( )( )( ) ⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
Y
D
i
X
D
i
Y
E
i
X
E
i
u
u
u
u
=
( )( )( )( ) ⎪⎪⎭
⎪⎪⎬
⎫
⎪⎪⎩
⎪⎪⎨
⎧
Y
D
i
X
D
i
Y
E
i
X
E
i
F
F
F
F
(2.34)
ou:
( )[ ] ( ){ } ( ){ }eieiei FUK = i=1,2,... (2.35)
onde:
( )[ ] =eiK
⎥⎥
⎥⎥
⎥
⎦
⎤
⎢⎢
⎢⎢
⎢
⎣
⎡
−
−−
−−
−−
iiiiii
iiiiii
iiiii
iiiiii
i
ii
sensencossensencos
sencoscossencoscos
sensencossensencos
sencoscossencoscos
L
AE
αααααα
αααααα
αααααα
αααααα
22
22
22
22
(2.36)
é a matriz de rigidez elementar do elemento de treliça.
Devemos observar que cada elemento possui quatro graus de liberdade (2
graus de liberdade por nó, correspondendo aos deslocamentos nas direções X e Y.
Uma vez desenvolvida a expressão geral da matriz de rigidez para um
elemento de barra genérico, podemos utilizar os procedimentos já detalhados na
Seção 2.1 para a montagem da matriz de rigidez global – traduzindo a
compatibilidade
de deslocamentos e equilíbrio de forças nos nós – e a imposição das
D.A.Rade Formulação do MEF pelo Processo Direto
29
condições de contorno. No que segue, estas etapas serão detalhadas a partir de um
exemplo numérico.
2.4 – Exemplo resolvido usando o MATLAB®.
Consideremos a treliça plana constituída de três barras, ilustrada na
Figura 2.7(a), cujas características físicas e geométricas são dadas na Tabela 2.1.
Interessa-nos determinas os deslocamentos nodais, as reações nos apoios e as
tensões em cada uma das barras.
A Figura 2.7(b) mostra a numeração escolhida para os nós e elementos, bem
como a indicação das componentes das forças e deslocamentos nodais no sistema de
referência global adotado.
(a)
(b)
Figura 2.7
45o 100 KN
Y
X
3
2
1
1
2
3 ( )Xu3
( )Yu3
( )Xu2
( )Yu2
( )XF3
( )YF3
( )YF2
( )XF2
( )Xu1
( )Yu1 ( )YF1
( )XF1
D.A.Rade Formulação do MEF pelo Processo Direto
30
Tabela 2.1
Elemento A(cm2) E(N/m2) L(m) Conectividade
(nós)
α (graus)
1
32,3
6,9×1010
2,54
2,3
90
2
38,7
20,7×1010
2,54
1,2
0
3
25,8
20,7×1010
3,59
1,3
135
O programa MATLAB® utilizado e os resultados obtidos são apresentados na
Figura 2.5.
Deve ser observado que, após o cálculo das forças e deslocamentos nodais, as
tensões nas barras foram calculadas de acordo com o seguinte procedimento:
Utilizando a equação (2.28), escrevemos:
( ) ( ) ( )Di D Ex ii i ix xi i
F E u u
A L
σ ⎡ ⎤= = −⎢ ⎥⎣ ⎦
Os valores de ( )Di xu e ( )Ei xu são calculados a partir dos deslocamentos nodais
expresso no sistema de referência global utilizando a equação (2.30).
Os resultados obtidos para este exemplo são os seguintes:
• Deslocamentos
( )1 0,6858mmXu = − ( )2 0,9100 mmXu = − ( )2 0,8059 mmYu = −
• Reações de apoio
( )2 70,711KNYF = ( )3 70,711KNXF = ( )3 0YF =
• Tensões nas barras
1 21,819MPaσ = 2 18,271MPaσ = 3 38,760MPaσ = −
D.A.Rade Formulação do MEF pelo Processo Direto
31
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Programa para análise estática por EF de treliças planas
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Domingos Alves Rade Última modificação:
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% ENTRADA DE DADOS
% propriedades físicas e geométricas das barras [A(m^2) E(N/m^2) L(m) alfa(grau)]
propriedades = [ 32.3e-4 6.9e10 2.54 90
38.7e-4 20.7e10 2.54 180
25.8e-4 20.7e10 2.59 135];
% número total de graus de liberdade
nb_gdl = 6;
% matriz de conectividade (graus de liberdade ordenados segundo:
% 1: nó 1, direção x
% 2: nó 1, direção y
% 3: nó 2, direção x
% 4: nó 2, direção y
% 5: nó 3, direção x
% 6: nó 3, direção y
mat_conect=[ 3 4 5 6
1 2 3 4
1 2 5 6];
% condições de contorno nos gdl impostos
cond_cont=[2 0
5 0
6 0];
% forças externas aplicadas nos gld livres(valores em Newtons)
forcas_aplic= [ 1 0
3 -100000*cos(pi/4)
4 -100000*sin(pi/4)];
[nb_ele,dummy]=size(propriedades);
% CONSTRUÇÃO DAS MATRIZES ELEMENTARES E MONTAGEM DA MATRIZ DE RIGIDEZ GLOBAL
K_global=zeros(nb_gdl);
for ii=1:nb_ele
coef=propriedades(ii,1)*propriedades(ii,2)/propriedades(ii,3);
alfa=propriedades(ii,4)*pi/180;
K_elementar=coef*[1 -1 ; -1 1];
mat_transf=[cos(alfa) sin(alfa) 0 0
0 0 cos(alfa) sin(alfa)];
K_elementar=mat_transf'*K_elementar*mat_transf;
mat_ident=eye(nb_gdl);
mat_transf1=[mat_ident(mat_conect(ii,1),:); mat_ident(mat_conect(ii,2),:); ...
mat_ident(mat_conect(ii,3),:); mat_ident(mat_conect(ii,4),:)];
K_global=K_global+mat_transf1'*K_elementar*mat_transf1;
end
% IMPOSIÇÃO DAS CONDIÇÕES DE CONTORNO PELO MÉTODO DO PARTICIONAMENTO
%identificacao dos gdl livres e gdl impostos
gdl_livres=forcas_aplic(:,1);
gdl_impostos = cond_cont(:,1);
% construção das submatrizes de rigidez
K_ll=K_global(gdl_livres,gdl_livres);
K_li=K_global(gdl_livres,gdl_impostos);
K_ii=K_global(gdl_impostos,gdl_impostos);
% construção dos vetores de forças nos gdl livres e de deslocamentos nos gdl
impostos
f_liv=forcas_aplic(:,2);
d_imp=cond_cont(:,2);
% CÁLCULO DOS DESLOCAMENTOS NOS GDL LIVRES E FORÇAS DE REAÇÃO NOS GDL IMPOSTOS
disp('Deslocamentos nos gdl livres (em milímetros)');
d_liv=inv(K_ll)*(f_liv-K_li*d_imp)*1000
disp('Reacoes de apoio (em kN) ');
f_imp=(K_li'*d_liv/1000+K_ii*d_imp)/1000
des_global=zeros(nb_gdl,1);
des_global(cond_cont(:,1))=cond_cont(:,2);
des_global(forcas_aplic(:,1))=d_liv/1000;
% CÁLCULO DAS TENSÕES NAS BARRAS
for ii=1:nb_ele
alfa=propriedades(ii,4)*pi/180;
D.A.Rade Formulação do MEF pelo Processo Direto
32
mat_transf=[cos(alfa) sin(alfa) 0 0
0 0 cos(alfa) sin(alfa)];
des_local=mat_transf*des_global(mat_conect(ii,:));
sigma(ii)=propriedades(ii,2)/propriedades(ii,3)*(des_local(2)-des_local(1));
end
disp('Tensões nas barras (em MN/m^2) ');
sigma/1e6
Figura 2.8
D.A.Rade Formulação do MEF pelo Processo Direto
33
2.5 - Aplicação do processo direto a outros tipos de problemas unidimensionais
A formulação do MEF pelo processo direto, ilustrado na seção 2.1 para
sistemas formados por associações de molas lineares pode ser utilizado para a
resolução de vários outros tipos de problemas unidimensionais em regime
permanente nas áreas de Mecânica dos Sólidos, Mecânica dos Fluidos e
Transferência de Calor. Para cada um destes tipos de problemas o processo
consiste em determinar, a partir da aplicação das leis Físicas pertinentes, as
equações equivalentes às equações de equilíbrio das molas. A partir daí, o
procedimento de montagem e resolução das equações globais segue o mesmo
procedimento estudado anteriormente. A obtenção das equações equivalentes é
ilustrada a seguir.
2.5.1 - Barras em solicitação axial
Consideramos aqui problemas formados por barras solicitadas por forças
axiais, como exemplificado na Figura 2.9, na qual, para cada segmento da barra,
, ,i i iE A l indicam, respectivamente, o módulo de elasticidade, área de seção
transversal e comprimento. ei iN U designam, respectivamente, as forças axiais e
os deslocamentos aplicados nos nós que delimitam cada segmento.
Figura 2.9
1 (1)
(2)
(3)
N1
U1
N2
U2
N3
U3
N4
U4
E1, A1, l1
E2, A2, l2
E3, A3, l3
2 3 4
(i)
NiE
Ui E
Ni D
Ui D
Ei, Ai, li
i i+1
ki= EiAi/li
fiE
ui E
fi D
ui D
D.A.Rade Formulação do MEF pelo
Processo Direto
34
Da Resistência dos Materiais sabemos que para cada segmento em
equilíbrio, considerado isoladamente, admitindo comportamento linear elástico,
valem as seguintes relações:
D D Ei ii i i i i i i i
i
E AN A E A U U
l
E D E Di ii i i i
i
E AN N U U
l
As duas equações acima por se postas na seguinte forma matricial:
1 1
1 1
E E
i ii i
D D
i i i
U NE A
l U N
(2.38)
Comparando as equações (2.1) e (2.38) concluímos que cada segmento da
barra solicitada axialmente pode ser considerado como uma mola linear cuja
constante de rigidez equivalente é dada por
i i
i
i
E Ak
l
Desta forma, o procedimento de montagem das equações globais e
imposição das condições de contorno e carregamentos é feito da mesma forma
adotada para os sistemas de molas considerados na Seção 2.1.
2.5.2 Eixos sujeitos a torção
Consideramos aqui sistemas formados por eixos solicitados por momentos
torsores (torques), conforme mostrado na Figura 2.10, na qual, para cada
segmento, , ,i i iG J l indicam, respectivamente, o módulo de cisalhamento, momento
polar de inércia da seção transversal e comprimento. Além disso, ei i
designam, respectivamente, os torques aplicados e ângulos de torção das seções
transversais nos nós que delimitam cada segmento.
D.A.Rade Formulação do MEF pelo Processo Direto
35
Figura 2.10
Da Resistência dos Materiais sabemos que para cada segmento em
equilíbrio, considerado isoladamente, admitindo comportamento linear elástico,
valem as seguintes relações:
D D Ei ii i i
i
G J
l
E D E Di ii i i i
i
G J
l
As duas equações acima por se postas na seguinte forma matricial:
1 1
1 1
E E
i ii i
D D
i i i
G J
l
(2.39)
Comparando as equações (2.1) e (2.39) concluímos que cada segmento da
barra solicitada em torção pode ser considerado como uma mola linear cuja
constante de rigidez equivalente é dada por
i i
i
i
G Jk
l
Para o caso de eixos de seção circular de raio iR , o momento polar de
inércia é dado por:
4
2
i
i
RJ
G1, J1, l1 G2, J2, l2
T1 ,1 T2 ,2
T3 , 3
(2) (1)
ki = GiJi,/li
TiE
iE
TiD
iD
Gi, Ji, li
TiE , i E TiD
(1)
, i D
D.A.Rade Formulação do MEF pelo Processo Direto
36
O procedimento de montagem das equações globais e imposição das
condições de contorno e carregamentos é feito da mesma forma adotada para os
sistemas de molas considerados na Seção 2.1.
2.5.3 Problemas de escoamento em regime permanente
Nesta seção tratamos o problema de escoamento unidimensional em
regime permanente em tubulações que podem ser consideradas constituída pela
associação de vários segmentos, como ilustrado na Figura 2.11(a) que mostra a
vazão total Q e as vazões que percorrem os segmentos, indicadas por qi , i=1,2, ...
O problema consiste em determinar os valores destas vazões, bem como das
pressões nos nós que delimitam os segmentos, os quais são também indicados na
Figura 2.11.
Figura 2.11(a)
A Figura 2.11(b) mostra um segmento arbitrário, delimitado por dois nós E
e D, nos quais estão indicadas as respectivas pressões Eip e Dip e vazões Eiq e Diq .
Figura 2.11(b)
Em virtude da hipótese de escoamento em regime permanente, para cada
segmento i deve-se ter:
A
E H
D Q Q
q1
q2 q2
q3
q4
F G
B Cq1
q1
q3
1/Ri
E
iq
D
iq
E
ip
D
ip
i
E
iq
E
ip D
ip
D
iq
D.A.Rade Formulação do MEF pelo Processo Direto
37
D E
i i iq q q (2.41)
Além disso, sabe-se da Mecânica dos Fluidos que:
D E
i i i ip p R q (2.42)
onde Ri é o coeficiente de resistência, que depende essencialmente do tipo de
fluido, da geometria da seção transversal e do acabamento da superfície interna
do duto.
Combinando as equações (2.41) e (2,42) escrevemos:
1E E Di i i
i
q p p
R
1D E Di i i
i
q p p
R
ou:
1 11
1 1
E E
i i
D D
ii i
q p
Rq p
(2.43)
Como nos casos anteriores, a comparação de (2.1) e (2.43) permite concluir
que o problema de escoamento unidimensional em regime permanente pode ser
modelado da mesma forma que o sistema de molas lineares considerado na Seção
2.1, sendo a rigidez de mola equivalente dada por 1/ Ri.
2.5.4 Transferência de calor unidimensional em regime permanente
Tratamos aqui problema de transferência de calor unidimensional em
regime permanente em uma barra constituída pela associação de vários
segmentos, como ilustrado na Figura 2.12(a), onde Ki, Ai e li designam,
respectivamente, o coeficiente de condutividade térmica, a área da seção
transversal e o comprimento dos segmentos que compõem a barra. O problema a
ser resolvido consiste em determinar a distribuição de temperatura T(x) para um
conjunto de condições de contorno impostas nas extremidades 1 e 4.
A Figura 2.12(b) mostra um segmento genérico i da barra e as grandezas a
ele associadas.
D.A.Rade Formulação do MEF pelo Processo Direto
38
Figura 2.12(a)
Figura 2.12(b)
Recordamos abaixo as equações pertinentes aos dois mecanismos de
transferência de calor (condução e convecção) que consideramos estar envolvidos
no problema em questão.
Condução (Lei de Fourier):
i i i
dTq K A
dx
(2.44.a)
Convecção:
i i i Pq h A T T (2.44.b)
qi
D
Ki, Ai, li EiT
E
i
D
iT
E
iq
D
iq
E D
i iT T
Ki
TiE
qiE
TiD
qiD
Q1
T1
Q4
K1, A1, l1
K2, A2, l2
K3, A3, l3
1 3 4
T2 T3 T4
2
x
T
TP
qi
D.A.Rade Formulação do MEF pelo Processo Direto
39
Análise em nível elementar
E
i iq q Di iq q
E Di ii i i
i
K Aq T T
l
Di iq q
Assim, para o elemento i:
1 1
1 1
E E
i ii i
D D
ii i
q TK A
lq T
(2.45)
Como nos casos anteriores, a comparação de (2.1) e (2.44) mostra que o
problema de transferência de calor unidimensional em regime permanente pode
ser modelado da mesma forma que o sistema de molas lineares considerado na
Seção 2.1, sendo a rigidez de mola equivalente dada por KiAi/li.
Uma especificidade a ser considerada neste problema é que são possíveis
combinações de três tipos de condição de contorno aplicadas nas extremidades da
barra:
1) valor do calor imposto;
2) valor da temperatura
imposta;
3) condição de contorno do tipo convecção com o valor da temperatura T
imposto;
A condição de contorno do tipo temperatura imposta deve ser tratada da
mesma forma que os deslocamentos impostos, apresentada na Seção 2.1.3 e a
condição de valor imposto é análoga à aplicação de uma força de valor conhecido
nos problemas estruturais.
O tratamento das condições de contorno do tipo convecção requer algumas
transformações adicionais das equações, uma vez que este tipo de condição de
contorno envolve as temperaturas das extremidades da barra, conforme pode ser
visto na Equação (2.44.b).
Considerando, para exemplificação, a barra discretizada com três
elementos finitos, as equações globais resultantes do procedimento de montagem,
tomam a forma:
D.A.Rade Formulação do MEF pelo Processo Direto
40
1 1 1 1
1 1
1 11 1 1 1 2 2 2 2
1 1 2 2 2
33 3 3 32 2 2 2
4 42 2 3 3
3 3 3 3
3 3
0 0
0
0
0
0
0 0
K A K A
l l
q TK A K A K A K A
l l l l T
TK A K AK A K A
q Tl l l l
K A K A
l l
(2.46)
Admitamos que na extremidade 1 tenhamos uma condição de contorno do
tipo temperatura imposta e na extremidade 4 tenhamos uma condição de
contorno de convecção, as quais são expressas segundo :
1 1T T
4 4 3 4q h A T T
Introduzindo estas expressões em (2.46), obtemos:
1 1 1 1
1 1
1 11 1 1 1 2 2 2 2
1 1 2 2 2
33 3 3 32 2 2 2
4 3 2 2 3 3 4
3 3 3 3
4 3
3 3
0 0
0
0
0
0
0 0
K A K A
l l
q TK A K A K A K A
l l l l T
TK A K AK A K A
h A T l l l l T
K A K A h A
l l
O sistema de equações acima pode ser decomposto em dois conjuntos de
equações cuja resolução fornece o valor das temperaturas desconhecidas T2, T3 e
T4 e o calor transferido através da seção transversal no nó 1, segundo:
1 1 2 2 2 2
1 2 2
2
3 3 3 32 2 2 2
3
2 2 3 3
4 3 4
3 3 3 3
4 3
3 3
0
0
0
0
K A K A K A
l l l
T
K A K AK A K A T
l l l l
h A T T
K A K A h A
l l
D.A.Rade Formulação do MEF pelo Processo Direto
41
1
21 1 1 1
1
1 1 3
4
0 0
T
TK A K Aq
l l T
T
D.A. Rade Fundamentos de Cálculo Variacional
42
CAPÍTULO 3
FUNDAMENTOS DE CÁLCULO VARIACIONAL
É possível, e muito conveniente, considerar o procedimento de obtenção de
modelos de elementos finitos como um processo variacional. Através deste
processo procura-se estabelecer as equações do movimento ou de equilíbrio em
nível elementar e em nível global mediante a busca dos pontos estacionários (ou
pontos críticos, ou pontos extremos) de certas funções, denominadas funcionais,
que são funções escalares de funções que definem a distribuição da variável de
campo sobre o domínio considerado. Neste contexto o Cálculo Variacional
constitui o conjunto de técnicas matemáticas destinadas à determinação dos
pontos estacionários de funcionais.
Os chamados Princípios Variacionais são leis que estabelecem que a
configuração de equilíbrio estático ou de equilíbrio dinâmico de um dado sistema
mecânico devem corresponder a pontos estacionários de determinados funcionais
associados ao problema. Em diversos problemas da Mecânica do Contínuo, os
funcionais representam a energia do sistema. As equações de equilíbrio ou do
movimento podem então ser obtidas desenvolvendo, matematicamente, as
condições de estacionaridade destes funcionais.
Neste Capítulo, apresentaremos os fundamentos do Cálculo Variacional e
ilustraremos os desenvolvimentos apresentados através de exemplos baseados no
uso de Princípios Variacionais da Mecânica, notadamente o Princípio da Energia
Potencial Estacionária e o Princípio de Hamilton. Tais princípios serão
apresentados detalhadamente no Capítulo 5.
3.1 Funcionais envolvendo uma variável dependente e uma variável
independente
Consideremos inicialmente o problema de determinação dos pontos
estacionários do seguinte funcional:
dxxu,xu,xFuI x
x 21 (3.1)
Admitiremos ainda que a função u deva satisfazer as seguintes condições
de contorno:
11 uxu ; 22 uxu (3.2)
Nas equações acima, x designa a variável independente, xu é a variável
dependente, que representa a variável de campo associada ao problema, u indica
a derivada de xu em relação à variável independente, xu,xu,xF indica
genericamente uma relação funcional envolvendo x, xu e xu e 1u e 2u
D.A. Rade Fundamentos de Cálculo Variacional
43
designam valores prescritos, assumidos pela variável dependente nos extremos
do intervalo em que o problema é definido.
Designando doravante por u=u(x) a função que corresponde à solução
estacionária exata do funcional definido em (3.1), o procedimento variacional
consiste em definir as chamadas funções de comparação, representadas sob a
forma:
xuxuxu~ 21 xxx (3.3.a)
com:
xxu (3.3.b)
onde 1 é um parâmetro escalar arbitrário e x é uma função contínua
arbitrária definida no intervalo 21 xxx , que assume valor nulo nos extremos
deste:
01 x 02 x (3.4)
Em decorrência de (3.4), as funções de comparação satisfazem as condições
de contorno do problema:
11 uxu~ ; 22 uxu~ (3.5)
Em (3.3.b), x é entendida como uma função arbitrária que representa
uma perturbação infinitesimal acrescentada à solução exata u(x) para cada valor
da variável independente x, conforme ilustra a Figura 3.1, onde ficam
evidenciadas as condições (3.4) e (3.5).
Figura 3.1
x 2x 1x
u, u
xy
xu
xu~
1u
2u
D.A. Rade Fundamentos de Cálculo Variacional
44
Na seqüência do processo variacional, devemos introduzir a solução
variada (3.3) no funcional dado por (3.1) e impor a condição de estacionariedade
do funcional, que determina que sua variação I , associada a qualquer variação
arbitrária u deve ser nula. Esquematicamente, este processo é ilustrado como
segue:
0 IIIFFuu
Buscaremos, primeiramente, expressar a variação da função
xu,xu,xF , designada por xu,xu,xF , resultante da substituição de xu
por xuxu . Assim, definimos:
xu,xu,xFuu,uu,xFxu,xu,xF (3.6)
Em decorrência de (3.3.b), temos:
xu (3.7)
Introduzindo (3.3.b) e (3.7) em (3.6), obtemos:
xu,xu,xFxu,xu,xFxu,xu,xF (3.8)
Negligenciando os termos de ordem igual e superior à segunda em ,
procedemos à expansão do primeiro termo do lado direito de (3.8) em série de
Taylor, obtendo:
x
u
Fx
u
Fxu,xu,xFxu,xu,xF (3.9)
Introduzindo (3.9) em (3.8), temos:
F FF x,u x ,u x x x
u u
(3.10)
A variação do funcional I correspondente à variação da função F é definida
segundo:
FIFFII = dxxu,xu,xFxu,xu,xFx
x 21
dxxu,xu,xFx
x 21 dxxu,xu,xFxx 21 (3.11)
Introduzindo (3.10) em (3.11), temos:
D.A. Rade Fundamentos de Cálculo Variacional
45
2
1
x
x
F FI x,u x ,u x x x dx
u u
(3.12)
A condição necessária para a estacionaridade do funcional I é que sua
variação deve anular-se, ou seja, 0 xu,xu,xI . Assim, com base em (3.12),
devemos ter:
2
1
0
x
x
F Fx x dx
u u
(3.13)
Procedemos, em seguida, à integração por partes da segunda parcela do
integrando:
2
1
2
1
2
1
x
x
x
x
x
x
dxx
u
F
dx
dx
u
Fdxx
u
F (3.14)
Introduzindo (3.14) em (3.13), temos:
22
1 1
0
xx
x x
F F d Fx x dx
u u dx u
(3.15)
Em decorrência de (3.4), a primeira parcela de (3.15) se anula:
02
1
x
x
x
u
F , (3.16)
de modo que:
2
1
0
x
x
F d F x dx
u dx u
(3.17)
Neste ponto, evocamos o seguinte lema:
Lema: Se 1x e 2x ( 1x < 2x ) são constantes e G(x) é uma função contínua de x
no intervalo 21; xx e se
02
1
dxxGx
x
x
para toda e qualquer função continuamente diferenciável x que satisfaz
D.A. Rade Fundamentos de Cálculo Variacional
46
01 x 02 x ,
conclui-se que 0xG identicamente no intervalo 21; xx .
Com base neste lema podemos afirmar que a equação (3.17) se verifica se:
0F d F
u dx u
(3.18)
A equação (3.18) é denominada Equação de Euler-Lagrange associada ao
funcional dado por (3.1). Por outro lado, a equação (3.16) indica todas as possíveis
combinações de condições de contorno impostas nos extremos do intervalo 21; xx .
No caso do funcional analisado, as condições de contorno que satisfazem (3.16)
são:
• Em x= 1x : 01 x (3.18.a)
0
1
xxu
F (3.18.b)
• Em x= 2x : 02 x (3.18.c)
0
2
xxu
F (3.18.d)
Neste ponto, introduzimos a seguinte classificação para as condições de
contorno:
a) condições de contorno geométricas ou de Dirichlet: são aquelas que
envolvem diretamente a função de comparação e suas derivadas.
b) condições de contorno naturais ou de Neumann: são aquelas
representadas por expressões que envolvem as derivadas da função
xu,xu,xF em relação a um de seus argumentos ou derivadas de
x .
De acordo com estas definições, as condições de contorno indicadas pelas
equações (3.18.a) e (3.18.c) são condições de contorno geométricas, ao passo que
aquelas indicadas pelas equações equações (3.18.b) e (3.18.d) são condições de
contorno naturais.
D.A. Rade Fundamentos de Cálculo Variacional
47
Observações:
• a equação de Euler-Lagrange (3.18) apresenta-se sob a forma de uma
equação diferencial, cuja resolução, complementada pelas condições de contorno
(3.2), fornece a função u(x) que corresponde ao ponto estacionário do funcional
(3.1).
• para um determinado ponto estacionário do funcional, pode-se
determinar sua natureza (máximo ou mínimo), a partir da análise do sinal da
segunda variação do funcional I. A segunda variação é obtida pela expansão de uuI em torno de uI utilizando uma série de Taylor limitada aos termos de
segunda ordem em u . Para o funcional em questão, escrevemos:
uu
uu
Fu
u
Fu
u
Fu
u
Fu
u
Fu,u,xFuu,uu,xF
2
2
2
2
2
2
2
2
1
2
1
ou:
u,u,xFu,u,xFu,u,xFuu,uu,xF 2 (3.19)
com:
u
u
Fu
u
Fu,u,xF
(3.20)
uu
uu
Fu
u
Fu
u
Fu,u,xF
2
2
2
2
2
2
2
2
2
1
2
1 (3.21)
Introduzimos agora a segunda variação do funcional I:
dxuu,uu,xFuuI x
x 21 (3.22)
Introduzindo (3.19) em (3.22), temos:
2
1
2
1
2
1
22
x
x
x
x
x
x
dxu,u,xFdxu,u,xFdxu,u,xFuuI (3.23)
Levando em conta (3.1) e (3.11), reescrevemos a equação acima sob a
forma:
IIuIuuI 2 (3.24)
onde:
D.A. Rade Fundamentos de Cálculo Variacional
48
2
1
222
x
x
dxu,u,xFI (3.25)
é a segunda variação do funcional I.
Introduzindo (3.21) em (3.25) e levando em conta (3.3.b), escrevemos:
2
1
2
2
2
2
2
2
2
22 2
2
1 x
x
dx
uu
F
u
F
u
FI (3.26)
Pode-se mostrar que um dado ponto estacionário do funcional I é um ponto
de mínimo se 02 I e um ponto de máximo se 02 I para qualquer função
arbitrária .
Exemplo 1: O problema de equilíbrio de uma corda tensionada
Consideremos a corda elástica de comprimento L ilustrada na Figura 3.2,
sujeita a uma tração axial constante T e a um carregamento transversal
distribuído (força por unidade de comprimento) q(x). Nesta mesma figura u=u(x)
indica o campo de deslocamento transversais correspondente à posição de
equilíbrio da corda.
Pode-se mostrar que a energia potencial total do sistema, incluindo a
energia de deformação da corda e o trabalho da carga distribuída é dada pelo
seguinte funcional:
dxquuTuI L 0
2
2
(a)
Figura 3.2
Confrontando a equação (a) com a equação (3.1), reconhecemos:
L
q(x)
x
y
u(x)
D.A. Rade Fundamentos de Cálculo Variacional
49
quuTu,uF 2
2
(b)
Após o cálculo das derivadas parciais necessárias, as equações (3.18) e
(3.16), nos fornecem, respectivamente, a equação de Euler-Lagrange e os termos
de contorno do problema:
T
xqu (c)
00 LuT (d)
Admitindo constanteqxq , integrando duas vezes (c) e introduzindo as
condições de contorno determinadas pelos vínculos mostrados na Figura 3.1:
00 u 0Lu ,
obtemos a seguinte expressão para o campo de deslocamentos transversais da
corda em equilíbrio:
xLx
T
qxu
2
Investiguemos agora a natureza do ponto estacionário, mediante a análise
do sinal da segunda variação do funcional (a), dada pela equação (3.26). Para
tanto, a partir de (b) avaliamos as seguintes derivadas:
02
2
u
F ; T
u
F
2
2
; 0
2
uu
F
Introduzindo as equações acima em (3.26), obtemos:
22
0
1 2
2
L
I T dx
Observamos que a segunda variação do funcional é sempre positiva,
independentemente da função x , o que indica que o ponto estacionário
determinado corresponde a um ponto de mínimo do funcional (a).
D.A. Rade Fundamentos de Cálculo Variacional
50
3.2 Funcionais envolvendo derivadas de ordem superior.
O procedimento desenvolvido na seção precedente pode ser estendido para
a determinação dos pontos extremos de um funcionais envolvendo derivadas de
ordem superior à primeira. Consideremos o funcional expresso dado por:
dxxu,,xu,xu,xFuI x
x
n 2
1
(3.27)
onde a função u x satisfaz às seguintes condições de contorno:
11 uxu ; 22 uxu
11 uxu ; 22 uxu (3.28)
1111 nn uxu ; 1221 nn uxu
Por procedimento similar àquele apresentado na seção anterior,
introduzimos a seguinte função de comparação:
xuxuxu~ 21 xxx
com:
xxu
onde 1 é um parâmetro escalar arbitrário e x é uma função contínua
arbitrária definida no intervalo 21 xxx , que satisfaz às seguintes condições:
01 x ; 02 x
01 x ; 02 x (3.29)
011 xn ; 021 xn
Em decorrência de (3.4), as funções de comparação satisfazem às condições
de contorno do problema dadas por (3.28).
A variação da função F resulta:
nnn uFuFuFu,,u,u,xF (3.30)
D.A. Rade Fundamentos de Cálculo Variacional
51
e a condição de estacionariedade fica:
02
1
dxu,u,u,xFI xx n (3.31)
Introduzindo (3.30) em (3.31), obtemos:
2
1
0
x
x
)n(
)n( dxu
F
u
F
u
F (3.32)
Desenvolvemos (3.32) fazendo integrações por partes sucessivas até a
eliminação, no integrando, das derivadas n,, . Para tanto, utilizamos o
seguinte desenvolvimento para a integração por partes:
2 2
1 1
1
x x k
kk
kk k
x x
F d Fdx dx
dxu u
2
1
2 1
0 1 2 11 2 3
2 11 1 1 1
xk
kk k k
kk k k k
x
F d F d F d F
dx dx dxu u u u
ou, de forma mais compacta:
22 2
1 1 1
1
1
0
1 1
xx x k jk
k jk k j
k jk k k
jx x x
F d F d Fdx dx
dx dxu u u
(3.33)
Introduzindo (3.33) em (3.32), a condição de estacionariedade do funcional
toma a forma:
2
1 0
1
x kn
k
k k
kx
d F dx
dx u
01
0
1
0
1
2
1
n
k
k
j
x
x
jk
kj
j
j
u
F
dx
d (3.34)
Em (3.34), a segunda parcela do lado esquerdo se anula em virtude das
condições (3.29), o que conduz a:
2
1 0
1 0
x kn
k
k k
kx
d F dx
dx u
(3.35)
D.A. Rade Fundamentos de Cálculo Variacional
52
Utilizando mais uma vez o Lema fundamental do Cálculo Variacional
apresentado na seção precedente, concluímos que a equação (3.35) se verifica
para qualquer função arbitrária x à condição que:
0
1 0
kn
k
k k
k
d F
dx u
(3.36)
A equação (3.36) é a Equação de Euler-Lagrange associada ao funcional
dado por (3.27).
Os termos de contorno que figuram em (3.34):
01
0
1
0
1
2
1
n
k
k
j
x
x
jk
kj
j
j
u
F
dx
d (3.37)
representam o conjunto de todas as possíveis condições de contorno do problema.
Exemplo 2: O problema de equilíbrio de uma viga de Euler-Bernoulli
Consideremos a viga de uniforme ilustrada na Figura 3.3, de módulo de
rigidez à flexão EI, submetida a um carregamento transversal distribuído q(x).
Neste caso, com base na teoria de vigas de Euler-Bernoulli, a energia potencial
total do sistema, incluindo a energia de deformação elástica e o trabalho da carga
distribuída, associada a um campo de deslocamentos transversais representado
por u(x), é dada pelo seguinte funcional:
dxxyxq
dx
udEII
L
0
2
2
2
2
1 (a)
Figura 3.3
L
q(x)
x
u
u(x)
D.A. Rade Fundamentos de Cálculo Variacional
53
Confrontando as equações (a) e (3.27), identificamos n=2 e
21
2
F u,u EI u q x u x (b)
Particularizando (3.36) e (3.37) para n=2, obtemos:
equação de Euler-Lagrange:
02
2
u
F
dx
d
u
F
dx
d
u
F (c)
termos de contorno:
0
0
L
u
F
u
F
dx
d
u
F (d)
Do desenvolvimento de (c) a partir de (b), resulta a seguinte equação de
Euler-Lagrange:
EI
xqu 4 , (e)
que é a bem conhecida equação diferencial da linha elástica de vigas de Euler-
Bernoulli.
Desenvolvendo os termos de contorno dados por (d), obtemos a expressão:
00 LuEIuEI (f)
Em (f), associamos a função ao deslocamento transversal da viga, à
rotação da seção transversal. Relembramos ainda as seguintes definições:
• uEIM : momento fletor
• uEIV : esforço cortante
Observamos que qualquer que seja o tipo de vínculo existente nas
extremidades da viga, em x=0 e x=L, a equação (f) será satisfeita. Com efeito:
• Para uma extremidade livre: M = 0 ; V = 0.
• Para uma extremidade engastada: = 0 ; = 0
D.A. Rade Fundamentos de Cálculo Variacional
54
• Para uma extremidade simplesmente apoiada: M = 0 ; = 0
3.3 Funcionais envolvendo várias variáveis independentes
Em diversos tipos de problemas de Engenharia, duas ou mais variáveis
independentes podem intervir. É o caso, por exemplo, de problemas de equilíbrio
bidimensionais nos quais a variável de campo são funções de duas coordenadas
espaciais.
Consideraremos, no desenvolvimento que segue, funcionais do tipo:
D
yx dydxu,u,y,x,uFuI (3.38)
onde x e y são as variáveis independentes, D é um domínio bidimensional
definido no plano x-y (ver Figura 3.4), y,xuu é a variável dependente e yx uu e
designam, respectivamente, as derivadas da função u, em relação a x e y.
Admitiremos, ainda, que a função u deva assumir um conjunto de valores
prescritos ao longo de uma curva fechada C que constitui a fronteira do domínio
D:
y,xuy,xu CC (3.39)
Figura 3.4
x
y
D
C
D.A. Rade Fundamentos de Cálculo Variacional
55
Para o estabelecimento das condições de estacionariedade do funcional
(3.38), procedemos de forma análoga àquela adotada na seção anterior,
introduzindo a função de comparação:
y,xuy,xuy,xu ,
com:
y,xy,xu (3.40)
A função de comparação deve satisfazer as seguintes condições de
fronteira:
y,xuy,xu~ CC (3.41)
De (3.39) e (3.41) decorre que a função y,x deve anular-se identicamente
na fronteira C:
0Cy,x (3.42)
A variação da função yx u,u,y,x,uF associada a y,xu é dada por:
y
y
x
x
y
y
x
x u
F
u
F
u
Fu
u
Fu
u
Fu
u
FF (3.43)
de modo que a condição necessária para que o funcional I tenha um ponto
estacionário é:
D
dydxFI 0
D
y
y
x
x u
F
u
F
u
F (3.44)
o que conduz a:
0
D
y
y
x
x u
F
u
F
u
F (3.45)
Para efetuar as integrações por partes das duas últimas parcelas do
integrando em (3.45), empregaremos o Teorema de Green em duas dimensões,
expresso sob a forma:
dxy,xFdyy,xGy,xdydx
v
y,xy,xF
x
y,xv,xG
CD
D.A. Rade Fundamentos de Cálculo Variacional
56
dydx
y
y,xF
x
y,xGy,x
D
(3.46)
Utilizando (3.46), procedemos à integração dos dois últimos termos de
(3.45) da seguinte forma:
D
y
y
x
x
dydx
u
F
u
F
=
dyuFdxuF yxC dydxu
F
dy
d
u
F
dx
d
D yx
(3.47)
Introduzindo (3.47) em (3.45), obtemos:
dyuFdxuF yxC dydxu
F
dy
d
u
F
dx
d
u
F
D yx
= 0 (3.48)
Em virtude da condição (3.42) e do Lema 1, de (3.48) resulta:
0
yx u
F
dy
d
u
F
dx
d
u
F (3.49)
dyuFdxuF yxC = 0 (3.50)
A equação (3.49) é a Equação de Euler-Lagrange associada ao funcional
dado por (3.38) e os termos de contorno (3.50) fornecem todas as combinações
possíveis de condições de contorno naturais e geométricas a serem satisfeitas na
fronteira do domínio de definição do problema.
D.A. Rade Fundamentos de Cálculo Variacional
57
Exemplo 3: O problema de equilíbrio de uma membrana tensionada
Consideremos a membrana mostrada na Figura 3.5, a qual se encontra em
equilíbrio estático, submetida a um carregamento transversal externo q(x,y) e a
uma tensão constante T em seu plano.
Figura 3.5
A energia potencial total, incluindo o trabalho da carga distribuída q(x,y) é
dada pelo funcional:
D yx dxdyquuuTuI 2221 (a)
onde u=u(x,y) designa o campo de deslocamentos transversais da membrana.
Confrontando a equação (a) com (3.38), identificamos:
yx u,u,u,y,xF quuuT yx 2221 (b)
Aplicando a equação (3.49), obtemos a seguinte equação de Euler-
Lagrange:
022 quuT yx (c)
Em particular, se D for um domínio circular de raio R, convém expressar
(c) em coordenadas polares, introduzindo as transformações (ver Figura 3.6):
y
x
z
q(x,y)
T
T
D.A. Rade Fundamentos de Cálculo Variacional
58
cosrx senrx 222 yxr
Figura 3.6
Assim procedendo, (c) assume a forma:
0
q
dr
dur
dr
d
r
T (d)
Admitindo que a carga distribuída seja constante, a integração da equação
(d), subsidiada pela condição de contorno:
0,ru para r = R,
conduz à seguinte expressão para o campo de deslocamentos transversais da
membrana:
22
4
rR
T
qru
Observemos que, devido à axissimetria geométrica e de carregamento, o
deslocamento transversal é independente do ângulo .
x
y
r
P(x,y)
R
D.A. Rade Fundamentos de Cálculo Variacional
59
3.4 Funcionais envolvendo várias variáveis dependentes
Nesta seção, consideraremos funcionais envolvendo mais de uma variável
dependente e uma única variável independente. No que segue, trataremos de
funcionais envolvendo apenas duas variáveis dependentes, sendo imediata a
extensão dos desenvolvimentos ao caso de funcionais envolvendo um número
maior de variáveis dependentes.
Buscaremos estabelecer as condições de estacionariedade do seguinte
funcional:
dxv,u,v,u,xFv,uI x
x xx 21 (3.51)
onde x é a variável independente, u=(x) e v=(x) são as variáveis dependentes e
xx vu e são, respectivamente, as derivadas de u e v em relação a x.
Admitiremos que as funções u e v satisfaçam as seguintes condições de
contorno:
11 uxu ; 22 uxu
(3.52) 11 vxv ; 22 vxv
Adotando procedimento análogo àquele empregado nas seções anteriores,
introduzimos as seguintes funções de comparação:
xuxuxu~ 21 xxx (3.53.a)
xvxvxv~ 21 xxx (3.53.b)
com:
xxu (3.53.a)
xxv (3.53.b)
onde 1 é um parâmetro escalar arbitrário e x e x são funções contínuas
arbitrárias definidas no intervalo 21 xxx , que devem assumir valores nulos
nos extremos deste intervalo:
01 x 02 x (3.55.a)
01 x 02 x (3.55.b)
D.A. Rade Fundamentos de Cálculo Variacional
60
Em decorrência de (3.55), as funções de comparação satisfazem as
condições de contorno do problema:
11 uxu~ ; 22 uxu~ (3.56.a)
11 uxv~ ; 22 uxv~ (3.56.b)
Empregando mais uma vez o procedimento baseado na expansão em série
de Taylor limitada aos termos de primeira ordem, obtemos a seguinte expressão
para a variação F :
x
x
x
x v
F
u
F
v
F
u
FF (3.57)
e a variação do funcional I correspondente à variação da função F resulta:
2
1
x
x
x
x
x
x
dx
v
F
u
F
v
F
u
FI (3.58)
De (3.58) resulta a seguinte condição de estacionariedade do funcional I:
00
2
1
x
x
x
x
x
x
dx
v
F
u
F
v
F
u
FI (3.59)
Procedemos, em seguida, à integração por partes das duas últimas parcelas
do integrando, conforme detalhado abaixo:
2
1
2
1
2
1
x
x x
x
xx
x
x
x
x
dx
u
F
dx
d
u
Fdx
u
F (3.60.a)
2
1
2
1
2
1
x
x x
x
xx
x
x
x
x
dx
v
F
dx
d
v
Fdx
v
F (3.60.b)
Introduzindo (3.60) em (3.50), obtemos:
2 2 2
1 1 1
0
x x x
x x x xx x x
F F F d F F d Fdx dx
u v u dx u v dx v
(3.61)
Em virtude de (3.55), os termos de contorno em (3.61) se anulam:
D.A. Rade Fundamentos de Cálculo Variacional
61
0
2
1
x
xxx v
F
u
F (3.62)
de modo que, de (3.61), resulta a condição:
2 2
1 1
0
x x
x xx x
F d F F d Fdx dx
u dx u v dx v
Mais uma vez, evocando o Lema Fundamental do Cálculo Variacional
empregado anteriormente, concluímos que a condição acima se verifica para
quaisquer funções arbitrárias x e x se, e somente se:
0
x
F d F
u dx u
(3.63.a)
0
x
F d F
v dx v
(3.63.b)
Os termos de contorno expressos por (3.62) representam todas as possíveis
condições de contorno geométricas e naturais do problema, enquanto (3.63)
constituem as Equações de Euler-Lagrange associadas ao funcional I.
O procedimento desenvolvido acima pode ser imediatamente estendido ao
caso de funcionais envolvendo um número qualquer de variáveis dependentes,
representados genericamente sob a forma:
dxz,,v,u,z,v,u,xFz,,v,uI x
x xxx 21 (3.64)
Neste caso, os termos de contorno resultam:
0
2
1
x
xxxx z
F
v
F
u
F (3.65)
e haverá tantas equações de Euler-Lagrange quantas forem as variáveis
dependentes:
0
x
F d F
u dx u
(3.66.a)
0
x
F d F
v dx v
(3.66.b)
D.A. Rade Fundamentos de Cálculo Variacional
62
0
x
F d F
z dx z
(3.66.c)
Exemplo 4: O problema de vibrações de um sistema mecânico discreto de dois
graus de liberdade.
Para o sistema mecânico ilustrado na Figura 3.7, dado um conjunto de
deslocamentos 1u e 2u e as respectivas velocidades 1u e 2u , são as seguintes as
expressões para a energia cinética e a energia de deformação do sistema:
22221121 umumT (a)
212122220211021 uukukukV (b)
Figura 3.6
Formemos o seguinte funcional:
2
1
t
t
dtVTI dtuukukukumumt
t
2
1
2
1221
2
220
2
110
2
22
2
112
1 (c)
Em (c), identificamos:
21221222021102222112121 21 uukukukumumu,u,u,uF (d)
Empregando as equações (3.66), obtemos as seguintes equações de Euler-
Lagrange:
10k 20k 12k
1m 2m
tu1 tu2
D.A. Rade Fundamentos de Cálculo Variacional
63
0121211011 uukukum (e.1)
0121222022 uukukum (e.2)
Podemos colocar as equações (e) sob a seguinte forma matricial:
10 12 121 1 1
12 20 122 2 2
0 0
0 0
k k km u u
k k km u u
(f.1)
ou ainda:
0 uKuM (f.2)
onde:
2
1
u
u
u (g.1)
2
1
0
0
m
m
M (matriz de inércia)
(g.2)
10 12 12
12 20 12
k k k
K
k k k
(matriz de rigidez) (g.3)
As equações (f) constituem as equações do movimento em regime livre (sem
forças externas) do sistema mecânico. Sua resolução, subsidiada por um conjunto
de condições iniciais, conduzem à determinação das funções tu1 e tu2 , que
descrevem o movimento das massas 1m e 2m em função do tempo.
D.A. Rade Princípios Variacionais da Mecânica
64
CAPÍTULO 4
PRINCÍPIOS VARIACIONAIS DA MECÂNICA
Neste Capítulo apresentaremos dois dos principais princípios variacionais
da Mecânica: o Princípio da Energia Potencial Estacionária, aplicável a corpos
deformáveis em equilíbrio estático e o Princípio de Hamilton, que rege os
problemas de dinâmica de sistemas rígidos ou flexíveis. A partir do princípio de
Hamilton, aplicando as técnicas do Cálculo Variacional, apresentadas no
Capítulo 3, obteremos as Equações de Lagrange do Movimento.
Será também apresentado o Método de Rayleigh-Ritz destinado à obtenção
de soluções aproximadas de problemas regidos por princípios variacionais. Nos
capítulos subseqüentes, o método de Rayleigh-Ritz, combinado com a noção de
aproximação por subdomínios, será utilizado para a formulação do método dos
elementos finitos, considerando alguns tipos particulares de elementos.
4.1 Fundamentos da Teoria da Elasticidade em pequenos deslocamentos
Nesta seção recapitularemos os fundamentos da Teoria da Elasticidade
com base na hipótese de pequenos deslocamentos. Esta hipótese admite que as
componentes do deslocamento de um ponto qualquer P do corpo, denotadas
doravante por u, v e w (ver Figura 4.1), sejam pequenas o bastante para justificar
a linearização das equações que governam o problema.
Figura 4.1
y
x
z
v
u
w
P’
P
D.A. Rade Princípios Variacionais da Mecânica
65
A Figura 4.2 ilustra o estado mais geral de tensões a que pode estar
submetido um elemento diferencial de volume de um corpo elástico, onde
zyx ,, são as componentes de tensões normais, zyyzzxxzyxxy ,,,,, são as
componentes de tensões cisalhantes e X , Y e Z são as componentes da força
externa
de volume (força por unidade de volume) atuante sobre o elemento.
Figura 4.2
• Tensor das tensões. O estado de tensões em um ponto qualquer do corpo
deformável é determinado por nove componentes formando o chamado tensor das
tensões, definido da seguinte forma:
zzyzx
yzyyx
xzxyx
(4.1)
y
x
z
x
xz
xy
z
zx
zy
y
yx
yz
Z
Y
X
dx
dy
dz
D.A. Rade Princípios Variacionais da Mecânica
66
• Equações de equilíbrio. Impondo as o equilíbrio de forças e momentos para o
elemento diferencial de volume mostrado na Figura 4.1, pode-se mostrar que as
componentes de tensão devem satisfazer as seguintes equações de equilíbrio:
yxxy (4.2a)
zxxz (4.2b)
zyyz (4.2c)
0
X
zyx
xzxyx (4.2d)
0
Y
zyx
yzyxy (4.2e)
0
Z
zyx
zyzxz (4.2f)
As equações (4.2a)-(4.2c) implicam a simetria do tensor das tensões que,
doravante, será expresso em termos de apenas seis componentes independentes
de tensões:
zyzxz
yzyxy
xzxyx
(4.3)
• Tensor das deformações. O estado de deformações em um ponto qualquer do
corpo deformável é determinado por nove componentes formando o chamado
tensor das deformações, definido da seguinte forma:
zzyzx
yzyyx
xzxyx
(4.4)
onde zyx ,, são as deformações normais e yzxzxy ,, são as deformações
cisalhantes.
D.A. Rade Princípios Variacionais da Mecânica
67
• Relações deformações-deslocamentos. Pode-se demonstrar que, admitindo a
hipótese de pequenos deslocamentos, são válidas as seguintes relações
envolvendo as componentes de deslocamento e as deformações:
x
u
x
(4.5a)
y
v
y
(4.5b)
z
w
z
(4.5c)
x
v
y
u
xy
(4.5d)
x
w
z
u
xz
(4.5e)
y
w
z
v
yz
(4.5f)
• Relações tensões-deformações. Para materiais obedecendo à Lei de Hooke
(proporcionalidade entre tensões e deformações), as chamadas equações
constitutivas, relacionando as componentes de tensões e as componentes de
deformações são expressas através da seguinte relação:
E (4.6)
onde:
Tyzxzxyzyx (4.7)
é o vetor contendo as componentes de tensões
Tyzxzxyzyx (4.8)
é o vetor contendo as componentes de deformações e
D.A. Rade Princípios Variacionais da Mecânica
68
666564636261
565554535251
464544434241
363534333231
262524222221
161514131211
eeeeee
eeeeee
eeeeee
eeeeee
eeeeee
eeeeee
E (4.9)
é o chamado tensor de rigidez ou tensor das constantes elásticas. As constantes
elásticas satisfazem às seguintes relações, que conferem simetria ao tensor de
rigidez.
jiij ee para i,j = 1,2, ..., 6 (4.10)
A relação (4.6) pode ser invertida, fornecendo:
S (4.11)
onde:
SE 1
666564636261
565554535251
464544434241
363534333231
262524232221
161514131211
ssssss
ssssss
ssssss
ssssss
ssssss
ssssss
(4.12)
é o chamado tensor de flexibilidade. Este é um tensor simétrico cujas
componentes satisfazem as relações:
jiij ss para i,j = 1,2, ..., 6 (4.13)
Para materiais isotrópicos, cujas características mecânicas independem da
direção, o número de constantes elásticas independentes fica reduzido a apenas
dois. Neste caso, as relações tensões-deformações assumem a forma:
zyxxx G 212 (4.14a)
zyxyy G 212 (4.14b)
zyxzz G 212 (4.14c)
D.A. Rade Princípios Variacionais da Mecânica
69
xyxy G (4.14d)
xzxz G (4.14e)
yzyz G (4.14f)
e, inversamente:
zyxx E 1 (4.15a)
zxyy E 1 (4.15b)
yxzz E 1 (4.15c)
com:
12
EG (4.16)
• Condições de contorno. No que diz respeito às condições de contorno, podemos
imaginar a superfície externa do corpo deformável dividida em duas regiões: a
região 1S , sobre a qual as condições de contorno são estabelecidas em termos de
forças externas aplicadas e a região 2S , sobre a qual as condições de contorno são
expressas em termos de deslocamentos impostos.
Denotando por Z,Y,X as componentes cartesianas as forças impostas
por unidade de área em 1S , as condições de contorno que traduzem o equilíbrio de
forças na fronteira são expressas segundo:
XX YY ZZ em 2S (4.17)
onde:
nmX xzxyx (4.18a)
nmY yzyxy (4.18b)
nm zyzxz (4.18c)
D.A. Rade Princípios Variacionais da Mecânica
70
As equações (4.18) são obtidas impondo o equilíbrio do tetraedro
infinitesimal mostrado na Figura 4.3, onde knjmin
é o vetor unitário na
direção da normal exterior à superfície 1S .
Figura 4.3
Denotando as componentes de deslocamentos impostos por w,v,u , as
condições de contorno ditas geométricas, aplicáveis à região 2S , são expressas
segundo:
uu (4.19a)
vv em 2S (4.19b)
ww (4.19c)
No problema geral de Elasticidade temos um total de 15 incógnitas: 3
componentes de deslocamento, 6 componentes de deformação e 6 componentes de
tensão. Estas incógnitas podem ser determinadas mediante a resolução de 15
equações diferenciais: 3 equações de equilíbrio (4.2), 6 relações deformações-
deslocamentos (4.5) e 6 relações tensões-deformações (4.6). A solução deve
também satisfazer as condições de contorno (4.17) e (4.19).
4.2 O Princípio do Trabalho Virtual em Elasticidade
Considerando um corpo deformável em equilíbrio sob a ação de forças
externas e um conjunto de restrições cinemáticas impostas, introduzimos o
conceito de deslocamentos virtuais:
kwjviur
(4.20)
x
y
z
n
z
zx
zy
y
yx
yz
x
xy
xz
dS1
D.A. Rade Princípios Variacionais da Mecânica
71
entendidos como um campo de deslocamentos infinitesimais arbitrários aplicados
aos pontos do corpo. Os deslocamentos virtuais são deslocamentos imaginários,
diferindo dos deslocamentos reais pelo fato que a eles não associamos nenhuma
variação do tempo. Deste modo, as forças atuantes sobre o corpo permanecem
inalteradas durante a ocorrência dos deslocamentos virtuais. Admite-se ainda
que os deslocamentos virtuais não violem as condições de contorno geométricas.
Desta forma, devemos ter:
0u 0v 0w em 2S (4.21)
Define-se o trabalho virtual de uma força F
como sendo o trabalho
realizado por esta força quando seu ponto de aplicação sobre um deslocamento
virtual:
rFW (4.22)
Em termos de componentes cartesianas:
wFvFuFkwjviukFjFiFW zyxzyx (4.23)
Com base nestas definições, é a seguinte a expressão para o trabalho
virtual das forças externas impostas sobre o contorno do corpo:
dSwZvYuXW
S
S
1
(4.24)
Considerando as condições de contorno na fronteira 1S , dadas por (4.17), a
equação (4.24) se reescreve:
dSwnmvnmunmW
S
zyzxzyzyxyxzxyx
S
1
(4.25)
Aplicando o Teorema de Gauss, expresso sob a forma:
dV
z
N
y
M
x
LdSNnMmL
S V
(4.26)
a equação (4.25) fica:
dVwvu
z
wvu
y
wvu
x
W
V
yzyzxzyzyxyxzxyx
S
Desenvolvendo as derivadas indicadas na equação acima, temos:
D.A. Rade Princípios Variacionais da Mecânica
72
dVw
zyx
v
zyx
u
zyx
W
V
zxyxzyzyxyxzxyxS
+
dV
z
w
z
v
y
w
y
v
z
u
x
w
y
u
x
v
x
u
V
zyzyxzxyx
(4.27)
Introduzimos agora as seguintes deformações associadas aos
deslocamentos virtuais:
z
w
y
v
x
u
zyx
(4.28)
z
v
y
w
z
u
x
w
y
u
x
v
yzxzxy
Utilizando as definições (4.28) e levando em conta as equações de equilíbrio
(5.2), reescrevemos (4.27) sob a forma:
dVwZvYuXW
V
S +
dV
V
yzyzxzxzxyxyzzyyxx (4.29)
Por fim, confrontando (4.24) e (4.29), escrevemos:
dSwZvYuX
S
1
dVwZvYuX
V
=
dV
V
yzyzxzxzxyxyzzyyxx (4.30)
ou ainda:
IVS WWW (4.31)
onde:
D.A. Rade Princípios Variacionais da Mecânica
73
• dSwZvYuXW
S
S
1
(4.32)
é o trabalho virtual das forças externas impostas na fronteira S1;
• dVwZvYuXW
V
V (4.33)
é o trabalho virtual das forças externas de volume;
• dVW
V
yzyzxzxzxyxyzzyyxx
I (4.34)
é o trabalho virtual das forças internas, associadas às tensões distribuídas sobre
o volume do corpo.
As equações (4.30) e (4.31) expressam o chamado Princípio do Trabalho
Virtual para corpos deformáveis em equilíbrio, sob a hipótese de pequenos
deslocamentos.
4.3 O Princípio da Energia Potencial Estacionária
Consideremos os casos em que seja possível definir uma função energia
potencial de deformação ,,,,,,UU yzxzxyzyx de modo que:
x
x
U
y
y
U
z
z
U
(4.35)
xy
xy
U
xz
xz
U
yz
yz
U
No caso de materiais obedecendo as relações constitutivas (4.6), podemos
facilmente verificar que a energia potencial de deformação é dada por:
EU TT
2
1
2
1 (4.36)
Mais particularmente no caso particular de materiais isotrópicos, regidos
pelas equações constitutivas (4.14), o desenvolvimento de (4.36) conduz à
expressão:
2222222 212121 yzxzxyzyxzyxU (4.37)
Com base em considerações físicas, pode-se mostrar que, em (4.36) U é uma
função definida positiva das componentes de deformação. Isto significa que U>0
para qualquer campo de deformações correspondente a um campo de
D.A. Rade Princípios Variacionais da Mecânica
74
deslocamentos cinematicamente admissível. Tal fato fica evidenciado na forma
particular (4.37).
Introduzindo (4.35) em (4.34), podemos expressar o trabalho virtual das
forças internas em termos da energia potencial de deformação da seguinte forma:
dVUUUUUUW
V
yz
yz
xz
xz
xy
xy
z
z
y
y
x
x
I
Na última equação acima observamos que o integrando representa o
diferencial total da função ,,,,,,UU yzxzxyzyx . Assim, podemos escrever:
dVUW
V
I (4.38)
com:
yz
yz
xz
xz
xy
xy
z
z
y
y
x
x
UUUUUUU
(4.39)
Podemos ainda permutar a ordem dos operadores integração e
diferenciação em (4.38), escrevendo:
dVUW
V
I (4.40)
Consideremos agora o trabalho virtual das forças externas (forças de
volume e forças impostas na fronteira S1. Admitindo que estas forças sejam
conservativas (o trabalho por elas realizado independe do caminho percorrido por
seu ponto de aplicação), suas componentes podem ser expressas em termos de
derivadas parciais de funções potenciais da seguinte forma:
u
X
v
Y
w
Z
(5.41)
u
X
vY
wZ
onde w,v,u e w,v,u são, respectivamente, as funções potenciais
associadas às forças externas de volume e às forças externas impostas na
fronteira S1.
Nestas condições, partindo de (4.32) e (4.33), o trabalho virtual das forças
externas pode ser expresso sob a forma:
D.A. Rade Princípios Variacionais da Mecânica
75
VS
IS dVw
w
v
v
u
u
dSw
w
v
v
u
u
WW
1
Os integrandos na equação acima correspondem aos diferenciais totais das
funções e , fato que nos permite escrever:
dVwvudSwvuWW
SS
IS
11
,,,, (4.42)
Assim, combinando as equações (4.38) e (4.42), o Princípio do Trabalho
Virtual, expresso por (5.31), assume a forma:
0
1
V V S
dSdVdVU (4.43a)
ou ainda:
0 (4.43b)
onde:
V V S
dSdVdVU
1
(4.44)
é a energia potencial total.
As equações (4.43) expressam o chamado Princípio da Energia Potencial
Estacionária, o qual estabelece que, dentre todos os campos de deslocamentos
(u,v,w) que satisfazem as condições de contorno geométricas, aquele que
corresponde à real configuração de equilíbrio é um ponto estacionário da energia
potencial total.
Do ponto de vista do Cálculo Variacional, podemos observar que as
equações de equilíbrio e as condições de contorno naturais em S1 são,
respectivamente, as equações de Euler-Lagrange e os termos de contorno do
problema de estacionariedade do funcional (4.44). Resta mostrar que o extremo
associado ao campo de deslocamentos da configuração de equilíbrio estático
corresponde a um ponto de mínimo do funcional (4.44). Para tanto, denotaremos,
a seguir, por (u,v,w) as componentes de deslocamento correspondentes à
configuração de equilíbrio e por wvu ,, um campo de deslocamentos
arbitrário, próximo a (u,v,w), que satisfazem as condições de contorno
geométricas. Relacionaremos ambos os campos de deslocamentos através das
relações:
uuu vvv www (4.45)
D.A. Rade Princípios Variacionais da Mecânica
76
Introduzindo (4.45) em (4.44) e empregando uma expansão em série de
Taylor limitada à segunda ordem, escrevemos:
wvuwvuwvuwvu ,,,,,,,, 2
Como (u,v,w) é um ponto estacionário do funcional , temos 0,, wvu
e a equação acima fica:
wvuwvuwvu ,,,,,, 2 (4.46)
A equação (5.46) nos mostra que a natureza do ponto crítico (mínimo ou
máximo) é determinada pelo sinal da segunda variação wvu ,,2 , de modo que:
• 0,,2 wvu ponto de mínimo
• 0,,2 wvu ponto de máximo
Podemos mostrar que:
dVwvuU
V
,,2 (4.47)
onde U é a energia potencial de deformação, dada por (4.36). Do fato de U ser
uma função definida positiva resulta que 02 para quaisquer valores
arbitrárias de wvu e, , o que demonstra que o campo de deslocamentos
correspondente à posição de equilíbrio do corpo deformável corresponde a um
ponto de mínimo do funcional que representa a energia potencial total. Por esta
razão, o Princípio da Energia Potencial Estacionária é às vezes denominado
Princípio da Energia Potencial Mínima.
4.4 – Princípio do Trabalho Virtual Aplicado à Dinâmica de Sistemas de
Partículas
Visando deduzir o Princípio Variacional de Hamilton, aplicável a sistemas
mecânicos em movimento, estabeleceremos primeiramente o Princípio do
Trabalho Virtual para um conjunto arbitrário de n partículas, ilustrado na
Figura 4.4, ficando subentendido que qualquer sistema mecânico discreto ou
contínuo poderá ser considerado como formas particulares deste tipo de conjunto.
Na Figura 4.4, sobre uma partícula genérica iP , de massa im , estão indicadas:
• iF
: a resultante das forças impostas sobre a partícula, representando as
ações mecânicas dos corpos a ela vizinhos, incluindo as forças exercidas pelas
outras partículas do sistema;
D.A. Rade Princípios Variacionais da Mecânica
77
• if
: a força de restrição exercida sobre pela superfície de restrição iS sobre
a qual iP é forçada a se mover. Negligenciaremos, por enquanto, o atrito, de modo
que a força de restrição terá sempre a direção normal à superfície de restrição,
indicada na Figura 7.1 por in .
• ir
: vetor posição da partícula em relação ao sistema de referência,
suposto fixo, OXYZ .
Figura 4.4
Aplicando a Segunda Lei de Newton para cada partícula do sistema,
considerada isoladamente, escrevemos:
iiii rmfF
i=1, 2, ..., n (4.48)
O Princípio de d’Alembert, apresentado na Seção 3.6, nos permite escrever
(4.49) sob a forma:
0 iiii rmfF
i=1, 2, ..., n (4.49)
com a interpretação de que, para um observador movimentando-se juntamente
com a partícula, esta é observada em estado de equilíbrio sob a ação das forças
Y
nP
1r
X
Z
ir
nr
iP
1P
O
1F
iF
nF
nf
if
1f
1r
ir
nr
it
in
1S
iS
nS
D.A. Rade Princípios Variacionais da Mecânica
78
iF
, if
e da força de inércia ii rm . Esta situação é denominada equilíbrio
dinâmico.
Introduzimos agora o conceito dos chamados deslocamentos virtuais, que
serão denotados por 1r
, ..., ir , ..., nr . Tais deslocamentos virtuais são
concebidos com as seguintes propriedades:
• são perturbações imaginárias, infinitesimais e arbitrárias das posições
das partículas, que não violam as restrições cinemáticas impostas pelas
superfícies de restrição. Esta última condição implica que os deslocamentos
virtuais devem ter a direção tangente à superfície de restrição, na posição
instantaneamente ocupada pela partícula. Assim, conforme ilustra a Figura 7.1,
os vetores if
e ir
são mutuamente perpendiculares.
• são admitidos ocorrer instantânea e simultaneamente, de modo que aos
deslocamentos virtuais não associamos nenhum lapso de tempo finito. Desta
forma, as forças aplicadas às partículas do sistema não variam em decorrência
dos deslocamentos virtuais.
Em termos de coordenadas cartesianas, os vetores posição das partículas
são dados por:
kZjYiXr iiii
, i=1, 2, ..., n
de modo que os deslocamentos virtuais podem ser expressos segundo:
kZjYiXr iiii
, i=1, 2, ..., n
onde iii ZYX ,, são entendidos como variações dadas às coordenadas que
representam a posição das partículas em relação ao sistema de referência
empregado.
Em virtude de (4.49), podemos afirmar que, na situação de equilíbrio
dinâmico, o trabalho realizado por todas as forças aplicadas sobre a partícula iP ,
incluindo a força de inércia, associado ao deslocamento virtual ir
é nulo, ou seja:
0 iiiii rrmfF i=1, 2, ..., n (4.50)
Do fato dos vetores if
e ir
serem mutuamente perpendiculares resulta
que as forças de restrição produzem trabalho virtual nulo ( 0 ii rf
) e (4.50) fica:
0 iiii rrmF i=1, 2, ..., n (4.51)
Adicionando as n equações (4.51), obtemos:
D.A. Rade Princípios Variacionais da Mecânica
79
0
1
n
i
iiii rrmF
(4.52)
ou ainda:
0 IF WW (4.53)
onde:
n
i
ii
F rFW
1
e IW
n
i
iii rrm
1
(4.54)
desingam o trabalho virtual das forças impostas e o trabalho virtual das forças de
inércia, respectivamente.
A equação (4.53) expressa o Princípio do Trabalho Virtual, o qual
estabelece que para qualquer posição de um sistema de partículas, o trabalho
virtual total realizado por todas as forças impostas e todas as forças de inércia
resulta nulo para todo e qualquer conjunto arbitrário de deslocamentos virtuais
introduzidos a partir daquela posição.
4.5 – Princípio de Hamilton
No desenvolvimento
que segue, ao operador , que designa uma variação
virtual, serão atribuídas as mesmas propriedades do operador diferencial do
Cálculo Diferencial.
A fim de desenvolver (4.52), introduzimos a identidade:
iiiiiiiiii rrrrrrrrdt rrd
2
1
donde:
iiiiii rrdt rrdrr
2
1 (4.55)
Introduzindo (4.55) em (4.52) e desenvolvendo a equação resultante,
escrevemos:
n
i
ii rF
1
0
2
1
11
n
i
iii
n
i
ii
i rrmdt
rrdm
(4.56)
Em (4.56), reconhecemos:
D.A. Rade Princípios Variacionais da Mecânica
80
n
i
iii rrmT
12
1 (energia cinética do sistema de partículas), (4.57)
de modo que, considerando (4.54), podemos reescrever (4.56) sob a forma:
TW F
n
i
ii
i dt
rrd
m
1
(4.57)
Observamos que, no caso geral, a energia cinética pode ser expressa como
uma função de várias variáveis. Em termos das componentes cartesianas dos
vetores posição e velocidade das partículas do sistema, tal função assume a
forma:
nnnnnn zyxzyxzyxzyxzyxzyxTT ,,,,,,,,,,,,,,,,,,, 222111222111 , (4.58)
de modo que T designa o diferencial total da função T, expresso segundo:
n
n
n
n
n
n
z
z
Ty
y
Tx
x
Tz
z
Ty
y
Tx
x
TT
1
1
1
1
1
1
n
n
n
n
n
n
z
z
Ty
y
Tx
x
Tz
z
Ty
y
Tx
x
T
1
1
1
1
1
1
(4.59)
Multiplicando (4.57) por dt e integrando entre dois instantes de tempo 1t e
2t , temos:
dtTWtt F21
2
1
2
1 11
t
t
n
i
iii
t
t
n
i
ii
i rrmdtdt
rrd
m
(4.60)
Neste ponto, admitiremos que os deslocamentos virtuais, embora
arbitrários, sejam nulos nos instantes 1t e 2t ., ou seja, 021 trtr ii ,
i=1,2,...,n, de modo que, da equação (4.60), resulta:
02
1
dtTWtt F (4.61)
Admitindo que as forças impostas sejam todas conservativas, podemos
escrever:
VW F (4.62)
onde V é a função potencial ou energia potencial associada às forças impostas,
que depende exclusivamente das posições ocupadas pelas partículas do sistema.
Em coordenadas cartesianas, esta função é expressa segundo:
D.A. Rade Princípios Variacionais da Mecânica
81
nnn z,y,x,,z,y,x,z,y,xVV 222111 (4.63)
de modo que V é interpretado como o diferencial total da função V, dado por:
n
n
n
n
n
n
z
z
Vy
y
Vx
x
Vz
z
Vy
y
Vx
x
VV
1
1
1
1
1
1
(4.64)
Assim, introduzindo (4.62) em (4.61), escrevemos:
2
1
0
t
t
dtVT
ou ainda:
0
2
1
t
t
dtL (4.65)
onde:
VTL (4.66)
é o chamado Lagrangeano do sistema.
Levando em conta (4.59) e (4.64), observamos que o Lagrangeano assume a
forma de uma função de várias variáveis, em termos das componentes dos vetores
posição e velocidade da partículas. Assim, escrevemos:
nnnnnn zyxzyxzyxzyxzyxzyxLL ,,,,,,,,,,,,,,,,,,, 222111222111 , (4.67)
e,
n
n
n
n
n
n
z
z
Ly
y
Lx
x
Lz
z
Ly
y
Lx
x
LL
1
1
1
1
1
1
n
n
n
n
n
n
z
z
Ly
y
Lx
x
Lz
z
Ly
y
Lx
x
L
1
1
1
1
1
1
(4.68)
Podemos permutar os operadores integração e diferenciação, expressando
então (4.65) sob a seguinte forma:
2
1
0
t
t
dtL (4.69)
D.A. Rade Princípios Variacionais da Mecânica
82
A equação (4.69) expressam o Princípio Variacional de Hamilton, que se
aplica aos casos em que todas as forças impostas são conservativas. Este princípio
estabelece que dentre todos os movimentos possíveis de ser realizados pelo
sistema mecânico, satisfazendo as restrições cinemáticas impostas, entre dois
instantes de tempo subseqüentes quaisquer t1 e t2, o movimento realmente
desenvolvido pelo sistema será aquele que corresponde a um ponto estacionário
do funcional:
2
1
t
t
dtLI (4.70)
Conforme ilustrado no exemplo a seguir, a imposição das condições
necessárias para a ocorrência de um ponto crítico, estudadas no âmbito do
Cálculo Variacional, abordado no Capítulo 3, conduz às equações do movimento
do sistema mecânico, que são determinadas pelas equações de Euler-Lagrange
associadas ao funcional expressso em (4.70).
4.6 O Princípio de Hamilton Estendido
Quando, além das forças impostas conservativas, houver forças impostas
não conservativas, podemos, em (4.61), separar os trabalhos realizados por estes
dois tipos de forças, escrevendo:
FncFcF WWW (4.71)
onde FncFc WW e designam, os trabalhos realizados pelas forças impostas
conservativas e não conservativas, respectivamente.
Fazendo uso da relação (4.62), para as forças impostas conservativas,
escrevemos:
VW Fc (4.72)
Associando (4.71) e (4.72) com (4.61), escrevemos:
02
1
2
1
dtWdtL tt Fnctt (4.73)
A equação (7.27) traduz o chamado Princípio de Hamilton Estendido,
aplicável aos casos em que existem forças impostas não conservativas.
É importante observar que, no lado esquerdo de (4.73), o termo FncW , que
representa o trabalho virtual das forças não conservativas não pode ser expresso
sob a forma de diferencial total de uma função de várias variáveis. Em
conseqüência (4.73) não pode ser interpretada como a condição de ocorrência de
um ponto crítico de uma função de várias variáveis, o que significa que,
estritamente, o Princípio de Hamilton Estendido não pode ser considerado como
D.A. Rade Princípios Variacionais da Mecânica
83
um princípio variacional. Contudo, (4.73) pode ser utilizada para obter as
equações do movimento, conforme mostraremos com o auxílio do exemplo
apresentado a seguir.
Uma outra observação importante a ser feita é que, embora tenham sido
desenvolvidos para sistemas de partículas, o Princípio de Hamilton e o Princípio
de Hamilton Estendido por ser imediatamente aplicáveis a sistemas mecânicos
formados por corpos rígidos e deformáveis, desde que sejam corretamente
computadas suas energias cinética e potencial e os trabalhos virtuais das forças e
momentos aplicadas a estes corpos.
Exemplo 4.1: Sabendo que no sistema ilustrado a mola está indeformada quando
r=ro, obteremos as equações diferenciais do movimento. Na figura abaixo, t
representa um torque externo, m e mo designam, respectivamente, a massa do
cursor e da barra uniforme OA e k é a constante de rigidez da mola. Obtenhamos
as equações do movimento do sistema em termos das grandezas r
e definidas
na figura.
Considerando o torque externo t como um esforço não conservativo,
utilizaremos o Princípio de Hamilton Estendido, dado pela equação (7.27),
repetida abaixo:
02
1
2
1
dtWdtL tt Fnctt (a)
mo
m
L
O
A
r
k
r
t Ref.
D.A. Rade Princípios Variacionais da Mecânica
84
Em termos das coordenadas radial-transversal indicadas na figura,
escrevemos as seguintes expressões para a energia cinética e potencial do
sistema:
barracursor TTT
onde:
2222 2121 rrmvvmT rcursor
222
6
1
2
1 LmJT oobarra (rotação não baricêntrica)
(elást.)(gravit.)(gravit.) molabarracursor VVVV
onde:
cosmgrVcursor
cos
2
LgmV obarra
2
2
1
omola rrkV
Com base nas equações acima, escrevemos o Lagrangeano sob a forma:
2221 rrmL 2261 Lmo cosmgr cos2Lgmo 221 orrk (b)
O trabalho virtual do torque não conservativo é dado por:
FncW (c)
Com o entendimento que o operador funciona como o diferencial total,
escrevemos:
LLr
r
Lr
r
LL (d)
D.A. Rade Princípios Variacionais da Mecânica
85
A partir de (b), avaliamos as derivadas parciais indicadas na equação
acima:
2mr
r
L
cosmgrrk o rmr
L
sen2sen
LgmmgrL o
22
3
1 LmmrL o
Introduzindo estas derivadas em (d), escrevemos:
rrmrmgrrkmrL o cos2
222
3
1sen
2
sen LmmrLgmmgr oo (e)
Introduzindo (c) e (e) em (a), obtemos:
rrmrmgrrkmrt
t
o
2
1
cos2
0
3
1
2
222
dtLmmrsenLgmsenmgr oo (f)
Desenvolvemos, em seguida, eliminar os diferenciais das variações
temporais das coordenadas que aparecem em (f). Para tanto, efetuamos as
seguintes integrações por partes:
• 2
1
2
1
2
1
t
t
t
t
t
t dtrrrrmdtrrm
•
2
1
2
1
2222
3
1
3
1t
t
t
t
oo LmmrdtLmmr
dtLmmrrmrt
t o
12 22 3
12
Introduzindo as duas últimas equações em (f), após rearranjo, escrevemos:
2
1
2
t
t
o rrmcosmgrrkmr
dtLmmrrmrsenLgmsenmgr oo 22 3122
D.A. Rade Princípios Variacionais da Mecânica
86
0
3
1 2
1
22
t
t
o Lmmrrrm (g)
Relembrando que e r devem ser nulas em t=t1 e t=t2, o último termo do
lado esquerdo da equação acima resulta nulo. Além disso, como a equação
resultante deve ser satisfeita para quaisquer incrementos arbitrários e r , (g) é
satisfeita se forem nulas as funções que multiplicam estes incrementos. Deste
fato, resultam as equações diferenciais do movimento do sistema:
02 cosmgrrkrrm o (h)
tgLmmrrmrLmmr oo
sen
2
2
3
1 22 (i)
Observemos que as equações do movimento são acopladas e não lineares. A
partir de um conjunto de condições iniciais e dada a função que define o torque
externo aplicado t , estas equações podem ser integradas numericamente (ver
Seção 3.8) para obtenção das funções que tr e t que definem a configuração
do sistema em função do tempo.
4.7 Número de graus de liberdade e coordenadas generalizadas
Considerando o sistema de n partículas mostrado na Figura 4.4, podemos
afirmar que, se todas as partículas puderem se mover livre e independentemente
uma das outras, a posição de cada uma delas fica determinada por três
coordenadas espaciais iii zyx ,, i=1,2,...,n, de modo que, para a completa
determinação da configuração espacial do sistema são necessárias
3n coordenadas. Contudo, se houver restrições cinemáticas que limitam o
movimento das partículas a uma dada região do espaço ou impedem algum tipo
de movimento relativo entre as partículas, algumas destas 3n coordenadas
tornam-se dependentes entre si. Neste caso, o número de coordenadas
independentes necessárias e suficientes para determinar a configuração espacial
do sistema é dado por:
N = 3n – p (4.74)
onde p é o número de equações de restrição impostas pelas restrições cinemáticas
e N é o denominado número de graus de liberdade do sistema.
No caso de sistemas formados por corpos rígidos, o número de coordenadas
necessárias para definir a posição espacial de cada corpo é 6 (3 coordenadas
lineares e 3 coordenadas angulares), de modo que o número de graus de liberdade
de um sistema formado por n corpos rígidos é dado por:
N = 6n – p (4.75)
D.A. Rade Princípios Variacionais da Mecânica
87
A escolha do conjunto de coordenadas independentes é arbitrária, embora
seu número permaneça invariável. Qualquer conjunto de N coordenadas
independentes constitui um conjunto de coordenadas generalizadas, que serão
aqui denotadas por Nqqq ,,, 21 . Ilustramos os conceitos de número de graus de
liberdade e coordenadas generalizadas com o auxílio do seguinte exemplo.
Exemplo 4.2: No sistema ilustrado na figura abaixo, as partículas P1 e P2 são
forçadas a se mover sobre o plano x-y, estando ligadas pela barra rígida de
comprimento L. Para este sistema, determinemos o número de graus de liberdade
e um conjunto de coordenadas generalizadas.
No caso em questão temos um sistema com duas partículas (n=2) e número
toal de coordenadas é 3n=6. Em termos das componentes de seus vetores posição
em relação ao sistema de eixos indicados na figura, estas seis coordenadas são
denotas por 222111 z,y,x,z,y,x . No entanto, estas seis coordenadas não são todas
indepentes, uma vez que as seguintes as equações de restrição devem ser
satisfeitas.
z1=0 z2=0 2212212 Lyyxx
P1
P2
L
y
x
z
D.A. Rade Princípios Variacionais da Mecânica
88
Temos então p = 3 (três equações de restrição), de modo que, de acordo com
a equação (4.74), o número de graus de liberdade do sistema resulta ser
N=323=3.
Alguns conjuntos possíveis de coordenadas generalizadas são:
( 11 xq , 12 yq , 23 xq )
( 11 xq , 12 yq , 13 yq )
( 11 xq , 12 yq , 3q )
4.8 Equações de Lagrange
Nesta seção, desenvolveremos a formulação que, partindo do Princípio de
Hamilton Estendido, conduz às chamadas Equações de Lagrange, as quais
constituem uma forma bastante elegante e expedita para a obtenção das equações
do movimento de sistemas dinâmicos.
As equações
de Lagrange são formuladas em termos de coordenadas
generalizadas, cujo conceito foi introduzido na Seção 4.7.
Para o conjunto de n partículas mostrado na Figura 4.4, podemos sempre
expressar os vetores posição das partículas em função de um conjunto
previamente escolhido de N coordenadas generalizadas através de relações do
tipo:
tqqqrr Nii ,,,, 21 i=1,2, ..., n (4.76)
Em termos de coordenadas cartesianas, a transformação (4.76) pode ser
detalhada como segue:
kzjyixr iiii
(4.77)
com:
tqqqxx Nii ,,, 21
tqqqyy Nii ,,, 21 (4.78)
tqqqzz Nii ,,, 21 i=1,2, ..., n
Buscaremos, em seguida, representar as velocidades das partículas em
termos das coordenadas generalizadas. Para tanto, empregando a regra da cadeia
da derivação, derivamos (4.77) em relação ao tempo, escrevendo:
kzjyixv iiii
(4.79)
D.A. Rade Princípios Variacionais da Mecânica
89
com:
t
x
q
q
x
t
x
q
q
x
q
q
x
q
q
x
dt
dx
x i
N
j
j
j
ii
N
N
iiii
i
1
2
2
1
1
t
y
q
q
y
q
q
y
q
q
y
dt
dy
y iN
N
iiii
i
2
2
1
1
t
yq
q
y iN
j
j
j
i
1
(4.80)
t
z
q
q
z
q
q
z
q
q
z
dt
dz
z iN
N
iiii
i
2
2
1
1
t
z
q
q
z iN
j
j
j
i
1
i=1,2, ..., n
Levando em conta (4.79), a expressão da energia cinética do sistema de
partículas, dada por (4.57), é desenvolvida da seguinte forma:
n
i
iiii
n
i
iii zyxmvvmT
1
222
1 2
1
2
1 (4.81)
Introduzindo (4.80) em (4.81), escrevemos:
2
1
2
1
2
112
1
t
zq
q
z
t
yq
q
y
t
xq
q
xmT i
N
j
j
j
ii
N
j
j
j
ii
N
j
j
j
i
n
i
i
Na equação acima podemos perceber que, de modo geral, a energia cinética
será função das coordenadas generalizadas nqqq ,, 21 , de suas derivadas
temporais nqqq ,, 21 e do tempo t. Assim, escrevemos:
tqqqqqqTT NN ,,,,,, 2121 (4.82)
Da mesma forma, a energia potencial, que é função exclusiva dos vetores
posição (e, eventualmente, do tempo), é expressa sob a forma:
tqqqqqqVV NN ,,,,,, 2121 (4.83)
Com base em (4.82) e (4.83), podemos expressar o Lagrangeano como uma
função das coordenadas generalizadas, das coordenadas generalizadas nqqq ,, 21 , de suas derivadas temporais nqqq ,, 21 e do tempo t:
tqqqqqqLVTL NN ,,,,,, 2121 (4.84)
A partir de (4.84), podemos expressar da seguinte forma a variação do
Lagrangeano associada a um conjunto arbitrário de deslocamentos virtuais:
D.A. Rade Princípios Variacionais da Mecânica
90
N
N
N
N
q
q
Lq
q
Lq
q
Lq
q
Lq
q
Lq
q
LL
2
2
1
1
2
2
1
1
ou:
N
j
j
j
j
j
q
q
Lq
q
LL
1
(4.85)
Notemos que, de acordo com as propriedades atribuídas aos deslocamentos
virtuais (ver Seção 4.4), a variação no tempo, representada pelo termo dttL não
é incluída no desenvolvimento que conduz a (4.85).
No desenvolvimento que segue, buscaremos expressar o trabalho virtual
das forças não conservativas em termos das coordenadas generalizadas. Para
tanto, escrevemos:
n
i
ii
F
nc rFW
1
(4.86)
Fazendo uso da equação (4.77), expressamos os deslocamentos virtuais sob
a forma:
ir N
N
iii q
q
r
q
q
r
q
q
r
2
2
1
1
N
j
j
j
i q
q
r
1
(4.87)
onde os termos jq , j=1,2, ..., N, são entendidos como variações arbitrárias e
independentes introduzidas nas coordenadas generalizadas.
Notemos, mais uma vez que, de acordo com as propriedades atribuídas aos
deslocamentos virtuais na Seção 4.4, a variação no tempo, representada pelo
termo dttri não é incluída no desenvolvimento (4.87).
Introduzindo (4.87) em (4.86), escrevemos:
FncW
N
j
j
n
i j
i
i qq
rF
1 1
(4.88)
ou ainda:
FncW
N
j
jj qQ
1
(4.89)
onde:
D.A. Rade Princípios Variacionais da Mecânica
91
n
i j
i
ij q
rFQ
1
(4.90)
são os denominados esforços generalizados.
Voltemos agora ao Princípio de Hamilton Estendido, expresso pela equação
(4.73), repetida abaixo:
02
1
2
1
dtWdtL tt Fnctt (4.91)
Introduzindo (7.85) e (4.89) em (4.91), escrevemos:
dtq
q
Lq
q
Lt
t
N
j
j
j
j
j
2
1 1
0
2
1 1
dtqQ
t
t
N
j
jj
Rearrajando:
0
2
1 1 1
dtq
q
LqQ
q
Lt
t
N
j
N
j
j
j
jj
j
(4.92)
Efetuamos em seguida as seguintes integrações por partes:
dtq
q
L
dt
dq
q
Ldtq
q
L t
t
N
j
j
j
t
t
j
j
t
t
N
j
j
j
2
1
2
1
2
1 11
(4.93)
Considerando a hipótese que os deslocamentos virtuais (e, portanto, as
variações correspondentes das coordenadas generalizadas) devem anular-se nos
instantes 1t e 2t , ou seja, 021 tqtq jj , o primeiro termo do lado direito de
(4.93) resulta nulo. Introduzindo então a equação resultante em (4.92),
escrevemos, após alguns rearranjos:
0
2
1 1
dtqQ
q
L
q
L
dt
dt
t
N
j
jj
jj
(4.94)
Uma vez que a equação acima deve ser satisfeita para todo e qualquer
conjunto de variações virtuais arbitrárias e independentes, jq , j=1,2, ..., N,
resulta que as funções que multiplicam cada uma das variações jq em (4.94)
devem anular-se identicamente. Assim, temos:
j
jj
Q
q
L
q
L
dt
d
j=1,2, ..., N (4.95)
D.A. Rade Princípios Variacionais da Mecânica
92
As equações (4.95) são as denominadas Equações de Lagrange do
movimento. Note-se que tais equações apresentam-se sob a forma de equações
diferenciais de segunda ordem na variável tempo, que constituem as equações do
movimento do sistema mecânico.
Devemos também observar que, embora tenham sido deduzidas para
sistemas de partículas,
as Equações de Lagrange na forma (4.95) podem ser
imediatamente aplicadas a sistemas constituídos de corpos rígidos e deformáveis,
bastando que se usem as expressões adequadas para o cálculo das energias
cinética e potencial que compõem o Lagrangeano.
O uso das equações de Lagrange para a obtenção das equações do
movimento é ilustrado no exemplo a seguir.
Exemplo 4.3: Considerando o sistema mecânico apresentado do Exemplo 4.1,
obtenhamos as equações diferenciais do movimento utilizando as Equações de
Lagrange.
Como as coordenadas ,r são independentes, podemos adotá-las como
coordenadas generalizadas e as duas equações de Lagrange do movimento,
expressas por (4.95) tomam as formas:
rQr
L
r
L
dt
d
(i)
Q
LL
dt
d
(ii)
Sendo o trabalho virtual do torque não conservativo, dado por:
FncW
e desenvolvendo (4.89) considerando a escolha feita para as coordenadas
generalizadas:
FncW QrQr ,
constatamos que:
0rQ Q (iii)
Efetuando as derivadas parciais indicadas em (i) e (ii), a partir da
expressão do Lagrangeano dada pela Eq. (a) do Exemplo 4.1, obtemos:
2221 rrmL 2261 Lmo cosmgr cos2Lgmo 221 orrk (b)
D.A. Rade Princípios Variacionais da Mecânica
93
r
L
rm rmr
L
dt
d
r
L orrkmgmr cos2
L 22
3
1 Lmmr o
L
dt
d 22
3
12 Lmmrrmr o
L senLgmmgrsen o 2 (iv)
Introduzindo as equações (iii) e (iv) em (i) e (ii), obtemos, as seguintes
equações do movimento, que são idênticas às equações (h) e (i) obtidas no
Exemplo 4.1:
0cos2 mgrrkrrm o (h)
tgsenLmmrrmrLmmr oo
2
2
3
1 22 (i)
4.9 Método de Rayleigh-Ritz
Conforme vimos nas seções anteriores deste capítulo, os métodos
variacionais conduzem às equações diferenciais de equilíbrio, no caso de
problemas estáticos, ou equações diferenciais do movimento, no caso de
problemas dinâmicos. Em diversas situações, notadamente no caso de sistemas
mecânicos complexos, as equações de equilíbrio ou do movimento podem ser
obtidas de forma mais simples e elegante, em comparação com os procedimentos
baseados nos princípios da Mecânica Newtoniana. Entretanto, a resolução do
problema completo requer a resolução destas equações, a qual pode ser feita
através de métodos analíticos apenas nos casos mais simples. Este fato justifica o
uso de métodos que forneçam soluções aproximadas dos problemas estudados,
conhecidos como métodos diretos do Cálculo Variacional.
O método de Rayleigh-Ritz, também conhecido como Método de Ritz é um
método destinado à obtenção de soluções aproximadas de problemas regidos por
princípios variacionais. Ilustraremos seus princípios considerando um problema
regido por um funcional do tipo:
dxxu,,xu,xu,xFuI x
x
n 2
1
(4.96)
com as condições de contorno:
11 uxu ; 22 uxu
D.A. Rade Princípios Variacionais da Mecânica
94
11 uxu ; 22 uxu (4.97)
1111 nn uxu ; 1221 nn uxu
Diante da impossibilidade (ou não conveniência) de se encontrar a solução
exata da função u x que extremiza o funcional (4.96), propõe-se uma solução
aproximada expressa da seguinte forma:
1
n
i i
i
u x c x
(4.98)
As funções , 1, 2, ,i x i n são escolhidas arbitrariamente e devem ser
tais que pelo menos as condições de contorno geométricas do problema sejam
satisfeitas por u x . Estas funções são denominadas funções de comparação (trial
functions). Note-se que, na aproximação (4.98), permanecem desconhecidas as
constantes , 1, 2, ,ic i n .
Introduzindo (4.98) em (4.96), e efetuando as operações de derivação e
integração indicadas, obtém-se o funcional sob a forma:
1 2, , , nI u I c c c
Observa-se que este funcional resulta expresso como uma função escalar de
n variáveis, para o qual as condições de estacionariedade são:
1 2
0; 0; ; 0
n
I u I u I u
c c c
(4.99)
As equações (4.99) representam um sistema de n equações algébricas com n
incógnitas que, uma vez resolvido, fornece os valores de , 1, 2, ,ic i n .
Introduzindo estes valores em (4.98) tem-se a solução aproximada do problema.
Conforme veremos no próximo capítulo, o Método de Rayleigh-Ritz é a base
para a formulação do Método de Elementos Finitos pelo processo variacional.
D.A. Rade Princípios Variacionais da Mecânica
95
Exemplo 4.4. Ilustraremos o método através do problema de flexão de uma viga
de Euler-Bernoulli, uniforme, de módulo de rigidez à flexão EI, submetida a um
carregamento transversal uniformemente distribuído q, conforme ilustrado na
Figura 4.5. Desejamos obter uma aproximação para o campo de deslocamentos
transversais da viga, indicado por y(x).
Figura 4.5
O problema de equilíbrio é regido pelo Princípio da Energia Potencial
Mínima, que, para o problema em questão, é formulado sob a seguinte forma
variacional:
22
2
0
1 0
2
L d yI y EI qy dx
dx
(a)
Vamos admitir uma solução aproximada do problema sob a forma:
1 2 2x xy x a sen a senL L
(b)
ou
1 1 2 2y x a x a x (c)
onde:
1 xx sen L
2 2 xx sen L
são denominadas funções de Ritz. Estas funções são escolhidas de modo que y x
satisfaça às condições de contorno geométricas do problema.
0 0y y L
y(x)
L
q
EI
x
D.A. Rade Princípios Variacionais da Mecânica
96
Visando formular a condição de estacionariedade (a), vamos primeiramente
avaliar:
22
2
0 2
L EI d yI y qy dx
dx
22 2
1 2 1 2
0
2 2 2
2
L EI x x x xa sen a sen q a sen a sen dx
L L L L L L
Desenvolvendo a expressão acima e efetuando as integrações, temos:
4 42 2 2 21 2
0
2 2
2
L EI x xI y a sen a sen
L L L L
2 2
1 2 1 2
2 2 22 x x x xa a sen sen qa sen qa sen dx
L L L L L L
Efetuando as integrações indicadas na equação acima, obtemos:
4 42 2 11 22 22 2 2 2
a LEI L EI LI y a a q
L L
(d)
Devemos agora impor a condição de estacionariedade do funcional I.
Observamos em (4.99), que ele é função de duas incógnitas a1 e a2. Assim, a
condição de estacionariedade conduz a:
1 2
1 2
0I II a a
a a
1 0
I
a
; 2
0I
a
1
0I
a
4
1 5
4 qla
EI
2
0I
a
2 0a
Assim, o campo de deflexões aproximado fica:
454 qL xy x senEI L
(4.100)
As soluções exata e aproximada para o deslocamento no ponto médio da
viga são:
D.A. Rade Princípios Variacionais da Mecânica
97
45
2 384
exato L qLy
EI
4
.
5
4
2
aprox L qLy
EI
O erro da aproximação é de apenas 0,39 %, o que mostra a boa precisão.
4.9 Bibliografia
FOWLES, G.R. CASSIDAY, G.L., Analytical Mechanics, 6th Edition, Harcourt
Brace College Publishers, 1998.
WASHIZU, K., Variational Methods in Elasticity and Plasticity (Monographs
in Aeronautics & Astronautics), 2nd Ed., Elsevier, 1975.
LEMOS, N., Mecânica Analítica, 2ª Edição, Editora Livraria da Física, 2007.
LANCZOS, C., The Variational Principles of Mechanics, 4th Edition,
University of Toronto Press, 1977.
D.A. Rade Formulação do MEF pelo Processo Variacional
98
CAPÍTULO 5
FORMULAÇÃO DO MEF PELO PROCESSO VARIACIONAL. APLICAÇÃO
À MODELAGEM DO COMPORTAMENTO DINÂMICO DE VIGAS DE
EULER –BERNOULLI
5.1 – Introdução
Neste capítulo, o método dos elementos finitos é formulado para obtenção das
equações do movimento para vigas de Euler-Bernoulli, utilizando o processo
variacional enfocado no Capítulo 3.
São admitidas as seguintes hipóteses:
1ª. Consideram-se exclusivamente os efeitos da flexão transversal, não sendo
considerados os efeitos associados ao movimento longitudinal;
2ª. São adotadas as simplificações da teoria de vigas de Euler-Bernoulli:
desprezam-se os efeitos do cisalhamento transversal e de inércia de rotação das
seções transversais;
3ª. Os deslocamentos e rotações são considerados pequenos e os materiais
apresentam comportamento linear elástico linear isotrópico.
4ª. Desprezam-se os efeitos dissipativos, isto é, considera-se o sistema sem
amortecimento.
Será considerada a viga mostrada na Figura 5.1 de comprimento e com
condições de contorno arbitrárias, para a qual se admitem as hipóteses da teoria de
Euler-Bernoulli, sendo m x (kg/m) e E(x) (N.m-2) a densidade linear e módulo de
elasticidade longitudinal do material que constitui a viga, respectivamente, I x
(m4) é o momento de inércia de área da seção transversal da viga em relação ao eixo
baricêntrico z, q(x,t) (N.m-1) é o carregamento transversal distribuído aplicado
externamente e v(x,t) (m) indica o campo de deslocamentos transversais da viga.
D.A. Rade Formulação do MEF pelo Processo Variacional
99
Figura 5.1
Serão a seguir detalhadas as etapas típicas de uma análise por elementos
finitos, a saber:
1ª. Obtenção das equações do movimento em nível elementar, empregando as
o Princípio Variaciona de Hamilton;
2ª. Obtenção das equações do movimento em nível global;
3ª. Imposição das condições de contorno;
4ª. Resolução das equações do movimento para obtenção de diferentes tipos de
respostas dinâmicas.
5.2 – Obtenção das equações do movimento em nível elementar
A viga mostrada na Figura 5.1 é admitida dividida em um certo número
arbitrário de elementos cada um deles delimitado por dois nós. A Figura 5.2 mostra
um destes elementos, identificado genericamente pelo índice i, contendo dois nós iL e
iR , cujas propriedades relevantes são indicadas por i , im x , iE x , iI x . Além disso
,iq x t e ,iv x t indicam o carregamento externo e o campo de deslocamentos
transversais no interior do elemento ( 0 ix ). As forças transversais ,L Ri iV t V t e
os momentos ,L Ri iM t M t representam as ações mecânicas exercidas pelos
elementos vizinhos, à esquerda e à direita, sobre o elemento considerado. Também
são indicados os graus de liberdade escolhidos para o elemento, que são constituídos
pelos valores dos deslocamentos transversais nos nós, indicados por Liv t e Riv t , e
v(x,t)
q(x,t)
y
x
L
D.A. Rade Formulação do MEF pelo Processo Variacional
100
pelas rotações das seções transversais nos nós, indicadas por Li t e Ri t . Estes
graus de liberdade são agrupados no vetor:
Te L L R Ri i i i iu t v t t v t t (5.1)
onde o superscrito (e) designa os graus de liberdade em nível elementar.
Deve-se observar que, de acordo com a teoria de Euler-Bernoulli, vale a
seguinte relação entre as rotações das seções transversais e os deslocamentos
transversais da viga:
,, ii v x tx t x
(5.2)
Figura 5.2
Os deslocamentos transversais são interpolados pelo seguinte polinômio
cúbico:
2 31 2 3 4, 0i iv x t c t c t x c t x c t x x (5.3)
Em conseqüência, levando em consideração a Eq. (5.2), as rotações das seções
transversais ficam aproximadas pela seguinte função quadrática:
22 3 4, 2 3 0i ix t c t c t x c t x x (5.4)
As equações (5.3) e (5.4) devem satisfazer as condições cinemáticas:
0, Ei iv t v t , Di i iv t v t (5.5.a)
y
x
Ei Di
Riv t
Ri t Li t
Liv t
Ei(x), Ii(x), i(x)
,i x t ,iv x t
qi(x)
RiV t LiV t
LiM t RiM t
D.A. Rade Formulação do MEF pelo Processo Variacional
101
0
Ei
i
x
v t
x
i
Di
i
x
v t
x
(5.5.b)
Impondo as condições (5.5) às aproximações (5.3) e (5.4), são obtidas quatro
equações algébricas lineares que, uma vez resolvidas, fornecem as expressões das
constantes ic t em função dos graus de liberdade que formam o vetor dado por (5.1)
e o campo de deslocamentos transversais resulta expresso da seguinte forma:
1 2 3 4, E E D Di i i i iv x t x v t x t x v t x t (5.6.a)
ou:
, ei iv x t x u t (5.6.b)
com:
1 2 3 4x x x x x (5.7)
2 31 1 3 2x xx l l
(5.8.a)
2 32 2 x xx x l ll l
(5.8.b)
2 33 3 2x xx l l
(5.8.c)
2 34 x xx l ll l
(5.8.d)
As funções , i = 1 a 4 são chamadas funções de forma ou funções de
interpolação. São mostradas graficamente na Figura 5.3.
D.A. Rade Formulação do MEF pelo Processo Variacional
102
Figura 5.3
As equações (5.6) mostram que o campo de deslocamentos transversais e
rotações das seções transversais são expressos como combinações lineares dos
deslocamentos e rotações nodais, sendo os coeficientes de combinação linear dados
pelas funções de forma.
Visando obter as equações do movimento em nível elementar através do
Princípio Variacional de Hamilton, devemos formular as energias cinética e
potencial, além do trabalho virtual do carregamento externo, conforme
desenvolvimento a seguir.
Energia Cinética
A energia cinética do elemento é dada pela expressão:
2
0
,1
2
i
i
i i
v x t
T m x dx
t
(5.9)
1(x)
x
1
4(x)
x
45o
45o
3(x)
x
1
2(x)
x
i i
i i
D.A. Rade Formulação do MEF pelo Processo Variacional
103
Introduzindo (5.6.b) em (5.9) após derivação de iv t em relação ao tempo,
escrevemos:
12 Te e ei i i iT u t M u t (5.10)
onde:
4 4
0
i
Te
i iM m x x x dx
(5.11)
é a chamada matriz de massa ou matriz de inércia do elemento de viga de Euler-
Bernoulli, cujo elemento geral é calculado segundo:
0
i
e
i i m nmn
M m x x x dx
, , 1a 4m n (5.12)
Para o caso em que im x é constante, as equações (5.11) e (5.12) conduzem à
seguinte matriz:
2 2
2
156 22 54 13
4 13 3
156 22420
4
i i
e i i ii i
i
i
i
mM
simétrica
(5.13)
Energia Potencial
Para a viga de Euler-Bernoulli, a energia potencial, incluindo a energia de
deformação e os trabalhos das forças e momentos externos aplicados, é expressa
segundo:
2
2
2
0
,1
2
i
i
i i i
v x t
V E x I x dx
x
0
, ,i i iq x t v x t dx
E E E E D D E Di i i i i i i iV t u t M t t V t u t M t t
D.A. Rade Formulação do MEF pelo Processo Variacional
104
Introduzindo (5.6.b) em (5.14), após derivação dupla de ,iv x t em relação a x,
escrevemos:
12 Te e ei i i iV u t K u t Te ei iu t Q t Te ei iu t Q t
(5.15)
onde:
4 4
0
i
Te
i i iK E x I x x x dx
, , 1a 4m n (5.16)
é a chamada matriz de rigidez do elemento de viga de Euler-Bernoulli, cujo
elemento geral é calculado segundo:
0
i
e
i i i m nmn
K E x I x x x dx
, , 1a 4m n , (5.17)
Além disso:
0 ,i Tei iQ t x q x t dx (5.18)
é o vetor de esforços generalizados em nível elementar, cujo elemento geral é dado
por:
0 ,iei m imQ t x q x t dx , 1a 4m , (5.19.a)
e
Te E E D Di i i i imQ t V t t V t t (5.19.b)
Para o caso em que o módulo de rigidez à flexão e o carregamento distribuído
são independentes de x no interior do elemento: i i i iE x I x E I , ,i iq x t q t , as
equações (5.16) e (5.18) conduzem às seguintes expressões:
D.A. Rade Formulação do MEF pelo Processo Variacional
105
2 2
3
2
12 6 12 6
4 6 2
12 6
4
i i
e i i ii i
i
i i
i
E IK
simétrica
(5.20)
2
2
2
12
2
12
i
i
i
i
e
i
i
i
i
i
q t
q t
Q t
q t
q t
(5.21)
O Princípio Variacional de Hamilton estabelece:
2
1
0
t
i it
T V dt (5.22)
Combinando as equações (5.10), (5.15) e (5.22), temos:
2
1
1 1 0
2 2
T T T Tt e e e e e e e e e e
i i i i i i i i i it
u t M u t u t K u t u t Q t u t Q t dt
Na equação acima reconhecemos um funcional do tipo dado pela equação
(3.51), no qual a variável independente é o tempo, as variáveis dependentes são os
graus de liberdade elementares que formam o vetor eiu t e as derivadas de mais
alta ordem são as derivadas de primeira ordem em relação ao tempo.
As equações de Euler-Lagrange associadas ao funcional, dadas por (3.63),
podem ser expressas, de forma alternativa, segundo:
0i ie ei iL Lddtu u
(5.23)
sendo i i iL T V o Lagrangeano.
Efetuando as operações indicadas em (5.23), obtemos as seguintes equações
do movimento em nível elementar:
D.A. Rade Formulação do MEF pelo Processo Variacional
106
( ) ( )e e e ee ei i i i i iM u t K u t Q t Q t (5.24)
5.2 Montagem das equações do movimento em nível global
Na seção anterior, as equações do movimento foram obtidas para um
elemento genérico considerado isolado. Entretanto, deve ser considerado que a viga
é discretizada em vários elementos de modo que elementos vizinhos compartilham
nós, conforme ilustrado na Figura 5.4 para dois elementos contíguos i e j.
Te L L R Ri i i i iu t v t t v t t Te L L R Rj j j j ju t v t t v t t
1 1 2 2 3 3 Tgu t v t t v t t v t t
Figura 5.4
graus de liberdade elementares
Li Ri
Riv t
Ri t Li t
Liv t
Lj Rj
Rjv t
Rj t Lj t
Ljv t
3
3v t
3 t 1 2
2v t
2 t 1 t
1v t
graus de liberdade globais
D.A. Rade Formulação do MEF pelo Processo Variacional
107
As equações do movimento dos dois elementos devem ser combinadas de
modo que seja assegurado o equilíbrio dinâmico e a continuidade de deslocamentos e
rotações nos nós que são compartilhados dois elementos vizinhos. Este
procedimento será ilustrado a seguir para os dois elementos mostrados na Figura
5.4.
Estabelecem-se, para cada elemento, transformações lineares relacionando os
graus de liberdade elementares e os graus de liberdade globais, sob a forma:
4 64 1 6 1e gi iu t L u t (5.25.a)
4 64 1 6 1e gj ju t L u t (5.25.b)
O esquema da Figura 5.4 permite facilmente construir as matrizes iL e
jL , conforme mostrado abaixo. Nota-se que são matrizes booleanas (formadas por
zeros e uns) nas quais as posições dos uns permitem localizar os graus de liberdade
elementares dentre os graus de liberdade globais.
Elemento i
1
1
2
2
3
3
1 0 0 0 0 0
0 1 0 0 0 0
0 0 1 0 0
0
0 0 0 1 0 0
L
i
L
i
R
i
R
i
v t
v t t
t v t
v t t
t v t
t
(5.26.a)
Elemento j
1
1
2
2
3
3
0 0 1 0 0 0
0 0 0 1 0 0
0 0 0 0 1 0
0 0 0 0 0 1
L
j
L
j
R
j
R
j
v t
v t t
t v t
v t t
t v t
t
(5.26.b)
Para o sistema formado pelos elementos i e j escreve-se a energia cinética
total sob a forma:
D.A. Rade Formulação do MEF pelo Processo Variacional
108
i jT T T (5.27)
Utilizando a Eq. (5.10) para ambos os elementos e as equações (5.25), obtém-
se a energia cinética total expressa sob a seguinte forma, em termos dos graus de
liberdade globais:
12 Tg g gT u t M u t (5.28)
onde a matriz de massa global é expressa segundo:
6 6
TTg e e
i i i j j jM L M L L M L
(5.39)
De modo semelhante, opera-se sobre a energia potencial total e o trabalho
virtual total do carregamento externo:
i jV V V = 12 T Tg g g g gu t K u t u t Q t (5.30)
com:
6 6
TTg e e
i i i j j jK L K L L K L
(5.31)
e:
6 1 TTg e ei i j jQ L Q L Q (5.32)
Vale observar que, de acordo com a terceira Lei de Newton, as forças e
momentos exercidos pelo elemento 1 sobre o elemento 2 são iguais aos
correspondentes exercidos pelo elemento 2 sobre o elemento 1. Isso faz com que as
operações que conduzem a (5.30) eliminem tais esforços correspondentes aos nós
intermediários da formulação, permanecendo apenas aqueles correspondentes aos
nós das extremidades da viga. Estes últimos são tratados quando são impostas as
condições de contorno, conforme detalhado na Seção 4.3.
Aplicando novamente o princípio de Hamilton expresso segundo:
2
1
0
t
i it
T V dt (5.33)
obtemos as seguintes equações de Euler-Lagrange associadas:
D.A. Rade Formulação do MEF pelo Processo Variacional
109
gg gi
L d L Q
dtu u
com L T V
(5.37)
Associando as equações (5.28), (5.30) e (5.37), obtêm-se as equações do
movimento em nível global, formado pelos dois elementos finitos conectados:
g g g g gM u t K u t Q t (5.38)
As equações (5.39) constituem um sistema de 6 equações diferenciais
ordinárias no tempo que possuem exatamente a mesma forma das equações do
movimento de sistemas discretos não amortecidos. Desta forma, os métodos
apresentados naquele capítulo podem ser utilizados para calcular os diferentes tipos
de respostas dinâmicas do sistema sujeito a diferentes formas de excitação, que
podem ser completadas por um conjunto de condições iniciais:
00g gu u 00g gu u (5.39)
O procedimento aqui ilustrado para o acoplamento de apenas dois elementos pode
ser aplicado sucessivamente para o acoplamento de um número qualquer de
elementos.
4.3 Imposição das condições do contorno
Ao contrário do que ocorre nos problemas de equilíbrio (problemas estáticos)
os problemas de análise dinâmica podem ser resolvidos sem dificuldades
particulares quando o sistema estrutural estiver livre no espaço, sem nenhuma
vinculação ao sistema de referência fixo. Este é o caso de estruturas aeronáuticas ou
espaciais. Nestes casos, o sistema de equações diferenciais (5.38) pode ser
processado diretamente. Por outro lado, caso haja vinculações cinemáticas
representadas por deslocamentos impostos em alguns pontos da estrutura, as
equações do movimento devem ser devidamente modificadas para levar em conta
tais vinculações. Apresenta-se, a seguir um procedimento para esta finalidade.
As matrizes e vetores que figuram em (5.38) são particionados segundo os
graus de liberdades ditos “livres”, que não são afetados pelos vínculos cinemáticos e
os graus de liberdade denominados “impostos” que são aqueles cujos valores são
prescritos. Os primeiros são doravante indicados por subscritos “L” e os últimos por
subscritos “I” nas equações do movimento, que se expressam segundo:
g g gg g g g
L L LLL LI LL LI
g g g gg g g
IL II IL III I I
u t u t Q tM M K K
M M K Ku t u t Q t
(5.40)
D.A. Rade Formulação do MEF pelo Processo Variacional
110
Os dois blocos de equações resultantes do particionamento indicado em (5.40)
podem ser postos sob a forma:
g g g g g g g g gLL L LL L L LI I LI IM u t K u t Q t M u t K u t (5.41)
g g g g g g g g gI IL L II I IL L II IQ t M u t M u t K u t K u t (5.42)
A resolução das equações (5.41) e (5.42) permite determinar a resposta
dinâmica nas coordenadas livres e, subseqüentemente, os esforços de reação
associados aos graus de liberdade impostos.
4.4 Realização de análise dinâmicas
A partir das equações do movimento globais resolução das equações (5.41) e
(5.42), podem-se computar diferentes tipos resposta dinâmica, cuja formulação é
apresentada a seguir:
Análise Modal
Para o problema homogêneo (respostas em regime livre):
0g g g gLL L LL LM u t K u t (5.43)
procuram-se soluções harmônicas do tipo:
g i tLu t U e (5.44)
Introduzindo (5.44) em (5.43), chega-se ao seguinte problema de autovalor:
2 0g gLL i LL iK M U (5.45)
cujas soluções ,i iU , i=1,2, ...são, respectivamente, as freqüências naturais e os
modos naturais de vibração.
Análise Harmônica
Considera-se o problema não-homogêneo (respostas em regime forçado):
g g g g gLL L LL L LM u t K u t Q t (5.43)
com as forças excitadoras e respostas harmônicas dadas, respectivamente, por:
D.A. Rade Formulação do MEF pelo Processo Variacional
111
g i tL LQ t Q e (5.44)
g i tL Lu t U e (5.45)
Introduzindo (5.44) e (5.45) em (5.43), obtêm-se a seguinte relação entre as
amplitudes de excitação e amplitudes das respostas harmônicas:
2g gLL LL L LK M U Q (5.46)
ou:
12g gL LL LL LU K M Q (5.47)
Análise Transiente
O problema de análise transiente consiste na determinação das respostas
temporais da estrutura sob a ação de forças externas variáveis no tempo e/ou a um
conjunto de condições iniciais (deslocamentos e/ou velocidades) impostas no instante
inicial. Alternativamente, a excitação pode ser provocada por deslocamentos
impostos variáveis no tempo.
Nestes casos, a resolução é feita integrando as equações do movimento (5.41),
utilizando algoritmos de integração numérica, tais como os métodos de Newmark,
Runge-Kutta, Huboldt, dentre outros existentes para esta finalidade. Ao final da
integração, são conhecidos os valores dos graus de liberdade do modelo de elementos
finitos em um conjunto discreto de pontos dentro do intervalo de tempo escolhido.
4.5. Bibliografia
CRAIG JR., R. R. e KURDILA, A. J. Fundamentals of Structural Dynamics. 2.ed.,
New Jersey – USA: John Wiley & Sons, 2006. 744p.
MEIROVITCH, L. Methods of Analytical Dynamics. New York – USA: McGraw-Hill,
1970. 524p.