From ce2bbf7aac55aabd536849c2b04a270f24b8a62d Mon Sep 17 00:00:00 2001 From: Luca Date: Tue, 9 Apr 2019 16:30:15 +0200 Subject: [PATCH] zero arguments and zeta function --- gpl_module.f90 | 12 +++++++++++- test.f90 | 31 ++++++++++++++----------------- 2 files changed, 25 insertions(+), 18 deletions(-) diff --git a/gpl_module.f90 b/gpl_module.f90 index 5d71265..c5fd612 100644 --- a/gpl_module.f90 +++ b/gpl_module.f90 @@ -12,6 +12,15 @@ CONTAINS res = merge(1,n*factorial(n-1),n==1) END FUNCTION factorial + FUNCTION zeta(n) + real(kind=prec) :: values(9), zeta + integer :: n + values = (/1.6449340668482262, 1.2020569031595942, 1.0823232337111381, & + 1.03692775514337, 1.0173430619844488, 1.008349277381923, & + 1.0040773561979441, 1.0020083928260821, 1.000994575127818/) + zeta = values(n-1) + END FUNCTION zeta + FUNCTION GPL_has_convergent_series(m,z,y,k) ! tests if GPL has a convergent series representation integer :: m(:), k @@ -34,6 +43,7 @@ CONTAINS ! used to compute the value of GPL when all zi are zero integer :: l complex(kind=prec) :: y, GPL_zero_zi + GPL_zero_zi = 1.0d0/factorial(l) * log(y) ** l END FUNCTION GPL_zero_zi @@ -48,7 +58,7 @@ CONTAINS ! are all z_i = 0 ? if(k == 1 .and. z(1) == 0) then ! for that we assume that only one argument was passed, the rest through m1 - GPL = GPL_zero_zi(m(1)-1,y) + GPL = GPL_zero_zi(m(1),y) return end if diff --git a/test.f90 b/test.f90 index 8545c96..6a96ecb 100644 --- a/test.f90 +++ b/test.f90 @@ -6,21 +6,22 @@ PROGRAM TEST use mpl_module use gpl_module implicit none - + complex(kind=prec) :: res real, parameter :: tol = 1.0e-14 logical :: tests_successful = .true. + - call do_MPL_tests() + ! call do_MPL_tests() call do_GPL_tests() - + if(tests_successful) then print*, 'All tests passed. ' else print*, 'Some tests failed. ' stop 1 end if - + CONTAINS subroutine check(res, ref) @@ -40,22 +41,19 @@ CONTAINS integer :: m(:) complex(kind=prec) :: x(:), ref character(len=*) :: test_id - + print*, ' ', 'testing MPL ', test_id, ' ...' res = MPL(m,x) call check(res,ref) end subroutine test_one_MPL subroutine do_MPL_tests() - integer :: m2(2), m3(3) - integer :: bla - complex(kind=prec) :: x2(2), x3(3) - complex(kind=prec) :: res, ref + complex(kind=prec) :: ref print*, 'doing MPL tests...' ref = dcmplx(0.022696600480693277651633) call test_one_MPL((/ 1,1 /),cmplx((/ 0.3156498673740053, 0.3431255827785649 /)),ref, '1.1') - + ref = dcmplx(0.00023134615630308335448329926098409) call test_one_MPL((/ 1,1 /),cmplx((/ 0.03, 0.5012562893380046 /)),ref, '1.2') @@ -75,19 +73,18 @@ CONTAINS end subroutine test_one_GPL subroutine do_GPL_tests() - integer :: m2(2), m1(1), k - complex(kind=prec) :: z2(2), z1(1), y, res, ref + complex(kind=prec) :: ref print*, 'doing GPL tests...' - + ref = dcmplx(0.0819393734128676) call test_one_GPL((/ 1,1 /),cmplx((/ 1.3d0, 1.1d0 /)),cmplx(0.4),2,ref,'2.1') - + ref = dcmplx(0.01592795952537145) call test_one_GPL((/ 3,2 /),cmplx((/ 1.3d0, 1.1d0 /)),cmplx(0.4),2,ref,'2.2') - - ref = dcmplx(0.0173042341866201179) + + ref = dcmplx(0.0020332632172573974) call test_one_GPL((/ 4 /),cmplx((/ 0 /)),cmplx(1.6),1,ref,'2.3') - + end subroutine do_GPL_tests END PROGRAM TEST -- GitLab