adding GiNaC method to simplify_rational

[sagetrac-bgoodri]

Currently simplify_rational() only offers 3 Maxima methods. GiNaC offers another possibility via its normal() method. This issue is discussed here

http://groups.google.com/group/sage-devel/browse_thread/thread/843c17dcbd9c2958

I have a patch and a benchmark but need to redownload sage because I am getting unrelated doctest failures with or without the patch.

EDIT: All tests pass now with the attached patch, as they should because the default behavior is not changed. Also, I am attaching a benchmark script using random rational expressions that simplify to 1. In this benchmark, the GiNaC option is about 10 times faster than the default option (Maxima's fullratsimp, without utilizing libraryness).

One limitation of this patch is that it does not support Maxima's map option. GiNaC has a map function, but utilizing it from sage would require a bit more effort.

CC: @williamstein @eviatarbach

Component: symbolics

Author: Ben Goodrich

Issue created by migration from https://trac.sagemath.org/ticket/10268

[sagetrac-bgoodri]

Description changed:

--- 
+++ 
@@ -3,3 +3,7 @@
 http://groups.google.com/group/sage-devel/browse_thread/thread/843c17dcbd9c2958
 
 I have a patch and a benchmark but need to redownload sage because I am getting unrelated doctest failures with or without the patch.
+
+EDIT: All tests pass now with the attached patch, as they should because the default behavior is not changed. Also, I am attaching a benchmark script using random rational expressions that simplify to 1. In this benchmark, the GiNaC option is about 10 times faster than the default option (Maxima's fullratsimp, without utilizing libraryness).
+
+One limitation of this patch is that it does not support Maxima's map option. GiNaC has a map function, but utilizing it from sage would require a bit more effort.

[sagetrac-bgoodri]

comment:1

Attachment: trac_10268_enhance_simplify_rational.patch.gz

[sagetrac-bgoodri]

Attachment: test.ginsh.gz

Attachment: test.sage.gz

[sagetrac-bgoodri]

comment:2

Okay, actually my patch does not work so well for the example that motivated it. I am now attaching a test.ginsh file that defines a rational expression, substitutes a variable with another big rational expression, calls normal and quits. The GiNaC shell finishes in about 1 minute on my laptop (in a shell execute: "time ginsh test.ginsh"). But after applying my patch to sage and then trying to do the equivalent thing via test.sage, it can go for hours without finishing. So, something is badly wrong, possibly my patch.

[kcrisman]

comment:3

I'm not a Cython expert, but maybe should you only use the sig_on business inside the else? It would make sense to start off with method='normal' and then do Ginac, otherwise do Maxima; there is no need to create the GEx unless you are actually going to use it, I would think. I don't know if this would have anything to do with the hang, but it's worth a shot.

[sagetrac-bgoodri]

comment:4

Replying to @kcrisman:

I am just learning sage, but it seems that the compiler does not like the GEx to be declared inside a conditional statement, which makes sense. The _sig_on and _sig_off thing I think is for catching segfaults, which doesn't seem to be a problem and when I comment those out, the behavior is the same. A slight possibility is the fact that when I use the GiNaC shell directly it is the most recent version, whereas Pynac forked off an older version, but the normal function has been in GiNaC for a long, long time.

More interesting is that when I interrupt sage, I get this traceback

KeyboardInterrupt                         Traceback (most recent call last)

/media/disk30/sage-4.6/<ipython console> in <module>()

/media/disk30/sage-4.6/local/lib/python2.6/site-packages/sage/symbolic/expression.so in sage.symbolic.expression.Expression.simplify_rational (sage/symbolic/expression.cpp:23989)()

/media/disk30/sage-4.6/local/lib/python2.6/site-packages/sage/symbolic/pynac.so in sage.symbolic.pynac.py_gcd (sage/symbolic/pynac.cpp:6440)()

/media/disk30/sage-4.6/local/lib/python2.6/site-packages/sage/rings/arith.pyc in gcd(a, b, **kwargs)
   1363 sigma = Sigma() 
   1364 
-> 1365 def gcd(a, b=None, **kwargs):
   1366     r"""
   1367     The greatest common divisor of a and b, or if a is a list and b is

/media/disk30/sage-4.6/local/lib/python2.6/site-packages/sage/interfaces/get_sigs.pyc in my_sigint(x, n)
      7 
      8 def my_sigint(x, n):
----> 9     raise KeyboardInterrupt
     10 
     11 def my_sigfpe(x, n):

KeyboardInterrupt: 

Why would it fall into the gcd function from /media/disk30/sage-4.6/local/lib/python2.6/site-packages/sage/rings/arith.pyc? The patch does not call it directly, and it is a waste because normal in GiNaC already cancels the greatest common factor from the numerator and the denominator. And then a related question is why does gcd seem to hang?

[burcin]

comment:5

Attachment: bench.sage.gz

I'm really happy to see some effort to use pynac/ginac to replace functionality we normally use maxima for. Unfortunately this is a really busy period for me so I can't help much. Thanks a lot for your effort Ben.

Replying to @sagetrac-bgoodri:

Why would it fall into the gcd function from /media/disk30/sage-4.6/local/lib/python2.6/site-packages/sage/rings/arith.pyc? The patch does not call it directly, and it is a waste because normal in GiNaC already cancels the greatest common factor from the numerator and the denominator. And then a related question is why does gcd seem to hang?

The numeric::gcd() method calls sage.symbolic.pynac.py_gcd(). See here:

http://pynac.sagemath.org/hg/file/b233d9dadcfa/ginac/numeric.cpp#l2526

It could be that our gcd() function doesn't work exactly like CLN's gcd() function which is used originally by ginac. This would effect the termination criteria used in the multivariate gcd code in ginac/pynac.

I haven't looked into the functionality in normal.cpp much, but one of William's goals was to make it call Singular (or the Factory library) to factor multivariate polynomials instead of the code in ginac. This library generally performs much better and it is actively being developed.

---

BTW, kcrisman was right about his comment on the use of sig_on/sig_off. These functions are used to override the signal handlers so we can catch CTRL-C in long running library code. You don't need them around the declaration of GEX, but you should use them in the call to normal(). AFAIR, Jeroen Demeyer recently wrote a section on this in the developers guide, but I don't have a link handy right now. The new call signature for these might be sig_on(), like a function call.

[burcin]

Author: Ben Goodrich

[sagetrac-bgoodri]

comment:6

Replying to @burcin: cc'ing William for clarification

I haven't looked into the functionality in normal.cpp much, but one of William's goals was to make it call Singular (or the Factory library) to factor multivariate polynomials instead of the code in ginac. This library generally performs much better and it is actively being developed.

Do you really want to do full factorization in simplify_rational(). I think neither Maxima nor GiNaC do that, only square-free factorization and gcd cancellation. We could add an option to do full factorization of the numerator and denominator before returning. If so, would it make more sense to first backport the functionality in factor.cpp from GiNaC 1.5.x to the pynac fork than to code a pynac-libSingular link?

Also, I think this might be a bit separate from the issue I was hitting. When I ran test.sage last night under trace() with the enter key wedged down, by the morning it had called gcd() over 30,000 times and hadn't even passed the rational expression to GiNaC yet. This is a waste because GiNaC's normal() function was going to do 1 gcd cancellation anyway. So, it seems what we need is an option to prevent sage from trying to find the gcd of every subexpression.