Colabore! Saiba mais screen_rotation

Ao navegar por este site , você concorda com a Política de Uso de Dados.

| | | | |

1 Problemas Unidimensionais 1.2 Problema Modelo 1.4 Malhas Auto-Adaptativas

1.3 Condições de Contorno

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

Em revisão

Nesta seção, vamos discutir sobre soluções de elementos finitos para a equações diferencial

\color[rgb]{0,0,1}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,1}% \pgfsys@color@rgb@stroke{0}{0}{1}\pgfsys@color@rgb@fill{0}{0}{1}-u^{\prime% \prime}=f,\quad x\in I=[0,L],

(1.105)

com diferentes condições de contorno.

1.3.1 Condições de Dirichlet

Em revisão

Consideramos o seguinte problema com condições de contorno de Dirichlet⁵⁵endnote: ⁵Johann Peter Gustav Lejeune Dirichlet, 1805 - 1859, matemático alemão. Fonte: Wikipédia: Johann Peter Gustav Lejeune Dirichlet.: encontrar $u$ tal que

	$\displaystyle-u^{\prime\prime}=f,\quad\forall x\in I=[0,L],$		(1.106)
	$\displaystyle u(0)=u_{0},\quad u(L)=u_{L},$		(1.107)

com $u_{0}$ , $u_{L}$ e $f$ dados.

Tomando uma função teste $v\in\color[rgb]{0,0,1}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,1}% \pgfsys@color@rgb@stroke{0}{0}{1}\pgfsys@color@rgb@fill{0}{0}{1}V_{0}:=H^{1}_{% 0}(I):=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=v(L)=0\}$ e multiplicando-a em (1.106), obtemos

-\int_{I}u^{\prime\prime}v\,dx=\int_{I}fv\,dx.

(1.108)

Aplicando a integração por partes, temos

\int_{I}u^{\prime}v^{\prime}\,dx=\int_{I}fv\,dx.

(1.109)

Desta forma, definimos o seguinte problema fraco associado: encontrar $u\in\color[rgb]{0,0,1}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,1}% \pgfsys@color@rgb@stroke{0}{0}{1}\pgfsys@color@rgb@fill{0}{0}{1}V:=\{v\in H^{1% }(I);\leavevmode\nobreak\ v(0)=u_{0},\leavevmode\nobreak\ v(L)=v_{L}\}$ tal que

\displaystyle\color[rgb]{0,0,1}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,1}% \pgfsys@color@rgb@stroke{0}{0}{1}\pgfsys@color@rgb@fill{0}{0}{1}a(u,v)=L(v),% \quad\forall v\in V_{0},

(1.110)

onde $a(u,v)$ é a forma bilinear

a(u,v)=\int_{I}u^{\prime}v^{\prime}\,dx

(1.111)

e $L(v)$ é a forma linear

L(v)=\int_{I}fv\,dx.

(1.112)

Exemplo 1.3.1.

Consideramos o problema

	$\displaystyle-u^{\prime\prime}=1,\quad x\in I=[0,1],$		(1.113)
	$\displaystyle u(0)=1/2,\quad u(1)=1.$		(1.114)

Sua solução analítica é $u(x)=-x^{2}/2+x+1/2$ .

Para obtermos uma aproximação de elementos finitos, consideramos o seguinte problema fraco: encontrar $u\in V:=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=1/2,\leavevmode\nobreak\ v(1% )=1\}$ tal que

a(u,v)=L(v),

(1.115)

para todo $v\in V_{0}=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=v(1)=0\}$ , onde

	$\displaystyle a(u,v)=\int_{I}u^{\prime}v^{\prime}\,dx,$		(1.116)
	$\displaystyle L(v)=\int_{I}fv\,dx.$		(1.117)

Então, o problema de elementos finitos no espaço das funções lineares por partes lê-se: encontrar $u_{h}\in V_{h}=\{v\in P_{1}(I);\leavevmode\nobreak\ v(0)=1/2,\leavevmode% \nobreak\ v(1)=1\}$ tal que

a(u_{h},v_{h})=L(v_{h}),

(1.118)

para todo $v_{h}\in V_{h,0}=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=v(1)=0\}$ .

Código 4: ex_mef1d_dirichlet.py

⬇

1from mpi4py import MPI

3# malha

4from dolfinx import mesh

5domain = mesh.create_unit_interval(MPI.COMM_WORLD,

6 nx = 5)

7# espaço

8from dolfinx import fem

9V = fem.functionspace(domain, ('P', 1))

11# condição de contorno

12import numpy as np

13uD = fem.Function(V)

14def dirichlet_bc(x):

15 y = np.full(x.shape[1], 0.5)

16 y[x[0,:] > 0.5] = 1.

17 return y

18uD.interpolate(dirichlet_bc)

20tdim = domain.topology.dim

21fdim = tdim - 1

22domain.topology.create_connectivity(fdim, tdim)

23boundary_facets = mesh.exterior_facet_indices(domain.topology)

24boundary_dofs = fem.locate_dofs_topological(V, fdim,

25 boundary_facets)

26bc = fem.dirichletbc(uD, boundary_dofs)

28# problema mef

29import ufl

30from dolfinx import default_scalar_type

31from dolfinx.fem.petsc import LinearProblem

32u = ufl.TrialFunction(V)

33v = ufl.TestFunction(V)

35f = fem.Constant(domain, default_scalar_type(1.))

37a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx

38L = f * v * ufl.dx

40problem = LinearProblem(a, L, bcs=[bc])

41uh = problem.solve()

43# armazena para visualização (paraview)

44from dolfinx import io

45from pathlib import Path

46results_folder = Path("results")

47results_folder.mkdir(exist_ok=True, parents=True)

48filename = results_folder / "u"

49with io.XDMFFile(domain.comm, filename.with_suffix(".xdmf"), "w") as xdmf:

50 xdmf.write_mesh(domain)

51 xdmf.write_function(uh)

1.3.2 Condições de Neumann

Em revisão

Consideramos o seguinte problema com condições de contorno de Neumann⁶⁶endnote: ⁶Carl Gottfried Neumann, 1832 - 1925, matemático alemão. Fonte: Wikipédia: Carl Neumann. homogênea em $x=L$ : encontrar $u$ tal que

	$\displaystyle-u^{\prime\prime}=f,\quad\forall x\in I=[0,L],$		(1.119)
	$\displaystyle u(0)=u_{0},\quad u^{\prime}(L)=0,$		(1.120)

com $u_{0}$ e $f$ dados. Trata-se de um problema com condição de contorno de Dirichlet à esquerda e condição de contorno de Neumann⁷⁷endnote: ⁷Carl Gottfried Neumann, 1832 - 1925, matemático alemão. Fonte: Wikipédia: Carl Neumann. homogênea à direita.

Tomando uma função teste $v\in V:=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=0\}$ e multiplicando-a em (1.119), obtemos

-\int_{I}u^{\prime\prime}v\,dx=\int_{I}fv\,dx.

(1.121)

Aplicando a integração por partes, temos

\int_{I}u^{\prime}v^{\prime}\,dx-\underbrace{u^{\prime}(L)v(L)}_{u^{\prime}(L)% =0}+\underbrace{u^{\prime}(0)v(0)}_{v(0)=0}=\int_{I}fv\,dx.

(1.122)

Desta forma, definimos o seguinte problema fraco associado: encontrar $u\in\tilde{V}:=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=u_{0}\}$ tal que

\displaystyle a(u,v)=L(v),\quad\forall v\in V,

(1.123)

onde $a(u,v)$ é a forma bilinear

a(u,v)=\int_{I}u^{\prime}v^{\prime}\,dx

(1.124)

e $L(v)$ é a forma linear

L(v)=\int_{I}fv\,dx.

(1.125)

Exemplo 1.3.2.

Consideramos o problema

	$\displaystyle-u^{\prime\prime}=1,\quad x\in I=[0,1],$		(1.126)
	$\displaystyle u(0)=0,\quad u^{\prime}(1)=0.$		(1.127)

Sua solução analítica é $u(x)=-x^{2}/2+x$ .

Podemos construir uma aproximação por elementos finitos do seguinte problema fraco associado: encontrar $u\in V=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=0\}$ tal que

a(u,v)=L(v),

(1.128)

para todo $v\in V$ , com as formas bilinear $a(\cdot,\cdot)$ e linear $L(\cdot)$ dadas em (1.124) e (1.125).

Então, considerando elementos lineares por partes, temos o seguinte problema de elementos finitos: encontrar $u_{h}\in V_{h}=\{v_{h}\in P_{1}(I);\leavevmode\nobreak\ v_{h}(0)=0\}$ tal que

a(u_{h},v_{h})=L(v_{h}),\quad\forall v_{h}\in V_{h}.

(1.129)

Código 5: ex_mef1d_neumann.py

⬇

1from mpi4py import MPI

3# malha

4from dolfinx import mesh

5domain = mesh.create_unit_interval(MPI.COMM_WORLD,

6 nx = 5)

7# espaço

8from dolfinx import fem

9V = fem.functionspace(domain, ('P', 1))

11# c.c. dirichlet

12import numpy as np

13from dolfinx.fem import dirichletbc, locate_dofs_geometrical

14uD = fem.Function(V)

15uD.interpolate(lambda x: np.full(x.shape[1], 0.))

17def boundary_D(x):

18 return np.isclose(x[0], 0.)

20dofs_D = locate_dofs_geometrical(V, boundary_D)

21bc = dirichletbc(uD, dofs_D)

23# problema mef

24import ufl

25from dolfinx import default_scalar_type

26from dolfinx.fem.petsc import LinearProblem

27u = ufl.TrialFunction(V)

28v = ufl.TestFunction(V)

30f = fem.Constant(domain, default_scalar_type(1.))

32a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx

33L = f * v * ufl.dx

35problem = LinearProblem(a, L, bcs=[bc])

36uh = problem.solve()

38# armazena para visualização (paraview)

39from dolfinx import io

40from pathlib import Path

41results_folder = Path("results")

42results_folder.mkdir(exist_ok=True, parents=True)

43filename = results_folder / "u"

44with io.XDMFFile(domain.comm, filename.with_suffix(".xdmf"), "w") as xdmf:

45 xdmf.write_mesh(domain)

46 xdmf.write_function(uh)

Agora, consideramos o seguinte problema com condições de Neumann não-homogênea em $x=L$ : encontrar $u$ tal que

	$\displaystyle-u^{\prime\prime}=f,\quad\forall x\in I=[0,L],$		(1.130)
	$\displaystyle u(0)=u_{0},\quad u^{\prime}(L)=\alpha,$		(1.131)

com $u_{0}$ , $\alpha$ e $f$ dados.

Tomando uma função teste $v\in V:=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=0\}$ e multiplicando-a em (1.130), obtemos

-\int_{I}u^{\prime\prime}v\,dx=\int_{I}fv\,dx.

(1.132)

Aplicando a integração por partes, temos

\int_{I}u^{\prime}v^{\prime}\,dx-\alpha v(L)=\int_{I}fc\,dx.

(1.133)

Desta forma, definimos o seguinte problema fraco associado: encontrar $u\in\tilde{V}:=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=u_{0}\}$ tal que

\displaystyle a(u,v)=L(v),\quad\forall v\in V,

(1.134)

onde $a(u,v)$ é a forma bilinear

a(u,v)=\int_{I}u^{\prime}v^{\prime}\,dx

(1.135)

e $L(v)$ é a forma linear

L(v)=\int_{I}fv\,dx+\alpha v(L).

(1.136)

Exemplo 1.3.3.

Consideramos o problema

	$\displaystyle-u^{\prime\prime}=1,\quad x\in I=[0,1],$		(1.137)
	$\displaystyle u(0)=0,\quad u^{\prime}(1)=1.$		(1.138)

Sua solução analítica é $u(x)=-x^{2}/2+2x$ .

Agora, consideramos o seguinte problema fraco associado: encontrar $u\in V=\{v\in H^{1}(I);\leavevmode\nobreak\ v(0)=0\}$ tal que

a(u,v)=L(v),\quad\forall v\in V,

(1.139)

com

a(u,v)=\int_{I}u^{\prime}v^{\prime}\,dx

(1.140)

L(v)=\int_{I}fv\,dx+1\cdot v(1).

(1.141)

Então, consideramos o seguinte problema de elementos finitos associado: encontrar $u_{h}\in V_{h}=\{v_{h}\in P_{1}(I);\leavevmode\nobreak\ v_{h}(0)=0\}$ tal que

a(u_{h},v_{h})=L(v_{h}),\quad\forall v_{h}\in V_{h}.

(1.142)

Código 6: ex_mef1d_neumann_nh.py

⬇

1from mpi4py import MPI

3# malha

4from dolfinx import mesh

5domain = mesh.create_unit_interval(MPI.COMM_WORLD,

6 nx = 5)

7# espaço

8from dolfinx import fem

9V = fem.functionspace(domain, ('P', 1))

11# c.c. dirichlet

12import numpy as np

13from dolfinx.fem import dirichletbc, locate_dofs_geometrical

14uD = fem.Function(V)

15uD.interpolate(lambda x: np.full(x.shape[1], 0.))

17def boundary_D(x):

18 return np.isclose(x[0], 0.)

20dofs_D = locate_dofs_geometrical(V, boundary_D)

21bc = dirichletbc(uD, dofs_D)

23# problema mef

24import ufl

25from dolfinx import default_scalar_type

26from dolfinx.fem.petsc import LinearProblem

27u = ufl.TrialFunction(V)

28v = ufl.TestFunction(V)

30f = fem.Constant(domain, default_scalar_type(1.))

31g = fem.Constant(domain, default_scalar_type(1.))

33a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx

34L = f * v * ufl.dx

35L += g * v * ufl.ds

37problem = LinearProblem(a, L, bcs=[bc])

38uh = problem.solve()

40# armazena para visualização (paraview)

41from dolfinx import io

42from pathlib import Path

43results_folder = Path("results")

44results_folder.mkdir(exist_ok=True, parents=True)

45filename = results_folder / "u"

46with io.XDMFFile(domain.comm, filename.with_suffix(".xdmf"), "w") as xdmf:

47 xdmf.write_mesh(domain)

48 xdmf.write_function(uh)

1.3.3 Condições de Robin

Em revisão

Consideramos o seguinte problema com condições de contorno de Robin⁸⁸endnote: ⁸Victor Gustave Robin, 1855 - 1897, matemático francês. Fonte: Wikipedia: Victor Gustave Robin.: encontrar $u$ tal que

	$\displaystyle-u^{\prime\prime}=f,\quad\forall x\in I=[0,L],$		(1.143)
	$\displaystyle u^{\prime}(0)=r_{0}(u(0)-s_{0}),\leavevmode\nobreak\ -u^{\prime}% (L)=r_{L}(u(L)-s_{L}),$		(1.144)

com $r_{0}$ , $r_{L}$ , $s_{0}$ , $s_{L}$ e $f$ dados.

Tomando uma função teste $v\in V=H^{1}(I)$ e multiplicando-a em (1.143), obtemos

-\int_{I}u^{\prime\prime}v\,dx=\int_{I}fv\,dx.

(1.145)

Aplicando a integração por partes, temos

\int_{I}u^{\prime}v^{\prime}\,dx-\underbrace{u^{\prime}(L)v(L)}_{-u^{\prime}(L% )=r_{L}(u(L)-s_{L})}+\underbrace{u^{\prime}(0)v(0)}_{u^{\prime}(0)=r_{0}(u(0)-% s_{0})}=\int_{I}fc\,dx.

(1.146)

ou, mais adequadamente,

\int_{I}u^{\prime}v^{\prime}\,dx+r_{L}u(L)v(L)+r_{0}u(0)v(0)=\int_{I}fc\,dx+r_% {L}s_{L}v(L)+r_{0}s_{0}v(0).

(1.147)

Desta forma, definimos o seguinte problema fraco associado: encontrar $u\in H^{1}(I)$ tal que

\displaystyle a(u,v)=L(v),\quad\forall v\in V,

(1.148)

onde $a(u,v)$ é a forma bilinear

a(u,v)=\int_{I}u^{\prime}v^{\prime}\,dx+r_{L}u(L)v(L)+r_{0}u(0)v(0)

(1.149)

e $L(v)$ é a forma linear

L(v)=\int_{I}fv\,dx+r_{L}s_{L}v(L)+r_{0}s_{0}v(0).

(1.150)

Exemplo 1.3.4.

Consideramos o problema

	$\displaystyle-u^{\prime\prime}=1,\quad x\in I=[0,1],$		(1.151)
	$\displaystyle u^{\prime}(0)=u(0),\quad-u^{\prime}(1)=u(1)-1.$		(1.152)

Sua solução analítica é $u(x)=-x^{2}/2+5x/6+5/6$ .

Aqui, tomamos o seguinte problema fraco: encontrar $u\in V=H^{1}(I)$ tal que

a(u,v)=L(v),\quad\forall v\in V,

(1.153)

onde

a(u,v)=\int_{I}u^{\prime}v^{\prime}\,dx+u(1)v(1)+u(0)v(0)

(1.154)

L(v)=\int_{I}fv\,dx+1\cdot v(1).

(1.155)

Então, uma aproximação por elementos finitos lineares por partes pode ser obtida resolvendo o seguinte problema: encontrar $u_{h}\in V_{h}=P_{1}(I)$ tal que

a(u_{h},v_{h})=L(v_{h}),\quad\forall v_{h}\in V_{h}.

(1.156)

⬇

1from mpi4py import MPI

3# malha

4from dolfinx import mesh

5domain = mesh.create_unit_interval(MPI.COMM_WORLD,

6 nx = 5)

7# espaço

8from dolfinx import fem

9V = fem.functionspace(domain, ('P', 1))

11# boundary colors

12from dolfinx.mesh import locate_entities

13from dolfinx.mesh import meshtags

14boundaries = [(0, lambda x: np.isclose(x[0], 0.)),

15 (1, lambda x: np.isclose(x[0], 1.))]

16facet_indices, facet_markers = [], []

17fdim = domain.topology.dim - 1

18for (marker, locator) in boundaries:

19 facets = locate_entities(domain, fdim, locator)

20 facet_indices.append(facets)

21 facet_markers.append(np.full_like(facets, marker))

22facet_indices = np.hstack(facet_indices).astype(np.int32)

23facet_markers = np.hstack(facet_markers).astype(np.int32)

24sorted_facets = np.argsort(facet_indices)

25facet_tag = meshtags(domain, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets])

27# problema mef

28import ufl

29from dolfinx import default_scalar_type

30from dolfinx.fem.petsc import LinearProblem

31u = ufl.TrialFunction(V)

32v = ufl.TestFunction(V)

34f = fem.Constant(domain, default_scalar_type(1.))

35g = fem.Constant(domain, default_scalar_type(1.))

37ds = ufl.Measure('ds', domain=domain, subdomain_data=facet_tag)

39a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx

40a += u * v * ds(1) + u * v * ds(0)

41L = f * v * ufl.dx

42L += g * v * ds(1)

44problem = LinearProblem(a, L, bcs=[])

45uh = problem.solve()

47# armazena para visualização (paraview)

48from dolfinx import io

49from pathlib import Path

50results_folder = Path("results")

51results_folder.mkdir(exist_ok=True, parents=True)

52filename = results_folder / "u"

53with io.XDMFFile(domain.comm, filename.with_suffix(".xdmf"), "w") as xdmf:

54 xdmf.write_mesh(domain)

55 xdmf.write_function(uh)

1.3.4 Exercícios

Em revisão

E. 1.3.1.

Considere o problema

	$\displaystyle-u^{\prime\prime}+u^{\prime}+2u=-\cos(x),\quad x\in(0,\pi/2),$		(1.157)
	$\displaystyle u(0)=-0,3,\quad u(\pi/2)=-0,1.$		(1.158)

Obtenha uma aproximação por elementos finitos para a solução deste problema, empregando o espaço de elementos finitos linear sobre uma malha uniforme com $10$ células. Então, compare a aproximação computada com sua solução analítica $u(x)=0,1(\operatorname{sen}(x)+3\cos(x))$ , bem como, compute o erro $\|u-u_{h}\|_{L^{2}}$ .

Resposta.

Código.

Envie seu comentário

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. Aproveito para agradecer a todas/os que de forma assídua ou esporádica contribuem enviando correções, sugestões e críticas!

| | | |