Note that I've written these in a form that makes them useful as
transformations.
I think you may want to try writing a custom simplification routine
that uses these. It would go something like this:
First, expand the original expression and get a list of the monomails.
You can use poly() to do this. For example:
In [100]: a = expressions.smallExpr
If you don't know before hand which trig functions are in your
expression, you can use something like this:
In [101]: from sympy.functions.elementary.trigonometric import
TrigonometricFunction
In [102]: a.atoms(TrigonometricFunction)
Out[102]: set([cos(q4), sin(q1), sin(q4), cos(q1), cos(q2), sin(q2),
sin(q3), cos(q3)])
In [116]: p = Poly(a, sin(q1), sin(q2), sin(q3), sin(q4), cos(q1),
cos(q2), cos(q3), cos(q4), sin(2*q1), sin(2*q2), sin(2*q3), sin(2*q4),
cos(2*q1), cos(2*q2), cos(2*q3), cos(2*q4))
In [117]: p.monoms()
Out[117]:
[(4, 4, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 3, 2, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 3, 1, 1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 2, 3, 0, 0, 2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 2, 2, 1, 0, 2, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 2, 2, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 2, 1, 0, 0, 2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 2, 1, 0, 0, 0, 2, 1, 0, 0, 0, 0, 0, 0, 0, 0),
(4, 2, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0),...
Note that I only send the generators to poly() that are needed. dq1,
etc. are treated as constants. I include the double argument terms
because these will be needed after the transformations.
This list means that you have a terms with the given monomials. For
example, (4, 4, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0) means that
there is a term sin(q1)**4*sin(q2)**4*sin(q3)*cos(q4). You can obtain
the coefficient of this term by using p.coeffs() (they are in the same
order as the monomials, you can put them together with zip).
The idea would be to work on the monomials, which is much more space
efficient than working on the expressions directly. The above rules
would then be something like (sorry, I've switched from numbering
starting at 1 to numbering starting at 0):
== a*sin(x)**2 + a*cos(x)**2 => a ==
- If (2 + s0, s1, s2, s3, c0, c1, c2, c3, ds0, ds1, ds2, ds3, dc0,
dc1, dc2, dc3) and (s0, s1, s2, s3, 2 + c0, c1, c2, c3, ds0, ds1, ds2,
ds3, dc0, dc1, dc2, dc3) are both monomails with the same coefficient
a, then they can be replaced with a single monomial (s0, s1, s2, s3,
c0, c1, c2, c3, ds0, ds1, ds2, ds3, dc0, dc1, dc2, dc3) with
coefficient a.
- If (s0, 2 + s1, s2, s3, c0, c1, c2, c3, ds0, ds1, ds2, ds3, dc0,
dc1, dc2, dc3) and (s0, s1, s2, s3, c0, 2 + c1, c2, c3, ds0, ds1, ds2,
ds3, dc0, dc1, dc2, dc3) are both monomails with the same coefficient
a, then they can be replaced with a single monomial (s0, s1, s2, s3,
c0, c1, c2, c3, ds0, ds1, ds2, ds3, dc0, dc1, dc2, dc3) with
coefficient a.
...
These rules can probably be generalized in to a single rule using some
more advanced logic. Note how I allow for things like
sin(q1)**4*cos(q2)**2 + cos(q1)**2*sin(q1)**2*cos(q2)**2 =>
sin(q1)**2*cos(q2)**2. 2 + s0 would be implemented as "if s0>= 2:".
== a*sin(x)*cos(x) => a/2*sin(2*x) ==
- If (n, s1, s2, s3, n, c1, c2, c3, ds0, ds1, ds2, ds3, dc0, dc1, dc2,
dc3) is a monomial with coefficient a, it can be replaced with (0, s1,
s2, s3, 0, c1, c2, c3, ds0 + n, ds1, ds2, ds3, dc0, dc1, dc2, dc3)
with coefficient a/2**n.
- If (s0, n, s2, s3, c0, n, c2, c3, ds0, ds1, ds2, ds3, dc0, dc1, dc2,
dc3) is a monomial with coefficient a, it can be replaced with (0, s1,
s2, s3, 0, c1, c2, c3, ds0, ds1 + n, ds2, ds3, dc0, dc1, dc2, dc3)
with coefficient a/2**n.
...
Note that you have to use a/2**n (not a/2!).
== a*cos(x)**2 - sin(x)**2 => a*cos(2*x) ==
- If (s0, s1, s2, s3, 2 + c0, c1, c2, c3, ds0, ds1, ds2, ds3, dc0,
dc1, dc2, dc3) is a monomial with coefficient a (2 + s0, s1, s2, s3,
c0, c1, c2, c3, ds0, ds1, ds2, ds3, dc0, dc1, dc2, dc3) is a monomial
with coefficient -a, then they can be replaced with a single monomial
(s0, s1, s2, s3, c0, c1, c2, c3, ds0, ds1, ds2, ds3, dc0 + 1, dc1,
dc2, dc3) with coefficient a.
...
The algorithm would then be to iterate these rules until no changes
can be made. The resulting expression is the fully simplified one.
You can probably be smart about the order that the rules are applied
in order to minimize the number of passes.
Be aware that I may have made typos or small logical errors above, but
you get the idea. This algorithm could probably be generalized and
put into trigsimp(), but even a slightly more custom version that fits
your needs should not be too hard to write, and would be much faster
(the only thing that will bog it down memory wise will be the
expansion of the original expression).
It also could be generalized to handle 2*n angle arguments by
repeating the double angle rules over the functions that already have
double angles, if you find that is necessary. And as you can imagine,
it can also be generalized in a bunch of other directions as well.
This logic could probably also be applied by using some fancy replace
calls, but the monomial version will be more precise, less bug prone
(because replace logic doesn't always work right for things like x**6
=> x**2*x**4), will be faster, and will use much less memory.
Sorry if this is more of a description of a solution rather than a
solution itself.
