Buscar

06_Integrales_definidas

Prévia do material em texto

NOTAS DE CLASE
Tema: Evaluación de Integrales definidas
 - Métodos de Newton-Cotes
 - Trapecios
 - Simpson
 - Cuadraturas Gaussianas
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
program integ_splines1d
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
implicit none
integer :: i,ndat,ncalc
real(8), allocatable :: x(:),y(:),y2(:)
real(8) :: yp1,ypn
real(8) :: x_int,y_int
real(8) :: func,integral
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
open(10,file='integ_splines.in',status='old')
read(10,*)ndat
allocate(x(ndat),y(ndat),y2(ndat))
read(10,*)x(1),x(ndat)
read(10,*)yp1,ypn
close(10)
do i=1,ndat
 x(i)=x(1)+(i-1)*(x(ndat)-x(1))/(ndat-1)
 y(i)=func(x(i))
enddo
call spline(x, y, ndat, yp1, ypn, y2)
integral=0._8
do i=1,ndat-1
 integral=integral+ &
 (x(i+1)-x(i))/2._8* &
 (y(i)+y(i+1)- &
 (x(i+1)-x(i))**3/12._8*(y2(i)+y2(i+1)))
enddo
write(*,*)'Integral=',integral
deallocate(x,y,y2)
end program integ_splines1d
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
real(8) function func(x)
implicit none
real(8) :: x
 func=1._8/(1._8+12._8*x*x)
end function func
 10 <- ndat
-1. 1. <-- x(1), x(ndat)
 1.e30 1.e30 <-- yp1, ypn (if > 0.99e+30 -> Natural splines)
File integ_splines.in
subroutine spline (x,y,n,yp1,ypn,y2) : VER CAPÍTULO: AJUSTE DE DATOS 
gfortran integ_splines.f90 spline.f90 -o integ_splines.x
Resultado
Ejercicio: Seleccionar varias integrales y comparar 
los resultados obtenidos con los métodos de Trapecios, 
Simpson y Splines Cúbicos, en similares condiciones 
(usar una grilla de número impar de datos) 
Integral= 0.74325
 
Carl Friedrich Gauss: 
30/04/1777, Brunswick, 
Alemania – 23/02/1855, 
Gottingen, Alemania
Matemático, astrónomo, 
geodesta y físico alemán 
que contribuyó 
significativamente en 
muchos campos, incluida 
la teoría de números, el 
análisis matemático, la 
geometría ...
 
 
 
 
 
 
 
 
 
Ejemplo: Área de un círculo con Simpson ⅓ y MonteCarlo (MC)
module commonvars
implicit none
integer(2) :: n=111 ! n debe ser impar y mayor que 1 !!!!
end module commonvars
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
program area_circulo_simpson
use commonvars
implicit none
real(8) :: a,b,simpson
a=0._8
b=1._8
write(*,*)'El area es:',simpson(a,b)
end program area_circulo_simpson
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
real(8) function simpson(a,b)
use commonvars
implicit none
integer(2) :: i
real(8), intent(in) :: a,b
real(8) :: h
real(8) :: fun
h=(b-a)/(n-1) ; simpson=(fun(a)+fun(b)+4*fun(a+h))*h/3
do i=2,n-2,2
 simpson=simpson+(2._8*fun(a+i*h)+4._8*fun(a+(i+1)*h))*h/3._8
enddo
end function simpson
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
real(8) function fun(x)
use commonvars
implicit none
integer(2) :: i
real(8) :: x,a,b,h
a=0.5_8-sqrt(0.25_8-(x-0.5_8)**2) ; b=-a+1._8 ; h=(b-a)/(n-1)
fun=(1._8+1._8+4._8)*h/3._8
do i=2,n-2,2
! y=a+i*h
 fun=fun+(2._8*1._8+4._8*1._8)*h/3._8
enddo
end function fun
1/2
1/2
x
y
program area_circulo_MC
implicit none
integer(2) :: i,n
integer :: idum=0
real(8) :: x,y,ran0,d,cont
n=1000
cont=0._8
x=ran0(idum)
y=ran0(idum)
d=sqrt((x-0.5)**2+(y-0.5)**2)
do i=1,n
 if ( d <= 0.5_8 ) cont=cont+1._8
enddo
cont=cont/n
write(*,*)'Area del círculo:',cont
end program area_circulo_MC
real(8) function ran0(idum)
implicit none
integer :: idum
integer,parameter :: IA=16807
integer,parameter :: IM=2147483647
integer,parameter :: IQ=127773
integer,parameter :: IR=2836
integer,parameter :: MASK=123459876
real(8),parameter :: AM=1._8/IM
integer :: k
idum=ieor(idum,MASK)
k=idum/IQ
idum=IA*(idum-k*IQ)-IR*k
if (idum.lt.0) idum=idum+IM
ran0=AM*idum
idum=ieor(idum,MASK)
end function ran0
 
FIN NOTAS DE CLASE
Tema: Evaluación de Integrales definidas
	Slide 1
	Slide 2
	Slide 3
	Slide 4
	Slide 5
	Slide 6
	Slide 7
	Slide 8
	Slide 9
	Slide 10
	Slide 11
	Slide 12
	Slide 13
	Slide 14
	Slide 15
	Slide 16
	Slide 17
	Slide 18
	Slide 19
	Slide 20
	Slide 21
	Slide 22
	Slide 23
	Slide 24
	Slide 25
	Slide 26
	Slide 27
	Slide 28
	Slide 29

Continue navegando