Disciplina de MECÂNICA E MODELAÇÃO COMPUTACIONAL

Mestrado Integrado em ENGENHARIA BIOMÉDICA

4º Ano, 1º Semestre 2010/11

ANÁLISE DE UMA ESTRUTURA DE BARRAS 2D ATRAVÉS DO

MÉTODO DOS ELEMENTOS FINITOS
Pinto, J.
63772
e-mail: joanasequeirapinto@gmail.com

Semedo, J.
62839
e-mail: joao-semedo@hotmail.com

Palavras-chave: Método dos Elementos Finitos, Barras, Deformação, Tensão, MATLAB

Resumo. Com este trabalho computacional pretende-se desenvolver um procedimento para

resolver problemas de estruturas de barras 2D através do Método dos Elementos Finitos. Com
este propósito, desenvolveu-se um programa no software MATLAB que permite obter os
deslocamentos generalizados, as reacções nos apoios e a tensão em cada barra. O programa
desenvolvido, além de apresentar uma interacção simples e intuitiva com o utilizador, consegue
resolver os problemas propostos em tempos computacionais muitíssimo reduzidos, provando a
enorme relevância deste método na mecânica estrutural.

1
 Pinto, J. e Semedo, J.

1. INTRODUÇÃO

1.1.Apresentação do problema
O objectivo deste trabalho computacional é implementar e aplicar o método dos elementos
finitos (MEF), de forma a resolver o problema apresentado pela estrutura de barras em 2D da
Figura 1. Pretende-se determinar a configuração deformada, as reacções nos apoios, o valor da
deformação no ponto e a tensão em cada uma das barras.
A implementação do MEF foi efectuada com recurso ao MATLAB, tendo como objectivo
criar uma aplicação geral, simples e intuitiva.

Figura 1. Estrutura de barras a utilizar.

2
 Pinto, J. e Semedo, J.

1.2.Método dos elementos finitos

Grande parte dos sistemas estudados pela engenharia são descritos por sistemas de equações
diferenciais parciais, sendo a resolução exacta das mesmas geralmente impossível. Uma das
primeiras abordagens à resolução deste tipo de problemas foi a utilização de métodos
variacionais tradicionais. No entanto estes apresentam sérias limitações, principalmente na
construção das funções interpoladoras, que é realizada de forma não sistemática,
impossibilitando uma implementação computacional eficiente dos mesmos [1].
É precisamente esta a vantagem do MEF, a forma sistemática como são construídas as
funções de interpolação. O MEF baseia-se na representação de um domínio por domínios mais
simples, elementos finitos (EF). Cada EF é constituído por nós, onde se consideram
explicitamente as variáveis do problema. Assim, as funções interpoladoras são construídas para
cada elemento de forma semelhante, efectuando-se em seguida a assemblagem da estrutura,
garantindo a continuidade da solução ao longo do domínio. As funções interpoladoras são
escolhidas à partida, dependendo tanto do tipo de problema em estudo como da sua geometria e
do número de nós considerados em cada elemento.
O MEF pode assim ser sistematizado em três passos principais:
1. Divisão do domínio em subdomínios, que permite a derivação sistemática das funções de
aproximação, assim como a representação de domínios complexos.
2. Derivação das funções de aproximação em cada elemento. Estas funções são geralmente
polinómios, obtidos recorrendo a métodos de interpolação.
3. Assemblagem dos elementos. Esta é baseada na continuidade da solução em todo o
domínio e nos fluxos internos.
Nos problemas de barras sujeitas a tensões axiais a equação diferencial regente é:

(1)

em que representa o deslocamento, a variável primária do problema, a

deformação e a força axial distribuída na barra, sendo o módulo de Young e a área da
secção de corte.
O objectivo é encontrar uma solução para cada EF. Para tal define-se a solução aproximada
(em cada elemento) como a combinação linear de um conjunto de funções que constituem uma
base:

(2)

onde são os valores da solução nos nós do EF e as funções de aproximação no mesmo.

3
 Pinto, J. e Semedo, J.

A solução aproximada não cumpre obviamente a equação (1) em todos os pontos, sendo por
isso utilizado um integral ponderado desta equação. Assim não é exigido o cumprimento da
equação em todo o seu domínio mas sim em média, ou seja:

(3)

sendo a função de ponderação utilizada.

Utiliza-se agora integração por partes para tornar e igualmente diferenciáveis,
diminuindo numa ordem a exigência da continuidade das funções de base (como se verá adiante,
estas também definirão ). O resultado desta integração é a formulação fraca da equação (1):

(4)

O segundo termo da equação (4), desconsiderando , representa a força de reacção das

extremidades do elemento, a variável secundária do sistema1. Num EF com nós interiores a
integração é efectuada por intervalos, , sendo o número de nós.
Assim, e efectuando a substituição , a equação (4) pode ser reescrita:

(5)

Uma vez que a equação (5) inclui não só a equação diferencial mas também as condições
impostas à variável secundária, resta exigir que a solução aproximada cumpra as condições de
fronteira. Para um EF com dois nós: e .
Regra geral são utilizados polinómios (completos) para aproximar a solução da equação
diferencial, não só pela facilidade na sua manipulação como também pela simplicidade com que
estas funções podem ser utilizadas na interpolação de dados. A solução aproximada, , deve ser
contínua e diferenciável ao longo do EF, tal como requerido pela formulação fraca, de forma a
desconsiderar soluções triviais. Para que as condições de fronteira sejam cumpridas há também
que garantir que interpola a variável primária pelo menos nos nós partilhados,
assegurando a continuidade da solução geral.
A solução aproximada é então um polinómio de grau nunca inferior a um.
(6)

1
Repare-se que tem dimensões de força. Para uma explicação completa consulte [1].

4
 Pinto, J. e Semedo, J.

Aplicando as condições de fronteira, mais uma vez para um EF com dois nós:

(7)

A equação (7) permite escrever as constantes polinomiais em função dos valores de fronteira.
Substituindo as expressões obtidas para e na equação (6) e comparando com a equação (2) é agora
possível identificar a base utilizada:

(8)

Estas funções fazem parte da família das funções interpoladoras de Lagrange.

Tal como referido anteriormente, as funções de interpolação são derivadas sistematicamente
no MEF: depois de escolhido o grau do polinómio interpolador derivam-se expressões para as
constantes polinomiais em função dos valores da variável principal nos nós, obtendo-se
finalmente a solução aproximada como combinação linear de polinómios de Lagrange. O
número de nós deve obviamente ser o suficiente para que se possa reescrever os coeficientes do
polinómio interpolador em função dos valores da variável primária nos nós, ou seja, o grau do
polinómio interpolador tem de ser igual a .
Procedendo de acordo com o método de Galerkin introduz-se a equação (2) na equação (5) e
substitui-se , obtendo-se assim equações
algébricas:

(9)

para , onde:

(10)

(11)

ou, em notação matricial2:

(12)

2
Por uma questão de simplicidade omitem-se as dependências, em , nas equações (9), (10) e (11).

5
 Pinto, J. e Semedo, J.

é a chamada matriz de rigidez, sendo o vector de forças distribuídas. A equação (12)

contém incógnitas: e , logo existe a necessidade de considerar,
pelo menos, condições. Algumas destas condições são simplesmente as condições de fronteira
do problema enquanto outras provêm da continuidade da solução ao longo do domínio, surgindo
com a assemblagem do sistema.
Obtidas as equações algébricas descritas pela equação (12) para cada EF há que realizar a
assemblagem, obtendo-se o sistema de equações geral:
(13)
Os vectores de incógnitas e referem-se agora ao número total de nós considerados no
domínio (tal como e ). Impondo finalmente as condições fronteira na equação (13)
estamos em condições de resolver o sistema.
O objectivo deste trabalho é a implementação e aplicação do MEF na resolução de problemas
de estruturas de barras bidimensionais, ou seja, os EF em estudo são barras, sendo por isso
constituídos por dois nós. Há então que derivar a matriz de rigidez para qualquer barra em 2D.
A matriz de rigidez de uma barra unidimensional é facilmente obtida avaliando a equação
(10) para as funções apresentadas na equação (8):

(14)

onde é o comprimento da barra.

A generalização para 2D, onde passa a haver dois graus de liberdade por cada nó, é imediata.
Considerando uma barra orientada segundo a horizontal, apenas sujeita a tensão axial, tem-se:

(15)

Para obter a matriz de rigidez de uma barra com qualquer orientação basta realizar a rotação
do sistema de coordenadas, através de uma matriz de transformação. A expressão geral para 2D
de uma barra é então:

(16)

em que representa o ângulo formado entre a barra e a direcção horizontal.

6
 Pinto, J. e Semedo, J.

2. IMPLEMENTAÇÃO COMPUTACIONAL

2.1.Organização de dados
Os dados do problema em análise são armazenados em duas matrizes principais:
 DATAnodes: é a matriz que armazena as posições iniciais de todos os nós do sistema,
além de informação sobre as suas características. Uma vez que se pretendem analisar
sistemas em 2D, a cada nó corresponderem duas coordenadas. Esta matriz tem por isso
linhas e duas colunas, em que é o número de nós utilizados na definição sistema.

(17)

onde são as coordenadas iniciais dos nós considerados, sendo que as

constantes assumem o valor ou conforme se trate de uma coordenada fixa
ou livre.
 DATAconec: nesta matriz são armazenadas as informações de cada EF considerado. Estas
incluem as informações de conectividade do sistema além das constantes características
de cada barra.
(18)
Na equação (18) exemplifica-se uma linha desta matriz, onde é o nó inicial da barra e
o nó final. e representam, respectivamente, a área da secção de corte e o módulo de
Young do material. Cada linha da matriz corresponde a um elemento.

2.2.Abordagem computacional
A implementação efectuada divide-se em três componentes distintas:
 Compute_K: função utilizada para construir a matriz de rigidez do sistema. Percorre
todos os elementos construindo a matriz de rigidez para cada um, assemblando
automaticamente a matriz de rigidez geral.
 FEM: função principal que, dadas as matrizes apresentadas na subsecção anterior, calcula
os deslocamentos generalizados, as reacções nos apoios e a tensão em cada barra.
 FEM_GUI: ficheiro de código MATLAB que executa a interface gráfica do programa.
Através da interface gráfica (abordada da próxima subsecção) são construídas, por esta
ordem, as matrizes e . As forças aplicadas são de seguida
adicionadas ao vector . A partir desta informação é gerada a matriz de rigidez geral, utilizada
para calcular os deslocamentos generalizados, as reacções nos apoios e a tensão em cada barra. A

7
 Pinto, J. e Semedo, J.

interface permite representar não só a configuração final do sistema, mas também informações
qualitativas, através de um sistema de cores, sobre a tensão em cada barra.

2.3.Guia do utilizador
O programa desenvolvido tira partido de uma interface gráfica para simplificar ao máximo a
interacção com o utilizador. Para iniciar a sua execução basta inserir FEM_GUI na linha de
comandos do MATLAB. Surge a janela apresentada na Figura 2.
Os dados do problema em análise devem ser inseridos da seguinte forma:
1. Adicionar todos os nós considerados no sistema, um a um, inserindo as suas coordenadas
através do campo “Add Node”. O tipo de apoio a que o nó adicionado é sujeito é
seleccionado através do pop-up menu presente (Figura 3). Depois de preenchidos
adequadamente todos os campos carregar “Add Node!”.
2. Adicionar todos as barras que compõem o sistema, uma a uma, utilizando o campo “Add
Bar”. Os campos a preencher para cada barra são os nós inicial e final assim como a área
da secção de corte e o módulo de Young do material de que é feita. Para adicionar a barra
pretendida ao sistema carregar “Add Bar!”.
3. Adicionar todas as força aplicadas ao sistema, uma a uma, inserindo as suas
componentes horizontal e vertical além do nó de aplicação, no campo “Add Force”.
Confirme a inserção carregando “Add Force!”.
O gráfico presente na janela do programa é automaticamente actualizado sempre que é
adicionado um novo dado ao sistema, sendo que completados os três passos descritos deve
encontrar-se desenhada a estrutura em análise, assim como representadas, nos nós correctos, as
forças aplicadas.
O sistema encontra-se agora completamente inserido no programa e resta efectuar a sua
resolução. Para tal basta carregar “Do it!”, por baixo do gráfico. É agora apresentada a
configuração final (utilizando um factor de escala, para que os deslocamentos sejam
perceptíveis) utilizando um código de cores para representar a tensão em cada barra.
Na linha de comandos do MATLAB são exportados os resultados os algoritmo, os seja, os
vectores , e (explicado em maior detalhe na subsecção seguinte). Os resultados
são apresentados desta forma de forma a facilitar a sua posterior utilização dentro do MATLAB
ou a exportação no formato que o utilizador considerar mais conveniente. No final o programa
verifica a correcção da solução encontrada executando o equilíbrio de forças e momentos. Os
resultados são apresentados no final do output na linha de comandos do MATLAB.

8
 Pinto, J. e Semedo, J.

Figura 2. Interface gráfica desenvolvida.

Figura 3. Pop-up menu onde estão disponíveis os vários tipos de apoio.

2.4.Exemplo
Mostra-se agora um exemplo de aplicação do programa. Na Figura 4 é visível a configuração
inicial da estrutura que serve de exemplo, apresentando-se na Figura 5 a configuração deformada
e na Figura 6 o output na linha de comandos.
Uma vez que o programa desenvolvido lida com problemas em 2D a cada nó corresponderam
duas coordenadas. Assim, nos vectores e o primeiro par de entradas correspondem ao nó
1, o segundo par ao nó 2 (indicados na Figura 4), etc.
O único nó livre deste sistema é o nó 2, sendo por isso o único que não conserva a posição
inicial. Na Figura 6 destaca-se o deslocamento deste nó e a força aplicada no mesmo.
De destacar que a força registada para o nó 2 é obviamente a força aplicada no input do
programa. As forças associadas aos restantes nós são as reacções nos apoios calculadas.

9
 Pinto, J. e Semedo, J.

Figura 4. Estrutura utilizada para a exemplificação do funcionamento do programa.

Figura 5. Configuração deformada (colorida) da estrutura de exemplo.

10
 Pinto, J. e Semedo, J.

Figura 6. Output na linha de comandos do MATLAB.

Quanto ao vector , cada valor de tensão está associado a uma das barras, sendo a
ordem aquela pela qual as barras foram inseridas (também indicada na representação da estrutura
deformada).

3. RESULTADOS E DISCUSSÃO

3.1.Resultados
O processamento do problema em mãos pelo programa implementado é de seguida
apresentado. É utilizado um factor de escala de 500 para destacar a deformação sofrida pela
estrutura.
No enunciado é pedida a configuração deformada, Figura 7, as reacções nos apoios, Tabela 1,
o valor da deformação em (que neste caso corresponde ao nó 6), Tabela 2, e a tensão em cada
barra, Tabela 3.

11
 Pinto, J. e Semedo, J.

Figura 7. Estrutura inicial, a tracejado, e configuração deformada (colorida).

Apoio
Nó 1 -0.8451 0.5000
Nó 5 -2.1549 4.5000
Tabela 1. Reacções nos apoios.

Nó
1 0 0
2 0.8356 -0.2562
3 0.7948 -0.1236
4 0.7540 -0.1837
5 0 0
6 0.4246 -0.1236
Tabela 2. Deformações obtidas.

12
 Pinto, J. e Semedo, J.

Barra
1 -11.773
2 -5.631
3 -5.633
4 -8.442
5 10.164
6 -25.896
7 0.008
8 10.161
9 -25.895
Tabela 3. Tensão aplicada a cada barra.

3.2.Discussão
Começando por observar a Figura 7, particularmente a configuração deformada, é possível
constatar que esta se encontra dentro do esperado. São aplicadas duas forças, uma na direcção
horizontal e outra na vertical mas no sentido negativo. Os apoios fixos permanecem na mesma
posição que ocupavam inicialmente além de que a deformação registada parece consistente com
as forças aplicadas. Os deslocamentos obtidos são pequenos, o que está de acordo com a hipótese
dos pequenos deslocamentos. Em particular o deslocamento obtido para o ponto , no programa
o nó 6, destacado a vermelho na Tabela 2, é precisamente positivo na horizontal e negativo na
vertical, apresentando deslocamentos na ordem dos decímetros.
A escala de cores apresentada também de encontra de acordo com o esperado, surgindo barras
comprimidas com cores mais frias e barras estendidas com cores quentes. Tal pode ser
confirmado comparando a Tabela 3 com a Figura 7.
O programa implementado realiza automaticamente uma verificação dos resultados obtidos,
através do equilíbrio estático de forças e de momentos. Tal como acontece na Figura 6, a
execução do problema sugerido pelo docente é aprovado nesta verificação. De referir que uma
vez que as forças obtidas são na ordem de grandeza dos kN, o cálculo dos equilíbrios,
principalmente o equilíbrio de momentos por incluir a multiplicação com os braços das forças,
não apresentará a precisão desejada. Este problema poderia ser atenuado aumentando o espaço
em memória para as variáveis do programa, aumentando a precisão dos resultado, mas tal resulta
num aumento do tempo de computação, o que pode ser drástico em estruturas mais complexas.
Ainda assim, tendo em conta a ordem das grandezas envolvidas o cálculo dos equilíbrios
continua a fornecer uma verificação importantíssima dos resultados. Os equilíbrios de forças e de
momentos para o problema sugerido são apresentados na Tabela 4. Como forma de validação do
programa foram ainda utilizadas outras estruturas tendo os resultados sido igualmente bons para
estas.

13
 Pinto, J. e Semedo, J.

3.2983E-012 9.0949E-013 -3.4095

Tabela 4. Equilíbrio de forças e momentos para o problema proposto.

4. CONCLUSÃO
Apesar de, como qualquer método numérico, ter um erro associado, o MEF apresenta
resultados bastante aproximados da realidade. Tal é de realçar uma vez que a versão
implementada é relativamente simples, com custos computacionais baixos. Se a isto se associar
uma interface simples e intuitiva com o utilizador, este torna-se no método de eleição na
resolução de estruturas de barras em 2D. Além disto é de referir que a mesma análise pode
facilmente ser estendida a 3D, aumentando ainda mais o potencial deste método.

REFERÊNCIAS
[1] J. N. Reddy, An Introduction to the Finite Element Method, 3ª edição, McGraw-Hill, 2006.

14