Buscar

Método das Diferenças Finitas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você viu 3, do total de 17 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você viu 6, do total de 17 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Você viu 9, do total de 17 páginas

Faça como milhares de estudantes: teste grátis o Passei Direto

Esse e outros conteúdos desbloqueados

16 milhões de materiais de várias disciplinas

Impressão de materiais

Agora você pode testar o

Passei Direto grátis

Prévia do material em texto

O Me´todo das Diferenc¸as Finitas
Equac¸o˜es na˜o lineares (Navier-Stokes)
Prof: Rafael Gabler Gontijo
November 29, 2014
Abstract
Esse documento apresenta um esquema utilizado para soluc¸a˜o de sis-
temas de equac¸o˜es na˜o lineares. O objetivo e´ que o aluno possa uti-
lizar o Me´todo de Newton para resolver um sistema de equac¸o˜es comum
na Mecaˆnica dos Fluidos, associado a` forma discretizada da equac¸a˜o de
Navier-Stokes em um problema de camada limite. O esquema se baseara´
em um modelo dida´tico de uma malha mais simples contendo 16 no´s a
t´ıtulo de ilustrac¸a˜o.
1 Introduc¸a˜o e motivac¸a˜o
Em muitos problemas da f´ısica precisamos resolver equac¸o˜es diferenciais or-
dina´rias ou parciais para entendermos como determinada grandeza se com-
porta. Muitas leis fundamentais da f´ısica cla´ssica sa˜o expressas em linguagem
matema´tica atrave´s de equac¸o˜es diferenciais. Essas equac¸o˜es geralmente esta˜o
associadas ao comportamento de func¸o˜es escalares, vetoriais ou tensoriais, no
espac¸o e no tempo. Nesse sentido o conceito de campo e´ muito utilizado para
descrever o comportamento de um sistema f´ısico.
Quando um determinado sistema f´ısico e´ descrito em linguagem matema´tica
atrave´s de equac¸o˜es diferenciais que possuem soluc¸a˜o anal´ıtica, esta e´ a te´cnica
utilizada para entender o comportamento deste. A vantagem de uma soluc¸a˜o
anal´ıtica de uma equac¸a˜o diferencial e´ que esta e´ exata e fornece uma expressa˜o
matema´tica como resposta. Esta expressa˜o por sua vez pode ser interpretada
para gerar uma teoria f´ısica. E´ um consenso entre os principais cientistas de
que a soluc¸a˜o anal´ıtica e´ sem du´vidas a busca primeira e primordial para o
desenvolvimento de uma teoria cient´ıfica. Entrentanto, frequentemente nos de-
paramos com problemas regidos por equac¸o˜es diferenciais com elevado grau
de complexidade que impedem o desenvolvimento de uma soluc¸a˜o anal´ıtica.
Nesses casos podemos abordar o problema com dois focos poss´ıveis: uma soluc¸a˜o
nume´rica (computacional) ou um estudo experimental para buscar correlac¸o˜es
emp´ıricas que representem o comportamento do sistema em estudo.
A vantagem de uma abordagem experimental consiste no fato de que se o
experimento for conduzido de forma meticulosa e com o menor grau poss´ıvel
de interfereˆncia do sistema de medic¸a˜o no sistema, ela representa uma exce-
lente visa˜o do que poder´ıamos chamar de realidade. As desvantagens de uma
abordagem experimental sa˜o: necessidade de ana´lises estat´ısticas custosas para
1
garantir a repetibilidade dos resultados, alto custo (dependendo da grandeza
f´ısica a ser mensurada), interfereˆncia do dispositivo de medida (senores, sondas)
no sistema, em muitos casos as medic¸o˜es sa˜o indiretas e dependem da inter-
pretac¸a˜o do experimentalista, dificuldade de controle de fatores externos que
podem influenciar no experimento (som, luz, temperatura, pressa˜o, ru´ıdos da
rede ele´trica, impurezas, entre outros).
A abordagem nume´rica muitas vezes e´ o caminho preferido por muitos ci-
entistas (ainda mais com a enorme evoluc¸a˜o tecnlo´gica recente que impulsiona
o poder de ca´lculo de nossos computadores modernos) por apresentar vanta-
gens como: facilidade no controle dos paraˆmetros f´ısicos do problema, acesso
facilitado a inu´meras ferramentas de pre´ e po´s processamento de dados, di-
versidade de linguagens de programac¸a˜o, acesso amplo a` inu´meras te´cnicas de
discretizac¸a˜o espacial para soluc¸a˜o de problemas de valor de contorno, entre
outros. Evidamentemente o me´todo nume´rico tambe´m apresenta desvantagens,
como: depedeˆncia dos resultados de acordo com a resoluc¸a˜o espacial e temporal
adotada, necessidade eventual de calibrac¸a˜o de paraˆmetros em aproximac¸o˜es de
func¸o˜es cont´ınuas por valores discretos e (provavelmente a principal delas) a
impossibilidade de obtenc¸a˜o de uma expressa˜o matema´tica geral como soluc¸a˜o
da equac¸a˜o diferencial a ser resolvida. Tudo que o me´todo nume´rico (como
o pro´prio nome diz) fornece sa˜o nu´meros. Nu´meros estes que podem ser or-
ganizados na forma de tabelas, gra´ficos e imagens de campo para estimar o
comportamento f´ısico aproximado de determinado sistema.
Figure 1: Exemplo do uso de um me´todo computacional para a soluc¸a˜o da
equac¸a˜o de Navier-Stokes para a determinac¸a˜o do campo de velocidades de um
escoamento na sa´ıda de uma turbina.
De um modo geral a lo´gica principal por tra´s das diversas estrate´gias de
soluc¸a˜o de equac¸o˜es diferenciais cont´ınuas atrave´s de um me´todo nume´rico con-
siste em discretizar um domı´nio espacial de ca´lculo cont´ınuo (com infinitos pon-
tos) em um domı´nio espacial discreto (com um nu´mero finito de pontos). Essa
transformac¸a˜o tambe´m se aplica ao tempo (em problemas transientes) que em
2
sistemas reais flui de forma cont´ınua, mas em sistemas nume´ricos deve ser dis-
cretizado em diferentes instantes sucessivos descont´ınuos.
Essa transformac¸a˜o cont´ınuo-discreto e´ tambe´m responsa´vel pela trans-
formac¸a˜o da necessidade de soluc¸a˜o de uma u´nica equac¸a˜o diferencial em N
equac¸o˜es alge´bricas para cada instante de tempo (aqui N representa o nu´mero
de no´s da malha que representa o espac¸o da soluc¸a˜o discretizado). Cada equac¸a˜o
alge´brica se refere a` maneira com a qual o valor de uma determinada varia´vel em
um no´ i qualquer da malha de ca´lculo se relaciona com seus no´s vizinhos. Essa
relac¸a˜o evidentemente e´ dada pela discretizac¸a˜o da equac¸a˜o diferencial gover-
nante do problema, que por sua vez e´ deduzida com base em um (ou um conjunto
de) princ´ıpio(s) f´ısico(s). Como um no´ de refereˆncia qualquer da malha possui
ao menos 4 vizinhos (em sistemas bidimensionais), este no´ em uma ana´lise pos-
terior sera´ o vizinho de um outro no´ de refereˆncia. Por este motivo a conexa˜o
entre os no´s da malha de ca´lculo e´ intricada, de modo que um no´ depende do
outro e vice-versa. Por esse motivo esse conjunto de N equac¸o˜es discretas (resul-
tado da transformac¸a˜o cont´ınuo-discreto da equac¸a˜o diferencial governante)
e´ comumente expresso atrave´s de um sistema linear de equac¸o˜es. Um exemplo
cla´ssico e´ o da equac¸a˜o da conduc¸a˜o de calor bidimensional permanente. Para
esta equac¸a˜o temos a seguinte expressa˜o cont´ınua:
∇2T = 0 →
∂2T
∂x2
+
∂2T
∂y2
= 0, eq. cont´ınua, (1)
em que ∇ e´ o operador nabla, T e´ o campo de temperatura, x e y representam
coordenadas espaciais Eulerianas. Para uma malha na qual o espac¸amento entre
os no´s na direc¸a˜o x, dado por ∆x e´ igual ao espac¸amento y dado por ∆y, a versa˜o
discretizada (via me´todo das diferenc¸as finitas) da equac¸a˜o (1) e´ dada por
Ti −
1
4
(Ti−1 + Ti+1 + Ti−n + Ti+n) = 0 (2)
aqui n se refere ao nu´mero de no´s na direc¸a˜o x. A equac¸a˜o (2) pode ser escrita
na forma de um sistema linear para os no´s internos como


· · ·
...
...
...
...
...
...
... · · ·
· · · −0.25 · · · −0.25 1 −0.25 · · · −0.25 · · ·
· · ·
...
...
...
...
...
...
... · · ·

 ·


...
Ti−n
...
Ti−1
Ti
Ti+1
...
Ti+n
...


=


...
...
...
...
0
...
...
...
...


(3)
O preenchimento das outras linhas da matriz A (primeiro termo do lado es-
querdo), do vetor x (segundo termo do lado esquerdo) e do vetor B (termo do
lado direito da equac¸a˜o) dependem das condic¸o˜es de contorno do problema.
3
De qualquer forma a soluc¸a˜o da equac¸a˜o cont´ınua (1) em um domı´nio dis-
cretizado com N no´s passa pela soluc¸a˜o de um sistema linear de ordem N do
tipo A·x = B. Existem diversos me´todos dispon´ıveis para a soluc¸a˜o de sistemas
lineares. Alguns desses me´todos demandam a inversa˜o da matriz A (proced-
imento muito caro do ponto de vista computacional), outrosexigem soluc¸o˜es
iterativas.
O problema da soluc¸a˜o de um sistema linear por me´todos convencionais
para soluc¸a˜o destes tipos de sistema e´ que estas sa˜o va´lidas para a soluc¸a˜o de
sistemas de equac¸o˜es discretizadas oriundos de equac¸o˜es diferenciais cont´ınuas
lineares. Em outras palavras: na˜o sa˜o va´lidas para sistemas na˜o lineares. Algu-
mas equac¸o˜es diferenciais importantes que regem uma gama enorme de sistemas
f´ısicos sa˜o na˜o lineares. Um exemplo cla´ssico e´ a equac¸a˜o de Navier-Stokes que,
para um fluido incompress´ıvel, e´ dada por
ρ
(
∂u
∂t
+ u · ∇u
)
= −∇p+ η∇2u+ ρg, (4)
em que ρ e η representam a massa espec´ıfica e a viscosidade dinaˆmica do fluido
respectivamente, p e u denotam os campos Eulerianos de pressa˜o e velocidade
do escoamento, t e´ a varia´vel tempo e g representa o vetor acelerac¸a˜o do campo
gravitacional. Note que o termo u ·∇u torna a equac¸a˜o (4) uma EDP na˜o linear
de segunda ordem na˜o homogeˆnea. O elevado grau de complexidade da equac¸a˜o
(4) nos impossibilita uma soluc¸a˜o anal´ıtica geral, de modo que esta so´ apresenta
soluc¸a˜o em casos assinto´ticos muito particulares. Dessa forma muitas vezes a
Mecaˆnica dos Fluidos avanc¸a em termos de compreensa˜o f´ısica atrave´s de sim-
ulac¸o˜es computacionais baseadas na soluc¸a˜o da equac¸a˜o (4) com condic¸o˜es ini-
ciais e de contorno bem definidas (problema bem posto) via me´todos nume´ricos
de discretizac¸a˜o temporal. O grande inconveniente de discretizarmos a equac¸a˜o
(4) cheia (com todos os termos) e´ o surgimento de um sistema de equac¸o˜es na˜o
lineares.
Este documento apresenta um esquema dida´tico mostrando como pode-
mos resolver uma equac¸a˜o diferencial parcial na˜o linear atrave´s do me´todo das
diferenc¸as finitas combinado com te´cnicas de soluc¸a˜o de sistemas de equac¸o˜es
na˜o lineares.
2 O me´todo de Newton - ra´ızes de func¸o˜es na˜o
lineares
Existem diversos me´todos utilizados para a soluc¸a˜o de equac¸o˜es na˜o lineares,
que va˜o desde uso de algoritmos gene´ticos, me´todo da descida ı´ngreme, me´todos
de ponto fixo, entre outros. Um me´todo de fa´cil compreensa˜o e ao mesmo tempo
muito eficaz consiste no me´todo de Newton. Geralmente alunos de cursos de
ca´lculo nume´rico esta˜o familiarizados com o me´todo de Newton para a deter-
minac¸a˜o de zeros de func¸o˜es e em certo sentido e´ essa a ideia por tra´s de seu uso
para resolver sistemas de equac¸o˜es na˜o lineares. A ideia do me´todo de Newton
e´ que uma func¸a˜o f(x) qualquer (linear ou na˜o) pode ser aproximada por uma
se´rie de Taylor no ponto x = x0 +∆x por seu valor e de sua derivada no ponto
x0 como
f(x) ≈ f(x0) +
df
dx
(x0)∆x (5)
4
aqui ∆x = x − x0. Caso desejemos encontrar alguma ra´ız da func¸a˜o f(x),
ou seja, determinar o valor (ou valores) de x que fazem com que f(x) = 0,
igualamos a equac¸a˜o (5) e isolamos x, este procedimento fornece
x = x0 −
f(x0)
f
′(x0)
(6)
aqui f
′
(x0) denota a primeira derivada de f com respeito a x avaliada em x0.
Como o me´todo se baseia na aproximac¸a˜o dada por (5) o uso da equac¸a˜o (6) na˜o
fornecera´ a resposta exata com apenas uma u´nica computac¸a˜o. Na realidade o
me´todo de Newton e´ um me´todo iterativo, em que
xi+1 = xi −
f(xi)
f
′(xi)
, i = 0, 1, 2, · · · (7)
A partir do momento em que a diferenc¸a entre dois valores consecutivos de x
for suficientemente pequena (dentro de uma toleraˆncia pre´-estabelecida) dizemos
que o me´todo convergiu. Considere por exemplo que desejemos determinar uma
ra´ız da func¸a˜o dada por
f(x) = x2 − 4, (8)
aplicando o me´todo descrito em (7) para a equac¸a˜o (8) obtemos
xi+1 = xi −
xi
2
− 4
2xi
, i = 0, 1, 2, · · · (9)
Sabemos que a equac¸a˜o (8) possui duas ra´ızes reais, dadas por 2 e -2. Con-
siderando x0 = 0.5 podemos montar a seguinte tabela:
i xi xi+1 |xi+1 − xi|
0 0.5 4.25 3.75
1 4.25 2.5956 1.6544
2 2.5956 2.0683 0.5273
3 2.0683 2.0011 0.0672
4 2.0011 2.0000003183 0.0011284441
5 2.0000003183 2.0 3.18 ×10−7
Table 1: Determinac¸a˜o de uma das ra´ızes de uma equac¸a˜o na˜o linear alge´brica
de segundo grau simples via Me´todo de Newton.
Note que a partir da sexta iterac¸a˜o a diferenc¸a entre dois valores consecutivos
da ra´ız obtida e´ da ordem de 10−7. Essa tabela mostra a efica´cia do me´todo e sua
ra´pida convergeˆncia. O uso do me´todo de Newton para a soluc¸a˜o de sistemas na˜o
lineares se baseia na generalizac¸a˜o deste para encontrar N ra´ızes de N equac¸o˜es
na˜o lineares organizadas na forma de um sistema na˜o linear. Um sistema na˜o
linear de N equac¸a˜o acopladas pode ser organizado em termos vetoriais como
uma func¸a˜o f(x) = 0. Essa func¸a˜o f e´ um vetor com N componentes em que
N e´ o nu´mero de equac¸o˜es a ser resolvida. Como exemplo considere o seguinte
sistema de duas equac¸o˜es na˜o lineares
f(x) =
(
f1
f2
)
=
(
x1 + 2x2 − 2
x21 + 4x
2
2
)
=
(
0
0
)
, (10)
5
considere ainda a confecc¸a˜o da seguinte matriz J
J =


∂f1
∂x1
∂f1
∂x2
∂f2
∂x1
∂f2
∂x2

 =

 1 2
2x1 8x2

 . (11)
Essa matriz e´ denominada matriz Jacobiana do sistema e e´ utilizada para escr-
ever a equac¸a˜o de recorreˆncia (9) em forma vetorial como
xi+1 = xi − J(xi)
−1
· f(xi). (12)
Note que na equac¸a˜o (12) surge a matriz inversa de J . Isso na˜o significa
necessariamente que precisamos inverter a matriz J para resolver o sistema.
Chamando o termo −J(xi)
−1
· f(xi) de si temos
si = −J(xi)
−1
· f(xi) (13)
multiplicando os dois lados da equac¸a˜o (13) por J(xi) obtemos
J(xi) · si = −J(xi) · J(xi)
−1
· f(xi) = −I · f(xi) = −f(xi), (14)
dessa forma temos
xi+1 = xi + si, i = 0, 1, 2, · · · (15)
com si dado pela soluc¸a˜o do sistema linear
J(xi) · si = −f(xi). (16)
Para o exemplo em questa˜o considere um vetor inicial x0 = [1 2]T , dessa forma
f(x0) = [3 13]T e
J(x0) =

 1 2
2 16

 (17)
Resolvendo o sistema 
 1 2
2 16

 · s0 =
(
−3
13
)
, (18)
obtemos s0 = [−1.83 −0.58]T e x1 = [−0.83 1.42]T . Para um pro´ximo ponto
temos f(x1) = [0 4.72]T e
J(x1) =

 1 2
−1.67 11.3

 (19)
Resolvendo o sistema
 1 2
−1.67 11.3

 · s1 =
(
0
4.72
)
, (20)
obtemos s1 = [0.64 − 0.32]T e x1 = [−0.19 1.1]T , repetindo o procedimento
diversas vezes a soluc¸a˜o convirgira´ para x = [0 1]T .
6
Note que o uso do Me´todo de Newton para a soluc¸a˜o de sistemas na˜o lineares
consiste na transformac¸a˜o de um u´nico sistema na˜o linear em diversos sitemas
lineares resolvidos de modo progressivo em um procedimento iterativo ate´ a
convergeˆncia. A ideia por tra´s deste me´todo e´ a de que uma func¸a˜o na˜o linear
pode ser aproximada como uma soma de diferentes func¸o˜es lineares em torno de
pequenos intervalos. A pro´xima sec¸a˜o mostra um exemplo do uso deste me´todo
na soluc¸a˜o da equac¸a˜o de Navier-Stokes para um problema de camada limite
laminar sobre placa plana de um fluido incompress´ıvel.
3 Exemplo - Soluc¸a˜o de uma EDP na˜o linear
Considere a seguinte equac¸a˜o diferencial parcial
u
∂u
∂x
+ v
∂u
∂x
=
1
Re
∂2u
∂y2
(21)
a equac¸a˜o (21) e´ conhecida como equac¸a˜o de Navier-Stokes (sua componente
apenas na direc¸a˜o x para um problema de camada limite). Considere a seguinte
imagem ilustrativa que representa o domı´nio de ca´lculo discretizado no qual
desejamos resolver a equac¸a˜o (21). A figura (3) ilustra um modelo esquema´tico
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
����������������������������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
�����������������
1 2 3 4
5 6 7 8
9 10 11 12
14 15 1613
u = 1.0
u = 0.0
u = 1.0
x
y
i i+1i−1
i+n
i−n
du/dx = 0
Figure 2: O domı´nio de ca´lculo discretizado consiste em uma malha com 16
no´s, sendo 4 internos e 12 externos sujeitos a`s condic¸o˜es de contorno descritas
na figura.
do problema de camda limite a ser abordado.
As func¸o˜es escalares u e v sa˜o func¸o˜es de x e y (coordenadas espaciais), de
modo que
u = u(x, y) e v = v(x, y) (22)
A versa˜o discretizada da equac¸a˜o (21) e´ dada por
ui
(
ui − ui−1
)
∆x
+ vi
(
ui − ui−n
)
∆x
=
1
Re
(
ui+n + ui−n − 2ui
)
∆y2
(23)
7
Figure 3: Modelo esquema´tico do problema de camada limite hidrodinaˆmica.
isolando ui temos
ui
(
ui
∆x
+
vi
∆y
+
2
Re∆y2
−
ui−1
∆x
)
−
viui−n
∆y
−
ui+n
Re∆y2
−
ui−n
Re∆y2
= 0 (24)
Note que a equac¸a˜o acima se refere a` um no´ i qualquer e a maneira com a qual
este se relaciona com seus vizinhos. Para fins de simplificac¸a˜o considere o caso
em que v = 0, de modo que
ui
(
ui
∆x
+
2
Re∆y2
−
ui−1
∆x
)
−
ui+n
Re∆y2
−
ui−n
Re∆y2
= 0 (25)
Podemos escrever esse sistema de equac¸o˜es na˜o lineares para a malha em
questa˜o (com 16 no´s) da forma
8
f(u) =


u1
u2
u3
u4
u5 − 1
u2
6
∆x
− u6u5
∆x
+ 2u6
Re∆y2
− u10
Re∆y2
− u2
Re∆y2
u2
7
∆x
− u7u6
∆x
+ 2u7
Re∆y2
− u11
Re∆y2
− u3
Re∆y2
u8 − u7
u9 − 1
u2
10
∆x
− u10u9
∆x
+ 2u10
Re∆y2
− u14
Re∆y2
− u6
Re∆y2
u2
11
∆x
− u11u10
∆x
+ 2u11
Re∆y2
− u15
Re∆y2
− u7
Re∆y2
u12 − u11
u13 − 1
u14 − 1
u15 − 1
u16 − 1


=


0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0


(26)
Considere a matriz Jacobiana J do sistema definida como
9
J =


∂f1
∂u1
∂f1
∂u2
∂f1
∂u3
∂f1
∂u4
∂f1
∂u5
∂f1
∂u6
∂f1
∂u7
∂f1
∂u8
∂f1
∂u9
∂f1
∂u10
∂f1
∂u11
∂f1
∂u12
∂f1
∂u13
∂f1
∂u14
∂f1
∂u15
∂f1
∂u16
∂f2
∂u1
∂f2
∂u2
∂f2
∂u3
∂f2
∂u4
∂f2
∂u5
∂f2
∂u6
∂f2
∂u7
∂f2
∂u8
∂f2
∂u9
∂f2
∂u10
∂f2
∂u11
∂f2
∂u12
∂f2
∂u13
∂f2
∂u14
∂f2
∂u15
∂f2
∂u16
∂f3
∂u1
∂f3
∂u2
∂f3
∂u3
∂f3
∂u4
∂f3
∂u5
∂f3
∂u6
∂f3
∂u7
∂f3
∂u8
∂f3
∂u9
∂f3
∂u10
∂f3
∂u11
∂f3
∂u12
∂f3
∂u13
∂f3
∂u14
∂f3
∂u15
∂f3
∂u16
∂f4
∂u1
∂f4
∂u2
∂f4
∂u3
∂f4
∂u4
∂f4
∂u5
∂f4
∂u6
∂f4
∂u7
∂f4
∂u8
∂f4
∂u9
∂f4
∂u10
∂f4
∂u11
∂f4
∂u12
∂f4
∂u13
∂f4
∂u14
∂f4
∂u15
∂f4
∂u16
∂f5
∂u1
∂f5
∂u2
∂f5
∂u3
∂f5
∂u4
∂f5
∂u5
∂f5
∂u6
∂f5
∂u7
∂f5
∂u8
∂f5
∂u9
∂f5
∂u10
∂f5
∂u11
∂f5
∂u12
∂f5
∂u13
∂f5
∂u14
∂f5
∂u15
∂f5
∂u16
∂f6
∂u1
∂f6
∂u2
∂f6
∂u3
∂f6
∂u4
∂f6
∂u5
∂f6
∂u6
∂f6
∂u7
∂f6
∂u8
∂f6
∂u9
∂f6
∂u10
∂f6
∂u11
∂f6
∂u12
∂f6
∂u13
∂f6
∂u14
∂f6
∂u15
∂f6
∂u16
∂f7
∂u1
∂f7
∂u2
∂f7
∂u3
∂f7
∂u4
∂f7
∂u5
∂f7
∂u6
∂f7
∂u7
∂f7
∂u8
∂f7
∂u9
∂f7
∂u10
∂f7
∂u11
∂f7
∂u12
∂f7
∂u13
∂f7
∂u14
∂f7
∂u15
∂f7
∂u16
∂f8
∂u1
∂f8
∂u2
∂f8
∂u3
∂f8
∂u4
∂f8
∂u5
∂f8
∂u6
∂f8
∂u7
∂f8
∂u8
∂f8
∂u9
∂f8
∂u10
∂f8
∂u11
∂f8
∂u12
∂f8
∂u13
∂f8
∂u14
∂f8
∂u15
∂f8
∂u16
∂f9
∂u1
∂f9
∂u2
∂f9
∂u3
∂f9
∂u4
∂f9
∂u5
∂f9
∂u6
∂f9
∂u7
∂f9
∂u8
∂f9
∂u9
∂f9
∂u10
∂f9
∂u11
∂f9
∂u12
∂f9
∂u13
∂f9
∂u14
∂f9
∂u15
∂f9
∂u16
∂f10
∂u1
∂f10
∂u2
∂f10
∂u3
∂f10
∂u4
∂f10
∂u5
∂f10
∂u6
∂f10
∂u7
∂f10
∂u8
∂f10
∂u9
∂f10
∂u10
∂f10
∂u11
∂f10
∂u12
∂f10
∂u13
∂f10
∂u14
∂f10
∂u15
∂f10
∂u16
∂f11
∂u1
∂f11
∂u2
∂f11
∂u3
∂f11
∂u4
∂f11
∂u5
∂f11
∂u6
∂f11
∂u7
∂f11
∂u8
∂f11
∂u9
∂f11
∂u10
∂f11
∂u11
∂f11
∂u12
∂f11
∂u13
∂f11
∂u14
∂f11
∂u15
∂f11
∂u16
∂f12
∂u1
∂f12
∂u2
∂f12
∂u3
∂f12
∂u4
∂f12
∂u5
∂f12
∂u6
∂f12
∂u7
∂f12
∂u8
∂f12
∂u9
∂f12
∂u10
∂f12
∂u11
∂f12
∂u12
∂f12
∂u13
∂f12
∂u14
∂f12
∂u15
∂f12
∂u16
∂f13
∂u1
∂f13
∂u2
∂f13
∂u3
∂f13
∂u4
∂f13
∂u5
∂f13
∂u6
∂f13
∂u7
∂f13
∂u8
∂f13
∂u9
∂f13
∂u10
∂f13
∂u11
∂f13
∂u12
∂f13
∂u13
∂f13
∂u14
∂f13
∂u15
∂f13
∂u16
∂f14
∂u1
∂f14
∂u2
∂f14
∂u3
∂f14
∂u4
∂f14
∂u5
∂f14
∂u6
∂f14
∂u7
∂f14
∂u8
∂f14
∂u9
∂f14
∂u10
∂f14
∂u11
∂f14
∂u12
∂f14
∂u13
∂f14
∂u14
∂f14
∂u15
∂f14
∂u16
∂f15
∂u1
∂f15
∂u2
∂f15
∂u3
∂f15
∂u4
∂f15
∂u5
∂f15
∂u6
∂f15
∂u7
∂f15
∂u8
∂f15
∂u9
∂f15
∂u10
∂f15
∂u11
∂f15
∂u12
∂f15
∂u13
∂f15
∂u14
∂f15
∂u15
∂f15
∂u16
∂f16
∂u1
∂f16
∂u2
∂f16
∂u3
∂f16
∂u4
∂f16
∂u5
∂f16
∂u6
∂f16
∂u7
∂f16
∂u8
∂f16
∂u9
∂f16
∂u10
∂f16
∂u11
∂f16
∂u12
∂f16
∂u13
∂f16
∂u14
∂f16
∂u15
∂f16
∂u16


(27)
Determinando as derivadas especificadas em (27) para o sistema definido em
(26) obtemos
10
J =


1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0
[
− 1
Re∆y2
]
0 0
[
− u6
∆x
] [
1
∆x
(2u6 − u5) +
2
Re∆y2
]
0 0 0
[
− 1
Re∆y2
]
0 0 0 0 0 0
0 0
[
− 1
Re∆y2
]
0 0
[
− u7
∆x
] [
1
∆x
(2u7 − u6) +
2
Re∆y2
]
0 0 0
[
− 1
Re∆y2
]
0 0 0 0 0
0 0 0 0 0 0 −1 1 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1


(28)
1
1
Considere o caso em que ∆x = ∆y = 0.1 e Re = 100, essa matriz e´ dada por
J =


1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 −1 0 0 −10u6 10 (2u6 − u5) + 2 0 0 0 −1 0 0 0 0 0 0
0 0 −1 0 0 −10u7 10 (2u7 − u6) + 2 0 0 0 −1 0 0 0 0 0
0 0 0 0 0 0 −1 1 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
...
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1


(29)
O me´todo de soluc¸a˜o consiste em supor uma condic¸a˜o inicial u0 para o vetor
u, substituir essa condic¸a˜o na matriz J(u0) e resolver o seguinte sistema linear
J(u0) · u
∗
0
= −f(u0) (30)
em seguida no´s atualizamos o novo valor do vetor u1 fazendo
u1 = u0 + u
∗
0
, (31)
com esse novo valor de u1 no´s resolvemos o seguinte sistema
J(u1) · u
∗
1
= −f(u1), (32)
em seguida fazemos
u2 = u1 + u
∗
1
, (33)
e assim sucessivamente ate´ que o me´todo convirja para um valor dentro de
umadeterminada faixa de toleraˆncia. Podemos resumir o me´todo nos seguintes
passos
1. Supor um valor inicial para o vetor u0
2. Determinar a matriz Jacobiana J(u0)
3. Determinar o vetor f (u0)
4. Resolver o sistema linear J(u0) · u
∗
0
= −f(u0)
12
5. Corrigir o pro´ximo valor de u1 com u1 = u0 + u
∗
0
6. Fazer para j de 1 a` N:
• Resolver o sistema linear J(uj) · u
∗
j = −f(uj)
• Corrigir o pro´ximo valor de uj com uj = uj + u
∗
j
7. Repetir o procedimento ate´ convergir;
Para fins de exemplificac¸a˜o considere o seguinte valor inicial para o vetor u0
u0 =


0
0
0
0
1
1
1
1
1
1
1
1
1
1
1
1


(34)
Para esses valores iniciais, com ∆x = ∆y = 0.1 e Re = 100, temos que J(u0) e´
dada por
J =


1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 −1 0 0 −10 12 0 0 0 −1 0 0 0 0 0 0
0 0 −1 0 0 −10 12 0 0 0 −1 0 0 0 0 0
0 0 0 0 0 0 −1 1 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
0 0 0 0 0 −1 0 0 −10 12 0 0 0 −1 0 0
0 0 0 0 0 0 −1 0 0 −10 12 0 0 0 −1 0
0 0 0 0 0 0 0 0 0 0 −1 1 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1


(35)
13
e o vetor f (u0) e´ dado por
f (u0) =


0
0
0
0
0
1
1
0
0
0
0
0
0
0
0
0


(36)
Dessa forma precisamos resolver inicialmente o seguinte sistema


1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0
0 −1 0 0 −10 12 0 0 0 −1 0 0 0 0 0 0
0 0 −1 0 0 −10 12 0 0 0 −1 0 0 0 0 0
0 0 0 0 0 0 −1 1 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0
0 0 0 0 0 −1 0 0 −10 12 0 0 0 −1 0 0
0 0 0 0 0 0 −1 0 0 −10 12 0 0 0 −1 0
0 0 0 0 0 0 0 0 0 0 −1 1 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1


·


u∗0
1
u∗0
2
u∗0
3
u∗0
4
u∗0
5
u∗0
6
u∗0
7
u∗0
8
u∗0
9
u∗0
10
u∗0
11
u∗0
12
u∗0
13
u∗0
14
u∗0
15
u∗0
16


=


0
0
0
0
0
−1
−1
0
0
0
0
0
0
0
0
0


(37)
A soluc¸a˜o desse sistema linear pelo me´todo de Jacobi (me´todo iterativo) com
14
500 iterac¸o˜es fornece
u∗
0
=


0
0
0
0
0
−0.083916068
−0.1548242
−0.1548242
0
−0.006993004
−0.018729515
−0.018729515
0
0
0
0


(38)
Fac¸amos agora a seguinte correc¸a˜o
u1 = u0 + u
∗
0
, (39)
de modo que
u1 =


0
0
0
0
1
0.9160839
0.8451759
0.8451759
1
0.9930070
0.9812705
0.9812705
1
1
1
1


(40)
Agora precisamos resolver o seguinte sistema
J(u1) · u
∗
1
= −f(u1), (41)
note que a matriz J muda a cada nova iterac¸a˜o. Montando o sistema e resol-
vendo para u∗
1
pelo mesmo me´todo (Jacobi com 500 iterac¸o˜es) obtemos o novo
15
valor de u2 dado por
u2 =


0
0
0
0
1
0.9092013
0.8277108
0.8277108
1
0.9923854
0.9790677
0.9790677
1
1
1
1


(42)
Repetindo esse procedimento 10 vezes obtemos u10, comparando o valor de u10
com o valor de u9 temos
u9 =


0
0
0
0
1
0.9091543
0.8274714
0.8274714
1
0.9923811
0.9790405
0.9790405
1
1
1
1


u10 =


0
0
0
0
1
0.9091544
0.8274716
0.8274716
1
0.9923812
0.9790406
0.9790406
1
1
1
1


(43)
16
Para ana´lise da convergeˆncia do processo repare no vetor dado por |u10 − u9|
|u10 − u9| =


0
0
0
0
0
1.1920929E− 07
1.1920929E− 07
1.1920929E− 07
0
5.9604645E− 08
5.9604645E− 08
5.9604645E− 08
0
0
0
0


(44)
O valor ma´ximo dessa diferenc¸a, nesse caso 1.1920929E− 07 pode ser utilizado
como o paraˆmetro que devera´ ser menor que determinada toleraˆncia pre´-definida
pelo usua´rio para garantir a convergeˆncia do me´todo.
4 Bibliografia
• Notas de aula do curso de fenoˆmenos de transporte - Prof. Rafael Gabler
Gontijo
• Micromecaˆnica e Microhidrodinaˆmica de Suspenso˜es Magne´ticas - Tese de
Doutorado - 2013 - Prof. Rafael Gabler Gontijo
17
	Introdução e motivação
	O método de Newton - raízes de funções não lineares
	Exemplo - Solução de uma EDP não linear
	Bibliografia

Outros materiais