attachment:poly_multiply_benchmark.sage of DavidHarvey/TimingExperiments

   1 """ This is a program for comparing times of NTL, MAGMA, PARI multiplication
   2 of polynomials over Z, over various degrees and coefficient sizes.
   3 
   4 AUTHOR: David Harvey (2006-08-21)
   5 """
   6 
   7 ## todo: docstring
   8 
   9 from sys import stdout
  10 
  11 # which systems to test:
  12 systems = ["magma", "pari", "ntl"]
  13 
  14 # initialise various subsystems
  15 if "magma" in systems:
  16     M = Magma()
  17     M.set("R<x>", "PolynomialRing(IntegerRing())")
  18 
  19 if "pari" in systems:
  20     pari.allocatemem()
  21     pari.allocatemem()
  22     pari.allocatemem()
  23     pari.allocatemem()
  24     pari.allocatemem()
  25     pari.allocatemem()
  26 
  27 
  28 # polynomials over Z
  29 
  30 
  31 def poly_to_string(coeffs):
  32     # returns string for polynomial with given coefficients
  33     return "+".join(["%s*x^%s" % (str(coeffs[i]), i) for i in range(len(coeffs))])
  34 
  35 
  36 def time_ntl(coeffs1, coeffs2, poly1_string, poly2_string, num_trials):
  37     poly1 = ntl.ZZX(coeffs1)
  38     poly2 = ntl.ZZX(coeffs2)
  39 
  40     t = cputime()
  41     for trial in range(num_trials):
  42         poly = poly1 * poly2
  43     return cputime(t)
  44 
  45 
  46 def time_magma(coeffs1, coeffs2, poly1_string, poly2_string, num_trials):
  47     M.set("poly1", poly1_string)
  48     M.set("poly2", poly2_string)
  49 
  50     t = M.cputime()
  51     for trial in range(num_trials):
  52         M.set("result", "poly1 * poly2")
  53     return M.cputime(t)
  54 
  55 
  56 def time_pari(coeffs1, coeffs2, poly1_string, poly2_string, num_trials):
  57     poly1 = pari(poly1_string)
  58     poly2 = pari(poly2_string)
  59     t = cputime()
  60     for trial in range(num_trials):
  61         poly3 = poly1 * poly2
  62     return cputime(t)
  63 
  64 
  65 
  66 # make a list of degrees to test
  67 degree_list = [10]
  68 while degree_list[-1] < 10000:
  69     degree_list.append(ceil(degree_list[-1] * 1.3))
  70 
  71 # make a list of coefficient sizes to test
  72 coeff_bits_list = [4]
  73 while coeff_bits_list[-1] < 10000:
  74     coeff_bits_list.append(ceil(coeff_bits_list[-1] * 1.3))
  75 
  76 for degree in degree_list:
  77     for coeff_bits in coeff_bits_list:
  78 
  79         # generate random polynomials to multiply
  80         max_coeff = 2**coeff_bits
  81         coeffs1 = [ZZ.random(max_coeff) for i in range(degree)]
  82         coeffs2 = [ZZ.random(max_coeff) for i in range(degree)]
  83 
  84         poly1_string = poly_to_string(coeffs1)
  85         poly2_string = poly_to_string(coeffs2)
  86 
  87         for system in systems:
  88             timing_function = eval("time_" + system)
  89 
  90             # determine an appropriate number of trials for each timing attempt
  91 
  92             num_trials = 0
  93             t = 0.0
  94             while t < 0.3:
  95                 if num_trials == 0:
  96                     num_trials = 1
  97                 else:
  98                     num_trials *= 2
  99 
 100                 t = timing_function(coeffs1, coeffs2, poly1_string, poly2_string, num_trials)
 101 
 102             # now take 5 samples with that number of trials each
 103 
 104             samples = [t / num_trials]
 105             for i in range(4):
 106                 t = timing_function(coeffs1, coeffs2, poly1_string, poly2_string, num_trials)
 107                 samples.append(t / num_trials)
 108 
 109             # get some stats on the sample timing
 110 
 111             average = sum(samples) / len(samples)
 112             value_range = max(samples) - min(samples)
 113 
 114             print "%d %d %s %.5f %.5f %d" % (degree, coeff_bits, system, average, value_range, num_trials)
 115             stdout.flush()