Calculo Numerico
267 pág.

Calculo Numerico


DisciplinaCálculo Numérico16.815 materiais273.848 seguidores
Pré-visualização35 páginas
! *** valores contido nesta tabela: ***
! 5., 49. / 10., 105. / 15., 172. / 20., 253. / 25., 352. / 30., 473. / 35., 619. / 40.,
793.
open(2,file=\u2019interpolados.dat\u2019,status=\u2019unknown\u2019)
do i = 1,n
read(1,1)(tab(i,j),j=1,2)
1 format(2f10.0)
enddo
do k=1,npi
y(k) = 0.
x(1) = 12.
x(2) = 22.
x(3) = 31.
do i = 1,n
254 Cap´\u131tulo 10 - Algoritmos implementados em FORTRAN 90
sum = 1.
do 20 j=1,n
if(i.eq.j) goto 20
sum = sum * ((x(k) - tab(j,1))/(tab(i,1) - tab(j,1)))
20 continue
y(k) = y(k) + sum * tab(i,2)
enddo
write(2,*)k,x(k),y(k)
enddo
stop
end program
10.3.5 Me´todo de Simpson
! ************************************************
! ****** Me´todo de Simpson ******
! ****** Calcula integral de func¸o\u2dces ******
! ************************************************
! *** Para\u2c6metros de entrada: ***
! n : Nu´mero de partic¸o\u2dces;
! h : espac¸amento;
! a : limite inferior do intervalo;
! b : limite superior do intervalo;
! *** Para\u2c6metros de Sa´\u131da: ***
Fund. de Ca´lculo Nume´rico para Engenheiros 255
! y : valor da integral.
program simpson
parameter (n=3,a=-0.5,b=0.5)
real, dimension(n+1) ::x
open(6,file=\u2019simpson.dat\u2019,status=\u2019unknown\u2019)
x(0)=a
y=0.
do i=1,n+1
h=(b-a)/n
x(i)=x(i-1)+h
enddo
do i=0,n/2-1
y=y+h/3.*(f(x(2*i),h)+4*f(x(2*i+1),h)+f(x(2*i+2),h))
write(6,*)i,y
enddo
write(*,*)\u2019Valor de h:\u2019,h
write(*,*)\u2019Solucao para Simpson:\u2019,y
end program
Real Function f (x,h)
if(x==0.)then
f = exp(-h)*log(h)
else
256 Cap´\u131tulo 10 - Algoritmos implementados em FORTRAN 90
f = exp(-x)*log(x)
endif
end
10.3.6 Me´todo de Runge-Kutta 2
! ****************************************************
! ****** Me´todo de Runge-Kutta Segunda Ordem ******
! ****** Resolve PVI de 1a Ordem ******
! ****************************************************
! *** Para\u2c6metros de entrada: ***
! ni : Nu´mero de partic¸o\u2dces;
! h : espac¸amento;
! a : limite inferior do intervalo;
! b : limite superior do intervalo;
! *** Para\u2c6metros de Sa´\u131da: ***
! x: valor das abcissas;
! y: valor da ordenadas.
program Rk2
parameter (ni=10, a=0, b=1)
real, dimension(ni+1) ::y,x,k1,k2
open(3, file =\u2019rk2.dat\u2019, status =\u2019unknown\u2019)
h=(b-a)/ni
x(0)=a
y(0)=1.
do i=1,ni+1
Fund. de Ca´lculo Nume´rico para Engenheiros 257
k1(i)=h*f(x(i-1),y(i-1))
k2(i)=h*f(x(i-1)+h,y(i-1)+h*k1(i))
y(i)=y(i-1)+(k1(i)+k2(i))/2.
x(i)= a+i*h
enddo
do i = 1,ni+1
write(3,20)x(i-1),y(i-1),k1(i),k2(i),f(x(i-1),y(i-1))
write(*,20)x(i-1),y(i-1),k1(i),k2(i),f(x(i-1),y(i-1))
enddo
20 format(4e16.7) end program
Real Function f (x,y)
f = (x**3)/2 + y/(x**2)
end
10.3.7 Me´todo de Runge-Kutta 4 para sistemas
Apresenta-se um programa computacional implementado na linguagem FORTRAN
90 para calcular um sistema ele´trico simples ( RCL ) utilizando o me´todo de Runge-
Kutta de quarta ordem para sistemas.
! ***************************************************
! ****** Me´todo de Runge-Kutta Quarta Ordem ******
! ****** Resolve Sistemas de PVI de 1a Ordem ******
! ***************************************************
258 Cap´\u131tulo 10 - Algoritmos implementados em FORTRAN 90
! *** Para\u2c6metros de entrada: ***
! ni : Nu´mero de partic¸o\u2dces;
! h : espac¸amento;
! a : limite inferior do intervalo
! b : limite superior do intervalo
! *** Para\u2c6metros de Sa´\u131da: ***
! x : valor das abcissas
! y : valor da ordenadas.
program sistemark42
external f,g
parameter(a=0,b=0.008,h=0.001)
integer:: i,ni
real:: k1,k2,k3,k4,L1,L2,L3,L4,x,y,z
open(2,file=\u2019sistema.dat\u2019,status=\u2019unknown\u2019)
x = 0.
y = 0.
z = 0.
ni = (b-a)/h
write(2,*) \u2019 tempo (t)\u2019,\u2019 Carga (q)\u2019,\u2019 Corrente (I)\u2019
2 format(3e16.7)
do i = 1,ni
k1 = h*f(x,y,z)
L1 = h*g(x,y,z)
k2 = h*f(x+0.5*h,y+0.5*k1,z+0.5*L1)
Fund. de Ca´lculo Nume´rico para Engenheiros 259
L2 = h*g(x+0.5*h,y+0.5*k1,z+0.5*L1)
k3 = h*f(x+0.5*h,y+0.5*k2,z+0.5*L2)
L3 = h*g(x+0.5*h,y+0.5*k2,z+0.5*L2)
k4 = h*f(x+h,y+k3,z+L3)
L4 = h*g(x+h,y+k3,z+L3)
y = y + 1./6.*(k1+2*(k2+k3)+k4)
z = z + 1./6.*(L1+2*(L2+L3)+L4)
x = x+h
write(2,2)x,y,z
write(*,*)x,y,z
enddo
end program
real function f(x,y,z)
f = z+0.*x +0.*y
return
end
real function g(x,y,z)
g = 5000. - 500.*z - 5000.*y + 0.*x
return
end
260 Cap´\u131tulo 10 - Algoritmos implementados em FORTRAN 90
10.3.8 Me´todo de Diferenc¸as Finitas
10.3.8.1 Equac¸a\u2dco do calor
A implementac¸a\u2dco a seguir resolve um problema de valor de contorno em diferenc¸as
finitas centrais para a varia´vel u(x,t) que pode representar a temperatura de um
corpo num meio homoge\u2c6neo.
! *******************************************************
! ****** Me´todo de Diferenc¸as Finitas Centrais ******
! ****** Resolve a equac¸a\u2dco do calor para uma varia´vel ******
! *******************************************************
! *** Para\u2c6metros de entrada: ***
! nx : Partic¸o\u2dces em x;
! nt : Passo de tempo;
! x : Espac¸amento em x;
! y : Espac¸amento em t;
! u : varia´vel dependente;
! *** Para\u2c6metros de sa´\u131da: ***
! u : Temperatura final.
program calor
integer :: nt, nx, t, i
real, dimension(:), allocatable :: u
real :: x, y
nx = 20
nt = 2000
open(2,file=\u2019eqcalor.dat\u2019,status=\u2019unknown\u2019)
y = 1.0/(nt+1)
Fund. de Ca´lculo Nume´rico para Engenheiros 261
x = 1.0/(nx+1)
allocate(u(0:nx+1))
u(0) = 0.0
u(1:nx) = 1.0
u(nx+1) = 0.0
do t = 1,nt u(1:nx) = 1.0-2.0*(y/(x**2))*u(1:nx) + (y/(x**2))*(u(2:nx+1) +
u(0:nx-1))
enddo
write(*,\u2019(/4(f16.14,1x))\u2019)(u(i),i=0,nx+1)
deallocate(u)
end
10.3.8.2 Equac¸a\u2dco da onda
! ******************************************************
! ****** Me´todo de Diferenc¸as Finitas Centrais ******
! ****** Resolve a equac¸a\u2dco do onda para uma varia´vel ******
! ******************************************************
! *** Para\u2c6metros de entrada: ***
! ni : Partic¸o\u2dces em x;
! nj : Partic¸o\u2dces em y;
! x : Espac¸amento nas abcissas;
! t : Espac¸amento nas ordenadas;
! w : varia´vel dependente;
! *** Para\u2c6metros de Sa´\u131da: ***
! w : Oscilac¸a\u2dco final.
262 Cap´\u131tulo 10 - Algoritmos implementados em FORTRAN 90
program onda
integer :: i,j
integer, parameter:: ni = 40 , nj = 60
real, dimension (ni+1,nj+1):: f,x,w,t
real:: h,z,k,pi
open(6, file = \u2019eqonda.dat\u2019, status = \u2019unknown\u2019)
pi= acos(-1.0)
h = 1./4
k = 1./4
z = 0.5/0.25
do i = 1,ni
do j = 1,nj
x(i,j) = i*h
enddo
enddo
do i = 1,ni
do j = 1,nj
f(i,j)= sin(pi*x(i,j))
enddo
enddo
do j = 1,nj
w(0,j) = 0.
w(ni,j) = 0.
enddo
w(0,0) = f(0,0)
Fund. de Ca´lculo Nume´rico para Engenheiros 263
w(ni,0) = f(1,0)
do i = 1, ni-1
w(i,0) = f(i*h,0)
w(i,1) = (1 - z**2)*f(i*h,0) + z**2/2*(f((i+1)*h,0) + f((i-1)*h,0))
enddo
do j = 1,nj-1
do i = 1,ni-1
w(i,j+1) = 2*(1-z**2)*w(i,j) + z**2*(w(i+1,j)+ w(i-1,j)) + w(i,j-1)
enddo
enddo
do j = 1,nj
do i = 1,ni
t(i,j) = j*k
x(i,j) = i*h
write(6,*)t(i,j),x(i,j),f(i,j)
enddo
enddo
12 format(3e15.5)
stop
end program
No que segue apresenta-se algumas refere\u2c6ncias bibliogra´ficas que serviram como
base no desenvolvimento do texto.
264 Refere\u2c6ncias Bibliogra´ficas
BIBLIOGRAFIA
[1] Anderson, D. A., Tannehill, J. C., and Pletcher, R. H. Compu-
tational Fluid Mechanics and Heat Transfer, Second ed. 1997.
[2] Barroso, L. C., Arau´jo, M. M., Campos, F. F., Carvalho, M.
L. B., and Maia, M. L. Ca´lculo Nume´rico com Aplicac¸o\u2dces, 2 ed. HAR-
BRA Ltda, 1990.
[3] Burden, R. L., and Faires, J. D. Ana´lise Nume´rica. Thompson Learn-
ing, 2003.
[4] Butcher, J. C. The Numerical Analysis of Ordinary Differential Equa-
tions - Runge-Kutta and General