program test_ze1
use fvn_fnlib
implicit none

complex(8) :: z,resz
complex(8), dimension(6) :: abramowitz_x,abramowitz_y
integer(4) :: i

abramowitz_x(1)=dcmplx(5.,5.)
abramowitz_y(1)=dcmplx(0.906058d0,0.070209d0)
abramowitz_x(2)=dcmplx(-5.,5.)
abramowitz_y(2)=dcmplx(1.079407d0,0.143879d0)
abramowitz_x(3)=dcmplx(11.,5.)
abramowitz_y(3)=dcmplx(0.932672d0,0.026361d0)
abramowitz_x(4)=dcmplx(-11.,5.)
abramowitz_y(4)=dcmplx(1.084526d0,0.049336d0)
abramowitz_x(5)=dcmplx(15.,18.)
abramowitz_y(5)=dcmplx(0.972359d0,0.029448d0)
abramowitz_x(6)=dcmplx(-15.,18.)
abramowitz_y(6)=dcmplx(1.026183d0,0.036552d0)

! Value of z*exp(z)*e1(z) from Abramowitz & Stegun
! to cover the different cases of the algo
!
! Case (5,5) modulus <10 with positive real part
! Case (-5,5) modulus <10 with negative real part
! Case (11,5) modulus between 10 and 20 with positive real part
! Case (-11,5) modulus between 10 and 20 with negative real part
! Case (15,18) modulus > 20 with positive real part
! Case (-15,18) modulus > 20 with negative real part

do i=1,6
  z=abramowitz_x(i)
  resz=z*exp(z)*e1(z)
  write(*,*) "Tabulated : ",abramowitz_y(i)
  write(*,*) "Calculated :",resz
  write(*,*)
end do
end program