test_ze1.f90 1.3 KB
program test_ze1
use fvn_fnlib
implicit none

complex(kind=dp_kind) :: z,resz
complex(kind=dp_kind), dimension(6) :: abramowitz_x,abramowitz_y
integer(kind=ip_kind) :: i

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

! 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