diff --git a/gpl_module.f90 b/gpl_module.f90 index 28855ec1d53f61e30e6f2e27dfef89c16c150896..e1eaf4f7e4b899e13f8e6f0c264a0cd715c1bdb0 100644 --- a/gpl_module.f90 +++ b/gpl_module.f90 @@ -44,12 +44,13 @@ CONTAINS ! used to compute the value of GPL when all zi are zero integer :: l complex(kind=prec) :: y, GPL_zero_zi - print*, 'computed value using zero' + print*, 'computed value using zi = 0' GPL_zero_zi = 1.0d0/factorial(l) * log(y) ** l END FUNCTION GPL_zero_zi FUNCTION G_with_flat_args(z_flat,y) result(res) + ! Calls G function with flat arguments, that is, zeroes not passed through the m's. complex(kind=prec) :: z_flat(:), y, res complex(kind=prec), allocatable :: z(:) integer :: m_prime(size(z_flat)), condensed_size @@ -73,29 +74,31 @@ CONTAINS ! computes the generalized polylogarithm G_{m1,..mk} (z1,...zk; y) ! assumes zero arguments expressed through the m's - integer :: m(:), k, i + integer :: m(:), k, i, kminusj complex(kind=prec) :: z(:), x(k), y, res, c(sum(m)+1,sum(m)+1), z_flat(sum(m)), a(sum(m)-1) ! print*, 'z = ', abs(get_flattened_z(m,z)) ! are all z_i = 0 ? - if(k == 1 .and. abs(z(1)) < zero) then + if(k == 1 .and. abs(z(1)) == 0) then ! assumes that the zeros at the beginning are passed through m1 res = GPL_zero_zi(m(1),y) return end if ! need to remove trailing zeros? - if(abs(z(k)) < zero ) then - print*, 'need to remove trailing zeros' + if(abs(z(k)) == 0 ) then + print*, 'need to remove trailing zeros' ! which we do in flat form ! flatten z z_flat = get_flattened_z(m,z) - a = z_flat(1: (size(z_flat)-1)) - c = shuffle_with_zero(a) - res = G_with_flat_args(a,y)*log(y) - do i = 2,k - res = res - G_with_flat_args(c(i,:),y) - end do - return + res = G_with_flat_args(z_flat,y) + + ! a = z_flat(1: (size(z_flat)-1)) + ! c = shuffle_with_zero(a) + ! res = G_with_flat_args(a,y)*log(y) + ! do i = 2,k + ! res = res - G_with_flat_args(c(i,:),y) + ! end do + ! return end if ! need make convergent? diff --git a/test.f90 b/test.f90 index 2f4a1570be0fd4fe658fabca61efbb7550a03e49..fdab131d5163071eaaf26b8c5e61cfe365ab6bbd 100644 --- a/test.f90 +++ b/test.f90 @@ -85,11 +85,9 @@ CONTAINS ref = dcmplx(0.0020332632172573974) call test_one_GPL((/ 4 /),cmplx((/ 0 /)),cmplx(1.6),1,ref,'2.3') - ! ref = dcmplx(0.0020332632172573974) - ! call test_one_GPL((/ 1,1 /),cmplx((/ 1.7,0.0 /)),cmplx(1.1),2,ref,'2.4') + ref = dcmplx(0.0020332632172573974) + call test_one_GPL((/ 1,1 /),cmplx((/ 1.7,0.0 /)),cmplx(1.1),2,ref,'2.4') - ! ref = dcmplx(0.0020332632172573974) - ! call test_one_GPL((/ 1 /),(/(0.7,epsilon)/),(1.1,1.4),1,ref,'2.4') end subroutine do_GPL_tests END PROGRAM TEST diff --git a/utils.f90 b/utils.f90 index 9f05acd7630ce46a2b327f0f4414213930706269..6d929e73f3416ca0105bf61177e1beaf3e7e3c61 100644 --- a/utils.f90 +++ b/utils.f90 @@ -12,7 +12,6 @@ MODULE utils implicit none integer, parameter :: prec = selected_real_kind(15,32) - real, parameter :: zero = 1.0e-50 ! logical :: print_enabled = .true. ! logical :: warnings_enabled = .true. @@ -61,6 +60,19 @@ CONTAINS end do END FUNCTION get_flattened_z + FUNCTION find_amount_trailing_zeros(z) result(res) + complex(kind=prec) :: z(:) + integer :: res, i + res = 0 + do i = size(z), 1, -1 + if( z(i) == 0 ) then + res = res + 1 + else + exit + end if + end do + END FUNCTION find_amount_trailing_zeros + FUNCTION find_first_zero(v) result(res) ! returns index of first zero, or -1 if there is no zero integer :: v(:), i, res @@ -119,33 +131,34 @@ END MODULE utils ! use utils ! implicit none -! complex(kind=prec), dimension(3) :: a = cmplx((/1,2,3/)) -! complex(kind=prec) :: z_flat(2) -! complex(kind=prec), allocatable :: z(:) -! integer :: m_prime(2), condensed_size -! integer, allocatable :: m(:) -! complex(kind=prec) :: b(size(a)+1,size(a)+1) - -! ! ! test shuffling -! ! b = 1 -! ! b = shuffle_with_zero(a) -! ! call print_as_matrix(b) - -! ! test condensing -! z_flat = cmplx((/4,0/)) -! m_prime = get_condensed_m(z_flat) -! if(find_first_zero(m_prime) == -1) then -! condensed_size = size(m_prime) -! else -! condensed_size = find_first_zero(m_prime)-1 -! end if -! allocate(m(condensed_size)) -! allocate(z(condensed_size)) -! m = m_prime(1:condensed_size) -! z = get_condensed_z(m,z_flat) -! z_flat = get_flattened_z(m,z) -! deallocate(m) -! deallocate(z) + +! ! complex(kind=prec), dimension(5) :: a = cmplx((/1,2,3/)) +! ! complex(kind=prec) :: z_flat(2) +! ! complex(kind=prec), allocatable :: z(:) +! ! integer :: m_prime(2), condensed_size +! ! integer, allocatable :: m(:) +! ! complex(kind=prec) :: b(size(a)+1,size(a)+1) + +! ! ! ! test shuffling +! ! ! b = 1 +! ! ! b = shuffle_with_zero(a) +! ! ! call print_as_matrix(b) + +! ! ! test condensing +! ! z_flat = cmplx((/4,0/)) +! ! m_prime = get_condensed_m(z_flat) +! ! if(find_first_zero(m_prime) == -1) then +! ! condensed_size = size(m_prime) +! ! else +! ! condensed_size = find_first_zero(m_prime)-1 +! ! end if +! ! allocate(m(condensed_size)) +! ! allocate(z(condensed_size)) +! ! m = m_prime(1:condensed_size) +! ! z = get_condensed_z(m,z_flat) +! ! z_flat = get_flattened_z(m,z) +! ! deallocate(m) +! ! deallocate(z) ! END PROGRAM test