program integ

   use fvn_integ

   implicit none
   real(8), external :: f1,f2,g,h
   real(8) :: a,b,epsabs,epsrel,abserr,res
   integer :: key,ier
   a=0.
   b=1.
   epsabs=1d-8
   epsrel=1d-8
   key=2
   call fvn_integ_1_gk(f1,a,b,epsabs,epsrel,key,res,abserr,ier)
   write(*,*) "Integration of x*x between 0 and 1 : "
   write(*,*) res
   call fvn_integ_2_gk(f2,a,b,g,h,epsabs,epsrel,key,res,abserr,ier)
   write(*,*) "Integration of x*y between 0 and 1 on both x and y : "
   write(*,*) res
  
  end program
  function f1(x)
   implicit none
      real(8) :: x,f1
      f1=x*x
  end function
  function f2(x,y)
   implicit none
      real(8) :: x,y,f2
      f2=x*y
  end function
  function g(x)
   implicit none
      real(8) :: x,g
      g=0.
  end function
  function h(x)
   implicit none
      real(8) :: x,h
      h=1.
  end function