84
Formulação de Problemas 2D e 3D Mecânica Estrutural (10371/10391/1411) 2011 Pedro V. Gamboa Departamento de Ciências Aeroespaciais

Formulação de Problemas 2D e 3D - webx.ubi.ptwebx.ubi.pt/~pgamboa/pessoal/10371/apontamentos/03-Problemas2De3D.pdf · Os passos principais são: 1. Descretização do domínio num

Embed Size (px)

Citation preview

Formulação de Problemas

2D e 3D

Mecânica Estrutural (10371/10391/1411)

2011

Pedro V. Gamboa Departamento de Ciências Aeroespaciais

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 2

1. Introdução

A análise de elementos finitos de problemas bidimensionais

envolve os mesmos passos básicos dos problemas

unidimensionais.

A análise é um pouco mais complicada pelo facto de os

problemas bidimensionais serem descritos por equações

diferenciais parciais sobre regiões geometricamente complexas.

A fronteira G de um domínio bidimensional W é, em geral, uma

curva. Logo, os elementos finitos são formas geométricas

bidimensionais simples que permitem aproximar um dado

domínio bidimensional bem como a solução sobre esse domínio.

Assim, num problema bidimensional não procuramos apenas a

solução aproximada de um dado problema mas também

aproximamos o domínio através de uma malha apropriada.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 3

1. Introdução

Consequentemente existirão erros de aproximação:

• devido à aproximação da solução;

• devido à descretização do domínio.

A malha de elementos finitos (descretização) consiste em

simples elementos bidimensionais, como triângulos, retângulos

e/ou quadriláteros que permitem a derivação única de funções

de interpolação.

Os elementos são ligados entre si nos pontos nodais nas

fronteiras dos elementos. A capacidade de representar

geometrias irregulares por uma coleção de elementos finitos

torna o método muito útil e prático para a solução de problemas

de valor de contorno, de valor inicial ou de valores próprios em

várias áreas da engenharia.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 4

1. Introdução

Neste capítulo os passos básicos dos problemas unidimensionais

vão ser estendidos a problemas bidimensionais com uma única

variável dependente.

Os passos básicos da análise de elementos finitos vão ser

descritos recorrendo a uma equação diferencial parcial de

segunda ordem que governa uma única variável.

Esta equação aparece em muitas áreas como mostra a tabela

seguinte.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 5

1. Introdução

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 6

2. Problemas de Valor de

Contorno

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 7

2.1. Equação do Modelo

Considere-se o problema de encontrar u(x,y) da equação

diferencial parcial de segunda ordem

dados aij (i,j=1,2), a00 e f, e condições de fronteira

especificadas.

As formas das condições de fronteira vão ficar claras através da

formulação fraca. Como caso especial, pode obter-se a equação

de Poissson a partir de (2.1) colocando a11=a22=k(x,y) e

a12=a21=a00=0:

onde é o operador gradiente.

00022211211

fua

y

ua

x

ua

yy

ua

x

ua

x (2.1)

W em, yxfuk (2.2)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 8

2.1. Equação do Modelo

Se i e j representarem os vetores unitários nas direções x e y,

respetivamente, o operador gradiente pode ser escrito

e (2.2) num sistema de coordenadas cartesiano fica

Seguidamente deriva-se o modelo de elementos finitos da

equação (2.1).

yx

ji

),( yxfy

uk

yx

uk

x

(2.3)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 9

2.1. Equação do Modelo

Os passos principais são:

1. Descretização do domínio num conjunto de elementos finitos;

2. Formulação fraca (ou integral ponderado) da equação

diferencial;

3. Derivação das funções de interpolação do elemento finito;

4. Desenvolvimento do modelo de elemento finito usando a forma

fraca;

5. Montagem dos elementos finitos para obter o sistema de

equações algébricas global;

6. Imposição das condições de fronteira;

7. Solução das equações;

8. Pós computação da solução e parâmetros de interesse.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 10

2.1. Equação do Modelo

Os passos 6 e 7 permanecem iguais ao problema de elmentos

finitos unidimensionais porque no final do passo 5 tem-se um

conjunto de equações algébricas cuja forma é independente da

dimensão do domínio ou da natureza do problema.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 11

2.2. Descretização em Elemento Finito

Em duas dimensões existe mais do que uma geometria que pode

ser usada como elemento finito:

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 12

2.2. Descretização em Elemento Finito

Como se irá ver, as funções de interpolação dependem não só do

número de nós do elemento e do número de incógnitas por nó,

mas também da forma do elemento.

A forma do elemento tem que ser tal que a sua geometria seja

definida por um conjunto de pontos, que são os nós do

elemento, no desenvolvimento das funções de interpolação. O

triângulo é a forma mais simples, seguida do retângulo.

A representação de uma dada região por um conjunto de

elementos (discretização ou geração da malha) é um passo

importante na análise de elementos finitos. A escolha do tipo de

elemento finito, número de elementos, e densidade de

elementos depende da geometria do domínio, do problema a ser

analisado e o grau de precisão desejado.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 13

2.2. Descretização em Elemento Finito

Obviamente, não existem fórmulas específicas para obter esta

informação. Em geral, o analista é guiado pelo seu

conhecimento técnico, compreensão do problema físico em

estudo, e experiência na modelação por elementos finitos.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 14

2.2. Descretização em Elemento Finito

As regras gerais para a geração da malha de uma formulação de

elementos finitos inclui:

1. Os elementos selecionados devem caracterizar as equações que

governam o problema;

2. O número, forma e tipo (linear ou quadrático) dos elementos deve ser

tal que a geometria do domínio seja representada com a precisão

desejada;

3. A densidade de elementos deve ser tal que as regiões com gradientes

elevados da solução sejam modelados adequadamente (mais elementos

ou elementos de ordem superior em regiões com gradientes elevados);

4. Os refinamentos da malha devem variar gradualmente das regiões com

mais densidade para as regiões com menos densidade. Se forem usados

elementos de transição, estes devem ser usados longe de regiões

críticas. Elementos de transição são aqueles que ligam elementos de

ordem inferior com elementos de ordem superior.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 15

2.3. Forma Fraca

Para o desenvolvimento da forma fraca basta considerar um

elemento típico.

Assumindo que We é um elemento típico, triangular ou

quadrilateral, da malha de elementos finitos, pode desenvolver-

se o modelo de elemento finito de (2.1) para We, seguindo o

procedimento de três-passos já apresentado.

O primeiro passo consiste em multiplicar (2.1) pela função de

ponderação w, que se assume ser diferenciável uma vez com

respeito a x e y, e depois integra-se a equação sobre o elemento

We:

dxdyfuaFy

Fx

w

e

W

00210 (2.4a)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 16

2.3. Forma Fraca

onde

No segundo passo distribui-se as derivadas igualmente por w e u.

Para conseguir isto integram-se os dois primeiros termos de

(2.4a) por partes. Primeiro, atendendo às igualdades

y

ua

x

uaF

y

ua

x

uaF

22212

12111

(2.4b)

2222

22

1111

11

ou

ou

wFy

Fy

w

y

Fw

y

FwF

y

wwF

y

wFx

Fx

w

x

Fw

x

FwF

x

wwF

x

(2.5a)

(2.5b)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 17

2.3. Forma Fraca

Depois, usando a forma dos componentes do teorema do

gradiente

onde nx e ny são as componentes (isto é, os cosenos diretores) do

vetor normal unitário

na fronteira Ge e ds é o comprimento de um elemento linha

infinitesimal ao longo da fronteira.

GW

GW

ee

ee

dsnwFdxdywFy

dsnwFdxdywFx

y

x

22

11 (2.6a)

(2.6b)

jsinicosjin yx nn (2.7)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 18

2.3. Forma Fraca

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 19

2.3. Forma Fraca

Usando (2.5a), (2.5b), (2.6a) e (2.6b) em (2.4a) obtém-se

Inspecionando o integral de fronteira em (2.8), pode ver-se que

a especificação de u constitui a condição de fronteira essencial

e, por isso, u é a variável principal.

A especificação do coeficiente da função de ponderação na

expressão da fronteira

(2.8)

(2.9)

G

W

e

e

dsy

ua

x

uan

y

ua

x

uanw

dxdywfwuay

ua

x

ua

y

w

y

ua

x

ua

x

w

yx 22211211

00222112110

y

ua

x

uan

y

ua

x

uanq yxn 22211211

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 20

2.3. Forma Fraca

constitui a condição de fronteira natural. Assim, qn é a variável

secundária da formulação. A função qn=qn(s) representa a

projeção do vetor a. u ao longo da normal unitária n.

Por definição, qn é positiva para fora da superfície quando nos

movemos no sentido anti-horário ao longo da fronteira Ge. Na

maior parte dos problemas, a variável secundária qn tem

interesse físico.

O terceiro e último passo da formulação é usar a definição (2.9)

em (2.8) e escrever a forma fraca de (2.1) como

(2.10)

G

W

e

e

dswq

dxdywfwuay

ua

x

ua

y

w

y

ua

x

ua

x

w

n

00222112110

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 21

2.3. Forma Fraca

ou

onde as formas bilinear Be(.,.) e linear le(.) são

A forma fraca (ou forma integral ponderada) em (2.10) ou

(2.11a) e (2.11b) é a base do modelo de elementos finito de

(2.1).

(2.11b)

GW

W

ee

e

dswqdxdywfwl

dxdywuay

ua

x

ua

y

w

y

ua

x

ua

x

wwuB

n

e

e

)(

),( 0022211211

(2.11a) )(),( wlwuB ee

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 22

2.3. Forma Fraca

Sempre que Be(w,u) é simétrico nos seus argumentos w e u, isto

é, Be(w,u)=Be(u,w), o funcional quadrático associado com o

problema variacional (2.11a) pode ser obtido com

A forma bilinear em (2.11b) é simétrica se e só se a12=a21. Então

o funcional é dado por

(2.12b)

GW

W

ee

e

dsuqdxdyuf

dxdyuay

ua

y

u

x

ua

x

uawI

n

e 2

00

2

2212

2

11 22

1)(

(2.12a) )(),()( wlwwBwI eee

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 23

2.3.1. Forma Vetorial do Problema

Variacional

É comum, especialmente na mecânica estrutural, exprimir as

formulações de elementos finitos na forma vetorial (em termos

de matrizes).

Embora a notação vetorial/matricial seja concisa, não é tão

transparente como a forma explícita.

A equação (2.11a) pode ser escrita na forma

onde, no presente caso, w é simplesmente w e u é u.

Agora é necessário exprimir Be(.,.) e le(.) em forma matricial.

(2.13) )(),( wuwee lB

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 24

2.3.1. Forma Vetorial do Problema

Variacional

Colocando

então Be e le de (2.11b) pode escreve-se na forma

(2.14)

T

yxa

aa

aa

1,

00

0

0

00

2221

1211

DC

(2.15a)

GW

W

ee

e

dsqwdxdyfwwl

dxdy

u

y

ux

u

a

aa

aa

w

y

wx

w

uwB

n

TTe

T

e

)(

00

0

0

),(

00

2221

1211

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 25

2.3.1. Forma Vetorial do Problema

Variacional

ou simplesmente

(2.15b)

GW

W

ee

e

dsdxdyl

dxdyB

TTe

Te

qwfww

DuCDwuw

)(

)(),(

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 26

2.4. Modelo de Elemento Finito

A forma fraca em (2.10) requer que a aproximação escolhida

para u seja no mínimo linear tanto em x como em y para que

nenhum termo em (2.10) seja zero. Uma vez que a variável

primária é a própria função, a família de funções de integração

de Lagrange é admissível.

Suponha-se que u é aproximado num elemento finito típico We

pela expressão

onde ue e Ye são vetores nx1

e uje é o valor de uh

e no nó j (xj,yj) do elemento;

(2.16a) eTee

h

e

j

n

j

e

j

e

h yxuyxuyxuyxu uY

),(ou),(),(),(1

(2.16b) Te

n

e

h

e

h

eeTe

n

e

h

e

h

ee uuuu 11 , Yu

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 27

2.4. Modelo de Elemento Finito

e je são as funções de interpolação de Lagrange com a seguinte

propriedade

Quando se derivam as equaçãos do elemento finito em termos

algébricos, não é necessário saber a forma do elemento We ou a

forma de ie. A forma específica de i

e é desenvolvida para

formas específicas dos elementos tal como triangular,

retangular, etc..

Substituindo a aproximação de elemento finito (2.16a) para u na

forma fraca (2.10) ou (2.13) obtém-se

(2.17) ),,2,1,),( njiyx ijjj

e

i

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 28

2.4. Modelo de Elemento Finito

ou

Esta equação deve ser válida para todas as escolhas admissíveis

da função de ponderação w.

Uma vez que são necessárias n equações algébricas

independentes para resolver as n incógnitas u1e, u2

e, …, une

escolhem-se n funções independentes lineares para w.

(2.18a)

G

W

e

e

dswqdxdywfuway

ua

xua

y

w

yua

xua

x

w

n

n

j

e

j

e

j

n

j

e

je

j

n

j

e

je

j

n

j

e

je

j

n

j

e

je

j

1

00

1

22

1

21

1

12

1

110

(2.18b) GWW

Y

eee

dsdxdydxdy TTeTTqwfwuDCDw )()(0

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 29

2.4. Modelo de Elemento Finito

w=1e, 2

e, …, ne (ou w={1

e 2e … n

e }=YT). Esta escolha

particular para a função de ponderação é natural quando a

função de ponderação é vista como uma variação virtual da

incógnita dependente, isto é

E o modelo de elemento finito resultante é conhecido como

modelo de elemento finito de forma fraca ou modelo de

elemento finito de Ritz.

Para cada escolha de w obtém-se uma relação algébrica entre

(u1e, u2

e, …, une).

n

i

iiuuw1

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 30

2.4. Modelo de Elemento Finito

A equação algébrica resultante da substituição de ie para w

dentro de (2.18a) é numerada como a primeira equação

algébrica, aquela resultante de w=2e é numerada como a

segunda equação, e assim consecutivamente.

Assim, a equação algébrica i é obtida substituindo w=ie em

(2.18a):

ou

GW

W

ee

e

dsqdxdyfdxdyja

ya

xa

yya

xa

x

e

in

e

i

e

i

e

i

n

j

e

j

e

je

i

e

j

e

je

i

00

1

222112110

(2.19a) ),,2,1(1

niQfuK e

j

e

j

n

j

e

j

e

ij

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 31

2.4. Modelo de Elemento Finito

onde

Na forma matricial (2.19a) toma a forma

G

W

W

e

e

e

dsqQ

dxdyff

dxdyja

ya

xa

yya

xa

xK

e

in

e

i

e

i

e

i

e

i

e

i

e

j

e

je

i

e

j

e

je

ie

ij

00

22211211

(2.19b)

eeeeeeee QfuK QfuK ou (2.20a)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 32

2.4. Modelo de Elemento Finito

onde

Note-se que Kij=Kji, isto é, [Ke] é uma matriz simétrica de ordem

nxn, apenas quando a12=a21.

Y

Y

G

W

W

e

n

ee

e

yn

e

y

e

y

e

xn

e

x

e

x

T

e

e

Te

e

e

e

ds

dxdy

dxdy

21

,,2,1

,,2,1

Ψ

DB

qQ

ff

CBBK

(2.20b)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 33

2.4. Modelo de Elemento Finito

As equações (2.20a) e (2.20b) representam o modelo do

elemento finito de (2.1).

Isto completa o desenvolvimento do modelo do elemento finito.

Antes de descrever a montagem das equações do elemento é

conveniente derivar as funções de interpolação ie para certos

elementos básicos e calcular as matrizes (2.19b) do elemento.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 34

2.5. Derivação da Função de

Interpolação

A aproximação de elemento finito uhe(x,y) no elemento We tem

que satisfazer as condições seguintes para que a solução

aproximada convirja para a solução verdadeira:

1. A representação de uhe tem que ser contínua como

necessário na forma fraca do problema (isto é, todos os

termos na forma fraca são representados como valores

diferentes de zero):

2. Os polinómios usados para representar uhe têm que ser

completos (isto é, todos os termos desde a constante até aos

termos da ordem superior têm que estar presentes em uhe);

3. Todos os termos do polinómio devem ser linearmente

independentes.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 35

2.5. Derivação da Função de

Interpolação

O número de termos linearmente independentes na

representação de uhe dita a forma e número de graus de

liberdade do elemento.

Em seguida serão vistos alguns polinómios básicos e elementos

relacionados para o problema do modelo com um único grau de

liberdade por nó.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 36

2.5.1. Elemento Triangular

Examinando a forma fraca (2.10) e as matrizes do elemento

finito (2.19b) pode ver-se que ie têm que ser pelo menos

funções lineares de x e y.

O polinómio linear completo em x e y de We tem a forma

onde cie são constantes.

O conjunto {1,x,y} é linearmente independente e completo. A

equação (2.21) define um plano único para valores fixos de cie.

Assim, se u(x,y) é uma superfície curva, uhe(x,y) aproxima a

superfície por um plano.

ycxccyxu eeee

h 321),( (2.21)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 37

2.5.1. Elemento Triangular

Em particular, uhe(x,y) é definido de forma única num triângulo

pelos três valores de uhe(x,y) nos vértices do triângulo.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 38

2.5.1. Elemento Triangular

Vamos escrever

onde (xi,yi) representa as coordenadas do vértice i do triângulo.

Notar que o triângulo é único definido pelos três pares de

coordenadas (xi,yi).

As três constantes cie (i=1,2,3) em (2.21) podem ser expressas

em termos dos três valores nodais uie (i=1,2,3). Assim, o

polinómio (2.21) é associado com um elemento triangular e

existem três nós identificados, nomeadamente os vértices do

triângulo.

ee

h

ee

h

ee

h uyxuuyxuuyxu 333222111 ),(,),(,),( (2.22)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 39

2.5.1. Elemento Triangular

As equações em (2.22) têm a forma explícita

onde o índice do elemento e é omitido por simplicidade. Esta

forma vai ser usada em seguida.

Em forma matricial tem-se

A solução de (2.23) para ci (i=1,2,3) precisa da inversa da matriz

de coeficientes A.

33321333

23221222

13121111

),(

),(

),(

ycxccyxuu

ycxccyxuu

ycxccyxuu

h

h

h

Acu

ou

c

c

c

yx

yx

yx

u

u

u

3

2

1

33

22

11

3

2

1

1

1

1

(2.23)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 40

2.5.1. Elemento Triangular

A inversa deixa de existir sempre que quaisquer duas linhas ou

colunas são iguais. Só quando todos os nós estão na mesma linha

é que se tem duas colunas ou duas linhas iguais na matriz de

coeficientes em (2.23).

Assim, em teoria, desde que os três vértices do triângulo sejam

distintos e não estejam na mesma linha, a matriz de coeficientes

pode ser invertida. No entanto, em computações reais, se dois

nós estão muito perto do terceiro nó ou os três nós formam

quase uma linha reta, a matriz de coeficientes pode ser quase

singular numericamente e não ter inversa.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 41

2.5.1. Elemento Triangular

Assim, os elementos com geometrias estreitas devem ser

evitados.

Invertendo a matriz de coeficientes em (2.23) obtém-se

onde 2A é o determinante da matriz A.

321

321

321

321

12,

2

1

AA

A

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 42

2.5.1. Elemento Triangular

A é a área do triângulo cujos vértices estão em (xi,yi) (i=1,2,3).

Resolvendo para ci em termos de ui obtêm-se

onde i, i e i são constantes que dependem apenas nas

coordenadas globais dos nós (xi,yi) do elemento

)(2

1

)(2

1

)(2

1

3322113

3322112

3322111

uuuA

c

uuuA

c

uuuA

c

(2.24a)

)naturalordemempermutam,,;(

)(

kjikji

xx

yy

yxyx

kji

kji

jkkji

(2.24b)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 43

2.5.1. Elemento Triangular

Substituindo para ci de (2.24a) em (2.21) obtém-se

onde ie são as funções de interpolação linear do elemento

e ie, i

e e ie são as constantes definidas em (2.24b).

3

1

332211

332211332211

),(

)(

)()(2

1),(

i

e

i

e

i

e

h

yxu

yuuu

xuuuuuuA

yxu

(2.25a)

)3,2,1(2

1 iyx

A

e

i

e

i

e

i

e

i (2.25b)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 44

2.5.1. Elemento Triangular

As funções de interpolação linear ie estão representadas na

figura.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 45

2.5.1. Elemento Triangular

As funções de interpolação ie têm as seguintes propriedades

Pode notar-se que (2.24a) determina uma superfície plana que

passa por u1, u2 e u3.

Assim, o uso de funções de interpolação lineares ie de um

triângulo resulta na aproximação da superfície curva u(x,y) por

uma função planar

como mostra a figura seguinte.

)3,2,1,(),( 2 jiyx ijj

e

j

e

i (2.26a)

(2.26b) 0,0,13

1

3

1

3

1

i

e

i

i

e

i

i

e

iyx

3

1i

e

i

e

i

e

h uu

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 46

2.5.1. Elemento Triangular

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 47

2.6. Avaliação das Matrizes e Vetores do

Elemento

O cálculo exato das matrizes [Ke] e {fe} do elemento em (2.19b)

regra geral não é fácil. As matrizes são avaliadas usando técnicas

de integração numérica.

No entanto, quando aij, a00 e f são constantes ao longo do

elemento, é possível avaliar os integrais exatamente ao longo

dos elementos triangulares ou retangulares lineares.

O integral de fronteira em {Qe} em (2.19b) pode ser avaliado

sempre que qn é conhecido. Para um elemento interior, a

contribuição do integral de fronteira cancela com contribuições

similares de elementos adjacentes da malha (semelhante ao Qie

nos problemas unidimensionais).

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 48

2.6. Avaliação das Matrizes e Vetores do

Elemento

Para simplificar a escrita, podemos escrever [Ke] em (2.19b)

como a soma de cinco matrizes básicas [S](,=0,1,2)

onde [.]T representa a transposta da matriz e

com i,∂i/∂x, x1=x e x2=y; i,0=i.

Todas as matrizes em (2.38) e funções de interpolação em (2.39)

são definidas sobre o elemento.

Agora é necessário calcular as matrizes em (2.39) e (2.19b)

usando as funções de interpolação linear.

(2.38) ][][][][][][ 22

22

12

21

12

12

11

11

00

00 SaSaSaSaSaK Te

(2.39) W

e

dxdyS jiij ,,

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 49

2.6.1. Matrizes de um Elemento

Triangular Linear

Primeiro, note-se que os integrais de polinómios sobre domínios

de triângulos arbitrários podem ser calculados com exatidão.

Para isso, assuma-se que Imn é o integral da expressão xmyn sobre

um triângulo arbitrário D

Então, sendo (xi,yi) as coordenadas dos vértices do triângulo,

pode ser mostrado que

(2.40) D dxdyyxI nm

mn

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 50

2.6.1. Matrizes de um Elemento

Triangular Linear

Podemos usar estes resultados para avaliar os integrais definidos

sobre elementos triangulares.

(2.41)

D

D

D

DD

DD

DD

23

1

22

02

23

1

22

20

3

1

11

3

1

10

01

3

1

01

10

00

00

ˆ912

ˆ912

ˆˆ912

3

1ˆ,ˆ

3

1ˆ,ˆ

triângulodoárea1

yyA

dxdyyI

xxA

dxdyxI

yxyxA

xydxdyI

yyyAydxdydxdyyxI

xxxAxdxdydxdyyxI

AdxdydxdyyxI

i

i

i

i

i

ii

i

i

i

i

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 51

2.6.1. Matrizes de um Elemento

Triangular Linear

Em seguida, calcula-se [Ke] e {fe} para um elemento triangular

linear assumindo que aij e f são constantes ao longo do

elemento. Tendo em conta, também, que

obtém-se

(2.42a) 0,0,23

1

3

1

3

1

i

e

i

i

e

ie

i

e

i A

(2.42b) ee

e

ie

e

i

e

i Ayx3

2ˆˆ

(2.43) e

e

ii

e

e

ii

AyAx 2,

2

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 52

2.6.1. Matrizes de um Elemento

Triangular Linear

Tendo em conta a identidade (2.42b) e para um valor constante

de f=fe ao longo do elemento, tem-se

(2.44)

jiijjiji

ijjiijjijiij

jiijjiijjiij

IIIA

yxA

S

AS

AS

AS

021120

00

221211

1

ˆˆ4

1

4

1,

4

1,

4

1

(2.45)

eee

e

ie

e

i

e

ie

ee

e

iee

e

ie

e

i

e

ee

i

e

i

e

i

e

e

e

i

e

i

e

i

e

ee

ie

e

i

Afyxf

yAxAAA

fIII

A

f

dxdyyxA

fdxdyyxff

ee

3

1ˆˆ

2

1

ˆˆ22

2),(

011000

DD

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 53

2.6.1. Matrizes de um Elemento

Triangular Linear

O resultado em (2.45) devia ser óbvio porque para uma fonte

constante fe a magnitude total da fonte (por exemplo calor) no

elemento é igual a feAe, que é distribuído igualmente por todos

os nós, resultando num valor nodal de feAe/3.

Uma vez conhecidas as coordenadas dos nós dos elementos,

pode calcular-se ie, i

e e ie a partir de (2.24b) e substituir na

equação (2.44) para se obter as matrizes do elemento que, por

sua vez, são usadas em (2.38) para se obter a matriz [Ke].

Em particular, quando a12, a21 e a00 são zero e a11 e a22 são

constantes no elemento, a equação (2.1) fica

efy

ua

x

ua

xW

em0

2

2

222

2

11 (2.46)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 54

2.6.1. Matrizes de um Elemento

Triangular Linear

e a matriz de coeficientes associados para um elemento

triangular linear é

(2.47) e

j

e

i

ee

j

e

i

e

e

e

ij aaA

K 22114

1

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 55

2.6.2. Avaliação dos Integrais de

Fronteira

Vamos considerar a avaliação do integral de fronteira do tipo

(ver equação (2.19b))

onde qne é uma função conhecida da distância s ao longo da

fronteira Ge.

Não é necessário calcular este integral quando uma porção de Ge

não coincide com a fronteira G do domínio total W. Em porções

de Ge que estão no interior do domínio W, qne no lado (i,j) do

elemento We cancela com qnf no lado (p,q) do elemento Wf

quando os lados (i,j) do elemento We e (p,q) do elemento Wf são

o mesmo (isto é no interface dos elementos We e Wf).

(2.56) G

e

dssqQ e

i

e

n

e

i )(

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 56

2.6.2. Avaliação dos Integrais de

Fronteira

Isto pode ser visto como o equilíbrio do “fluxo” interno.

Quando Ge coincide com a fronteira do domínio W, qn ou é

conhecido como função de s ou deve ser determinado no pós-

processamento. A variável primária deve ser especificada na

porção da fronteira onde qn não é especificado.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 57

2.6.2. Avaliação dos Integrais de

Fronteira

A fronteira Ge de um elemento bi-dimensional consiste em

segmentos de reta, que podem ser considerados como elementos

uni-dimensionais. Assim, a avaliação dos integrais de fronteira

em problemas bi-dimensionais consiste em avaliar integrais de

linha.

Por exemplo, considere-se um elemento triangular linear como

mostra a figura.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 58

2.6.2. Avaliação dos Integrais de

Fronteira

As funções de interpolação linear para este elemento são dadas

por (2.25).

Agora, escolhe-se o sistema de coordenadas locais (s,t) com a

origem no ponto 1 e a coordenada s paralela ao lado que liga os

nós 1 e 2. Os dois sistemas de eixos (x,y) e (s,t) são relacionados

com

As constantes a1, b1, c1, a2, b2 e c2 podem ser determinadas com

as seguintes condições

tcsbay

tcsbax

222

111

33

22

11

,,,quando

,,0,quando

,,0,0quando

yyxxbtcs

yyxxtas

yyxxts

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 59

2.6.2. Avaliação dos Integrais de

Fronteira

Obtém-se

As equações (2.57) permitem-nos exprimir i(x,y) como i(s,t)

que pode ser avaliada no lado que liga os nós 1 e 2 colocando

t=0 em i(s,t):

(2.57)

b

tyy

a

cy

a

c

a

syyytsy

b

txx

a

cx

a

c

a

sxxxtsx

321121

321121

1)(),(

1)(),(

a

syyysy

a

sxxxsx

sysxss iii

)()(,)()(

))0,(),0,(()0,()(

121121

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 60

2.6.2. Avaliação dos Integrais de

Fronteira

Por exemplo, tem-se

onde as definições de 1, 1 e 1 são usadas para reescrever a

expressão completa.

De forma idêntica tem-se

onde a=h12 é o comprimento do lado 1-2.

a

s

a

s

A

ya

sy

a

sx

a

sx

a

s

As

112

1

112

1)(

321

21121111

0)(,)( 32 sa

ss

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 61

2.6.2. Avaliação dos Integrais de

Fronteira

Note-se que 1(s) e 2(s) são exatamente as funções lineares e

unidimensionais associadas com o elemento linha qua liga os nós

1 e 2.

De forma similar, quando i(x,y) são avaliadas no lado 3-1 do

elemento, obtém-se

onde a coordenada s é medida ao longo do lado 3-1, com origem

no ponto 3, e h13 é o comprimento do lado 1-3. Assim, a

avaliação de Qie involve a utilização das funções de interpolação

unidimensionais e as incógnitas de qn na fronteira.

13

32

13

1 1)(,0)(,)(h

sss

h

ss

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 62

2.6.2. Avaliação dos Integrais de

Fronteira

Em geral, o integral (2.56) ao longo da fronteira do elemento

triangular linear pode ser expresso como

onde representa o integral ao longo da linha que liga o nó i ao

nó j, a coordenada s é medida do nó i para o nó j, com a origem

no nó i, e Qije é definido como sendo a contribuição de qn no

lado J do elemento We em Qie:

Onde i refere-se ao nó i do elemento e J refere-se ao lado J do

elemento.

(2.58a) e

i

e

i

e

i

ninini

e

i

QQQ

dssqsdssqsdssqsQ

321

133221

)()()()()()(

(2.58b) ladoJ

ni

e

iJ dsqQ

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 63

2.6.2. Avaliação dos Integrais de

Fronteira

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 64

2.6.2. Avaliação dos Integrais de

Fronteira

Por exemplo, tem-se

A contribuição do lado 2-3 é zero porque 1 é zero no lado 2-3

de um elemento triangular. Para um elemento retangular, Q1e

tem quatro partes mas apenas as contribuições dos lados 1-2 e 4-

1 são diferentes de zero porque 1 é zero nos lados 2-3 e 3-4.

G

13

113

21

12111 )()()( dsqdsqdssqQ nnn

e

e

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 65

2.7. Montagem das Equações do

Elemento

A montagem das equações do elemento finito é baseada nos

mesmos dois princípios que foram usados nos problemas

unidimensionais:

1. Continuidade das variáveis primárias

2. “Equilíbrio” das variáveis secundárias

Pode ilustrar-se o procedimento considerando uma malha de

elementos finitos constituída por elementos triangulares e um

elemento quadrilateral.

Vamos considerar que Kij1 (i,j=1,2,3) representa a matriz de

coeficientes do elemento triangular e que Kij2 (i,j=1,…,4)

representa a matriz de coeficientes do elemento quadrilateral.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 66

2.7. Montagem das Equações do

Elemento

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 67

2.7. Montagem das Equações do

Elemento

Da malha de elementos finitos da figura podemos notar a

correspondência seguinte (isto é, as relações de conectividade)

entre os nós globais e do elemento:

onde x indica que não existe valor.

A correspondência entre os valores nodais locais e globais é

que consiste em impor a continuidade das variáveis primárias

nos nós comuns dos elementos 1 e 2.

(2.59)

3542

321B

(2.60) 5

2

34

2

23

2

4

1

32

2

1

1

21

1

1 ,,,, UuUuUuuUuuUu

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 68

2.7. Montagem das Equações do

Elemento

Note-se que a continuidade das variáveis primárias nos nós inter-

elemento garante a continuidade da variável primária ao longo

de toda a fronteira inter-elemento. Para o caso da figura, o

requisito de u21=u1

2 e u31=u4

2 garante que uh1(s)=uh

2(s) no lado

que liga os nós globais 2 e 3.

Isto pode ser mostrado da seguinte forma. A solução uh1(s) ao

longo da linha que liga os nós globais 2 e 3 é linear e é dada por

onde s é a coordenada local com a sua origem no nó 2 e h é o

comprimento do lado 2-3 (ou lado 2).

h

su

h

susuh

1

3

1

2

1 )1()(

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 69

2.7. Montagem das Equações do

Elemento

De forma idêntica, a solução de elemento finito ao longo da

mesma linha mas no elemento 2 é

Uma vez que u21=u1

2 e u31=u4

2, então uh1(s)=uh

2(s) para todos os

valores de s ao longo do interface dos dois elementos.

Em seguida usa-se o equilíbrio das variáveis secundárias. No

interface entre os dois elementos, o fluxo dos dois elementos

deve ser igual em magnitude e opostos em sinal. Para os dois

elementos da figura, o interface é ao longo do lado que liga os

nós globais 2 e 3.

h

su

h

susuh

2

4

2

1

2 )1()(

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 70

2.7. Montagem das Equações do

Elemento

Assim, o fluxo interno qn1 no lado 2-3 do elemento 1 deve

equilibrar o fluxo qn2 no lado 4-1 do elemento 2:

No método dos elementos finitos impõe-se a relação acima na

forma de integrais ponderados:

onde hpqe representa o comprimento do lado que liga o nó p ao

nó q do elemento We.

As equações acima podem ser escritas na forma

41

2

32

1

14

2

32

1 ou

nnnn qqqq (2.61)

214

123

214

123

2

4

21

3

12

1

21

2

1 ,

h

n

h

n

h

n

h

n dsqdsqdsqdsq (2.62a)

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 71

2.7. Montagem das Equações do

Elemento

ou

onde QiJe representa a parte de Qi

e que vem do lado J do

elemento e:

Os lados de elementos triangulares ou quadrilaterais são

numerados como mostra a figura. Estas relações de equilíbrio

devem ser impostas na montagem das equações do elemento.

(2.62b) 0,0214

123

214

123

2

4

21

3

12

1

21

2

1 h

n

h

n

h

n

h

n dsqdsqdsqdsq

(2.62c) 0,0 1

44

1

32

1

14

1

22 QQQQ

ladoJ

e

i

e

n

e

iJ dsqQ

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 72

2.7. Montagem das Equações do

Elemento

Note-se que QiJe é apenas uma porção de Qi

e [ver (2.56) e

(2.58b)].

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 73

2.7. Montagem das Equações do

Elemento

As equações do elemento da malha de dois elementos da figura

são escritas em primeiro lugar. Para o problema em questão,

existe apenas um grau de liberdade primário por nó. Para o

elemento triangular, as equações do elemento têm a forma

(2.63a)

1

3

1

3

1

3

1

33

1

2

1

32

1

1

1

31

1

2

1

2

1

3

1

23

1

2

1

22

1

1

1

21

1

1

1

1

1

3

1

13

1

2

1

12

1

1

1

11

QfuKuKuK

QfuKuKuK

QfuKuKuK

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 74

2.7. Montagem das Equações do

Elemento

Para o elemento quadrilateral, as equações do elemento são

dadas por

Para se impor o equilíbrio das variáveis secundárias em (2.62c) é

necessário que se adicione a segunda equação do elemento 1 à

primeira equação do elemento 2 e, também, que se adicione a

terceira equação do elemento 1 à quarta equação do elemento

2:

(2.63b)

2

4

2

4

2

4

2

44

2

3

2

43

2

2

2

42

2

1

2

41

2

3

2

3

2

4

2

34

2

3

2

33

2

2

2

32

2

1

2

31

2

2

2

2

2

4

2

24

2

3

2

23

2

2

2

22

2

1

2

21

2

1

2

1

2

4

2

14

2

3

2

13

2

2

2

12

2

1

2

11

QfuKuKuKuK

QfuKuKuKuK

QfuKuKuKuK

QfuKuKuKuK

2

4

2

4

1

3

1

3

2

4

2

44

2

3

2

43

2

2

2

42

2

1

2

41

1

3

1

33

1

2

1

32

1

1

1

31

2

1

2

1

1

2

1

2

2

4

2

14

2

3

2

13

2

2

2

12

2

1

2

11

1

3

1

23

1

2

1

22

1

1

1

21

QfQfuKuKuKuKuKuKuK

QfQfuKuKuKuKuKuKuK

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 75

2.7. Montagem das Equações do

Elemento

Usando a notação global em (2.60), podemos escrever as

equações em cima como

Agora pode impor-se as condições em (2.62c), igualando a zero

as porções apropriadas das expressões dentro de parênteses do

lado direito das equações acima:

O termos sublinhados são zero pelas relações de equilíbrio

(2.62c).

2

4

1

3

2

4

1

35

2

434

2

423

2

44

1

232

2

41

1

321

1

31

2

1

1

2

2

1

1

25

2

134

2

123

2

14

1

232

2

11

1

221

1

21

QQffUKUKUKKUKKUK

QQffUKUKUKKUKKUK

2

43

2

42

2

41

2

44

1

32

1

33

1

31

2

44

2

43

2

42

2

41

1

33

1

32

1

31

2

4

1

3

2

13

2

12

2

11

2

14

1

22

1

23

1

21

2

14

2

13

2

12

2

11

1

23

1

22

1

21

2

1

1

2

QQQQQQQ

QQQQQQQQQ

QQQQQQQ

QQQQQQQQQ

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 76

2.7. Montagem das Equações do

Elemento

Os outros termos de cada equação ou são conhecidos porque qn é

conhecido na fronteira ou permanecem incógnitos porque a

variável primária é especificada na fronteira.

Em geral, quando vários elementos estão ligados, a montagem

dos elementos é feita através da colocação dos coeficientes dos

elementos Kije, fi

e e Qie nas posições adequadas da matriz e vetor

da direita globais. Isto é feito com as relações de conectividade,

isto é, com a as relações de correspondência do número local do

nó com o número global do nó.

Por exemplo, se o número global do nó 3 corresponder ao nó 3

do elemento 1 e ao nó 4 do elemento 2, então tem-se

2

44

1

3333

2

4

1

3

2

4

1

3

2

4

1

33 , KKKQQffFFF

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 77

2.7. Montagem das Equações do

Elemento

Se os números globais dos nós 2 e 3 corresponderem,

respetivamente, aos nós 2 e 3 do elemento 1 e aos nós 1 e 4 do

elemento 2, então os coeficientes globais K22, K23 e K33 são dados

por

De forma idêntica, os componentes fonte dos nós globais 2 e 3

são adicionados:

Para a malha de dois elementos da figura anterior, pode obter-

se as equações montadas.

2

44

1

3333

2

14

1

2323

2

11

1

2222 ,, KKKKKKKKK

2

4

1

33

2

1

1

22 , FFFFFF

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 78

2.7. Montagem das Equações do

Elemento

Assim

O procedimento de montagem descrito acima pode ser usado

para montar elementos de qualquer forma e tipo.

Por exemplo, considere-se a malha de elementos finitos

mostrada na figura anterior. A posição (4,4) da matriz global de

coeficientes contém K331+K11

2+K113. A posição 4 no vetor montado

contém F31+F1

2+F13.

(2.64)

2

3

2

2

2

4

1

3

2

1

1

2

1

1

5

4

3

2

1

2

33

2

32

2

34

2

31

2

23

2

22

2

24

2

21

2

43

2

42

2

44

1

33

2

41

1

32

1

31

2

13

2

12

2

14

1

23

2

11

1

22

1

21

1

13

1

12

1

11

0

0

00

F

F

FF

FF

F

U

U

U

U

U

KKKK

KKKK

KKKKKKK

KKKKKKK

KKK

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 79

2.7. Montagem das Equações do

Elemento

As posições (1,5), (1,6), (1,7), (2,5), (2,6), (2,7), (3,6), (3,7) e

(4,7) da matriz global contêm zeros porque KIJ=0 quando os nós

globais I e J não correspondem aos nós do mesmo elemento da

malha.

Estão completos os cinco primeiros passos na modelação por

elementos finitos da equação (2.1). Os dois passos seguintes da

análise, nomeadamente a imposição das condições de fronteira e

a solução das equações são iguais ao caso dos problemas

unidimensionais.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 80

2.8. Cálculos Posteriores

A solução do elemento finito em qualquer ponto (x,y) num

elemento We é dada por

e as suas derivadas são calculadas a partir de (2.65) como

As equações (2.65) e (2.66) podem ser usadas para calcular a

solução e as suas derivadas em qualquer ponto (x,y) dentro do

elemento.

(2.65)

n

j

e

j

e

j

e

h yxuyxu1

),(),(

(2.66)

n

j

e

je

j

e

h

n

j

e

je

j

e

h

yu

y

u

xu

x

u

1

1

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 81

2.8. Cálculos Posteriores

É útil gerar, pela interpolação de (2.65), a informação

necessária para representar graficamente uhe e as suas

derivadas.

As derivadas de uhe não serão contínuas nas fronteiras entre

elementos porque a continuidade das derivadas não é imposta

durante o procedimento de montagem. A forma fraca das

equações sugere que a variável primária é u, que é considerada

a variável nodal. Se outras variáveis, como derivadas de ordem

superior da incógnita dependente, forem consideradas como

variáveis nodais para as tornar contínuas na fronteira entre

elementos, o grau de interpolação (ou a ordem do elemento)

aumenta.

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 82

2.8. Cálculos Posteriores

Por outro lado, a continuidade das derivadas de ordem superior

que não são identificadas como variáveis primárias pode violar

os princípios físicos do problema.

Por exemplo, se ∂u/∂x for tornado contínuo viola-se o requisito

de que qx (=a11 ∂u/∂x) é contínuo no interface de dois materiais

diferentes porque a11 é diferente para os dois materiais no

interface.

Para o elemento triangular linear, as derivadas são constantes

dentro do elemento:

(2.67)

n

j e

j

e

je

hn

j e

j

e

je

h

j

e

e

j

j

e

e

j

jjj

e

e

j

A

u

y

u

A

u

x

u

AyAxyx

A

11 2,

2

2

1,

2

1,

2

1

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 83

2.8. Cálculos Posteriores

Para o elemento retangular linear, ∂Ue/∂x é linear em y e

∂uhe/∂y é linear em x:

onde x e y são as coordenadas locais.

Apesar de ∂uhe/∂x e ∂uh

e/∂y serem funções lineares de x e y,

respetivamente, em cada elemento, elas são descontínuas nas

fronteiras entre elementos. Consequentemente, as quantidades

calculadas usando as derivadas da solução do elemento finito uhe

são descontínuas nas fronteiras entre elementos.

(2.68)

n

j

je

j

je

hn

j

je

j

je

h

j

e

jj

e

j

a

xxu

by

u

b

yyu

ax

u

a

xx

byb

yy

ax

1

2

1

2 1)1(1

,1)1(1

11

,11

Faculdade de Engenharia

Universidade da Beira Interior

Mecânica Estrutural - 2011

Departamento de Ciências Aeroespaciais

Pedro V. Gamboa 84

2.8. Cálculos Posteriores

Por exemplo, se calcularmos qxe=a11

e ∂uhe/∂x num nó partilhado

por três elementos diferentes, são esperados três valores

diferentes de qxe. A diferença entre os três valores vai-se

reduzindo à medida que a malha é refinada. Alguns programas

comerciais de elementos finitos fornecem apenas um valor de qx

no nó através da média obtida dos valores dos vários elementos

ligados no nó.