Calculo Numerico
35 pág.

Calculo Numerico


DisciplinaCálculo Numérico12.094 materiais248.254 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