title: Why is pow(a, d, n) so much faster than a**d % n?
body: <p>I was trying to implement a <a href="http://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test">Miller-Rabin primality test</a>, and was puzzled why it was taking so long (> 20 seconds) for midsize numbers (~7 digits). I eventually found the following line of code to be the source of the problem:</p>
  
  <pre><code>x = a**d % n
  </code></pre>
  
  <p>(where <code>a</code>, <code>d</code>, and <code>n</code> are all similar, but unequal, midsize numbers, <code>**</code> is the exponentiation operator, and <code>%</code> is the modulo operator)</p>
  
  <p>I then I tried replacing it with the following:</p>
  
  <pre><code>x = pow(a, d, n)
  </code></pre>
  
  <p>and it by comparison it is almost instantaneous.</p>
  
  <p>For context, here is the original function:</p>
  
  <pre><code>from random import randint
  
  def primalityTest(n, k):
      if n &lt; 2:
          return False
      if n % 2 == 0:
          return False
      s = 0
      d = n - 1
      while d % 2 == 0:
          s += 1
          d &gt;&gt;= 1
      for i in range(k):
          rand = randint(2, n - 2)
          x = rand**d % n         # offending line
          if x == 1 or x == n - 1:
              continue
          for r in range(s):
              toReturn = True
              x = pow(x, 2, n)
              if x == 1:
                  return False
              if x == n - 1:
                  toReturn = False
                  break
          if toReturn:
              return False
      return True
  
  print(primalityTest(2700643,1))
  </code></pre>
  
  <p>An example timed calculation:</p>
  
  <pre><code>from timeit import timeit
  
  a = 2505626
  d = 1520321
  n = 2700643
  
  def testA():
      print(a**d % n)
  
  def testB():
      print(pow(a, d, n))
  
  print("time: %(time)fs" % {"time":timeit("testA()", setup="from __main__ import testA", number=1)})
  print("time: %(time)fs" % {"time":timeit("testB()", setup="from __main__ import testB", number=1)})
  </code></pre>
  
  <p>Output (run with PyPy 1.9.0):</p>
  
  <pre><code>2642565
  time: 23.785543s
  2642565
  time: 0.000030s
  </code></pre>
  
  <p>Output (run with Python 3.3.0, 2.7.2 returns very similar times):</p>
  
  <pre><code>2642565
  time: 14.426975s
  2642565
  time: 0.000021s
  </code></pre>
  
  <p>And a related question, why is this calculation almost twice as fast when run with Python 2 or 3 than with PyPy, when usually PyPy is <a href="http://speed.pypy.org/">much faster</a>?</p>