Peter Bingham Peter Bingham - 4 months ago 20
Python Question

Sympy never returns when processing complicated determinants

At least this it how it appears. The following code behaved correctly for 3x3 and 6x6 matrices.

deter = mat.det('method':'berkowitz')
#self.resultFileHandler.writeLogStr(str(deter))
a = sy_polys.Poly(deter, k)


For 3x3 it takes ~0.8s to execute this code, for 6x6 it takes ~288s (with only 650ms for the det function, the rest for the Poly). For 10x10, either the complexity has ramped at a colossal rate or some other reason is preventing it returning from the Poly call (I waited a week). No exceptions are thrown.

The elements of the determinants consist of large symbolic polynomials.

I was on 0.7.1 and just upgraded to 1.0 (problem in both versions).

I added the logging to try and get the determinant to file but it sticks again in the str(deter) function call. If I break my debugger can't display the deter (prob too large for the debugger).

Here is a stack:

MainThread - pid_135368_id_42197520
_print [printer.py:262]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
parenthesize [str.py:29]
_print_Mul [str.py:290]
_print [printer.py:257]
_print_Add [str.py:56]
_print [printer.py:257]
doprint [printer.py:233]
sstr [str.py:748]
__str__ [basic.py:396]
_getRoots_sympy_Poly_nroots [__init__.py:91]
getRoots [__init__.py:68]
findPolyRoots [__init__.py:697]
_getNroots [polefinder.py:97]
_doForN [polefinder.py:60]
_incN [polefinder.py:52]
__init__ [polefinder.py:39]
_doPoleFind [polefinderwrap.py:33]
_polesForPos [polefinderwrap.py:47]
<module> [polefinderwrap.py:60]
run [pydevd.py:937]
<module> [pydevd.py:1530]


OK, I've got an exception from the str function. Seems likely that the polynomial has become too large.

Traceback (most recent call last):
File "E:\Peter's Documents\PhD\Code\Git\ProtoQScat\multichannel\qscat\ratsmat.\polefinder.py", line 60, in _doForN
roots = self._getNroots(N)
File "E:\Peter's Documents\PhD\Code\Git\ProtoQScat\multichannel\qscat\ratsmat.\polefinder.py", line 97, in _getNroots
roots = ratSmat.findPolyRoots(False)
File "E:\Peter's Documents\PhD\Code\Git\ProtoQScat\multichannel\qscat\numerical/..\ratsmat\__init__.py", line 697, in findPolyRoots
roots = self.polyRootSolve.getRoots(mat, k)
File "E:\Peter's Documents\PhD\Code\Git\ProtoQScat\multichannel\qscat\numerical/..\ratsmat\__init__.py", line 68, in getRoots
ret = self._getRoots_sympy_Poly_nroots(mat, k)
File "E:\Peter's Documents\PhD\Code\Git\ProtoQScat\multichannel\qscat\numerical/..\ratsmat\__init__.py", line 91, in _getRoots_sympy_Poly_nroots
self.resultFileHandler.writeLogStr(str(deter))
File "C:\Python27\lib\site-packages\sympy\core\basic.py", line 396, in __str__
return sstr(self, order=None)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 748, in sstr
s = p.doprint(expr)
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 233, in doprint
return self._str(self._print(expr))
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 257, in _print
return getattr(self, printmethod)(expr, *args, **kwargs)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 56, in _print_Add
t = self._print(term)
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 257, in _print
return getattr(self, printmethod)(expr, *args, **kwargs)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 290, in _print_Mul
a_str = [self.parenthesize(x, prec) for x in a]
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 29, in parenthesize
return "(%s)" % self._print(item)
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 257, in _print
return getattr(self, printmethod)(expr, *args, **kwargs)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 56, in _print_Add
t = self._print(term)
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 257, in _print
return getattr(self, printmethod)(expr, *args, **kwargs)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 290, in _print_Mul
a_str = [self.parenthesize(x, prec) for x in a]
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 29, in parenthesize
return "(%s)" % self._print(item)
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 257, in _print
return getattr(self, printmethod)(expr, *args, **kwargs)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 56, in _print_Add
t = self._print(term)
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 257, in _print
return getattr(self, printmethod)(expr, *args, **kwargs)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 290, in _print_Mul
a_str = [self.parenthesize(x, prec) for x in a]
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 29, in parenthesize
return "(%s)" % self._print(item)
File "C:\Python27\lib\site-packages\sympy\printing\printer.py", line 257, in _print
return getattr(self, printmethod)(expr, *args, **kwargs)
File "C:\Python27\lib\site-packages\sympy\printing\str.py", line 69, in _print_Add
return sign + ' '.join(l)
MemoryError


EDIT:
Following from answer below here is a profile plot with the determinant size (channel). Ignore N (on y axis) it is another parameter of the calculation (governs the size of the polys in the elements).
enter image description here

Answer

OK, so I returned to this after reading literature and feeling (emphasis here) that the Berkowitz should perform between O(n^2) and O(n^3).

What I found was that the input to the det and Poly makes a huge difference to the performance (I admit that this was not alluded to in my question). Wrapping the original expression in a poly drastically improves performance.

Consider the following three codes

1: using MatrixSymbol. det takes 1.1s then still stuck at str after 20min

from sympy import MatrixSymbol, Matrix

X = MatrixSymbol('X', 10, 10)
Xmat = Matrix(X)

deter = Xmat.det(method='berkowitz')
str(deter)

2: Represents my original problem code. det takes 1.8s then still stuck at Poly after 20min

import sympy
from sympy import Matrix, I
from sympy import poly
from sympy.polys import Poly

matSz = 10

m = Matrix([[0.0]*matSz]*matSz)
x = sympy.symbols('x')
for i in range(matSz):
  for j in range(matSz):
    m[i,j] = 2.0*float(i+1)*x**2 + 2.0*float(j+1)*x - 5.0*float(i+1)

deter = m.det(method='berkowitz')
deter_poly = Poly(deter, x)  #Required or exception
roots = deter_poly.nroots()

3: Above but with m[i,j] = poly(. det takes 3.0s, Poly 0.04 and nroots 0.27

import sympy
from sympy import Matrix, I
from sympy import poly
from sympy.polys import Poly

matSz = 10

m = Matrix([[0.0]*matSz]*matSz)

x = sympy.symbols('x')
for i in range(matSz):
  for j in range(matSz):
    m[i,j] = poly(2.0*float(i+1)*x**2 + 2.0*float(j+1)*x*I - 5.0*float(i+1))

deter = m.det(method='berkowitz')
deter_poly = Poly(deter, x)  #Required or exception
roots = deter_poly.nroots()