On Wed, Apr 18, 2012 at 1:27 PM, Tom Bachmann <[email protected]> wrote: > Hi all, > > I managed to improve the situation quite a lot. The new code is this time > attached as a patch, or alternatively see my (ness01) branch trigsimp. > > Executive summary: run mytrigsimp, test the function mytrigsimp. Ignore the > debugging output, pass order=grlex or grevlex for sensible results, pass > n=2, n=3 etc to improve the results. > > > The function mytrigsimp > ----------------------- > > As outlined in the previous mail, mytrigsimp can be used to simplify a > fraction involving cos and sin terms. As currently designed it specifically > simplifes cos(n*x) and sin(n*x) stuff. (You can pass pretty much anything, > but nothing will happen to other terms except for slowing things down :D). > There is currently quite a lot of debug output which you can ignore. I also > added some short docstrings. Specifically, mytrigsimp is invoked as follows: > > mytrigsimp(expr, **opts) > > Opts is mostly passed through to polys. You probably shouldn't use the > gens=.. option. It is probably a good idea to pass order=grlex or > order=grevlex -- for reasons that are not clear to me the default order in > polys is lex. In any case, unless you pass a graded order (i.e. grlex or > grevlex), the result will *not* be of minimal total degree, in the naive > sense.
I guess just because lex is the default for printing. > > There is one non-polys option, which is currently called n. It can be used > to extend the search space. Namely, suppose you want so simplify > 2*sin(x)*cos(x). We want the answer to be sin(2*x). But currently mytrigsimp > is rather conservative, and since it sees no sin(2*x) it does not try to > find an expression involving that. Passing n=2 will double the search space, > in this case it will try sin(2*x) and cos(2*x). > > You will notice that performance degrades quickly as you increase either the > number of generators or n. Also, if n is "big" funny things can happen (e.g. > try n=4 for sin(x)*cos(x) :D). The reason is that ratsimpmodprime returns a > fraction of minimal degree, but does no other selection. > > > Changes to ratsimpmodprime > -------------------------- > > In the last mail I reported a bug in ratsimpmodprime. That was really more a > bug in my code and a failure to properly signal a problem in > ratsimpmodprime. (Specifically, there is a certain situation in which > ratsimpmodprime divides by zero, but this can only happen if the ideal is > not prime, or if the generating set is not a groebner basis.) There also > *is* a bug in ratsimpmodprime, but it is fairly minor (results correct up to > a rational number as prefactor -- fixed in my branch). > > I also implemented a (tiny) change which resulted in huge performance > improvements. I'm not 100% sure it is correct, so I'd be glad if someone > could verify the following: > > Ratsimpmodprime simplifies f/g mod I, essentially by writing out a general > fraction p/q of degree n in the generators (say (a + b*x + c*y)/(e + f*x + > g*y) for n=1, generators x,y), and then checking if there exists a choice of > coefficients such that f*q - p*g is in I. What the original code did was to > change the ground field from QQ to QQ(a, b, c, e, f, g) and then reduce > f*q-p*g modulo the groebner basis (one may easily see that changing the > ground field in this way does indeed preserve groebner bases). This will > result in an expression linear in the new coefficients, and finding a > solution is a matter of linear algebra. The problem is that computations > over this field of fractions are slow. > > My new method adds a, b, c, e, f, g to the generators of our polynomial > ring, that is we now work in k[x, y, a, .., g] instead of k[x, y]. Again > groebner bases are preserved under this operaton, and I believe a (the) > "completely reduced normal form" has exactly the properties we want, as > outlined above. I'm still trying to wrap my head around this stuff, so excuse me if this doesn't make sense, but doesn't this change the ordering, and the total degree of each term? > > > Discussion and Future > --------------------- > > It seems to me that this method is quite powerful. I believe the most > important thing now that this proof of concept code exists is to investigate > what it can and cannot do. So please test it thoroughly and report all > surprises :-). I believe we can get a fairly good trigsimp along the > following lines: > > 1. Replace all trig expressions by sin and cos equivalents. > 2. Run mytrigsimp, probably some clever improved version. > 3. Run some clever heuristics to get tan and similar back. > 4. Run the old trigsimp. Why do we still need the old trigsimp? > > I see the following problems/limitations with the code, at least in its > current form: > > 1. It only handles sin and cos. > 2. It quickly gets slow. > 3. Answers can be "funny". > > Let me look at these in the opposite order *g*. > > Regarding (3), I think the solutions ratsimpmodprime can be more complicated > than necessary (i.e. of minimal degree, but it could be possible using fewer > terms). > > For example, consider > > expr = (4*sin(x)*cos(x) + 12*sin(x) + 5*cos(x)**3 + 21*cos(x)**2 + 23*cos(x) > + 15)/(-sin(x)*cos(x)**2 + 2*sin(x)*cos(x) + 15*sin(x) + 7*cos(x)**3 + > 31*cos(x)**2 + 37*cos(x) + 21) I tried this for n=2 and got the same thing as for n=1. For n=3 I got AttributeError: 'NoneType' object has no attribute 'iteritems' (coming from line 231). > > Using n=1 you get something quite short. For n=2 you get something longer > (but of the same total degree). The reason here is that the normal form the > algorithm *starts* with depends on the ideal, and may be nasty. But it is of > degree (1,1), so the algorithm only tries (1, 0) and (0, 1) - neither of > which works. Even if we change it to also try (1,1), the result returned is > still not that nice. The reason is that in the algorithm, when presented > with choice, has currently no heuristics (I'll try to change that). > > [This expr is actually quite fun. Try simplifying various powers with > various choices of n ...] > > Regarding (2), I don't see much to do about that. We could try to precompute > groebner bases, or at least cache them intelligently. This would probably > speed up common cases. As the exrpession gets larger (in particular higher > degree or more variables), the reduction computations also quickly become > slow. Also we could add heuristics, to make things faster in the common case where it would otherwise be slow. Something like the fu et. al. algorithm could be used as a fast preprocessor, for example. What slows it down more, adding multiple independent generators (like sin(x) and sin(y)), or multiple dependent generators (like sin(2*x), ..., sin(10*x))? > > Regarding (1), there are actually some possibilities. ratsimpmodprime is > quite general, so various simplification strategies are possible. The > problem is that we have to produce a prime ideal of polynomial relations, > which is not quite trivial... Is the problem that the equations have to be independent of one another? I tried adding tan(x)*cos(x) - sin(x) to the ideal and (after fixing the above bug), I got a lot of cool results, like >>> mytrigsimp(sin(x).rewrite(tan).rewrite(cos), n=2) sin(x) >>> mytrigsimp(sin(x)/cos(x)) tan(x) >>> mytrigsimp(cos(x)/sin(x)) 1/tan(x) >>> mytrigsimp((cos(x)*sin(x) - cos(x))/sin(x)) (sin(x) + 1)/tan(x) I tried adding tan double angle identities to the generators list by adding tan(k*x).rewrite(cos).expand(trig=True), but I couldn't get just sin(x).rewrite(tan) to simplify (though I could get some other stuff to work, like mytrigsimp(tan(2*x).rewrite(cos).expand(trig=True), n=2)). I'm curious if we can extend the algorithm to work with sum and difference formula (like sin(x + y) = ...). > > > Anyhow, I have spent much longer than I intended already on writing this > mail. Let me know what you think! I think it's awesome! > > > > Best, > Tom By the way, why is mytrigsimp in a separate file? Why not just put it in simplify.py (or create a new trigsimp.py)? Aaron Meurer -- You received this message because you are subscribed to the Google Groups "sympy" group. To post to this group, send email to [email protected]. To unsubscribe from this group, send email to [email protected]. For more options, visit this group at http://groups.google.com/group/sympy?hl=en.
