title: Why is this iterative Collatz method 30% slower than its recursive version in Python?
body: <h3>Prelude</h3>
  
  <p>I have two implementations for a particular problem, one recursive and one iterative, and I want to know what causes the iterative solution to be ~30% slower than the recursive one.</p>
  
  <p>Given the recursive solution, I write an iterative solution making the stack explicit.<br>
  Clearly, I simply mimic what the recursion is doing, so of course the Python engine is better optimized to handle the bookkeeping. But can we write an iterative method with similar performance?</p>
  
  <p>My case study is <a href="http://projecteuler.net/problem=14" rel="nofollow">Problem #14</a> on Project Euler.  </p>
  
  <blockquote>
    <p>Find the longest Collatz chain with a starting number below one million.</p>
  </blockquote>
  
  <h3>Code</h3>
  
  <p>Here is a parsimonious recursive solution (credit due to <em>veritas</em> in the problem thread plus an optimization from <em>jJjjJ</em>):</p>
  
  <pre><code>def solve_PE14_recursive(ub=10**6):
      def collatz_r(n):
          if not n in table:
              if n % 2 == 0:
                  table[n] = collatz_r(n // 2) + 1
              elif n % 4 == 3:
                  table[n] = collatz_r((3 * n + 1) // 2) + 2
              else:
                  table[n] = collatz_r((3 * n + 1) // 4) + 3
          return table[n]
      table = {1: 1}
      return max(xrange(ub // 2 + 1, ub, 2), key=collatz_r)
  </code></pre>
  
  <p>Here's my iterative version:</p>
  
  <pre><code>def solve_PE14_iterative(ub=10**6):
      def collatz_i(n):
          stack = []
          while not n in table:
              if n % 2 == 0:
                  x, y = n // 2, 1
              elif n % 4 == 3:
                  x, y = (3 * n + 1) // 2, 2
              else:
                  x, y = (3 * n + 1) // 4, 3
              stack.append((n, y))
              n = x
          ysum = table[n]
          for x, y in reversed(stack):
              ysum += y
              table[x] = ysum
          return ysum
      table = {1: 1}
      return max(xrange(ub // 2 + 1, ub, 2), key=collatz_i)
  </code></pre>
  
  <p>And the timings on my machine (i7 machine with lots of memory) using IPython:</p>
  
  <pre><code>In [3]: %timeit solve_PE14_recursive()
  1 loops, best of 3: 942 ms per loop
  In [4]: %timeit solve_PE14_iterative()
  1 loops, best of 3: 1.35 s per loop
  </code></pre>
  
  <hr>
  
  <h3>Comments</h3>
  
  <p>The recursive solution is awesome:</p>
  
  <ul>
  <li>Optimized to skip a step or two depending on the two least significant bits.<br>
  My original solution didn't skip any Collatz steps and took ~1.86 s</li>
  <li>It is difficult to hit Python's default recursion limit of 1000.<br>
  <code>collatz_r(</code><a href="http://en.wikipedia.org/wiki/Collatz_conjecture#cite_ref-11" rel="nofollow"><code>9780657630</code></a><code>)</code> returns 1133 but requires less than 1000 recursive calls.</li>
  <li>Memoization avoids retracing</li>
  <li><code>collatz_r</code> length calculated on-demand for <code>max</code></li>
  </ul>
  
  <p>Playing around with it, timings seem to be precise to +/- 5 ms.<br>
  Languages with static typing like C and Haskell can get timings below 100 ms.<br>
  I put the initialization of the memoization <code>table</code> in the method by design for this question, so that timings would reflect the "re-discovery" of the table values on each invocation.</p>
  
  <p><code>collatz_r(2**1002)</code> raises <code>RuntimeError: maximum recursion depth exceeded</code>.<br>
  <code>collatz_i(2**1002)</code> happily returns with <code>1003</code>.</p>
  
  <p>I am familiar with generators, coroutines, and decorators.<br>
  I am using Python 2.7. I am also happy to use Numpy (1.8 on my machine).</p>
  
  <h3>What I am looking for</h3>
  
  <ul>
  <li>an iterative solution that closes the performance gap</li>
  <li>discussion on how Python handles recursion</li>
  <li>the finer details of the performance penalties associated with an explicit stack</li>
  </ul>
  
  <p>I'm looking mostly for the first, though the second and third are very important to this problem and would increase my understanding of Python.</p>