Performance improvements

Performance improvements

1 - The method for computing the root of the indicial polynomial now is cached.
This means that operators with the same indicial polynomial take less time to
compute their roots.
2 - Computing th roots of the indicial polynomial when there are parameters is done
using a modular approximation: we evaluate the parameters and when the polynomial
has the same degree we compute the roots of a univariate polynomial. That is a superset
of the roots we want, so we check and we keep those that were roots for the indicial polynomial.

Minor changes:
   matrix.py now use the method is_Matrix instead of sage.matrix.Matrix (deprecated in
SAGE: https://trac.sagemath.org/ticket/24096