New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Faster CSE #2355

Merged
merged 36 commits into from Sep 3, 2013

Conversation

Projects
None yet
7 participants
@cdsousa
Contributor

cdsousa commented Aug 5, 2013

Here is a proposal for a faster CSE.

For small expressions the CSE computation time is the same, but there is a huge improvement on large expressions, e.g., http://nbviewer.ipython.org/5986996.

I've let previous cse function, as prev_cse(), inside the file for easier compare.

I plan to do some comments on code in this PR tomorrow, explaining some decisions and asking for some help.

@moorepants

This comment has been minimized.

Member

moorepants commented Aug 6, 2013

Very cool, thanks for adding this!

@asmeurer

View changes

sympy/simplify/cse_main.py Outdated
@@ -210,7 +211,7 @@ def _remove_singletons(reps, exprs):
reps[:] = u_reps # change happens in-place
def cse(exprs, symbols=None, optimizations=None, postprocess=None):
def prev_cse(exprs, symbols=None, optimizations=None, postprocess=None):

This comment has been minimized.

@asmeurer

asmeurer Aug 6, 2013

Member

We need a better migration path than this.

@asmeurer

This comment has been minimized.

Member

asmeurer commented Aug 6, 2013

Can you list what the major differences are between the new and old cse are, from a user stand-point (I'm trying to figure out what the best thing to do with the old cse() is).

@jrioux

This comment has been minimized.

Member

jrioux commented Aug 6, 2013

SymPy Bot Summary: Passed after merging cdsousa/fast_cse (5af022aebc1f2bd79f371c98aa039a9551fc79c1) into master (7291221).
PyPy 2.0.0-beta-1; 2.7.3-final-42: pass
Python 2.7.2-final-0: pass
Python 3.2.1-final-0: pass
Sphinx 1.1.3: pass
Doc Coverage: decreased
36.047% of functions have doctests (compared to 36.065% in master)
41.531% of functions are imported into Sphinx (compared to 41.552% in master)

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Aug 6, 2013

Ok, I've removed prev_cse().
The results of the new CSE implementation have no major differences compared to the previous one.
Slight differences in the output, which imo is at least not worse than before, implied the chances in the unit tests.
I'll do some comments in the code diff to explain the differences.

The major difference in the new implementation is that subs() calls are completely avoided.
Instead, the whole expression tree is rebuild making replacement of repeated subexpressions where appropriate (function tree_cse()).
Those repetitions are spotted by previously parsing the expression tree.

Apart from the "optimizations" parameter, which remain the same, the optimizations on Adds, Muls and inverse powers were also reimplemented (function opt_cse()).
These optimizations rely on a opt_subs dictionary with pairs in the form
Pow(x, -y) -> Pow(Pow(x, y), -1, evaluate=False)
and
Add(x, y, z) -> Add(Add(x, y), z, evaluate=False)
Add(x, y, w) -> Add(Add(x, y), w, evaluate=False)
(...)
which are taken into account in repeated subexpression count and tree rebuild steps.

Add and Mul optimization is similar to previous one, having two major differences.

  • non-commutative arguments are not optimized, e.g., the non-commutative B*C part in x*y*A*B*C + x*y*B*C will not be replaced by an intermediate variable (although the commutative x*y part will be); I think this is not a big problem since even in previous cse() non-commutatives had some issues (see https://code.google.com/p/sympy/issues/detail?id=3618); for the given example, previous CSE will output [(x0, B*C)], [x*y*(x0 + x0*A)] where x0*A doesn't respect non-commutativeness;
  • common terms in Adds and Muls can be recursively matched, e.g.,
    • before: [w*x, w*x*y, w*x*y*z] results in [(x0, w*x)], [x0, x0*y, x0*y*z];
    • now: [w*x, w*x*y, w*x*y*z] results in [(x0, w*x), (x1, x0*y)], [x0, x1, x1*z].

Apart from this, there are some changes, which imo would improve CSE, that I've included before but then removed for inclusion in a future PR:

  • handle lists of Matrix's as input;
  • when input is a Basic, output Basic instead of [Basic].
@cdsousa

View changes

sympy/core/tests/test_expand.py Outdated
(x4, sin(x3)), (x5, atan2(0, x2)/4), (x6, sin(x5)), (x7, x4*x6),
(x8, cos(x3)), (x9, x6*x8), (x10, cos(x5)), (x11, x10*x4),
(x12, x10*x8)], [sqrt(2)*(x11 + I*x11 - x12 + I*x12 + x7 - I*x7 + x9 +
I*x9)/(8*pi**(3/2)*x2**(1/4))])''')

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Different, but equally valid, result from new cse implementation.

This comment has been minimized.

@asmeurer

asmeurer Aug 7, 2013

Member

That's fine. Even minor changes to cse or to the internal representation can change the output. It was never guaranteed to be what it was. OTOH, make sure that the output is always independent of hash values and arg orderings.

This comment has been minimized.

@cdsousa

cdsousa Aug 7, 2013

Contributor

OTOH, make sure that the output is always independent of hash values and arg orderings.

I do that for Mul and Add terms, see https://github.com/sympy/sympy/pull/2355/files#r5629536.
Is it necessary for any other Functions?

@@ -549,7 +549,7 @@ def test_issue_3460():
e = -log(-12*sqrt(2) + 17)/24 - log(-2*sqrt(2) + 3)/12 + sqrt(2)/3
# XXX modify cse so x1 is eliminated and x0 = -sqrt(2)?
assert cse(e) == (
[(x0, sqrt(2)), (x1, -x0)], [x0/3 - log(2*x1 + 3)/12 - log(12*x1 + 17)/24])
[(x0, sqrt(2))], [x0/3 - log(-12*x0 + 17)/24 - log(-2*x0 + 3)/12])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Different, but equally valid, result from new cse implementation.

@cdsousa

View changes

sympy/simplify/cse_main.py Outdated
c, nc = expr.args_cnc()
args = list(ordered(c)) + nc
elif Op is Add:
args = list(ordered(args))

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Parse Muls and Adds arguments by order so output order doesn't depend on machine dependent hash. See test_numbered_symbols() and mattpap@e72d3ac.

assert reduced == [x0*x1, x0, x1]
assert substs == [(x0, -z), (x1, x + x0), (x2, x0 + y)]
assert rsubsts == [(x0, -z), (x1, x0 + y), (x2, x + x0)]
assert reduced == [x1*x2, x1, x2]

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Two issues:

  • Results are now consistent: since cse([x - z, y - z]) returns [(x0, -z)], [x + x0, x0 + y] (both in previous and new CSE implementations), such replacement must also appear in cse([(x - z)*(y - z), x - z, y - z]); (Arguable not optimal transformation from [x - z, y - z] to [(x0, -z)], [x + x0, x0 + y] may be handled by future consistent optimizations.)
  • Intermediate variables order is now dependent on their appearance in the input expression list (this is a side effect of new implementation, not necessarily bad).
l = [w*y + w + x + y + z, w*x*y]
assert cse(l) == ([(x0, w*y)], [w + x + x0 + y + z, x*x0])
assert cse([x + y, x + y + z]) == ([(x0, x + y)], [x0, z + x0])
assert cse([x + y, x + z]) == ([], [x + y, x + z])
assert cse([x*y, z + x*y, x*y*z + 3]) == \
([(x0, x*y)], [x0, z + x0, 3 + x0*z])
@XFAIL
def test_non_commutative_cse():
A, B, C = symbols('A B C', commutative=False)
l = [A*B*C, A*C]
assert cse(l) == ([], l)
l = [A*B*C, A*B]
assert cse(l) == ([(x0, A*B)], [x0*C, x0])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Non-commutative Mul arguments not CSEed anymore...

This comment has been minimized.

@asmeurer

asmeurer Aug 7, 2013

Member

But if the previous results were wrong, this is a good thing. Tests should be added for the new behavior that is correct now.

This comment has been minimized.

@cdsousa

cdsousa Aug 7, 2013

Contributor

I've added a unit test to current behavior.
I've also discovered that replacements of non-commutative sub-expressions fails even on new implementation, see comment on diff of cdsousa@99819d6.

@cdsousa

View changes

sympy/simplify/tests/test_cse.py Outdated
A, B, x0 = symbols('A B, x0', commutative=False)
l = [A*A, B*A*A]
assert cse(l) == ([(x0, A**2)], [x0, b*x0])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

If non-commutative Mul arguments are to be CSEed, then intermediate variables for those must also be non-commutative.
This test is wrong anyway, I'll correct it latter.

assert cse(cos(1/x**2) + sin(1/x**2)) == \
([(x0, x**2)], [sin(1/x0) + cos(1/x0)])
([(x0, x**(-2))], [sin(x0) + cos(x0)])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Better CSE now.

@@ -201,7 +209,7 @@ def test_pow_invpow():
assert cse(exp(x**2) + x**2*cos(1/x**2)) == \
([(x0, x**2)], [x0*cos(1/x0) + exp(x0)])
assert cse((1 + 1/x**2)/x**2) == \
([(x0, x**2)], [(1 + 1/x0)/x0])
([(x0, x**(-2))], [x0*(x0 + 1)])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Better CSE now.

@@ -211,7 +219,7 @@ def test_postprocess():
assert cse([eq, Eq(x, z + 1), z - 2, (z + 1)*(x + 1)],
postprocess=cse_main.cse_separate) == \
[[(x1, y + 1), (x2, z + 1), (x, x2), (x0, x + 1)],
[x0 + exp(x0/x1) + cos(x1), x2 - 3, x0*x2]]
[x0 + exp(x0/x1) + cos(x1), z - 2, x0*x2]]

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Slightly better CSE now.

B(x0 + x1, x4)), (x6, G(b)), (x7, G(x3)), (x8, -x0), (x9,
(x4/2)**(x8 + 1)), (x10, x6*x7*x9*B(b - 1, x4)), (x11, x6*x7*x9*B(b,
x4)), (x12, B(x3, x4))], [(a, a + S(1)/2, x0, b, x3, x10*x5,
x11*x4*x5, x10*x12*x4, x11*x12, 1, 0, S(1)/2, z/2, x2, b + x8, x8)])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Different, yet equally valid, result from new cse implementation.

@cdsousa

View changes

sympy/simplify/cse_main.py Outdated
opt_subs[ops[i]] = Op(Op(*diff_i), com_op, evaluate=False)
diff_j = op_args[j].difference(com_args)
op_args[j] = diff_j | set([com_op])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

In previous Mul/Add optimization, when a common set of terms/symbols was found it was removed from expressions so it doesn't appear again. Previous code:

# remove this set of symbols so it doesn't appear again
            adds[i] = adds[i].difference(com)

Now, the common term operation is included in the result set for possible future matches:

           com_op = Op(*com_args)
           diff_i = op_args[i].difference(com_args)
           op_args[i] = diff_i | set([com_op])

, allowing cse([w*x, w*x*y, w*x*y*z]) to result in [(x0, w*x), (x1, x0*y)], [x0, x1, x1*z].

This comment has been minimized.

@asmeurer

asmeurer Aug 7, 2013

Member

You might add this sort of thing as a comment in the code.

This comment has been minimized.

@cdsousa
@cdsousa

View changes

sympy/simplify/cse_main.py Outdated
if isinstance(expr, Basic) and expr.is_Atom:
return
if isinstance(expr, bool):
return

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

if isinstance(expr, bool) was added here and in two more places as a dirty fix to test_issue_3164().
A better fix must be found.

This comment has been minimized.

@asmeurer

asmeurer Aug 7, 2013

Member

Maybe just return None on any non-Basic.

This comment has been minimized.

@cdsousa
else:
used = j - len(exprs)
if type(used) is int:

This comment has been minimized.

@smichr

smichr Aug 6, 2013

Member

Does cse(2_x+2_y+2*z) return the expression unchanged?

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Yes, indeed...

preprocess_for_cse with factor_terms firstly transforms 2*x+2*y+2*z into 2*(x+y+z), but this factorization is lost when tree is rebuilt due to new_expr = Op(*map(_recreate, args))...

A quick and possibly dirty fix would be putting
optimizations=[(factor_terms, factor_terms)] instead of optimizations=[(factor_terms, None)]

Nonetheless, current fast_cse branch code outputs
([(x0, x + y + z)], [2*x0, 3*x0])
for
cse([2*x+2*y+2*z, 3*x+3*y+3*z], optimizations=[(factor_terms, None)]),
so, if factor_terms is used just for increasing optimization opportunities, such opportunities are still taken into account.

Moreover, for master branch cse
cse([2*sin(x) + 2*y + 2*z, sin(x)])
outputs
[(x0, sin(x))], [2*x0 + 2*y + 2*z, x0],
rather than
[(x0, sin(x))], [2*(x0 + y + z), x0].

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

In fact there is a good way to make proposed cse to better match current cse output:
doing
new_expr = type(expr)(*map(_recreate, expr.args))
only if map(_recreate, expr.args) != expr.args.

Done and committed: cdsousa@be29279

So now, cse(2*x + 2*y + 2*z) returns 2*(x + y + z) when the default optimizations=[(factor_terms, None)] is used, just like master's cse implementation.

@cdsousa

View changes

sympy/simplify/tests/test_cse.py Outdated
A, B, C = symbols('A B C', commutative=False)
x0 = symbols('x0', commutative=False)
l = [B*C, A*B*C]
assert cse(l) == ([(x0, B*C)], [x0, A*x0])

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

If non-commutative Mul arguments are to be CSEed, then the intermediate variables for those must also be non-commutative.

@cdsousa

View changes

sympy/simplify/cse_main.py Outdated
if expr in opt_subs:
expr = opt_subs[expr]
Func, args = type(expr), expr.args

This comment has been minimized.

@cdsousa

cdsousa Aug 6, 2013

Contributor

Is this line ok, or it must have expr.func instead of type(expr)?

This comment has been minimized.

@asmeurer

asmeurer Aug 7, 2013

Member

Use .func. Currently, they are identical, but this might change in the future. The rebuild identity is expr.func(*expr.args).

This comment has been minimized.

@cdsousa
@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Aug 6, 2013

Updated large expression example at http://nbviewer.ipython.org/5986996.

@certik

This comment has been minimized.

Member

certik commented Aug 6, 2013

Thanks @cdsousa for this!! It will be very helpful to a lot of people. I'll try to go over it soon.

@jrioux

This comment has been minimized.

Member

jrioux commented Aug 7, 2013

SymPy Bot Summary: Passed after merging cdsousa/fast_cse (be292794c0efdcf138482d61f5d515797fabafea) into master (84ea466).
PyPy 2.0.0-beta-1; 2.7.3-final-42: pass
Python 2.7.2-final-0: pass
Python 3.2.1-final-0: pass
Sphinx 1.1.3: pass
Doc Coverage: decreased
36.053% of functions have doctests (compared to 36.065% in master)
41.538% of functions are imported into Sphinx (compared to 41.552% in master)

@jrioux

This comment has been minimized.

Member

jrioux commented Aug 7, 2013

SymPy Bot Summary: Passed after merging cdsousa/fast_cse (be292794c0efdcf138482d61f5d515797fabafea) into master (84ea466).
PyPy 2.0.0-beta-1; 2.7.3-final-42: pass
Python 2.7.2-final-0: pass
Python 3.2.1-final-0: pass
Sphinx 1.1.3: pass
Doc Coverage: decreased
36.053% of functions have doctests (compared to 36.065% in master)
41.538% of functions are imported into Sphinx (compared to 41.552% in master)

@hazelnusse

This comment has been minimized.

Contributor

hazelnusse commented Aug 7, 2013

@cdsousa this is awesome, thanks for the work on this. I am about to test on an example I have that currently never finishes (after 20+ hours using sympy master) to see how it does. I'll report back with results.

@asmeurer

This comment has been minimized.

Member

asmeurer commented Aug 7, 2013

In your notebook at the end, are larger or smaller numbers better?

@cdsousa cdsousa closed this Aug 7, 2013

@cdsousa cdsousa reopened this Aug 7, 2013

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Aug 7, 2013

Oops... Sorry, wrong click.

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Aug 7, 2013

I've added a verification step to the example notebook, and I've also added a test to the main/raw CSE routine (tree_cse()).

@asmeurer

In your notebook at the end, are larger or smaller numbers better?

Those counts give a fast and rough idea about the number of operations/functions per CSEed expression (it just counts characters...).
Small numbers 'mean' less operations, so small numbers are better...?
Anyway, even if for a given expression a CSE implementation gives less operations than another one, it doesn't necessarily mean that the same happens for other expressions.

@hazelnusse

I am about to test on an example I have that currently never finishes [...]

Let me suggest you to start doing a raw CSE with

from sympy.simplify import cse_main
raw_cse = cse_main.tree_cse
import sympy
from sympy.utilities import numbered_symbols

raw_cse(expressions, numbered_symbols())

and then a non-optimized CSE with

import sympy
sympy.cse(expressions, optimizations=[])

, before doing a full CSE.
This will almost surely give increasing times, and you can analyse the optimization/performance tradeoff.
It is also possible that for large expressions, optimizing CSE gives less optimized output, as it happens in my example notebook.

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Aug 28, 2013

Ok, I've made a master rebase and I've tried to improve documentation.

@jrioux

This comment has been minimized.

Member

jrioux commented Aug 28, 2013

SymPy Bot Summary: 🔴 Failed after merging cdsousa/fast_cse (3dc2845) into master (e968bc7).
@cdsousa: Please fix the test failures.
🔴 PyPy 2.0.0-beta-1; 2.7.3-final-42: fail
Python 2.7.2-final-0: pass
Python 3.2.1-final-0: pass
Sphinx 1.1.3: pass
Doc Coverage: increased
36.320% of functions have doctests (compared to 36.315% in master)
41.665% of functions are imported into Sphinx (compared to 41.645% in master)

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Aug 28, 2013

PyPy 2.0.0-beta-1; 2.7.3-final-42: fail @#2355 (comment)

I hadn't test with Pypy before, but now I installed it and fixed the bug.

@asmeurer

This comment has been minimized.

Member

asmeurer commented Aug 30, 2013

Yes, never use is to compare unless you really do want identity comparison, or unless you are comparing against a type that is guaranteed to be singletonized (True, False, None in Python, S.anything in SymPy).

@asmeurer

This comment has been minimized.

Member

asmeurer commented Aug 30, 2013

What remains to be done here?

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Aug 30, 2013

I think it is quite ok now.
Though, it would be best if you guys do a quick review, at least ...

@certik

This comment has been minimized.

Member

certik commented Sep 3, 2013

This is +1 from me. Let me review it some more before I merge it though.

@certik

This comment has been minimized.

Member

certik commented Sep 3, 2013

It looks good. Tests pass, so I am merging this. We can fix further problems in a new PR if needed. Thanks for all the hard work @cdsousa !

certik added a commit that referenced this pull request Sep 3, 2013

@certik certik merged commit df897c3 into sympy:master Sep 3, 2013

@asmeurer

This comment has been minimized.

Member

asmeurer commented Sep 3, 2013

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Sep 3, 2013

@certik Me, I thank all you guys for SymPy ;)

@asmeurer Sure :)

@cdsousa cdsousa deleted the cdsousa:fast_cse branch Sep 3, 2013

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Sep 4, 2013

@asmeurer

@cdsousa can you add a bit to https://github.com/sympy/sympy/wiki/release-notes-for-0.7.4 about this?

Done. (Maybe it got too verbose...?)

@asmeurer

This comment has been minimized.

Member

asmeurer commented Sep 4, 2013

I think that's fine. Thanks.

@moorepants

This comment has been minimized.

Member

moorepants commented Sep 4, 2013

@cdsousa Thanks for all the work on this. It is a great improvement that will help lots of folks using SymPy.

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Sep 4, 2013

I'm really happy for contributing to SymPy.

I think CSE can be used in SymPy's code generation modules.

One day I would like to study the inclusion of MatrixExprs in CSE.
For that, one needs Symbols which can somehow clone expressions assumptions...

@cdsousa

This comment has been minimized.

Contributor

cdsousa commented Sep 8, 2013

@asmeurer

This comment has been minimized.

Member

asmeurer commented Sep 22, 2013

SymPy Bot Summary: Could not fetch the branch cdsousa/fast_cse.
@cdsousa: Please make sure that cdsousa/fast_cse has been pushed to GitHub and run the sympy-bot tests again.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment