title: How to accelerate slow matrix multiplication in SymPy?
body: <p>I was writing a tool to solve specific recurrence equations with SymPy, and found that one of the steps that involved a matrix multiplication was taking an extraordinarily long time. For instance, if I try the following in the iPython console,</p>
  
  <pre><code>In [1]: from sympy import *
  
  In [2]: A = Matrix(500, 500, lambda i,j: 2 + abs(i-j) if i-j in [-1, 0, 1] else 0)
  
  In [3]: A[0:7, 0:7]
  Out[3]: 
  Matrix([
  [2, 3, 0, 0, 0, 0, 0],
  [3, 2, 3, 0, 0, 0, 0],
  [0, 3, 2, 3, 0, 0, 0],
  [0, 0, 3, 2, 3, 0, 0],
  [0, 0, 0, 3, 2, 3, 0],
  [0, 0, 0, 0, 3, 2, 3],
  [0, 0, 0, 0, 0, 3, 2]])
  
  In [4]: A * A
  ...
  </code></pre>
  
  <p>I've never waited long enough for it to finish.</p>
  
  <p>I appreciate that symbolic manipulation considerably slower than numerical evaluation, but this seems pretty ridiculous. Using Octave or other linear algebra packages can perform this calculation in a fraction of a second.</p>
  
  <p>Does anyone have experience using SymPy's matrix class for rows &amp; columns ~ 1000, with Rational entries? </p>