| | | |

Método dos elementos finitos

Ajude a manter o site livre, gratuito e sem propagandas. Colabore!

2.2 Problema modelo

Vamos estudar a aplicação do método de elementos finitos a equação de Poisson777Siméon Denis Poisson, 1781 - 1840, matemático francês. Fonte: Wikipédia: Siméon Denis Poisson. com condições de Dirichlet888Johann Peter Gustav Lejeune Dirichlet, 1805 - 1859, matemático alemão. Fonte: Wikipédia: Johann Peter Gustav Lejeune Dirichlet., que usaremos como um problema modelo. Mais precisamente, definimos o chamdo problema forte: encontrar u tal que

−Δ⁢u=f,(x,y)∈Ω:=[0,1]2, (2.31)
u=0,(x,y)∈∂Ω, (2.32)

onde Δ=∂2/∂x2+∂2/∂y2 é o operador de Laplace999Pierre-Simon Laplace, 1749 - 1827, matemático francês. Fonte: Wikipédia: Pierre-Simon Laplace. e f é uma função dada.

2.2.1 Formulação Fraca

A aplicação do método de elementos finitos é construída sobre a formulação fraca do problema (2.31)-(2.32). Para a obtermos, multiplicamos (2.31) por uma função teste v em um espaço adequado V0 e integramos no domínio Ω, obtendo

−∫ΩΔ⁢u⁢v⁢𝑑x=∫Ωf⁢v⁢𝑑x. (2.33)

Então, no lado esquerdo, aplicamos a fórmula de Green101010George Green, 1793 - 1841, matemático britânico. Fonte: Wikipédia: George Green.

∫ΩΔ⁢u⁢v⁢𝑑x=−∫Ω∇u⋅∇v⁢d⁢x+∫∂Ω(𝒏⋅∇u)⁢v⁢𝑑s. (2.34)

donde temos

∫Ω∇u⋅∇v⁢d⁢x−∫∂Ω(𝒏⋅∇u)⁢v⁢𝑑s=∫Ωf⁢v⁢𝑑x. (2.35)

Então, observando critérios de regularidade e a condição de contorno (2.32), escolhemos o espaço teste

V0:={v∈H1⁢(Ω):v|∂Ω=0}. (2.36)

Lembramos que H1⁢(Ω)={v:‖v‖L2⁢(Ω)+‖∇v‖L2⁢(Ω)<∞}.

Com isso, temos o seguinte problema fraco associado a (2.31)-(2.32): encontrar u∈V0 tal que

a⁢(u,v)=L⁢(v),∀v∈V0, (2.37)

onde a⁢(u,v) é chamada de forma bilinear e definida por

a⁢(u,v):=∫Ω∇u⋅∇v⁢d⁢x (2.38)

e L⁢(v) é chamada de forma linear e definida por

L⁢(v):=∫Ωf⁢v⁢𝑑x. (2.39)

2.2.2 Formulação de Elementos Finitos

A formulação de elementos finitos é obtida da formulação fraca (2.37) pela aproximação do espaço teste V0 por uma espaço de dimensão finita. Tomando uma triangulação 𝒦⊂Ω e considerando o espaço contínuo dos polinômios afins por partes

Vh:={v:v∈C0⁢(Ω),v|K∈P1⁢(K)⁢∀K∈𝒦}, (2.40)

assumimos o espaço de elementos finitos

Vh⁢,0:={v∈Vh:v|∂Ω=0}. (2.41)

Com isso, temos o seguinte problema de elementos finitos associado (2.37): encontrar uh∈Vh⁢,0 tal que

a⁢(uh,vh)=L⁢(vh),∀vh∈Vh⁢,0. (2.42)

Observemos que (2.42) é equivalente ao problema de encontrar uh∈Vh⁢,0 tal que

a⁢(uh,φi)=L⁢(φi), (2.43)

com i=0,1,⋯,nin−1, onde {φi}i=0nin−1 é a base nodal de Vh⁢,0 e nin é o número de funções bases (igual ao número de nodos internos da triangulação 𝒦). Ainda, como

uh=∑j=0nin−1ξj⁢φj, (2.44)

temos

a⁢(uh,φi)=a⁢(∑j=0nin−1ξj⁢φj,φi) (2.45)
=∑j=0nin−1ξja(φj,φi). (2.46)

Com isso, o problema de elementos finitos é equivalente a resolver o seguinte sistema linear

∑j=0ni−1ξj⁢a⁢(φj,φi)=L⁢(φi),i=0,1,⋯,ni−1, (2.47)

para as incógnitas ξj, j=0,1,⋯,nin−1. Ou, equivalentemente, temos sua forma matricial

A⁢𝝃=𝒃, (2.48)

onde A=[ai,j]i,j=0ni−1 é chamada de matriz de rigidez com

ai,j=a⁢(φj,φi) (2.49)

e 𝒃=(b0,b1,⋯,bni−1) é o vetor de carga com

bi=L⁢(φi). (2.50)
Exemplo 2.2.1.

Consideremos o seguinte problema de Poisson

−Δ⁢u=2⁢π2⁢sen⁡(π⁢x)⁢sen⁡(π⁢y),x∈Ω:=(0,1)2, (2.51)
u=0,x∈∂Ω. (2.52)

A solução exata deste problema é dada por u⁢(x,y)=sen⁡(π⁢x)⁢sen⁡(π⁢y). Na Figura 2.4 temos um esboço da aproximação de elementos finitos obtida em uma malha uniforme com 16×16 nodos. As isolinhas correspondem à solução exata u.

Refer to caption
Figura 2.4: Solução de elementos finitos versus solução analítica do problema do Exemplo 2.2.1.

O seguinte Código 16 computa a solução de elementos finitos deste problema. Verifique!

Código 16: mef2d_pm.py
1from mpi4py import MPI
2
3# malha
4from dolfinx import mesh
5domain = mesh.create_unit_square(MPI.COMM_WORLD,
6 nx = 16, ny=16)
7# espaço
8from dolfinx import fem
9V = fem.functionspace(domain, ('P', 1))
10
11# condição de contorno
12import numpy as np
13uD = fem.Function(V)
14uD.interpolate(lambda x: np.full(x.shape[1], 0.))
15
16tdim = domain.topology.dim
17fdim = tdim - 1
18domain.topology.create_connectivity(fdim, tdim)
19boundary_facets = mesh.exterior_facet_indices(domain.topology)
20boundary_dofs = fem.locate_dofs_topological(V, fdim,
21 boundary_facets)
22bc = fem.dirichletbc(uD, boundary_dofs)
23
24# problema mef
25import ufl
26from dolfinx import default_scalar_type
27from dolfinx.fem.petsc import LinearProblem
28u = ufl.TrialFunction(V)
29v = ufl.TestFunction(V)
30
31# termo de fonte
32x = ufl.SpatialCoordinate(domain)
33f = 2 * ufl.pi**2 * ufl.sin(ufl.pi * x[0]) * ufl.sin(ufl.pi * x[1])
34
35# forma bilinear
36a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx
37
38# forma linear
39L = f * v * ufl.dx
40
41# problema linear
42problem = LinearProblem(a, L, bcs=[bc],
43 petsc_options_prefix="ex_mef2d_modelo")
44uh = problem.solve()

2.2.3 Análise de elementos finitos

Vamos estudar os fundamentos da análise de elementos finitos para o problema modelo (2.31)-(2.32). Vamos estabelecer a existência e unicidade da solução do problema de elementos finitos (2.42) e, posteriormente, obter uma estimativas a priori e a posteriori do erro ‖|u−uh‖|L2⁢(Ω).

Existência e unicidade

Para mostrarmos a existência e unicidade da solução de elementos finitos, precisamos do seguinte resultado sobre a matriz de rigidez.

Lema 2.2.1.(Matriz positiva definida)

A matriz de rigidez é positiva definida.

Demonstração.

A matriz de rigidez A=[a⁢(φj,φi)]i⁢j=0nin−1 é obviamente simétrica. Além disso, para todo 𝝃∈ℝnin, 𝝃≠0, temos

𝝃T⁢A⁢𝝃=∑i,j=0nin−1ξj⁢a⁢(φj,φi)⁢ξi (2.53)
=∑i,j=0nin−1ξj∫Ω∇φj⋅∇φidxξi (2.54)
=∫Ω∇(∑j=0nin−1ξjφj)⋅∇(∑i=0nin−1ξiφi)dx (2.55)
=∥∇(∑j=0nin−1ξjφj)∥L2⁢(Ω)2. (2.56)

Portanto, 𝝃T⁢A⁢𝝃≥0 e é nulo se, e somente se, v=∑j=0nin−1ξj⁢φj for constante. Como v∈Vh⁢,0, temos que v constante implica v≡0, mas então 𝝃=0, o que é uma contradição. Logo, 𝝃T⁢A⁢𝝃>0 para todo 𝝃∈ℝnin, 𝝃≠0. ∎

Teorema 2.2.1.(Existência e unicidade)

O problema de elementos finitos (2.42) tem solução única.

Demonstração.

O problema de elementos finitos (2.42) se resume a resolver o sistema linear A⁢𝝃=𝒃. Do Lema 2.2.1, temos que A é uma matriz definida positiva e, portanto, invertível. Daí segue, imediatamente, que o problema (2.42) tem solução única. ∎

Estimativa a priori do erro

Teorema 2.2.2.(Ortogonalidade de Galerkin)

A solução uh do problema de elementos finitos (2.42) satisfaz

a⁢(u−uh,vh)=0,∀vh∈Vh⁢,0, (2.57)

onde u é a solução do problema fraco (2.37).

Demonstração.

Segue, imediatamente, do fato de que Vh⁢,0⊂V0 e, portanto,

a⁢(u,vh)=L⁢(vh),∀vh∈Vh⁢,0, (2.58)

bem como

a⁢(uh,vh)=L⁢(vh),∀vh∈Vh⁢,0. (2.59)

∎

No decorrer do texto, vamos usar a norma da energia, definida por

‖v‖E:=(∫Ω∇v⋅∇v⁢d⁢x)1/2=‖∇v‖L2⁢(Ω), (2.60)

para todo v∈V0. Nesta norma, a solução de elementos finitos uh é a melhor aproximação da solução exata u no espaço Vh⁢,0, como mostra o seguinte resultado.

Teorema 2.2.3.(Melhor aproximação)

A solução uh do problema de elementos finitos satisfaz

‖u−uh‖E≤‖u−vh‖E,∀vh∈Vh⁢,0. (2.61)
Demonstração.

Observando que u−uh=u−vh+vh−uh e usando a ortogonalidade de Galerkin (Teorema 2.2.2), temos:

‖u−uh‖E2=∫Ω∇(u−uh)⋅∇(u−uh)⁡d⁢x (2.62)
=∫Ω∇(u−uh)⋅∇(u−vh)dx
+∫Ω∇(u−uh)⋅∇(vh−uh)⁡d⁢x (2.63)
=∫Ω∇(u−uh)⋅∇(u−vh)dx (2.64)
≤∥∇(u−uh)∥L2⁢(Ω)∥∇(u−vh)∥L2⁢(Ω) (2.65)
=∥u−uh∥E∥u−vh∥E. (2.66)

∎

Teorema 2.2.4.(Estimativa a priori)

A solução uh do problema de elementos finitos (2.42) satisfaz

‖u−uh‖E2≤C⁢∑K∈𝒦hK2⁢‖D2⁢u‖L2⁢(K)2. (2.67)
Demonstração.

Do resultado da melhor aproximação (Teorema 2.2.3) e da estimativa do erro de interpolação (Proposição 2.1.2), segue que

‖u−uh‖E2≤‖u−π⁢u‖E2 (2.68)
=∥D(u−πu)∥L2⁢(Ω)2 (2.69)
≤C∑K∈𝒦hK2∥D2u∥L2⁢(K)2. (2.70)

∎

Para obtermos uma estimativa na norma L2⁢(Ω), podemos usar a desigualdade de Poincaré.

Teorema 2.2.5.(Desigualdade de Poincaré)

Seja Ω⊂ℝ2 um domínio limitado. Então, existe uma constante C=C⁢(Ω), tal que

‖v‖L2⁢(Ω)≤C⁢‖∇v‖L2⁢(Ω),∀v∈V0. (2.71)
Demonstração.

Se Ω tem contorno suficientemente suave, então existe ϕ tal que −Δ⁢ϕ=1 em Ω com supx∈Ω|∇ϕ|<C. Com isso, temos

‖v‖L2⁢(Ω)2=∫Ωv2⁢𝑑x (2.72)
=−∫Ωv2Δϕdx. (2.73)

Agora, usando o Teorema de Green111111George Green, 1793 - 1841, matemático britânico. Fonte: Wikipédia: George Green. e a desigualdade de Cauchy121212Augustin-Louis Cauchy, 1789-1857, matemático francês. Fonte: Wikipédia: Augustin-Louis Cauchy.-Schwarz131313Karl Hermann Amandus Schwarz, 1843-1921, matemático alemão. Fonte: Wikipédia: Hermann Amandus Schwarz., obtemos

‖v‖L2⁢(Ω)2=−∫∂Ωv2⁢(𝒏⋅∇ϕ)⁢𝑑s+∫Ω∇v2⋅∇ϕ⁢d⁢x (2.74)
=∫Ω2v∇v⋅∇ϕdx (2.75)
≤2supx∈Ω|∇ϕ|∥v∥L2⁢(Ω)∥∇v∥L2⁢(Ω). (2.76)

∎

Com a desigualdade de Poincaré e da estimativa a priori do erro (Teorema 2.2.4), temos

‖u−uh‖L2⁢(Ω)≤C⁢h⁢‖D2⁢u‖L2⁢(Ω), (2.77)

onde h é o tamanho global da malha . Entretanto, esta estimativa pode ser melhorada.

Teorema 2.2.6.(Estimativa ótima a priori do erro)

A solução uh do problema de elementos finitos (2.42) satisfaz

‖u−uh‖L2⁢(Ω)≤C⁢h2⁢‖D2⁢u‖L2⁢(Ω). (2.78)
Demonstração.

Seja e=u−uh o erro e ϕ a solução do problema dual (ou problema adjunto)

−Δ⁢ϕ=e,∀x∈Ω (2.79)
ϕ=0,∀x∈∂Ω. (2.80)

Então, usando a fórmula de Green141414George Green, 1793 - 1841, matemático britânico. Fonte: Wikipédia: George Green., a ortogonalidade de Galerkin151515Boris Galerkin, 1871 - 1945, engenheiro e matemático soviético. Fonte: Wikipédia: Boris Galerkin. e, então, a desigualdade de Cauchy161616Augustin-Louis Cauchy, 1789-1857, matemático francês. Fonte: Wikipédia: Augustin-Louis Cauchy.-Schwarz171717Karl Hermann Amandus Schwarz, 1843-1921, matemático alemão. Fonte: Wikipédia: Hermann Amandus Schwarz., temos

‖e‖L2⁢(Ω)2=−∫Ωe⁢Δ⁢ϕ⁢𝑑x (2.81)
=∫Ω∇e⋅∇ϕdx−∫∂Ωe(𝒏⋅∇ϕ)ds (2.82)
=∫Ω∇e⋅∇(ϕ−πϕ)dx (2.83)
≤∥∇e∥L2⁢(Ω)∥∇(ϕ−πϕ)∥L2⁢(Ω). (2.84)

Da estimativa a priori (2.77), temos

‖∇e‖L2⁢(Ω)≤C⁢h⁢‖D2⁢u‖L2⁢(Ω). (2.85)

Agora, da regularidade elíptica ‖D2⁢ϕ‖L2⁢(Ω)≤C⁢‖Δ⁢ϕ‖L2⁢(Ω) [6, Seção 6.3] e da estimativa do erro de interpolação (Proposição 2.1.2), temos

‖∇(ϕ−π⁢ϕ)‖L2⁢(Ω)≤C⁢h⁢‖D2⁢ϕ‖L2⁢(Ω) (2.86)
≤Ch∥Δϕ∥L2⁢(Ω) (2.87)
≤Ch∥e∥L2⁢(Ω). (2.88)

Então, temos

‖e‖L2⁢(Ω)2≤C⁢h⁢‖D2⁢u‖L2⁢(Ω)⁢C⁢h⁢‖e‖L2⁢(Ω). (2.89)

∎

Exemplo 2.2.2.

A Tabela 2.1 apresenta o estudo de convergência de malha para o problema do Exemplo 2.2.1. Verifique!

Tabela 2.1: Erros de aproximações por elementos finitos referente ao problema dado no Exemplo 2.2.2.
#nodos h ‖u−uh‖L2⁢(Ω) Taxa
4×4 3.54⁢e−1 7.91⁢e−2 –
8×8 1.77⁢e−1 2.11⁢e−2 1.90
16×16 8.84⁢e−2 5.38⁢e−3 1.97
32×32 4.42⁢e−2 1.35⁢e−3 1.99
64×64 2.21⁢e−2 3.38⁢e−4 2.00

2.2.4 Exercícios

E. 2.2.1.

Compute uma aproximação de elementos finitos para a solução do seguinte problema de Poisson com condições de Dirichlet

−Δ⁢u=2⁢π2⁢sen⁡(π⁢x)⁢cos⁡(π⁢y),x∈Ω:=(0,1)2, (2.90)
u=0,(x,y)∈{0,1}×(0,1), (2.91)
u=sen⁡(π⁢x),(x,y)∈(0,1)×{0}, (2.92)
u=−sen⁡(π⁢x),(x,y)∈(0,1)×{1}. (2.93)

A solução exata deste problema é dada por u⁢(x,y)=sen⁡(π⁢x)⁢cos⁡(π⁢y). Faça gráfico comparativo entre a solução de elementos finitos e a solução exata, bem como, faça um estudo de convergência de malha.

E. 2.2.2.

Compute uma aproximação de elementos finitos para a solução do seguinte problema de Poisson

−Δ⁢u=π2⁢(x2+y2)⁢cos⁡(π⁢x⁢y),x∈Ω:=(0,1)2, (2.94)
∂u∂x=0,(x,y)∈{0}×(0,1), (2.95)
∂u∂x=−π⁢sen⁡(π⁢x),(x,y)∈{1}×(0,1), (2.96)
u=1,(x,y)∈(0,1)×{0}, (2.97)
u=−cos⁡(π⁢x),(x,y)∈(0,1)×{1}. (2.98)

A solução exata deste problema é dada por u⁢(x,y)=cos⁡(π⁢x⁢y). Faça gráfico comparativo entre a solução de elementos finitos e a solução exata, bem como, faça um estudo de convergência de malha.

E. 2.2.3.

Compute uma aproximação de elementos finitos para a solução do seguinte problema de Poisson com condições de Robin

−Δ⁢u=2⁢π2⁢sen⁡(π⁢x)⁢sen⁡(π⁢y),x∈Ω:=(0,1)2, (2.99)
∂u∂x+u=sen⁡(π⁢x),(x,y)∈{0}×(0,1), (2.100)
∂u∂x+u=−sen⁡(π⁢x),(x,y)∈{1}×(0,1), (2.101)
∂u∂y+u=sen⁡(π⁢x),(x,y)∈(0,1)×{0}, (2.102)
∂u∂y+u=−sen⁡(π⁢x),(x,y)∈(0,1)×{1}. (2.103)

A solução exata deste problema é dada por u⁢(x,y)=sen⁡(π⁢x)⁢sen⁡(π⁢y). Faça gráfico comparativo entre a solução de elementos finitos e a solução exata, bem como, faça um estudo de convergência de malha.

E. 2.2.4.

Compute uma aproximação de elementos finitos para a solução do seguinte problema do calor

∂u∂t−Δ⁢u=2⁢π2⁢(1−e−t)⁢sen⁡(π⁢x)⁢sen⁡(π⁢y)
+e−t⁢sen⁡(π⁢x)⁢cos⁡(π⁢y),(x,y)∈Ω:=(0,1)2,t>0, (2.104)
u=0,(x,y)∈∂Ω,t>0, (2.105)
u⁢(0,x,y)=0,(x,y)∈Ω. (2.106)

A solução exata deste problema é dada por u⁢(x,y,t)=(1−e−t)⁢sen⁡(π⁢x)⁢sen⁡(π⁢y). Faça gráfico comparativo entre a solução de elementos finitos e a solução exata, bem como, faça um estudo de convergência de malha e de tempo.

E. 2.2.5.

Compute uma aproximação de elementos finitos para a solução do seguinte problema

∂u∂t+𝒗⋅∇u=Δ⁢u,(x,y)∈Ω:=(−1,1)2,t>0, (2.107)
u⁢(0,x,y)={1,(x,y)∈(−0.1,0.1)20,caso contrário,,(x,y)∈Ω,t>0, (2.108)
u⁢(x,y,t)=0,(x,y)∈∂Ω. (2.109)

Teste a aproximação de elementos finitos para diferentes campos de velocidade 𝒗, como 𝒗=(1,0), 𝒗=(0,1) e 𝒗=(1,1)/2.


Envie seu comentário

Aproveito para agradecer a todas/os que de forma assídua ou esporádica contribuem enviando correções, sugestões e críticas!

Opcional. Preencha seu nome para que eu possa lhe contatar.
Opcional. Preencha seu e-mail para que eu possa lhe contatar.
As informações preenchidas são enviadas por e-mail para o desenvolvedor do site e tratadas de forma privada. Consulte a política de uso de dados para mais informações.

Licença Creative Commons
Este texto é disponibilizado nos termos da Licença Creative Commons Atribuição-CompartilhaIgual 4.0 Internacional. Ícones e elementos gráficos podem estar sujeitos a condições adicionais.

Pedro H A Konzen
| | | |