Commit a528d95f authored by Luca's avatar Luca

shuffle algebra

parent ce2bbf7a
......@@ -51,10 +51,10 @@ CONTAINS
FUNCTION GPL(m,z,y,k)
! computes the generalized polylogarithm G_{m1,..mk} (z1,...zk; y)
! assumes zero arguments expressed through the m's
integer :: m(:), k, i
complex(kind=prec) :: z(:), x(k), y, GPL
! 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
......@@ -62,9 +62,14 @@ CONTAINS
return
end if
! do they have convergent series rep?
! need make convergent?
if(.not. GPL_has_convergent_series(m,z,y,k)) then
print*, ' ', 'does not have convergent series representation'
print*, 'need to make convergent'
if(k == 1 .and. m(1) == 1) then
print*, 'now we use the easy case. ha. ha.'
else
print*, ' ', 'does not have convergent series representation'
end if
end if
do i = 1,k
......@@ -76,19 +81,4 @@ CONTAINS
END MODULE gpl_module
! PROGRAM test
! ! used to test this module
! use gpl_module
! integer :: m(2) = (/ 1,1 /)
! complex(kind=prec) :: z(2) = dcmplx((/ 1.3d0, 1.1d0 /))
! complex(kind=prec) :: y = 0.4
! complex(kind=prec) :: res, ref
! res = GPL(m,z,y,2)
! ref = dcmplx(0.0819393734128676)
! print*, 'res=',res
! print*, 'ref=',ref
! END PROGRAM test
! This is currently a stand alone program which will merely be used as a
! guide for the implementation of the shuffle algebra for GPL functions
PROGRAM shuffle_algebra
implicit none
! Currently words can be no longer than the following values
! Might need to be adjusted
integer, parameter :: max_word_size = 6
integer, parameter :: max_word_sum_size = 1000
type word
character(len=max_word_size) :: letters
integer :: length
endtype word
type word_sum
type(word), dimension(max_word_sum_size) :: words
integer :: length
end type word_sum
type(word) :: v1 = word("abcd",4)
type(word) :: v2 = word("1234",4)
type(word) :: w1
type(word_sum) :: ws, ws1, ws2
ws = shuffle_product(v1,v2)
call print_word_sum(ws)
CONTAINS
RECURSIVE FUNCTION shuffle_product(v1, v2) result(res)
! takes two words and returns shuffle product as a word sum
type(word_sum) :: res, p1, p2, s1, s2
type(word) :: v1,v2,w1,w2
type(character) :: alpha, beta
! print*, '----------------------'
! print*, 'v1 = ', v1
! print*, 'v2 = ', v2
if(v1%length == 0) then
res = word_sum((/ v2 /),1)
else if(v2%length == 0) then
res = word_sum((/ v1 /),1)
else
alpha = v1%letters(1:1)
beta = v2%letters(1:1)
w1 = word(v1%letters(2:),v1%length-1)
w2 = word(v2%letters(2:),v2%length-1)
p1 = shuffle_product(w1,v2)
p2 = shuffle_product(v1,w2)
s1 = times(alpha, p1)
s2 = times(beta, p2)
res = combined_word_sum(s1,s2)
end if
END FUNCTION shuffle_product
FUNCTION combined_word_sum(s1, s2)
! combines word sums s1 and s2 into s
type(word_sum) :: s1, s2, combined_word_sum
type(word), dimension(s1%length + s2%length) :: combined_words
integer :: length
length = s1%length + s2%length
combined_words(1:s1%length) = s1%words(1:s1%length)
combined_words(s1%length+1:length) = s2%words(1:s2%length)
combined_word_sum = word_sum(combined_words,length)
END FUNCTION combined_word_sum
FUNCTION times(l, ws)
! computes word sum from letter times word sum, e.g. a(bc + cd) = abc + acd
character :: l
type(word_sum) :: ws
type(word_sum) :: times
integer :: i
times%length = ws%length
do i=1,ws%length
times%words(i)%letters = l // ws%words(i)%letters
end do
END FUNCTION times
SUBROUTINE print_word_sum(ws)
! prints a sum of word with plusses for easy readibility
integer :: i
type(word_sum) :: ws
do i = 1,ws%length
write(*, fmt="(1xai0)", advance="no") ws%words(i)%letters
if(i /= ws%length) then
write(*, fmt="(1xai0)", advance="no") " + "
end if
end do
END SUBROUTINE print_word_sum
END PROGRAM shuffle_algebra
......@@ -10,9 +10,8 @@ PROGRAM TEST
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
......@@ -58,8 +57,7 @@ CONTAINS
call test_one_MPL((/ 1,1 /),cmplx((/ 0.03, 0.5012562893380046 /)),ref, '1.2')
ref = dcmplx(0.000023446106415452030937059124671151)
call test_one_MPL((/ 2,1,2 /),cmplx((/ 0.03, 0.5012562893380046, 55.3832 /)),ref, '1.3')
call test_one_MPL((/ 2,1,2 /),cmplx((/ 0.03, 0.5012562893380046, 55.3832 /)),ref, '1.3')
end subroutine do_MPL_tests
subroutine test_one_GPL(m,z,y,k,ref,test_id)
......@@ -74,6 +72,7 @@ CONTAINS
subroutine do_GPL_tests()
complex(kind=prec) :: ref
complex(kind=prec), parameter :: epsilon = 1E-14
print*, 'doing GPL tests...'
ref = dcmplx(0.0819393734128676)
......@@ -84,7 +83,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 /),(/(0.7,epsilon)/),(1.1,1.4),1,ref,'2.4')
end subroutine do_GPL_tests
END PROGRAM TEST
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment