title: Make numpy.sum() return a sum of matrices instead of a single number
body: <p>I am doing a fairly complicated summation using a matrix with numpy. 
  The shape of the matrix is <code>matrix.shape = (500, 500)</code> and the shape of the array is <code>arr.shape = (25,)</code>. The operation is as follows:</p>
  
  <pre><code>totalsum = np.sum([i * matrix for i in arr])
  </code></pre>
  
  <p>Here is what I don't understand:</p>
  
  <p><code>np.sum()</code> is very slow and returns a single float, <code>float64</code>. Doing the same operation with Python's <code>sum.()</code>, i.e. </p>
  
  <pre><code>totalsum2 = sum([i*matrix for i in arr])
  </code></pre>
  
  <p>Preserves the shape of the matrix. That is, the resulting shape is <code>totalsum2.shape() = (500, 500)</code>. Huh?</p>
  
  <p>I also think it is strange that <code>np.sum()</code> takes longer than <code>sum()</code>, particularly when we are working with numpy ndarrays. </p>
  
  <p>What exactly is going on here? 
  How is <code>np.sum()</code> summing the above values in comparison to <code>sum()</code>?</p>
  
  <p>I would like <code>np.sum()</code> to preserve the matrix shape. How can I set the dimensions such that <code>np.sum()</code> preserves the matrix size and does not return a single float? </p>