Repository: sympy/sympy

Issue:
rref() should return a Matrix
Currently calling `.rref()` on a matrix produces a tuple with a matrix and a list of pivots.  I think this is unexpected for new users and the default should be that `rref()` returns the row-reduced matrix only.  An `include_pivots` flag could be added to get the current behavior where a tuple is returned.

Discussion:

Patch:
diff --git a/sympy/concrete/summations.py b/sympy/concrete/summations.py
index 0a6fb13bdf..5b6012e1cb 100644
--- a/sympy/concrete/summations.py
+++ b/sympy/concrete/summations.py
@@ -265,41 +265,9 @@ def _eval_difference_delta(self, n, step):
 
         return Sum(f, (k, upper + 1, new_upper)).doit()
 
-    def _eval_simplify(self, ratio=1.7, measure=None):
-        from sympy.simplify.simplify import factor_sum, sum_combine
-        from sympy.core.function import expand
-        from sympy.core.mul import Mul
-
-        #split the function into adds
-        terms = Add.make_args(expand(self.function))
-        s_t = [] # Sum Terms
-        o_t = [] # Other Terms
-
-        for term in terms:
-            if term.has(Sum):
-                #if there is an embedded sum here
-                #it is of the form x * (Sum(whatever))
-                #hence we make a Mul out of it, and simplify all interior sum terms
-                subterms = Mul.make_args(expand(term))
-                out_terms = []
-                for subterm in subterms:
-                    #go through each term
-                    if isinstance(subterm, Sum):
-                        #if it's a sum, simpify it
-                        out_terms.append(subterm._eval_simplify())
-                    else:
-                        #otherwise, add it as is
-                        out_terms.append(subterm)
-
-                #turn it back into a Mul
-                s_t.append(Mul(*out_terms))
-            else:
-                o_t.append(term)
-
-        #next try to combine any interior sums for further simplification
-        result = Add(sum_combine(s_t), *o_t)
-
-        return factor_sum(result, limits=self.limits)
+    def _eval_simplify(self, ratio, measure):
+        from sympy.simplify.simplify import sum_simplify
+        return sum_simplify(self)
 
     def _eval_summation(self, f, x):
         return None
diff --git a/sympy/core/exprtools.py b/sympy/core/exprtools.py
index 09ac9f62fb..bff811ece9 100644
--- a/sympy/core/exprtools.py
+++ b/sympy/core/exprtools.py
@@ -1144,8 +1144,6 @@ def factor_terms(expr, radical=False, clear=False, fraction=False, sign=True):
 
     """
     def do(expr):
-        from sympy.concrete.summations import Sum
-        from sympy.simplify.simplify import factor_sum
         is_iterable = iterable(expr)
 
         if not isinstance(expr, Basic) or expr.is_Atom:
@@ -1161,9 +1159,6 @@ def do(expr):
                 return expr
             return expr.func(*newargs)
 
-        if isinstance(expr, Sum):
-            return factor_sum(expr, radical=radical, clear=clear, fraction=fraction, sign=sign)
-
         cont, p = expr.as_content_primitive(radical=radical, clear=clear)
         if p.is_Add:
             list_args = [do(a) for a in Add.make_args(p)]
diff --git a/sympy/matrices/matrices.py b/sympy/matrices/matrices.py
index f458f7d1e7..2fd4efa8ca 100644
--- a/sympy/matrices/matrices.py
+++ b/sympy/matrices/matrices.py
@@ -4083,6 +4083,13 @@ def rref(self, iszerofunc=_iszero, simplify=False):
         (Matrix([
         [1, 0],
         [0, 1]]), [0, 1])
+        >>> rref_matrix, rref_pivots = m.rref()
+        >>> rref_matrix
+        Matrix([
+        [1, 0],
+        [0, 1]])
+        >>> rref_pivots
+        [0, 1]
         """
         simpfunc = simplify if isinstance(
             simplify, FunctionType) else _simplify
diff --git a/sympy/simplify/simplify.py b/sympy/simplify/simplify.py
index fc803f8994..a8413332a6 100644
--- a/sympy/simplify/simplify.py
+++ b/sympy/simplify/simplify.py
@@ -615,55 +615,37 @@ def shorter(*choices):
 def sum_simplify(s):
     """Main function for Sum simplification"""
     from sympy.concrete.summations import Sum
-    from sympy.core.function import expand
 
-    terms = Add.make_args(expand(s))
+    terms = Add.make_args(s)
     s_t = [] # Sum Terms
     o_t = [] # Other Terms
 
     for term in terms:
         if isinstance(term, Mul):
+            constant = 1
             other = 1
-            sum_terms = []
-
-            if not term.has(Sum):
-                o_t.append(term)
-                continue
-
-            mul_terms = Mul.make_args(term)
-            for mul_term in mul_terms:
-                if isinstance(mul_term, Sum):
-                    r = mul_term._eval_simplify()
-                    sum_terms.extend(Add.make_args(r))
+            s = 0
+            n_sum_terms = 0
+            for j in range(len(term.args)):
+                if isinstance(term.args[j], Sum):
+                    s = term.args[j]
+                    n_sum_terms = n_sum_terms + 1
+                elif term.args[j].is_number == True:
+                    constant = constant * term.args[j]
                 else:
-                    other = other * mul_term
-            if len(sum_terms):
-                #some simplification may have happened
-                #use if so
-                s_t.append(Mul(*sum_terms) * other)
+                    other = other * term.args[j]
+            if other == 1 and n_sum_terms == 1:
+                # Insert the constant inside the Sum
+                s_t.append(Sum(constant * s.function, *s.limits))
+            elif other != 1 and n_sum_terms == 1:
+                o_t.append(other * Sum(constant * s.function, *s.limits))
             else:
-                o_t.append(other)
+                o_t.append(term)
         elif isinstance(term, Sum):
-            #as above, we need to turn this into an add list
-            r = term._eval_simplify()
-            s_t.extend(Add.make_args(r))
+            s_t.append(term)
         else:
             o_t.append(term)
 
-
-    result = Add(sum_combine(s_t), *o_t)
-
-    return result
-
-def sum_combine(s_t):
-    """Helper function for Sum simplification
-
-       Attempts to simplify a list of sums, by combining limits / sum function's
-       returns the simplified sum
-    """
-    from sympy.concrete.summations import Sum
-
-
     used = [False] * len(s_t)
 
     for method in range(2):
@@ -672,93 +654,43 @@ def sum_combine(s_t):
                 for j, s_term2 in enumerate(s_t):
                     if not used[j] and i != j:
                         temp = sum_add(s_term1, s_term2, method)
-                        if isinstance(temp, Sum) or isinstance(temp, Mul):
+                        if isinstance(temp, Sum):
                             s_t[i] = temp
                             s_term1 = s_t[i]
                             used[j] = True
 
-    result = S.Zero
+    result = Add(*o_t)
+
     for i, s_term in enumerate(s_t):
         if not used[i]:
             result = Add(result, s_term)
 
     return result
 
-def factor_sum(self, limits=None, radical=False, clear=False, fraction=False, sign=True):
-    """Helper function for Sum simplification
-
-       if limits is specified, "self" is the inner part of a sum
-
-       Returns the sum with constant factors brought outside
-    """
-    from sympy.core.exprtools import factor_terms
-    from sympy.concrete.summations import Sum
-
-    result = self.function if limits is None else self
-    limits = self.limits if limits is None else limits
-    #avoid any confusion w/ as_independent
-    if result == 0:
-        return S.Zero
-
-    #get the summation variables
-    sum_vars = set([limit.args[0] for limit in limits])
-
-    #finally we try to factor out any common terms
-    #and remove the from the sum if independent
-    retv = factor_terms(result, radical=radical, clear=clear, fraction=fraction, sign=sign)
-    #avoid doing anything bad
-    if not result.is_commutative:
-        return Sum(result, *limits)
-
-    i, d = retv.as_independent(*sum_vars)
-    if isinstance(retv, Add):
-        return i * Sum(1, *limits) + Sum(d, *limits)
-    else:
-        return i * Sum(d, *limits)
 
 def sum_add(self, other, method=0):
     """Helper function for Sum simplification"""
     from sympy.concrete.summations import Sum
-    from sympy import Mul
-
-    #we know this is something in terms of a constant * a sum
-    #so we temporarily put the constants inside for simplification
-    #then simplify the result
-    def __refactor(val):
-        args = Mul.make_args(val)
-        sumv = next(x for x in args if isinstance(x, Sum))
-        constant = Mul(*[x for x in args if x != sumv])
-        return Sum(constant * sumv.function, *sumv.limits)
-
-    if isinstance(self, Mul):
-        rself = __refactor(self)
-    else:
-        rself = self
 
-    if isinstance(other, Mul):
-        rother = __refactor(other)
-    else:
-        rother = other
-
-    if type(rself) == type(rother):
+    if type(self) == type(other):
         if method == 0:
-            if rself.limits == rother.limits:
-                return factor_sum(Sum(rself.function + rother.function, *rself.limits))
+            if self.limits == other.limits:
+                return Sum(self.function + other.function, *self.limits)
         elif method == 1:
-            if simplify(rself.function - rother.function) == 0:
-                if len(rself.limits) == len(rother.limits) == 1:
-                    i = rself.limits[0][0]
-                    x1 = rself.limits[0][1]
-                    y1 = rself.limits[0][2]
-                    j = rother.limits[0][0]
-                    x2 = rother.limits[0][1]
-                    y2 = rother.limits[0][2]
+            if simplify(self.function - other.function) == 0:
+                if len(self.limits) == len(other.limits) == 1:
+                    i = self.limits[0][0]
+                    x1 = self.limits[0][1]
+                    y1 = self.limits[0][2]
+                    j = other.limits[0][0]
+                    x2 = other.limits[0][1]
+                    y2 = other.limits[0][2]
 
                     if i == j:
                         if x2 == y1 + 1:
-                            return factor_sum(Sum(rself.function, (i, x1, y2)))
+                            return Sum(self.function, (i, x1, y2))
                         elif x1 == y2 + 1:
-                            return factor_sum(Sum(rself.function, (i, x2, y1)))
+                            return Sum(self.function, (i, x2, y1))
 
     return Add(self, other)
 


Test patch:
diff --git a/sympy/concrete/tests/test_sums_products.py b/sympy/concrete/tests/test_sums_products.py
index 32007dc680..5d277cd414 100644
--- a/sympy/concrete/tests/test_sums_products.py
+++ b/sympy/concrete/tests/test_sums_products.py
@@ -702,51 +702,34 @@ def test_simplify():
     y, t, v = symbols('y, t, v')
 
     assert simplify(Sum(x*y, (x, n, m), (y, a, k)) + \
-        Sum(y, (x, n, m), (y, a, k))) == Sum(y * (x + 1), (x, n, m), (y, a, k))
+        Sum(y, (x, n, m), (y, a, k))) == Sum(x*y + y, (x, n, m), (y, a, k))
     assert simplify(Sum(x, (x, n, m)) + Sum(x, (x, m + 1, a))) == \
         Sum(x, (x, n, a))
     assert simplify(Sum(x, (x, k + 1, a)) + Sum(x, (x, n, k))) == \
         Sum(x, (x, n, a))
     assert simplify(Sum(x, (x, k + 1, a)) + Sum(x + 1, (x, n, k))) == \
-        Sum(x, (x, n, a)) + Sum(1, (x, n, k))
+        Sum(x, (x, k + 1, a)) + Sum(x + 1, (x, n, k))
     assert simplify(Sum(x, (x, 0, 3)) * 3 + 3 * Sum(x, (x, 4, 6)) + \
-        4 * Sum(z, (z, 0, 1))) == 4*Sum(z, (z, 0, 1)) + 3*Sum(x, (x, 0, 6))
+        4 * Sum(z, (z, 0, 1))) == Sum(4*z, (z, 0, 1)) + Sum(3*x, (x, 0, 6))
     assert simplify(3*Sum(x**2, (x, a, b)) + Sum(x, (x, a, b))) == \
-        Sum(x*(3*x + 1), (x, a, b))
+        Sum(3*x**2 + x, (x, a, b))
     assert simplify(Sum(x**3, (x, n, k)) * 3 + 3 * Sum(x, (x, n, k)) + \
         4 * y * Sum(z, (z, n, k))) + 1 == \
-            4*y*Sum(z, (z, n, k)) + 3*Sum(x**3 + x, (x, n, k)) + 1
+            y*Sum(4*z, (z, n, k)) + Sum(3*x**3 + 3*x, (x, n, k)) + 1
     assert simplify(Sum(x, (x, a, b)) + 1 + Sum(x, (x, b + 1, c))) == \
         1 + Sum(x, (x, a, c))
     assert simplify(Sum(x, (t, a, b)) + Sum(y, (t, a, b)) + \
-        Sum(x, (t, b+1, c))) == x * Sum(1, (t, a, c)) + y * Sum(1, (t, a, b))
+        Sum(x, (t, b+1, c))) == Sum(x + y, (t, a, b)) + Sum(x, (t, b+1, c))
     assert simplify(Sum(x, (t, a, b)) + Sum(x, (t, b+1, c)) + \
-        Sum(y, (t, a, b))) == x * Sum(1, (t, a, c)) + y * Sum(1, (t, a, b))
+        Sum(y, (t, a, b))) == Sum(x + y, (t, a, b)) + Sum(x, (t, b+1, c))
     assert simplify(Sum(x, (t, a, b)) + 2 * Sum(x, (t, b+1, c))) == \
         simplify(Sum(x, (t, a, b)) + Sum(x, (t, b+1, c)) + Sum(x, (t, b+1, c)))
     assert simplify(Sum(x, (x, a, b))*Sum(x**2, (x, a, b))) == \
         Sum(x, (x, a, b)) * Sum(x**2, (x, a, b))
     assert simplify(Sum(x, (t, a, b)) + Sum(y, (t, a, b)) + Sum(z, (t, a, b))) \
-        == (x + y + z) * Sum(1, (t, a, b))          # issue 8596
+        == Sum(x + y + z, (t, a, b))          # issue 8596
     assert simplify(Sum(x, (t, a, b)) + Sum(y, (t, a, b)) + Sum(z, (t, a, b)) + \
-        Sum(v, (t, a, b))) == (x + y + z + v) * Sum(1, (t, a, b))  # issue 8596
-    assert simplify(Sum(x * y, (x, a, b)) / (3 * y)) == \
-        (Sum(x, (x, a, b)) / 3)
-    assert simplify(Sum(Function('f')(x) * y * z, (x, a, b)) / (y * z)) \
-        == Sum(Function('f')(x), (x, a, b))
-    assert simplify(Sum(c * x, (x, a, b)) - c * Sum(x, (x, a, b))) == 0
-    assert simplify(c * (Sum(x, (x, a, b))  + y)) == c * (y + Sum(x, (x, a, b)))
-    assert simplify(c * (Sum(x, (x, a, b)) + y * Sum(x, (x, a, b)))) == \
-        c * (y + 1) * Sum(x, (x, a, b))
-    assert simplify(Sum(Sum(c * x, (x, a, b)), (y, a, b))) == \
-                c * Sum(x, (x, a, b), (y, a, b))
-    assert simplify(Sum((3 + y) * Sum(c * x, (x, a, b)), (y, a, b))) == \
-                c * Sum((3 + y), (y, a, b)) * Sum(x, (x, a, b))
-    assert simplify(Sum((3 + t) * Sum(c * t, (x, a, b)), (y, a, b))) == \
-                c*t*(t + 3)*Sum(1, (x, a, b))*Sum(1, (y, a, b))
-    assert simplify(Sum(Sum(d * t, (x, a, b - 1)) + \
-                Sum(d * t, (x, b, c)), (t, a, b))) == \
-                    d * Sum(1, (x, a, c)) * Sum(t, (t, a, b))
+        Sum(v, (t, a, b))) == Sum(x + y + z + v, (t, a, b))  # issue 8596
 
 
 def test_change_index():
@@ -856,7 +839,7 @@ def test_issue_2787():
     res = s.doit().simplify()
     assert res == Piecewise(
         (n*p, p/Abs(p - 1) <= 1),
-        ((-p + 1)**n*Sum(k*p**k*(-p + 1)**(-k)*binomial(n, k), (k, 0, n)),
+        (Sum(k*p**k*(-p + 1)**(-k)*(-p + 1)**n*binomial(n, k), (k, 0, n)),
         True))
 
 
diff --git a/sympy/core/tests/test_exprtools.py b/sympy/core/tests/test_exprtools.py
index 77541b8630..b41ff8f892 100644
--- a/sympy/core/tests/test_exprtools.py
+++ b/sympy/core/tests/test_exprtools.py
@@ -288,11 +288,6 @@ def test_factor_terms():
     assert factor_terms(e, sign=False) == e
     assert factor_terms(exp(-4*x - 2) - x) == -x + exp(Mul(-2, 2*x + 1, evaluate=False))
 
-    # sum tests
-    assert factor_terms(Sum(x, (y, 1, 10))) == x * Sum(1, (y, 1, 10))
-    assert factor_terms(Sum(x, (y, 1, 10)) + x) == x * (1 + Sum(1, (y, 1, 10)))
-    assert factor_terms(Sum(x*y + x*y**2, (y, 1, 10))) == x*Sum(y*(y + 1), (y, 1, 10))
-
 
 def test_xreplace():
     e = Mul(2, 1 + x, evaluate=False)