SHARE
TWEET

Untitled

a guest Dec 13th, 2012 18 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. #!/usr/local64.hg/bin/sage -python
  2.  
  3. # Computes the Plancherel integral of the trace of the r-fold symmetric representation of SU(3,CC)
  4.  
  5. import sys
  6. from sage.all import *
  7.  
  8. maxdim = 1
  9. maxprime = 10
  10.  
  11. # specify polynomials
  12. s = SymmetricFunctions(ZZ).schur()
  13. representations = [s([r+1]).expand(3) for r in range(maxdim)]
  14.  
  15. # specify primes p
  16. #primes = [2,3,5,7] # works
  17. primes = list(primes(2, maxprime)) # doesn't work, out of memory.
  18.  
  19. # function setup
  20. def e(x):
  21.     return exp(2*pi*i*x)
  22.  
  23. # gamma(s) * gamma(-s) for a prime p and an angle theta
  24. def gamma_sq(p,theta):
  25.     return ( (1-1/p*e(-theta)) * (1-1/p*e(theta)) ) \
  26.             / ( (1-e(-theta)) * (1-e(theta)) )
  27.  
  28. # |c(s)|^2 for all roots of sl3
  29. def c_sq(p,a,b):
  30.     return gamma_sq(p,a) * gamma_sq(p,b) * gamma_sq(p,-a-b)
  31.  
  32. # integrand = tr(repr'n) * factors for d mu(omega_s) -> ds
  33. #           = tr(repr'n) * factors for ds->dtheta1,2 / |c(theta1,theta2)|^2
  34. def integrand(p,poly,theta1,theta2):
  35.     return  poly(e(theta1),e(theta2),e(-theta1-theta2))\
  36.             * e(theta1+theta2)\
  37.             / c_sq(p,theta1,theta2)
  38.  
  39. def calc_integrals(reps,primes):
  40.     results = []
  41.     theta1,theta2 = var('theta1,theta2')
  42.     for rep in reps:
  43.         print "Representation Polynomial: %s" % rep
  44.         result_by_prime = []
  45.         for prime in primes:
  46.             print "Prime number %d: " % prime
  47.             integral_by_theta1 = -4 *(pi**2) * integral(integrand(prime,rep,theta1,theta2),theta1,0,1)
  48.             print "First integral: %s " % integral_by_theta1
  49.             integral_by_theta2 = integral(integral_by_theta1, theta2,0,1)
  50.             print "Second Integral: %s " % integral_by_theta2
  51.             result_by_prime.append(integral_by_theta2)
  52.         results.append(result_by_prime)
  53.     return results
  54.  
  55. # annnd..... integrate!
  56. results = calc_integrals(representations,primes)
  57. o = open('calc_res_2.txt','w')
  58. o.write(str(results))
  59. o.close()
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
Top