! Metodos Numericos y Simulaciones en Astrofisica
! FCEFN - UNSJ, 2010
!
! ejemplo 5b
!
! Metodo de Euler para problemas de valor inicial
!
! Calculo de la Ec. diferencial ordinaria:
!
! dy/dt = 2 t + 1 para y(0) = 1
!
! cuya solucion exacta es:
!
! y = t^2 + t + 1
!
!
program eje5b
implicit none
  integer, parameter ::p=selected_real_kind(16,200)
  real(p) :: hh,time,yy,yytrue
  
! ingrese el paso
!
  write(*,*)'Ingrese h: '
  read(*,*)hh
  
! valor inicial para y y t
!
  yy=1._p
  yytrue=1._P
  time=0._p
  open(7,file='eje5b.out')
  write(*,'(f8.3,2f23.18,g17.7)')time,yy,yytrue,0.
  write(7,*)time,yy,yytrue,0.
  
! inicia bucle de calculo
!
  do while(time < 10._p)
    time=time+hh
	
! calcula valor real de y
!
	yytrue=time**2+time+1._p

! calcula aproximacion por euler
!
   yy=yy+hh*(2._p*time+1._p)	
   
   write(*,'(f8.3,2f23.18,g17.7)')time,yy,yytrue,abs(yytrue-yy)
   write(7,'(f8.3,2f23.18,g17.7)')time,yy,yytrue,abs(yytrue-yy)
  enddo
  close(7)
  
end program eje5b
