title: N-D version of itertools.combinations in numpy
body: <p>I would like to implement <a href="http://docs.python.org/2/library/itertools.html#itertools.combinations" rel="noreferrer">itertools.combinations</a> for numpy.  Based on <a href="http://numpy-discussion.10968.n7.nabble.com/itertools-combinations-to-numpy-td16635.html" rel="noreferrer">this discussion</a>, I have a function that works for 1D input:</p>
  
  <pre><code>def combs(a, r):
      """
      Return successive r-length combinations of elements in the array a.
      Should produce the same output as array(list(combinations(a, r))), but 
      faster.
      """
      a = asarray(a)
      dt = dtype([('', a.dtype)]*r)
      b = fromiter(combinations(a, r), dt)
      return b.view(a.dtype).reshape(-1, r)
  </code></pre>
  
  <p>and the output makes sense:</p>
  
  <pre><code>In [1]: list(combinations([1,2,3], 2))
  Out[1]: [(1, 2), (1, 3), (2, 3)]
  
  In [2]: array(list(combinations([1,2,3], 2)))
  Out[2]: 
  array([[1, 2],
         [1, 3],
         [2, 3]])
  
  In [3]: combs([1,2,3], 2)
  Out[3]: 
  array([[1, 2],
         [1, 3],
         [2, 3]])
  </code></pre>
  
  <p>However, it would be best if I could expand it to N-D inputs, where additional dimensions simply allow you to speedily do multiple calls at once.  So, conceptually, if <code>combs([1, 2, 3], 2)</code> produces <code>[1, 2], [1, 3], [2, 3]</code>, and  <code>combs([4, 5, 6], 2)</code> produces <code>[4, 5], [4, 6], [5, 6]</code>, then <code>combs((1,2,3) and (4,5,6), 2)</code> should produce <code>[1, 2], [1, 3], [2, 3] and [4, 5], [4, 6], [5, 6]</code> where "and" just represents parallel rows or columns (whichever makes sense).  (and likewise for additional dimensions)</p>
  
  <p>I'm not sure:</p>
  
  <ol>
  <li>How to make the dimensions work in a logical way that's consistent with the way other functions work (like how some numpy functions have an <code>axis=</code> parameter, and a default of axis 0.  So probably axis 0 should be the one I am combining along, and all other axes just represent parallel calculations?)</li>
  <li>How to get the above code to work with ND (right now I get <code>ValueError: setting an array element with a sequence.</code>)</li>
  <li>Is there a better way to do <code>dt = dtype([('', a.dtype)]*r)</code>?</li>
  </ol>