Trabalho teorico e computacional

Prévia do material em texto

UNIVERSIDADE FEDERAL DA BAHIA 
INSTITUTO DE CIÊNCIAS AMBIENTAIS E 
DESENVOLVIMENTO SUSTENTÁVEL- (ICADS) 
COLEGIADO DE FÍSICA 
DISCIPLINA: CALCULO NUMÉRICO 
 
 
 
 
 
 
 
 
 
 
 Dérick	
  Gabriel	
  F.	
  Borges	
  	
  
 
 
 Calculo	
  numérico	
  abordagem	
  teórica	
  e	
  computacional	
  
Professor: Kennedy Morais Fernandes 
 
 
 
 
 
 
 
 
 
 
 
Barreiras-ba, 8 de outobro de 2012 
	
  
	
  
	
  
	
  
 
	
  
	
  
	
  
Sumário	
  	
  	
  
1- Introdução 
-Conceito de calculo numérico 
-Aplicações 
-erros 
 
2- Raizes 
-Método da Bisseção 
-Método de Newton 
-Método da Secante 
 
3- Sistemas Lineares 
-Método da Eliminação Gaussiana 
-Método da Gauss Seidel 
 
4- Ajuste de curvas 
-Métodos do Mínimo Quadrado 
 
5- Interpolação 
-Polinômio de Lagrange 
-Spline Cubico 
 
6- Integração Numérica 
-Regra do Trapézio 
-Método de Simpson 
 
7- Resolução de EDO 
-Método Euler 
-Métodos de Runge Kutta 
 
8- Conclusão 
 
9- Referências 	
  	
  	
  	
  	
  	
  	
  	
  	
  
1. Introdução	
  	
   Cálculo	
   Numérico	
   e	
   o	
   ramo	
   da	
   matemática	
   que	
   desenvolve	
   algoritmos	
  para	
   simular	
   processos	
   matemáticos	
   complexos	
   provenientes	
   do	
   mundo	
   real.	
  Esses	
  algoritmos	
  são	
  seqüência	
  de	
  instruções	
  ordenadas,	
  de	
  modo	
  a	
  fornecer	
  em	
  seu	
   decurso,	
   a	
   solução	
   para	
   um	
   problema	
   específico.	
   No	
   contexto	
   de	
   Cálculo	
  Numérico,	
   os	
   algoritmos	
   são	
   denominados	
   de	
   Métodos	
   Numéricos.[1]	
   Esses	
  métodos	
   se	
   aplicam	
   principalmente	
   a	
   problemas	
   que	
   não	
   apresentam	
   uma	
  solução	
  exata(situação	
  real),	
  portanto	
  precisam	
  ser	
  resolvidos	
  numericamente.	
  	
   Calculo	
   numerico	
   aplica-­‐se	
   as	
   ciências	
   exatas	
   e	
   as	
   engenharias,	
  geralmente	
  em	
  problemas	
  fisicos	
  do	
  contexto	
  real,	
  como:	
  -­‐ elétricidade	
  -­‐ petróleo	
  -­‐ cores	
  	
  -­‐ pressão	
  	
  -­‐ resistência	
  	
  A	
  utilização	
  de	
  simuladores	
  numéricos	
  para	
  determinação	
  da	
  solução	
  de	
  um	
  problema	
  requer	
  a	
  execução	
  da	
  seguinte	
  seqüência	
  de	
  etapas:	
  Etapa	
  1:	
  Definir	
  o	
  problema	
  real	
  a	
  ser	
  resolvido	
  	
  Etapa	
  2:	
  Observar	
  fenômenos,	
   levantar	
  efeitos	
  dominantes	
  e	
  fazer	
  referências	
  a	
  conhecimentos	
  prévios	
  físicos	
  e	
  matemáticos	
  Etapa	
  3:	
  Modelagem,	
  fase	
  de	
  obtenção	
  de	
  um	
  modelo	
  matemático	
  que	
  descreve	
  um	
  problema	
  físico	
  em	
  questão.	
  Etapa	
  4:	
  Resolução,	
  fase	
  de	
  obtenção	
  da	
  solução	
  através	
  da	
  aplicação	
  de	
  métodos	
  numéricos.	
  	
   𝑝𝑟𝑜𝑏𝑙𝑒𝑚𝑎   ⇒𝒎𝒐𝒅𝒆𝒍𝒂𝒈𝒆𝒎   ⇒ 𝑚𝑜𝑑𝑒𝑙𝑜  𝑚𝑎𝑡𝑒𝑚á𝑡𝑖𝑐𝑜   ⇒ 𝒓𝒆𝒔𝒐𝒍𝒖çã𝒐   ⇒ 𝑠𝑜𝑙𝑢çã𝑜	
  	
   Ao	
   se	
   tentar	
   representar	
   um	
   fenômeno	
   físico	
   por	
   meio	
   de	
   um	
   método	
  matemático,	
   raramente	
   se	
   tem	
   uma	
   descrição	
   correta	
   deste	
   fenômeno,	
  normalmente,	
  são	
  necessárias	
  várias	
  simplificações	
  do	
  mundo	
  físico	
  para	
  que	
  se	
  tenha	
  um	
  modelo.	
  Tais	
  simplificações	
  resultam	
  em	
  erros.	
  A	
  resolução	
  de	
  modelos	
  matemáticos	
   exige,	
   muitas	
   vezes,	
   o	
   uso	
   de	
   instrumentos	
   de	
   cálculo	
   que	
  necessitam	
  de	
  aproximações,	
  tais	
  aproximações	
  podem	
  gerar	
  erros,	
  tais	
  como:	
  [1]	
  -­‐ Erros	
  inerentes.	
  -­‐ Erros	
  de	
  arredondamento.	
  -­‐ Erros	
  de	
  truncamento.	
  	
   Erro	
   inerentes,	
   são	
   erros	
   que	
   existem	
   nos	
   dados	
   e	
   são	
   causados	
   pelos	
  erros	
  dos	
  próprios	
  equipamentos	
  utilizados	
  na	
  captaçao	
  dos	
  dados.	
  Erros	
  de	
  arredondamento,	
  são	
  os	
  erros	
  originados	
  pela	
  representação	
  do	
  números	
  reais	
  utilizando-­‐se	
  apenas	
  um	
  numero	
  finito	
  de	
  casas	
  decimais.	
  Erros	
   de	
   truncamento,	
   são	
   os	
   erros	
   causados	
   quando	
   utilizamos	
   num	
  processo	
  algoritmo	
  infinito	
  apenas	
  uma	
  parte	
  finita	
  do	
  processo.	
  A	
  noção	
  de	
  erro	
  sempre	
  estará	
  presente	
  em	
  todos	
  os	
  campos	
  do	
  Cálculo	
  Numérico,	
   pois	
   sempre	
   temos	
   dados	
   nem	
   sempre	
   exatos,	
   as	
   operações	
   sobre	
  valores	
   não	
   exatos	
   propagam	
  os	
   erros	
   e	
   os	
   	
  Métodos	
  Numéricos	
   são,	
   em	
  geral	
  métodos	
  de	
  aproximação.	
  [1]	
  
2. Raizes	
  	
  	
   Nas	
  mais	
  diversas	
  áreas	
  das	
  ciências	
  exatas	
  ocorrem,	
  frequentemente,	
  situações	
  que	
  envolvem	
  a	
  resolução	
  de	
  uma	
  equação	
  do	
  tipo	
  𝑓(𝑥) = 0	
  	
  A	
  idéia	
  central	
  dos	
  métodos	
  numéricos	
  para	
  encontrar	
  as	
  raizes	
  de	
  uma	
  equação,	
   parti	
   de	
   uma	
   aproximação	
   inicial	
   para	
   a	
   raiz	
   (um	
   intervalo	
   onde	
  imaginamos	
  a	
  raiz	
  estar	
  contida)	
  e	
  em	
  seguida	
  refinar	
  essa	
  aproximação	
  através	
  de	
  um	
  processo	
  iterativo.	
  [2]	
  A	
   forma	
   como	
   se	
   efetua	
   o	
   refinamento	
   é	
   que	
   diferencia	
   os	
   métodos.	
  Temos	
  varios	
  método,	
  serão	
  apontados	
  nesse	
  trabalho	
  o	
  	
  da	
  Bisseçao,	
  de	
  Newton,	
  da	
  Secante.	
  	
  
Método	
   da	
   Bisseção,	
   parte	
   de	
   um	
   	
   intervalo	
   inicial	
   [a,b]	
   que	
   se	
   sabe	
  conter	
   o	
   zero	
   de  𝑓,	
   suposto	
   único.	
   Em	
   cada	
   iteração	
   é	
   produzido	
   um	
   intervalo	
  com	
  metade	
  do	
  comprimento	
  do	
  intervalo	
  actual.	
  Para	
  tal,	
  divide-­‐se	
  o	
  intervalo	
  actual	
  a	
  meio	
  e	
  escolhe-­‐se	
  o	
  subintervalo	
  esquerdo	
  ou	
  direito	
  de	
  forma	
  a	
  que	
  a	
  função	
  tenha	
  sinais	
  diferentes	
  nos	
  extremos	
  do	
  subintervalo	
  escolhido.	
  Ou	
  seja,	
  sendo	
   [𝑎!,  𝑏!]	
  o	
   intervalo	
  na	
   iteração	
  n,	
   calcula-­‐se	
  𝑥!!! = !!!!!!  .	
  O	
  valor	
  𝑥𝑛 + 1	
  substitui	
  𝑎! 	
  ou	
  𝑏! 	
  consoante	
  𝑓 𝑥!!!    𝑓(𝑏!  ) < 0 	
  ou	
  𝑓 𝑥!!!    𝑓(𝑎!  ) < 0 .	
   Desta	
  forma,	
   assegura-­‐se	
   que	
   s	
  ∈	
   [𝑎!,	
  𝑏!]	
   em	
   qualquer	
   iteração.	
   O	
   erro	
   absoluto	
   da	
  estimativa	
  de	
  𝑥!	
  e	
  dado	
  por	
  	
  !!!!! 	
  .[3]	
  	
  Calculando	
  a	
  raiz	
  da	
  função	
  𝑓 𝑥 = 1+ 𝑥 + 𝑒!	
  pelo	
  método	
  da	
  bisseção,	
  com	
  erro	
  de	
  10!!.	
  	
  
Fortran	
  1	
  (implementação)	
  	
  
	
  
program bissecao 
implicit none 
 
real xm,xe,xd,f,e 
integer result,i 
result=21 
open(unit=result,file='result.txt') 
xe=-2. 
xd=-1. 
e=10**(-5.) 
xm=(xe+xd)/2. 
 
do while (abs(xe-xd)>e) 
if (f(xe)*f(xm)>0.) then 
xe=xm 
else 
xd=xm 
end if 
xm=(xe+xd)/2. 
write(*,*)'raiz',xm 
end do 
 
end 
 
function f(x) 
real f,x 
f=1+x+exp(x) 
return 
end 
	
  	
  
Resultados	
  obtidos	
  	
  	
  	
   n   𝑎!   f(𝑎!)   𝑏!   f(𝑏!)   𝑥!!!	
   f(𝑥!!!)  0-­‐2.000	
   -­‐0.865	
   -­‐1.000	
   +0.368	
   -­‐1.500	
   -­‐0.277	
  1	
   -­‐1.500	
   -­‐0.277	
   -­‐1.000	
   +0.368	
   -­‐1.250	
   +0.037	
  2	
   -­‐1.500	
   -­‐0.277	
   -­‐1.250	
   +0.037	
   -­‐1.375	
   -­‐0.122	
  3	
   -­‐1.375	
   -­‐0.122	
   -­‐1.250	
   +0.037	
   -­‐1.313	
   -­‐0.043	
  4	
   -­‐1.313	
   -­‐0.043	
   -­‐1.250	
   +0.037	
   -­‐1.281	
   -­‐0.004	
  5	
   -­‐1.281	
   -­‐0.004	
   -­‐1.250	
   +0.037	
   -­‐1.266	
   +0.016	
  6	
   -­‐1.281	
   -­‐0.004	
   -­‐1.266	
   +0.016	
   -­‐1.276	
   +0.006	
  7	
   -­‐1.281	
   -­‐0.004	
   -­‐1.273	
   +0.006	
   -­‐1.277	
   +0.001	
  Tabela	
  1	
  	
  A	
  solução	
  da	
  equação	
  será	
  s	
  =	
  −1.277  ±  4  ×  10!!,	
  ou	
  seja,	
  s	
  ∈	
  [−1.281,	
  −1.273].	
  	
  
	
  
Método	
   da	
   Newton,	
   e	
   um	
   dos	
   métodos	
   mais	
   poderosos	
   para	
   resolver	
  equações	
   do	
   tipo	
  𝑓(𝑥)  =  0.	
   Este	
   método	
   parte	
   de	
   uma	
   estimativa	
   inicial	
   x0	
   e	
  gera	
  uma	
  sucessão	
  {𝑥!}	
  de	
  uma	
  forma	
  recorrente.	
  Cada	
  novo	
  valor	
  da	
  sucessão,	
  𝑥!!!,	
  e	
  determinado	
  como	
  sendo	
  a	
  abcissa	
  do	
   ponto	
   de	
   interseção	
   com	
  o	
   eixo	
   dos	
   xx	
   da	
   reta	
   tangente	
   a	
   função	
   no	
   ponto	
  (𝑥!, (𝑓(𝑥!)),	
  ou	
  seja,	
  no	
  ponto	
  correspondente	
  ao	
  valor	
  anterior	
  da	
  sucessão.	
  A	
  expressão	
  de	
   interação	
  que	
  permite	
  determinar	
  𝑥!!!	
  em	
   função	
  de	
  𝑥!	
  e	
  dado	
  por	
  :	
  [3]	
  	
   y   =  f(𝑥!)  +  f′(𝑥!)   ·  (𝑥  −  𝑥!)	
  	
   𝑥! + 1 = 𝑥! −   𝑓(𝑥!)𝑓′(𝑥!)	
  	
  	
  Calculando	
  a	
  raiz	
  da	
  função	
  𝑓 𝑥 = 1+ 𝑥 + 𝑒!	
  pelo	
  método	
  de	
  Newton,	
  erro	
  de	
  10!!.	
  	
  
Fortran	
  2	
  (implementação)	
  	
  
program newton 
implicit none 
real xo,xn,f,df,e,h 
integer i 
xo=0. 
 
e=10**(-6.) 
xn=xo-(f(xo)/df(xo)) 
do while (abs(f(xn))>e) 
xo=xn 
xn=xo-(f(xo)/df(xo)) 
print*,xn 
end do 
 
end 
 
function f(x) 
real f,x 
f=1+x+exp(x) 
return 
end 
 
function df(x) 
real df,x,f,h 
h=0.01 
df=(f(x+h)-f(x))/h 
return 
end 	
  
Resultados	
  obtidos	
  	
  
	
   n   𝑥!   f(𝑥!)   f’(𝑥!)   𝑥!!!   erro  0   -­‐1.0000   +3.68  10!!   +1.368   -­‐1.26894   1.2  10!!  1   -­‐1.26894   +1.22  10!!   +1.281   -­‐1.27845   1.5  10!!  2   -­‐1.27845   +1.27  10!!   +1.278   -­‐1.27846   1.6  10!!!  Tabela	
  2	
  	
  A	
  solução	
  aproximada	
  será	
  s	
  ≃	
  −1.27846	
  	
  	
  
Método	
  da	
  Secante,	
  e	
  semelhante	
  ao	
  método	
  de	
  Newton,	
  com	
  a	
  diferença	
  de	
  que	
  a	
  reta	
  tangente	
  a	
  função	
  é	
  substituída	
  pela	
  reta	
  secante	
  nos	
  dois	
  últimos	
  pontos.	
   Este	
   método	
   obriga	
   a	
   que	
   em	
   cada	
   iteração	
   sejam	
   guardadas	
   as	
   duas	
  últimas	
  estimativas	
  da	
  solução	
  a	
  determinar.	
  Essa	
  reta	
  e	
  descrita	
  pela	
  equação:	
  [3]	
  
𝑦 = 𝑓 𝑥!!! + 𝑓 𝑥𝑛 − 𝑓 𝑥!!!𝑥𝑛 − 𝑥!!! (𝑥 − 𝑥!!!)	
  	
  
𝑥!!! = 𝑥𝑛 + 1  𝑓 𝑥𝑛 − 𝑥𝑛  𝑓(𝑥!!!)𝑓 𝑥𝑛 − 𝑓(𝑥!!!) 	
  	
  Calculando	
  a	
  raiz	
  da	
  função	
  𝑓 𝑥 = 1+ 𝑥 + 𝑒!	
  pelo	
  método	
  da	
  Secante,	
  com	
  erro	
  de	
  10!!.	
  	
  
	
  
Fortran	
  3	
  (implementação)	
  	
  
program newton 
implicit none 
real a,b,xn,f,e,h 
integer i 
a=-1. 
b=0. 
e=10**(-5.) 
xn=a-(f(a)/(f(a)-f(b)))*(a-b) 
do while (abs(f(xn))>e) 
a=xn 
xn=a-(f(a)/(f(a)-f(b)))*(a-b) 
print*,xn 
end do 
 
end 
 
function f(x) 
real f,x 
f=1+x+exp(x) 
return 
end 
	
  	
  
Resultados	
  obtidos	
  	
  	
  	
   n   𝑥!!!   𝑥!   𝑥!!!   f(𝑥!!!)   erro  0   -­‐1.00000   -­‐1.10000   -­‐1.27249   7.65  10!!   7.6  10!!  1   -­‐1.10000   -­‐1.27249   -­‐1.27834   1.55  10!!   1.7  10!!  2   -­‐1.27249   -­‐1.27834   -­‐1.27846   1.01  10!!   1.2  10!!  Tabela	
  3	
  
 A	
  estimativa	
  obtida	
  	
  é	
  s	
  ≃	
  −1.27846.	
  	
  	
  	
  
3. Sistemas	
  lineares	
  	
  
	
   A	
   resolução	
   de	
   sistemas	
   de	
   equações	
   lineares	
   	
   é	
   um	
   dos	
   problemas	
  numéricos	
   mais	
   comuns	
   em	
   aplicações	
   científicas.	
   Tais	
   sistemas	
   surgem,	
   por	
  exemplo,	
   em	
   conexão	
   com	
   a	
   solução	
   de	
   equações	
   diferenciais	
   parciais,	
  determinação	
  de	
  caminhos	
  ótimos	
  em	
  redes	
  e	
   interpolação	
  o	
  de	
  pontos,	
  dentre	
  outros.	
   Vamos	
   abordar	
   nesse	
   capitulo	
   dois	
   métodos	
   para	
   solução	
   de	
   sistemas	
  lineares	
  o	
  método	
  de	
  Eliminação	
  Gaussiana	
  e	
  Gauss	
  Seidel.[4]	
  
	
  	
  
Eliminação	
  Gaussiana,	
  e	
  um	
  método	
  exato	
  também	
  chamado	
  de	
  método	
  de	
   Gauss	
   Simples,	
   consiste	
   em	
   transformar	
   o	
   sistema	
   dado	
   num	
   sistema	
  triangular	
   equivalente	
   através	
   de	
   uma	
   sequência	
   de	
   operações	
   elementares	
  sobre	
  as	
  linhas	
  do	
  sistema	
  original,	
  isto	
  é,	
  o	
  sistema	
  equivalente	
  é	
  obtido	
  através	
  
da	
  aplicação	
  repetida	
  da	
  operação,	
  substituir	
  uma	
  equação	
  pela	
  diferença	
  entre	
  essa	
   mesma	
   equação	
   e	
   uma	
   outra	
   equação	
   multiplicada	
   por	
   uma	
   constante	
  diferente	
  de	
  zero.	
  É	
  claro	
  que	
  tal	
  operação	
  não	
  altera	
  a	
  solução	
  do	
  sistema,	
  isto	
  	
  é,	
  	
  obtem-­‐se	
  com	
  ela	
  outro	
  sistema	
  equivalente	
  ao	
  original.	
  O	
  objetivo	
  é	
  organizar	
  essa	
   sequência	
   de	
   operações	
   de	
   tal	
   forma	
   que	
   o	
   sistema	
   linear	
   resultante	
   seja	
  triangular	
  superior.	
  Utilizando	
   o	
  método	
   de	
   Eliminação	
   Gaussiana	
   para	
   encontrar	
   a	
   solução	
  para	
  	
  𝑥!,  𝑥!  𝑒  𝑥!.	
   3 1 11 4 20 2 5   𝑥!𝑥!𝑥! = 745    ⇒  𝑝𝑎𝑟𝑎  𝑥! = 000 	
  
resultados	
  obtidos:	
  Solução	
  exata	
   2.0000	
   0.0000	
   1.0000	
  Tabela	
  4	
  	
  	
   Os	
  resultados	
  mostrados	
  na	
  tabela	
  4	
  corresponde	
  a	
  solução	
  exata	
  para	
  𝑥!,  𝑥!  𝑒  𝑥!	
  nessa	
  ordem.	
  
	
  
Fortran	
  4	
  (implementação)	
  	
  
program elininacaogauss 
implicit none 
 
real A(100,100), B(100), X(100), max, aux 
integer i, j, k, n, entrada,saida 
 
entrada = 12 
saida=11 
n = 3 
 
open (unit = entrada, file = 'entrada.txt') 
open (unit = saida, file = 'saida.txt') 
do i=1,n,1 
read (entrada,*) (A(i,j), j=1,n,1) 
end do 
 
do i=1,n,1 
read (entrada,*) B(i) 
end do 
 
do i=1,n-1,1 
 
k = i 
max = abs(A(i,i)) 
 
do j=i+1,n,1 
if(abs(A(j,i)) > max) then 
max = abs(A(j,i)) 
k = j 
end if 
end do 
 
if (i .ne. k) then 
 
do j=1,n,1 
 
aux = A(i,j) 
A(i,j) = A(k,j) 
A(k,j) = aux 
 
end do 
 
aux = B(i) 
B(i) = B(k) 
B(k) = aux 
 
end if 
 
do j=i+1,n,1 
 
aux = A(j,i)/A(i,i) 
A(j,i) = 0 
 
do k=i+1,n,1 
A(j,k) = A(j,k) - aux*A(i,k) 
end do 
 
B(j) = B(j) - aux*B(i) 
 
end do 
end do 
 
do i=n, 1, -1 
 
X(i) = B(i) 
 
do j=i+1,n,1 
X(i) = X(i) - A(i,j)*X(j) 
end do 
 
X(i) = X(i)/A(i,i) 
 
end do 
 
do i=1,n,1 
write(saida,*) X(i) 
end do 
 
close (entrada) 
close (saida) 
end 
 
 Método	
  de	
  Gauss	
  Seidel,	
  sendo	
  um	
  métodointerativo,	
  tem	
  melhor	
  comportamento	
  do	
  que	
  os	
  métodos	
  exatos	
  quando	
  se	
  trabalha	
  com	
  matrizes	
  esparsas.	
  É	
  um	
  dos	
  métodos	
  iterativos	
  mais	
  comuns	
  e	
  simples	
  para	
  ser	
  programado	
  em	
  computador.	
  A	
  idéia	
  principal	
  é,	
  cada	
  coordenada	
  do	
  vetor	
  correspondente	
  à	
  nova	
  aproximação	
  é	
  calculada	
  a	
  partir	
  da	
  respectiva	
  equação	
  do	
  sistema,	
  utilizando-­‐se	
  as	
  coordenadas	
  do	
  vetor	
  aproximação	
  da	
  iteração	
  anterior,	
  quando	
  essas	
  ainda	
  não	
  foram	
  calculadas	
  na	
  iteração	
  corrente,	
  e	
  as	
  coordenadas	
  do	
  vetor	
  aproximação	
  da	
  iteração	
  corrente,	
  no	
  caso	
  contrário.	
  [5] A	
  função	
  interação	
  do	
  método	
  de	
  Gauss	
  Seidel	
  pode	
  ser	
  expressada	
  por:	
  
𝑥!! = 𝑏! − 𝑎!"𝑥!!!!!!!! − 𝑎!"𝑥!!!!!!!!!!𝑎!!  𝑝𝑎𝑟𝑎  𝑖 = 1…𝑛	
  
	
  
	
  
Fortran	
  5	
  (implementação)	
  	
  
program gausseidel 
implicit none 
 
real A(100,100), B(100), X(100) 
integer n, entrada, i, j, k,saida 
character * 15 texto 
 
entrada = 12 
 
open (unit = entrada, file = 'dados.txt') 
open (unit = saida, file = 'saida.txt') 
read (entrada,*) texto 
 
read (entrada,*) n 
 
read (entrada,*) texto 
 
do i=1,n,1 
read (entrada,*) (A(i,j), j=1,n,1) 
end do 
 
read (entrada,*) texto 
 
do i=1,n,1 
read (entrada,*) B(i) 
end do 
 
do i=1,n,1 
X(i) = 0 
end do 
 
do k=1,100,1 
 
do i=1,n,1 
 
X(i) = B(i) 
 
do j=1,i-1,1 
X(i) = X(i) - A(i,j)*X(j) 
end do 
 
do j=i+1,n,1 
X(i) = X(i) - A(i,j)*X(j) 
end do 
 
X(i) = X(i)/A(i,i) 
 
end do 
end do 
 
do i=1,n,1 
write(saida,*) X(i) 
end do 
end 	
  Aplicando-­‐se	
  Gaus	
  Seidel,	
  para	
  calcular	
  𝑥!,  𝑥!  𝑒  𝑥!	
  .	
  	
   3 1 11 4 20 2 5   𝑥!𝑥!𝑥! = 745    ⇒  𝑝𝑎𝑟𝑎  𝑥! = 000 	
  	
  
Resultados	
  obtidos	
  	
  
	
   Interações	
  	
   𝑥!	
   𝑥!	
   𝑥!	
  1	
   2.3333	
   0.4167	
   0.8333	
  2	
   1.9167	
   0.1042	
   0.9583	
  3	
   1.9792	
   0.0260	
   0.9896	
  4	
   1.9942	
   0.0065	
   0.9974	
  5	
   1.9987	
   0.0016	
   0.9993	
  6	
   1.9997	
   0.0004	
   0.9998	
  7	
   1.9999	
   0.0001	
   1.0000	
  …	
   …	
   …	
   …	
  Solução	
  exata	
   2.0000	
   0.0000	
   1.0000	
  Tabela	
  5	
  	
   Observa-­‐se,	
   que	
   com	
   poucas	
   interações	
   se	
   obteve	
   o	
   solução	
   exata	
  𝑥! = 2.000,  𝑥! = 0.000, 𝑥! = 1.000	
  	
  
4. Ajuste	
  de	
  curvas	
  	
  	
  	
  	
  	
   Seja	
   um	
   conjunto	
   de	
   dados	
   contendo	
   n	
  pares	
   de	
   valores	
   (x,	
   y)	
   obtidos	
  numérica	
  ou	
  experimentalmente.	
  De	
  modo	
  a	
  calcular	
  qualquer	
  valor	
  de	
  y	
  distinto	
  dos	
  valores	
  tabelados,	
  ajustamos	
  uma	
  função	
  y	
  =	
  f(x)	
  através	
  do	
  chamado	
  Método	
  dos	
  Mínimos	
  Quadrados.	
  
Considere	
   uma	
   equação	
   relacionando	
   a	
   variável	
   y	
   com	
   a	
   variável	
  independente	
  x,	
  como	
  y	
  = f	
  (x)	
  ,	
  onde	
  y	
  indica	
  que	
  este	
  é	
  o	
  valor	
  aproximado	
  de	
  y.	
  Queremos	
  encontrar	
  a	
  função	
  y	
  = f	
  (x)	
  ,	
  cujo	
  desvio	
  em	
  relação	
  aos	
  valores	
  y	
  seja	
  expresso	
  como	
  𝛿! = 𝑦! − 𝑦! .	
  Por	
  uma	
  questão	
  de	
  conveniência	
  trabalharemos	
  com	
  o	
  desvio	
  quadrático	
  𝛿!! = (𝑦! − 𝑦!)!	
  A	
   função	
   y	
   = f	
   (	
   x	
   )	
   que	
   melhor	
   ajusta	
   os	
   pontos	
   (x,	
   y)	
   dados	
   é	
   aquela	
   que	
  minimiza	
  o	
  somatório	
  dos	
  desvios	
  quadráticos	
  S:	
  
	
   𝑆 = 𝑦!!!!!! = (𝑦! − 𝑦!)!!!!! 	
  	
   A	
  condição	
  de	
  minimização	
  da	
  função	
  S	
  é	
  satisfeita	
  fazendo-­‐se	
  dS	
  =	
  0,	
  ou	
  seja,	
  necessitamos	
  calcular	
  a	
  derivada	
  da	
  função	
  S	
  em	
  relação	
  aos	
  parâmetros	
  de	
  ajuste	
   da	
   função	
   y	
  =	
   f(x)	
   para	
   que	
   possamos	
   encontrar	
   o	
   sistema	
   de	
   equações	
  denominado	
   equações	
   normais	
   que	
   conduz	
   ao	
  melhor	
   ajuste	
   dos	
   pontos	
   (x,y)	
  pela	
   função	
   y	
   =	
   f(x)	
   escolhida.	
   Para	
   cada	
   tipo	
   de	
   função	
   de	
   ajuste	
   existe	
   um	
  sistema	
  de	
  equações	
  normais	
  que	
  minimiza	
  a	
  soma	
  dos	
  desvios	
  quadráticos	
  S.	
  	
   Temos	
   varios	
   ajustes	
   que	
   se	
   fundamentam	
   nos	
   métodos	
   dos	
   Mínimos	
  Quadrados,	
  nesse	
   trabalho	
   será	
  apresentado	
  ajuste	
   linear,	
   ajuste	
  exponencial	
   e	
  ajuste	
  potencial.	
  	
  	
  Dados	
  os	
  pontos	
  abaixo,	
  faremos	
  3	
  ajustes:	
  linear,	
  exponencial	
  e	
  potencial.	
  	
  x	
   1.0	
   1.2	
   1.4	
   1.6	
   1.8	
   2.0	
   2.2	
   2.4	
   2.6	
   2.8	
   3.0	
  y	
   0.525	
   0.8448	
   1.2807	
   1.8634	
   2.6326	
   3.6386	
   4.944	
   6.6258	
   8.7768	
   11.5076	
   14.9484	
  Tabela	
  6	
  
	
  
Ajuste	
   linear,	
   função	
  na	
   forma	
  𝑦 = 𝑎! + 𝑎!𝑥	
  .	
   A	
   determinação	
  dos	
   coeficientes	
  𝑎!  𝑒  𝑎!	
  pelo	
  método	
  dos	
  mínimos	
  quadrados	
  e	
  dado	
  pelas	
  expressões:	
  	
   𝑎! = 𝑥! 𝑦!!!! − 𝑥!!!! 𝑥𝑦!!!!!!!!𝑛 𝑥! − ( 𝑥!!!! )!!!!! 	
  	
   𝑎! = 𝑛 𝑥𝑦!!!! − 𝑥!!!! 𝑦!!!!𝑛 𝑥! − ( 𝑥!!!!!!!! )! 	
  	
  
Fortran	
  6	
  (implementação)	
  
	
  program	
  ajuste	
  implicit	
  none	
  	
  real	
  X(100),	
  Y(100),	
  A(100,100),	
  B(100),	
  max,	
  aux	
  integer	
  i,	
  j,	
  k,	
  entrada	
  ,	
  saida,	
  n,	
  pontos	
  	
  entrada	
  =	
  12	
  
saida	
  =	
  14	
  	
  open	
  (unit	
  =	
  entrada,	
  file	
  =	
  'entrada.txt',status='old')	
  open	
  (unit	
  =	
  saida,	
  file	
  =	
  'saida.txt')	
  	
  read(entrada,*)	
  n	
  	
  	
  !	
  n	
  determina	
  o	
  tipo	
  de	
  ajuste	
  polinomial	
  	
  read(entrada,*)	
  pontos	
  	
  	
  	
  !quantidade	
  de	
  pontos	
  	
  do	
  i	
  =	
  1,	
  pontos,	
  1	
  read(entrada,*)	
  X(i),	
  Y(i)	
  end	
  do	
  	
  close(entrada)	
  	
  do	
  i=1,n,1	
  do	
  j=1,n,1	
  A(i,j)	
  =	
  0	
  end	
  do	
  end	
  do	
  	
  A(1,1)	
  =	
  pontos	
  	
  do	
  i=1,n,1	
  do	
  j=	
  1,pontos,1	
  A(i,n)	
  =	
  A(i,n)	
  +	
  X(j)**(n	
  +	
  i	
  -­‐	
  2)	
  end	
  do	
  end	
  do	
  	
  do	
  j=2,n-­‐1,1	
  do	
  i=1,pontos,1	
  A(1,j)	
  =	
  A(1,j)	
  +	
  X(i)**(j	
  -­‐	
  1)	
  end	
  do	
  end	
  do	
  	
  do	
  i=2,n,1	
  do	
  j=1,n-­‐1,1	
  A(i,j)	
  =	
  A(i	
  -­‐	
  1,	
  j	
  +	
  1)	
  end	
  do	
  end	
  do	
  	
  do	
  i=1,n,1	
  B(i)	
  =	
  0	
  end	
  do	
  	
  do	
  i=1,n,1	
  do	
  j=1,pontos,1	
  B(i)	
  =	
  B(i)	
  +	
  (X(j)**(i	
  -­‐	
  1))*Y(j)	
  end	
  do	
  
end	
  do	
  	
  !	
  inicio	
  da	
  resolução	
  do	
  sistema	
  	
  do	
  i=1,	
  n-­‐1,	
  1	
  	
  k=i	
  	
  max=abs(A(i,i))	
  	
  do	
  j=i+1,n,1	
  if	
  (abs(A(j,i))	
  >	
  max)	
  then	
  max	
  =	
  abs(A(j,i))	
  k	
  =	
  j	
  end	
  if	
  end	
  do	
  	
  if	
  (i	
  .ne.	
  k)	
  then	
  	
  do	
  j=1,	
  n,	
  1	
  aux	
  =	
  A(i,j)	
  A(i,j)	
  =	
  A(k,j)	
  A(k,j)	
  =	
  aux	
  end	
  do	
  	
  aux	
  =	
  B(i)	
  B(i)	
  =	
  B(k)B(k)	
  =	
  aux	
  end	
  if	
  	
  do	
  j=i+1,n,1	
  aux=A(j,i)/A(i,i)	
  A(j,i)=0	
  	
  do	
  k=i+1,n,1	
  A(j,k)=A(j,k)	
  -­‐	
  A(i,k)*aux	
  end	
  do	
  	
  B(j)=B(j)	
  -­‐	
  B(i)*aux	
  end	
  do	
  end	
  do	
  	
  do	
  i=n,1,-­‐1	
  	
  X(i)	
  =	
  B(i)	
  	
  do	
  j=i+1,n,1	
  X(i)=X(i)	
  -­‐	
  A(i,j)*X(j)	
  end	
  do	
  
	
  X(i)=X(i)/A(i,i)	
  	
  end	
  do	
  	
  write	
  (saida,*)	
  'Coeficientes	
  do	
  polinomio'	
  	
  do	
  i=1,n,1	
  write(saida,*)	
  X(i)	
  end	
  do	
  	
  close(saida)	
  	
  end	
  	
  
soluções	
  	
  	
  logo	
  o	
  ajuste	
  Linear	
  é:	
  𝑦 = −8.3187+ 6.770𝑥	
  	
  	
  O	
  gráfico	
  1	
  apresenta	
  o	
  ajuste	
  Linear	
  com	
  os	
  pontos	
  da	
  tabela	
  6.	
  
	
  
	
  
	
  
	
  
	
  
Ajuste	
  Exponencial,	
  função	
  na	
  forma	
  𝑦 = 𝑎!𝑎!!	
  ,	
  precisa	
  ser	
  linearizada.	
  Linearização	
  y	
  =	
  ln(y)	
  ;	
  	
   ln 𝑦 = 𝑙𝑛(𝑎!)+ 𝑥𝑙𝑛(𝑎!)	
  	
  Logo	
  𝑦!  =  𝑎!! + 𝑎!!𝑥	
  	
  
Solução	
  do	
  ajuste	
  exponencial	
  é:	
  𝑦 = 0.12349  𝑒!.!"#"!	
  	
  O	
  gráfico	
  2	
  apresenta	
  o	
  ajuste	
  exponencial	
  com	
  os	
  pontos	
  da	
  tabela	
  6.	
  	
  
	
  	
   Ajuste	
  potencial,	
  função	
  na	
  forma	
  𝑦 = 𝑎!𝑥!! ,	
  precisa	
  ser	
  linearizada.	
  Linearização	
  y	
  =	
  ln(y)	
  ;	
  x	
  =	
  ln(x).	
  	
   ln 𝑦 = ln 𝑎! + 𝑎!ln  (𝑥)	
  logo	
  𝑦!  =  𝑎!! + 𝑥!𝑎!	
  	
  solução	
  do	
  ajuste	
  potencial	
  é:	
  𝑦 = 0.47227  𝑥!.!"#$	
  	
  O	
  gráfico	
  3	
  apresenta	
  o	
  ajuste	
  potencial	
  com	
  os	
  pontos	
  da	
  tabela	
  6.	
  
	
  
O	
   gráfico	
   4	
   apresenta	
   todos	
   os	
   ajustes	
   com	
   os	
   pontos	
   da	
   tabela	
   6.	
  Utilizando	
  –	
  se	
  conceitos	
  de	
  estatística,	
  regra	
  de	
  Pearson	
  obtemos	
  que	
  o	
  melhor	
  ajuste	
  e	
  o	
  exponencial.	
  Visualmente	
  nota	
  –	
  se	
  também	
  que	
  o	
  ajuste	
  exponencial	
  aproxima	
  –	
  se	
  melhor	
  dos	
  pontos	
  do	
  que	
  os	
  outros	
  ajustes.	
  	
  
	
  	
  	
  	
  
5. Interpolação	
  	
  	
   O	
   problema	
   de	
   interpolação	
   consiste	
   em,	
   dado	
   um	
   conjunto	
   de	
   pares	
  ordenados	
   (x0,	
   y0),	
   (x1,	
   y1),	
   .	
   .	
   .	
   ,	
   (xn,	
   yn),	
   determinar	
   uma	
   função	
   g,	
   designada	
  função	
  interpoladora,	
  tal	
  que	
  :	
  g xi =  yi  , i = 0,1, . . . ,n.	
  	
   Os	
   valores	
   x0,x1,...,xn	
   designam-­‐se	
   por	
   nós	
   de	
   interpolação	
   e	
   devem	
  satisfazer	
   a	
   condição	
   i≠ j	
   ⇒	
   xi	
  ≠   xj,	
   ou	
   seja,	
   serem	
   todos	
   diferentes.	
   Os	
  correspondentes	
  valores	
  y0,	
  y1,	
  .	
  .	
  .	
  ,	
  yn	
  designam-­‐se	
  por	
  valores	
  nodais.	
  Perante	
   um	
   dado	
   problema	
   de	
   interpolação	
   será	
   necessário	
   ter	
   em	
  consideração	
   diversas	
   questões,	
   das	
   quais	
   se	
   destacam	
   a	
   escolha	
   da	
   classe	
   de	
  funções	
   interpoladoras	
   a	
   utilizar	
   e	
   a	
   forma	
   de	
   determinar	
   concretamente	
   a	
  função	
  interpoladora.	
  	
  O	
  problema	
  de	
  interpolação	
  tem	
  aplicações	
  em	
  diversas	
  situações	
  como:	
  -­‐ O	
   cálculo	
   de	
   funções	
   fornecidas	
   por	
   tabelas	
   quando	
   se	
   pretende	
   avaliar	
  função	
  em	
  pontos	
  não	
  tabelados.	
  -­‐ Quando	
  apenas	
  se	
  conhecem	
  os	
  valores	
  de	
  uma	
  função	
  em	
  certos	
  pontos,	
  por	
   exemplo	
   resultantes	
   de	
   medidas	
   experimentais,	
   e	
   se	
   pretende	
   avaliar	
   a	
  função	
  em	
  novos	
  pontos.	
  -­‐ Aproximação	
   de	
   funções	
   cujo	
   cálculo	
   seja	
   complexo	
   ou	
   exija	
   grande	
  esforço.	
  
-­‐ A	
  base	
  de	
  muitos	
  métodos	
  numéricos.	
  	
   Dentre	
   os	
   tipos	
   de	
   interpolação	
   se	
   destacam	
   interpolação	
   linear,	
  interpolação	
   polinomial	
   e	
   interpolação	
   trigonométrica.	
   Neste	
   trabalho	
   iremos	
  apresentar	
  a	
   interpolação	
  polinomial	
  pela	
  forma	
  de	
  Lagrange	
  que	
  será	
  definido	
  mas	
  adiante.	
  Mas	
  muitas	
  vezes	
  porem	
  ocorre	
  o	
   inconveniente	
  de	
  uma	
   função	
  g	
  qualquer	
  apresentar	
  um	
  valor	
  muito	
  discrepante	
  em	
  relação	
  a	
  um	
  dado	
  valor	
  de	
  f,	
   valor	
   real,	
   é	
   ai	
   que	
   surge	
   a	
   necessidade	
   de	
   se	
   utilizar	
   um	
   outro	
   tipo	
   de	
  interpolação	
   o	
   Spline,	
   que	
   interpola	
   f	
   em	
   grupos	
   de	
   poucos	
   pontos.	
   Nesse	
  trabalho	
  usaremos	
  o	
  Spline	
  Cúbico	
  que	
  também	
  será	
  definido	
  mais	
  adiante.	
  	
  
	
  
	
  
Forma	
  Lagrange	
  	
  
	
   Num	
   conjunto	
   de	
   pontos	
   distintos	
  (𝑥!)!!!! .	
   Os	
   polinómios	
   (de	
   grau	
   n)	
  definidos	
  pela	
  expressão:	
  	
   𝐿!  (𝑥) =   𝑥 − 𝑥!𝑥! − 𝑥!      𝑘 = 0,1…!!!!,!!! 	
  	
  	
  designam-­‐se	
  por	
  polinómios	
  de	
  Lagrange,	
  relativos	
  aos	
  pontos	
  x0,	
  x1,	
  .	
  .	
  .	
  ,	
  xn.	
  O	
   polinómio	
   interpolador	
   na	
   forma	
   de	
   Lagrange	
   é	
   obtido	
   como	
   uma	
  combinação	
  linear	
  dos	
  polinómios	
  de	
  Lagrange	
  relativos	
  aos	
  pontos	
  em	
  questão.	
  Os	
  coeficientes	
  desta	
  combinação	
  linear	
  serão	
  os	
  valores	
  de	
  y	
  a	
  interpolar.	
  O	
  polinómio	
  p,	
  de	
  grau	
  menor	
  ou	
   igual	
   a	
  n,	
  que	
   interpola	
  o	
   conjunto	
  de	
  valores	
  y0,	
  y1,	
  ...	
  ,	
  yn	
  nós	
  pontos	
  distintos	
  x0,	
  x1,	
  ...	
  ,	
  xn	
  é	
  dado	
  por:	
  [3]	
  	
  	
   𝑝! 𝑥 = 𝑦!𝐿!(𝑥)!!!! =  𝑦! 	
  	
  	
  Seja	
  a	
  função	
  𝑓 𝑥 =   𝑥 ∗ sen 𝑥	
  conhecida	
  nos	
  pontos	
  da	
  tabela	
  7.	
  	
  	
  	
  
x	
   y	
  0	
   0,031570	
  2	
   1.278745	
  4	
   -­‐1.521074	
  6	
   -­‐0.591799	
  8	
   2.793026	
  Tabela	
  7	
  	
   Temos	
   no	
   gráfico	
   5	
   a	
   curva	
   gerada	
   pela	
   interpolação	
   de	
   Lagrange,	
  comparado	
  com	
  a	
  curva	
  exata.	
  
	
  
	
  
	
  
	
  
Fortran	
  7	
  (implementação)	
  	
  
program Lagrange 
implicit none 
 
real a, h, x(100), f(100), Paux(100), P 
integer i, j, entrada, saida, n 
 
entrada = 40 
saida = 20 
h = 0.1 !intervalo com o qual os pontos ser∆o imprimidos na 
sa°da 
 
open (unit=entrada,file='lagrangeentrada.txt') 
open (unit=saida,file='lagrangesaida.txt') 
 
read (entrada,*) n !n£mero de pontos 
 
do i=1,n,1 
read(entrada,*) x(i), f(i) 
end do 
 
a = x(1) 
 
do while (a < x(n)) 
 
do i=1,n,1 
Paux(i) = 1 
end do 
 
do j=1,n,1 
 
do i=1,j-1,1 
Paux(j) = Paux(j)*(a - x(i))/(x(j) - x(i)) 
end do 
 
do i=j+1,n,1 
Paux(j) = Paux(j)*(a - x(i))/(x(j) - x(i)) 
end do 
 
Paux(j) = Paux(j)*f(j) 
 
end do 
 
P=0 
 
do i=1,n,1 
P = P + Paux(i) 
end do 
 
write(saida,*) a,P 
 
a = a + h 
 
end do 
 
close (entrada) 
close (saida) 
 
end 
 	
  	
  
Forma	
  Spline	
  Cúbica	
  
 Suponha	
   que	
   f(x)	
   esteja	
   tabelada	
   nos	
   pontos	
  𝑥! ,	
  𝑖 = 0,1,2,… ,𝑛	
  afunção	
  𝑆!(𝑥)	
  é	
   chamada	
   spline	
   cúbica	
   interpolante	
   de	
   f(x)	
   nos	
   nós	
  𝑥! ,	
  𝑖 = 0,1,2,… ,𝑛	
  se	
  existem	
  n	
  polinômios	
  de	
  grau	
  3,	
  𝑆!(𝑥),	
  𝑘 = 1,… ,𝑛  tal	
  que:	
  	
  𝑆! 𝑥 = 𝑎!(𝑥 − 𝑥!)! + 𝑏!(𝑥 − 𝑥!)! + 𝑐!(𝑥 − 𝑥!)! + 𝑑!,𝑘 = 1,2,… ,𝑛.  	
  
	
  𝑆!(𝑥)	
  Tem	
  derivadas	
   de	
   primeira	
   e	
   segunda	
   ordem	
   contínuas,	
   o	
   que	
   faz	
  com	
  que	
  a	
  curva	
  𝑆!  (𝑥)	
  não	
  tenha	
  picos	
  e	
  nem	
  troque	
  abruptamente	
  de	
  curvatura	
  nos	
  nós.	
  [6]	
  	
  Seja	
   a	
   mesma	
   função	
  𝑓 𝑥 =   𝑥 ∗ sen 𝑥 	
  do	
   exemplo	
   anterior,	
   conhecida	
   nos	
  pontos	
  da	
  tabela	
  7.	
  Temos	
  o	
  gráfico	
  da	
  curva	
  gerada	
  pela	
  interpolação	
  de	
  Spline	
  Cúbica.	
  	
   	
  Temos	
  no	
  gráfico	
  6	
  a	
  curva	
  gerada	
  por	
  Spline	
  Cúbica,	
   comparado	
  com	
  a	
  curva	
  da	
  solução	
  exata.	
  	
  
	
  
	
  
	
  
Fortran	
  8	
  (implementação)	
  	
  	
  
 
program splinecubico 
implicit none 
 
real X(0:100), Y(0:100), A(100,100), B(100), G(0:100), h, t, z 
real ak, bk, ck, dk, s 
integer n, i, j, k, entrada, saida 
 
entrada = 100 
saida = 200 
 
open (unit=entrada,file='splineentrada.txt') 
open (unit=saida,file='splinesaida.txt') 
 
read(entrada,*) n !numero de pontos 
read(entrada,*) h !intervalo entre os pontos 
 
!pontos 
 
do i=0,n-1,1 
read(entrada,*) X(i), Y(i) 
end do 
 
close (entrada) 
 
do i=1, n-2, 1 
do j=1, n-2,1 
A(i,j) = 0 
end do 
end do 
 
do i=1, n-2,1 
A(i,i) = 4*h 
end do 
 
do i=1,n-3,1 
A(i,i+1)=h 
A(i+1,i)=h 
end do 
 
!imprimindo matriz 
!do i=1, n-2, 1 
! write(saida,*) (A(i,j), j=1, n-2, 1) 
!end do 
 
do i=1,n-2,1 
B(i) = (6/h)*(Y(i+1) -2*Y(i) + Y(i-1)) 
end do 
 
!encontrando os gi com gausseidel 
 
do i=0,n-1,1 
G(i) = 0 
end do 
 
do k=1,100,1 
 
do i=1,n-2,1 
 
G(i) = B(i) 
 
do j=1,i-1,1 
G(i) = G(i) - A(i,j)*G(j) 
end do 
 
do j=i+1,n,1 
G(i) = G(i) - A(i,j)*G(j) 
end do 
 
G(i) = G(i)/A(i,i) 
 
end do 
end do 
 
!encontrando coeficientes para cada Si e imprimindo os pontos 
 
z = x(0) 
t = 0.1 !intervalo com o qual os pontos serao mostrados na 
saida 
 
do k=1, n-1,1 
 
ak = (G(k) - G(k-1))/(6*h) 
bk = G(k)/2 
ck = (Y(k) - Y(k-1))/h + (2*h*G(k) + G(k-1)*h)/6 
dk = Y(k) 
 
do while (z <= x(k)) 
 
s = ak*(z - x(k))**3 + bk*(z - x(k))**2 + ck*(z - x(k)) + dk 
 
write (saida,*) z,s 
 
z = z + t 
 
end do 
 
end do 
 
close (saida) 
 
end 
 No	
  gráfico	
  7	
  foi	
  exposto	
  as	
  curvas	
  de	
  Lagrange,	
  Spline	
  cúbica	
  e	
  a	
  curva	
  da	
  solução	
  exata.	
  	
  
	
  	
   Pode-­‐se	
  observar	
  visualmente	
  que	
  o	
  método	
  Spline	
  cubica	
  se	
  aproximou	
  melhor	
  da	
  curva	
  exata.	
  Sendo	
  assim	
  e	
  o	
  melhor	
  método	
  para	
  interpolar	
  os	
  pontos	
  da	
  tabela	
  7.	
  	
  	
  	
  	
  
6. Integração	
  numérica	
  	
  	
   A	
   integração	
   numérica	
   é	
   uma	
   técnica	
   comumente	
   empregada	
   na	
  determinação	
  de	
  uma	
   integral	
  definida,	
  cuja	
   função	
  ou	
  não	
  é	
  disponível	
  ou	
  não	
  
possui	
   uma	
   solução	
   analítica.	
   Ela	
   consiste	
   na	
   aproximação	
   de	
   uma	
   integral	
  definida	
  do	
  tipo:	
  	
   𝐼 = 𝑓 𝑥 𝑑𝑥!! 	
  por	
  uma	
  soma	
  do	
  tipo:	
   𝐼 = 𝑓(𝑥)𝑑𝑥!! ≅ 𝑤!𝑓(𝑥!)∆𝑥!!!! 	
  	
  na	
   qual	
  𝑓 𝑥! 	
  são	
   os	
   valores	
   da	
   função	
  𝑓(𝑥), 	
  ∆𝑥 = 𝑥!!! − 𝑥!  e	
  𝑤!  é	
   um	
   valor	
  numérico	
  de	
  ponderação	
  que	
  também	
  é	
  conhecido	
  por	
  função	
  peso.	
  Nesse	
  trabalho,	
  nos	
  restringiremos	
  aos	
  métodos	
  de	
   integração	
  numérica	
  para	
   intervalo	
   ∆𝑥 	
  constante.	
   A	
   determinação	
   desses	
   métodos	
   consiste	
  basicamente	
  em	
  avaliar	
  o	
  valor	
  da	
  função	
  peso	
  𝑤! .	
  São	
  os	
  métodos	
  	
  do	
  Trapézio	
  e	
  Simpson.	
  	
  
Regra	
  do	
  Trapézio,	
  consiste	
  em	
  se	
  aproximar	
  o	
  valor	
  da	
  função	
  contínua	
  de	
   𝑓(𝑥) 	
  no	
   intervalo	
   [a,b]	
   por	
   uma	
   função	
   de	
   primeira	
   ordem;	
   isto,	
  geometricamente,	
  é	
  equivalente	
  a	
  aproximar	
  uma	
  curva	
  qualquer	
  por	
  uma	
  reta.	
  Desta	
  forma,	
  a	
  área	
  sob	
  a	
  função	
  𝑓(𝑥),	
  que	
  é	
  equivalente	
  à	
  integral	
  dessa	
  função,	
  é	
  aproximada	
  pela	
  área	
  do	
  trapézio	
  cuja	
  largura	
  é	
  igual	
  a	
  (b	
  –	
  a)	
  e	
  a	
  altura	
  média	
  igual	
  a	
  [𝑓(𝑎)  +  𝑓(𝑏)]/𝑛.	
  Fazendo-­‐se	
  ∆𝑥 =  𝑏  –  𝑎,	
  a	
  fórmula	
  para	
  a	
  integral	
  pode	
  ser	
  escrita	
  como:	
  [5]	
   𝑓 𝑥 𝑑𝑥!! ≅ ℎ2 [𝑓 𝑎 + 𝑓 𝑏 ]	
  	
  Calcular	
  	
   𝑒!𝑑𝑥  !! empregando	
  a	
  regra	
  do	
  trapézios.	
  Solução	
  Exata	
   Solução	
  pelo	
  método	
  do	
  trapézio	
   Erro	
  gerado	
  1.71828	
   172722	
   0.0089	
  Tabela	
  8	
  Como	
   ilustrado	
   na	
   tabela	
   8	
   a	
   solução	
   exata	
   da	
   equação	
   𝑒!𝑑𝑥  !! é  1.71828 .	
  Utilizando	
  n=4	
  obtemos	
  com	
  o	
  método	
  do	
  trapézio	
  uma	
  solução	
  para	
  a	
  mesma	
  de	
  1.72722.	
   Calculamos	
   o	
   	
   erro	
   pelo	
   modulo	
   da	
   diferença	
   entre	
   as	
   soluções	
  𝑠!"#$# − 𝑠!"#$%&'( ,	
  tivemos	
  um	
  erro	
  de	
  0.0089.	
  
Fortran	
  9	
  (implementação)	
  
program trapézio 
implicit none 
real h,s,a,b,f ! a e b sao limite de integrcao 
integer n,solucao 
solucao=21 
open(unit=solucao,file='solucao.txt') 
a=0. 
b=1. 
n=4 
h=(b-a)/n 
s=(h/2)*(f(0*h)+(2*f(1*h))+(2*f(2*h))+(2*f(3*h))+f(4*h)) 
write(solucao,*)'regra trapezio' ,s 
end 
 
function f(x) 
real f,x 
f = exp(x) 
return 
end 
	
  
Regra	
  de	
  Simpson	
  consiste	
  na	
  aproximação	
  da	
  função	
  contínua	
  f(x)	
  no	
  intervalo	
  [a,b]	
  por	
  uma	
  função	
  de	
  segunda	
  ordem,	
  ou	
  seja,	
  na	
  aproximação	
  de	
  uma	
  curva	
  por	
  uma	
  parábola.	
  A	
  fórmula	
  para	
  a	
  integral	
  tem	
  a	
  forma: 	
   𝑓 𝑥 𝑑𝑥!! ≅ ℎ3 [𝑓 𝑥! + 4𝑓 𝑥! + 𝑓(𝑥!)	
  	
  Esta	
  fórmula	
  também	
  é	
  conhecida	
  como	
  regra	
  de	
  Simpson	
  1/3	
  por	
  causa	
  do	
  fator	
  que	
  multiplica	
  h.	
  Observar	
  que	
  são	
  necessários,	
  no	
  mínimo,	
  três	
  valores	
  de	
  f(𝑥!)	
  para	
  se	
  calcular	
  a	
  integral	
  pela	
  regra	
  de	
  Simpson.	
  Na	
  notação	
  usada	
  aqui,	
  𝑥!	
  =	
  a,	
  𝑥!	
  =	
   b	
   e	
  𝑥!	
  é	
   um	
   ponto	
   equidistante	
   de	
  𝑥!	
  e	
  𝑥!.	
   Para	
   n	
   intervalos	
  ∆𝑥, a	
   fórmula	
  pode	
  ser	
  escrita	
  como:	
  	
   𝑓 𝑥 𝑑𝑥!! ≅ ℎ3 [𝑓 𝑥! + 4𝑓 𝑥! + 2𝑓 𝑥! +⋯+ 2𝑓 𝑥!!! + 4𝑓 𝑥!!! + 𝑓(𝑥!)	
  	
  na	
  qual	
  n	
  é	
  par	
  (correspondendo	
  a	
  um	
  número	
  par	
  de	
  intervalos	
  de	
  integração)	
  ou,	
   equivalentemente,	
   a	
   regra	
   de	
   Simpson	
   1/3	
   só	
   pode	
   ser	
   aplicada	
   para	
   um	
  número	
  ímpar	
  de	
  pontos	
  𝑥! ,	
  f(𝑥!).	
  [6]	
  	
  Calcular	
  	
   𝑒!𝑑𝑥  !! empregando	
  a	
  regra	
  de	
  Simpson	
  1/3.	
  Solução	
  Exata	
   Solução	
  pelométodo	
  do	
  Simpson	
   Erro	
  gerado	
  1.71828	
   1.71832	
   0.00004	
  Tabela	
  9	
  	
   Como	
   ilustrado	
   na	
   tabela	
   9,	
   temos	
   uma	
   solução	
   exata	
   para	
   a	
   equação	
  𝑒!𝑑𝑥  !! de	
   1.71828.	
   Utilizando	
   n=4	
   pela	
   regra	
   do	
   Simpson	
   obtivemos	
   uma	
  solução	
   de	
   1.7832.	
   Calculamos	
   o	
   	
   erro	
   pelo	
   modulo	
   da	
   diferença	
   entre	
   as	
  soluções	
   𝑠!"#$# − 𝑠!"#$!%& ,	
  tivemos	
  um	
  erro	
  de	
  0.00004.	
  	
  
Fortran	
  10	
  (implementação)	
  	
  
program simpson 
implicit none 
real h,s,a,b,f ! a e b sao limite de integrcao 
integer n,solucao 
solucao=21 
open(unit=solucao,file='solucao.txt') 
a=0. 
b=1. 
n=4 
h=(b-a)/n 
s=(h/2)*(f(0*h)+(4*f(1*h))+(2*f(2*h))+(4*f(3*h))+f(4*h)) 
write(solucao,*)'regra simpson' ,s 
end 
 
function f(x) 
real f,x 
 
f = exp(x) 
return 
end 	
  	
  	
  
	
  
7. Resolução	
  de	
  EDO	
  
	
   	
  Muitos	
   problemas	
   encontrados	
   em	
   engenharia	
   e	
   outras	
   ciências	
   podem	
  ser	
   formulados	
   em	
   termos	
   de	
   equações	
   diferenciais.	
   Por	
   exemplo,	
   trajetórias	
  balísticas,	
  teoria	
  dos	
  satélites	
  artificiais,	
  estudo	
  de	
  redes	
  elétricas,	
  curvaturas	
  de	
  vigas,	
   estabilidade	
   de	
   aviões,	
   teoria	
   das	
   vibrações,	
   reações	
   químicas	
   e	
   entre	
  	
  outras.	
   Esse	
   capitulo	
   7	
   tem	
   como	
   objetivo	
   apresentar	
  métodos	
   numericos	
   que	
  resolvam	
   equações	
   diferenciais	
   ordinárias,	
   vamos	
   abordar	
   nesse	
   trabalho	
   o	
  método	
  de	
  Euler	
  e	
  métodos	
  de	
  Runge	
  Kutta.	
  	
  
O	
  método	
  de	
  Euler,	
  consiste	
  em	
  tomar	
  a	
  aproximação	
  de	
  primeira	
  ordem	
  de  𝑥(𝑡 + ℎ):	
  	
   𝑥 𝑡 + ℎ = 𝑥 𝑡 + 𝑥! 𝑡 ℎ + 𝑜(ℎ)	
  	
   	
  O	
  resto	
  o(h)	
  indica	
  que	
  o	
  erro	
  em	
  tomar	
  essa	
  aproximação	
  será,	
  em	
  geral	
  da	
   ordem	
   de	
  ℎ! 	
  ou	
   mais.	
   Como	
  𝑥! = 𝑓 ,	
   o	
   método	
   sugere	
   que,	
   na	
   pratica,	
  obtenhamos	
  𝑥!!! = 𝑥 𝑡!!! = 𝑥(𝑡! + ℎ)	
  como:	
  	
  	
   𝑥!!! = 𝑥! + 𝑓 𝑡! , 𝑥! ℎ	
  	
   	
  A	
   interação	
   a	
   partir	
   do	
  𝑥!acumulará	
   um	
   erro	
   da	
   ordem	
   de	
  ℎ!	
  em	
   cada	
  etapa	
   podendo	
   acumular	
   ao	
   final	
   um	
   erro	
   de	
  𝑛ℎ!.	
   Como	
  𝑛 = !!!! 	
  ,	
   então	
   o	
   erro	
  acumulado	
  será	
  da	
  ordem	
  de	
   𝑏 − 𝑎 ℎ	
  .	
  [3]	
  
Considerando	
   ,	
  𝑥!  𝑒  𝑣!igual	
   a	
   zero.	
  Calcule	
   a	
  distancia	
   “x”	
  percorrida	
  em	
  10s,	
  estabelecendo	
  uma	
  aceleração	
  “a”	
  constante	
  de	
  valor	
  5𝑚𝑠!!. Sendo a função 
𝑥(𝑡)	
  =2,5𝑡. Através	
   do	
   método	
   de	
   Euler,	
   encontramos	
   a	
   solução	
   aproximada	
   da	
  equação	
  Diferencial	
  Ordinária	
  do	
  problema	
  acima. 
 𝑥!  =  𝑥!!!  +  ℎ𝑣(𝑡)  
 
Fortran	
  11	
  (implementação)	
  	
  
Program mruv " 
implicit none " 
real xo,vo,h,t,to,x,tf,v,ex 
integer i,n,dado " 
dado=20 
open(unit=dado,file='exct.txt') 
tf= 10 " 
xo=0 " 
vo=0 " 
to=0 
h=0.1 " 
i=0 
n=(tf‐to)/h 
x=xo 
t=to 
write(dado,*)i,t,x,ex(t) 
doi=1,n,1 " 
x=x+h*v(t) " 
t=t+h 
write(dado,*)i,t,x,ex(t) 
end do " 
call system('pause') 
end " 
!.............. " 
function v(t) " 
real v,t,vo " 
v=vo+5*t " 
return " 
end " 
!...................... " 
function ex(t) " 
real ex,t " 
ex=2.5*t**2 " 
return " 
end 	
  	
  Utilizado	
  varios	
  valores	
  para	
  “h”	
  ,	
  observamos	
  no	
  gráfico	
  abaixo	
  visualmente	
  que	
  o	
  h=1	
  e	
  o	
  que	
  se	
  aproxima	
  mais	
  da	
  solução	
  exata.	
  
 
	
  	
  	
   Métodos	
   de	
   Runge-­‐Kutta,	
   são	
   métodos	
   de	
   passo	
   simples	
   requerem	
  apenas	
  derivadas	
  de	
  primeira	
  ordem	
  e	
  pode	
  fornecer	
  aproximações	
  precisas	
  com	
  erros	
  de	
  truncamento	
  da	
  ordem	
  de	
  ℎ!, ℎ!, ℎ!,	
  etc.	
  Todas	
  os	
  métodos	
  de	
  Runge-­‐Kutta	
  têm	
  a	
  seguinte	
  forma	
  geral:	
  	
  
	
   𝑦!!! = 𝑦! + ℎ𝜙(𝑥! ,𝑦! , ℎ)	
  	
  onde	
  𝜙, chamado	
   de	
   função	
   incremento,	
   é	
   uma	
   aproximação	
   conveniente	
   para	
  
f(x,y)	
  no	
  intervalo	
  𝑥! ≤ 𝑥 ≤ 𝑥!!!.[5]	
  	
   O	
   método	
   de	
   Runge-­‐Kutta	
   de	
   4ª	
   ordem	
   é	
   o	
   mais	
   usado na	
   solução	
  numérica	
   de	
   problemas	
   com	
   equações	
   diferenciais	
   ordinárias.	
   Nesse	
   trabalho	
  será	
  apresentado	
  um	
  algoritmo	
  para	
  Runge	
  Kutta	
  de	
  4ª	
  ordem.	
  	
  
Algoritmo	
  	
  	
  Início	
  Se	
  n	
  =	
  0	
  Então	
  Escreva(‘Divisão	
  imprópria’)	
  Senão	
  Se	
  n	
  <	
  0	
  Então	
  Escreva(‘Intervalo	
  inválido’)	
  Senão	
  Início	
  h	
  ← (xn	
  –	
  x0)/n	
  	
  y	
  ← y0	
  x	
  ← x0	
  	
  Para	
  i	
  =	
  1	
  até	
  n	
  faça	
  Início	
  K1	
  ← f(x,y)	
  	
  K2	
  ← f(x	
  +	
  h	
  /	
  2,	
  y	
  +	
  h	
  /	
  2	
  *	
  K1)	
  	
  
K3	
  ← f(x	
  +	
  h	
  /	
  2,	
  y	
  +	
  h	
  /	
  2	
  *	
  K2)	
  	
  K4	
  ← f(x	
  +	
  h,	
  y	
  +	
  h	
  *	
  K3)	
  	
  y←y+h	
  /	
  6*(K1+2*K2+2*K3+K4)	
  x←x+h	
  Fim-­‐Para	
  Escreva(‘Valor	
  da	
  função	
  f	
  no	
  ponto	
  k	
  é:	
  ’,	
  y)	
  	
  Fim	
  Fim	
  	
  	
  Variáveis	
  utilizadas	
  no	
  algoritmo:	
  	
  Reais:	
  h,	
  y,	
  x,	
  y0,	
  x0,	
  K1,	
  K2,	
  K3,	
  K4	
  Inteiras:	
  i,	
  n	
  	
  	
  
8. Conclusão Foi	
   demostrado	
   com	
   esse	
   trabalho	
   que	
   calculo	
   numérico,	
   é	
   usado	
   para	
  resolver	
  qualquer	
  tipo	
  de	
  problema	
  matemático	
  complexo.	
  	
  No	
   capítulo	
   2	
   o	
   trabalho	
   apresenta	
   resolução	
   de	
   raizes	
   através	
   de	
  métodos,	
  entre	
  os	
  três	
  métodos	
  apresentado	
  Bisseção,	
  Newton	
  e	
  Secante,	
  nota	
  –	
  se	
  que	
  a	
  melhor	
  aproximação	
  feita	
  pelos	
  métodos	
  foi	
  a	
  Secante	
  e	
  Newton,	
  sendo	
  que	
  o	
  método	
  da	
  secante	
  e	
  menos	
  complexo	
  para	
  se	
   feito	
  do	
  que	
  Newton,	
  pois	
  não	
  precisamos	
  da	
  derivada	
  da	
  função.	
  	
  No	
  capítulo	
  3	
  o	
  trabalho	
  apresento	
  resolução	
  para	
  sistemas	
  lineares	
  onde	
  trabalhamos	
   com	
   os	
   métodos	
   de	
   Eliminação	
   Gaussiana	
   e	
   o	
   de	
   Gauss	
   Seidel,	
  podemos	
  observar	
  com	
  as	
  tabelas	
  4	
  e	
  5	
  que	
  os	
  dois	
  método	
  são	
  eficientes.	
  	
  No	
  capítulo	
  4	
  foi	
  apresentado	
  ajuste	
  de	
  curvas	
  pelo	
  método	
  dos	
  mínimos	
  quadrados,	
   em	
   geral	
   todo	
   ajuste	
   fundamenta-­‐se	
   nesse	
   método,	
   o	
   Fortran	
   6	
  apresentado	
  faz	
  qualquer	
  ajuste;	
   linear,	
  parabólico	
  entre	
  outros,	
  só	
  depende	
  de	
  seus	
  dados	
  de	
  entrado	
  para	
  a	
  leitura	
  do	
  programa.	
  	
  No	
  capítulo	
  5	
  foi	
  apresentado	
  dois	
  métodos	
  para	
  interpolação,	
  Lagrange	
  e	
  Spline	
   cubica	
   bastantes	
   eficientes,	
   no	
   exemplo	
   apresentado	
   o	
   método	
   Splinecubica	
  foi	
  mas	
  eficiente,	
  mostrado	
  no	
  gráfico	
  7.	
  	
  No	
   capítulo	
   6	
   foi	
   apresentado	
   dois	
   métodos	
   para	
   se	
   fazer	
   integrações,	
  Trapézio	
  e	
  Simpson.	
  Olhando	
  –	
  se	
  as	
  tabelas	
  8	
  e	
  9	
  tem	
  –	
  se	
  que	
  para	
  o	
  exemplo	
  apresentado	
   o	
  método	
  de	
   Simpson	
   teve	
  melhor	
   solução	
   por	
   apresentar	
  menor	
  erro.	
   No	
  capítulo	
  7	
  apresentamos	
  métodos	
  para	
  resolução	
  de	
  EDO,	
  observa	
  –	
  se	
  com	
   o	
   gráfico	
   8	
   que	
   quanto	
  menor	
   o	
   valor	
   “h”	
   escolhido	
  melhor	
   aproximação	
  pelo	
  método	
  de	
  Euler.	
  Conceituo	
  –	
  se	
  o	
  método	
  de	
  Runge	
  Kutta	
  e	
  foi	
  apresentado	
  um	
  algoritmo	
  para	
  resoluções	
  de	
  4ª	
  ordem.	
  Enfim	
   calculo	
   numérico	
   facilita	
   a	
   resolução	
   de	
   problemas	
   físicos	
  matemáticos,	
  com	
  soluções	
  próximas	
  da	
  realidade	
  que	
  em	
  geral	
  seria	
  infinito	
  os	
  cálculos	
  na	
  mão	
  .	
  	
  
 	
  
9. Referências	
  	
  	
  [1]	
  Brito,	
  N.	
  S.	
  D.	
  Calculo	
  Numérico	
  [2]	
  Pilling,	
  S.	
  Calculo	
  Numérico	
  [3]	
  Matos,	
  A.	
  C.	
  C.	
  de.	
  Apontamentos	
  de	
  Análise	
  Numérica	
  [4]	
  A.L.	
  de	
  Bortoli,	
  C.	
  Cardoso,	
  M.P.G.	
  Fachin	
  &	
  R.D.	
  da	
  Cunha.	
  Calculo	
  numérico	
  com	
  aplicações.	
  [5]	
  Lira,	
  W.	
  W.	
  M.	
  L.	
  Apostila	
  de	
  Calculo	
  Numérico	
  [6]	
  Ruggiero,	
  M.	
  A.	
  G.	
  &	
  Lopes,	
  V.	
  L.	
  da	
  R.	
  Calculo	
  numérico	
  aspectos	
  teóricos	
  e	
  computacionais