9

If I ask SymPy to row-reduce the singular matrix

nu = Symbol('nu')
lamb = Symbol('lambda')
A3 = Matrix([[-3*nu, 1, 0, 0],
             [3*nu, -2*nu-1, 2, 0],
             [0, 2*nu, (-1 * nu) - lamb - 2, 3],
             [0, 0, nu + lamb, -3]])
print A3.rref()

then it returns the identity matrix

(Matrix([
[1, 0, 0, 0],
[0, 1, 0, 0],
[0, 0, 1, 0],
[0, 0, 0, 1]]), [0, 1, 2, 3])

which it shouldn't do, since the matrix is singular. Why is SymPy giving me the wrong answer and how can I get it to give me the right answer?

I know SymPy knows the matrix is singular, because when I ask for A3.inv(), it gives

raise ValueError("Matrix det == 0; not invertible.")

Furthermore, when I remove lamb from the matrix (equivalent to setting lamb = 0), SymPy gives the correct answer:

(Matrix([
[1, 0, 0, -1/nu**3],
[0, 1, 0, -3/nu**2],
[0, 0, 1,    -3/nu],
[0, 0, 0,        0]]), [0, 1, 2])

which leads me to believe that this problem only happens with more than one variable.

EDIT: Interestingly, I just got the correct answer when I pass rref() the argument "simplify=True". I still have no idea why that is though.

Jessica
  • 2,335
  • 2
  • 23
  • 36

1 Answers1

3

The rref algorithm fundamentally requires the ability to tell if the elements of the matrix are identically zero. In SymPy, the simplify=True option instructs SymPy to simplify the entries first at the relevant stage of the algorithm. With symbolic entries, this is necessary, as you can easily have symbolic expressions that are identically zero but which don't simplify to such automatically, like x*(x - 1) - x**2 + x. The option is off by default because in general such simplification can be expensive, through this can be controlled by passing in a less general simplify function than simplify (for rational functions, use cancel). The defaults here could probably be smarter.

asmeurer
  • 86,894
  • 26
  • 169
  • 240