Problem doing symbolic computations (bug in Pynac ?)

36 views
Skip to first unread message

Jean-Pierre Flori

unread,
Sep 1, 2010, 11:25:07 AM9/1/10
to sage-support
Hello,

I'm using Sage to do some symbolic computations.
What I am basically trying to do is to compute sums made of terms of
the form :
cst_{i_0,...,i_d)*b_0^{i_0}*...*b_d^{i_d}*2^{-b_0-...-b_d}
where the b_i are variables
(i.e. a product of an exponential with b_i in the exponent by a
multivariate polynomial in b_i)
using functions (especially symbolic_sum) from sage.calculus.calculus.
I am regularly calling expand() and simplify() functions as well as
some custom functions using substitute() to get a more compact
expression and avoid redundant terms.
Everything works fine as long as the sums are not too "complicated",
but Sage segfaults when they become too big.
Running Sage with the -gdb flag, gives different errors according to
which of my functions I am calling, but both errors seem to be related
to a "compare" function in Pynac (SIGSEGV or SIGBUS):

First error message :
Program received signal SIGSEGV, Segmentation fault.
GiNaC::power::compare (this=0x3281c20, other=...) at power.cpp:899
899 power.cpp: Aucun fichier ou dossier de ce type.
in power.cpp

Second error message :
Program received signal SIGBUS, Bus error.
GiNaC::mul::compare (this=0x83490e0, other=...) at mul.cpp:1174
1174 mul.cpp: Aucun fichier ou dossier de ce type.
in mul.cpp


My box is running a Debian unstable with packages from the
experimental repository.
My kernel version is 2.6.35, compiled by myself.
My cpu is a Intel(R) Core(TM)2 Quad CPU Q6600.
The error is happening with Sage 4.5.2 and 4.4.4.
I compiled both of them myself with gcc (Debian 4.4.4-11) 4.4.5
20100824 (prerelease).
I also tried the precompiled sage-4.5.2-linux-64bit-ubuntu_10.04_lts-
x86_64-Linux which gives me the same result.

I didn't manage to reproduce the bug with a short and simple piece of
code.
However I can mail the source code and instructions to reproduce the
bug to anyone who can help.
I will also provide any information needed to investigate this
problem.
I am currently trying to compile Sage on another box to do another
test.
Any help will be appreciated, and thanks for a such a great piece of
software.

Best regards,

Jean-Pierre Flori

Burcin Erocal

unread,
Sep 1, 2010, 12:01:46 PM9/1/10
to sage-s...@googlegroups.com
Hi Jean-Pierre,

<snip>


> I didn't manage to reproduce the bug with a short and simple piece of
> code.
> However I can mail the source code and instructions to reproduce the
> bug to anyone who can help.
> I will also provide any information needed to investigate this
> problem.

Thanks for the report.

If you send me your example, I can take a look at it tomorrow.

Cheers,
Burcin

Burcin Erocal

unread,
Sep 7, 2010, 5:26:39 PM9/7/10
to sage-s...@googlegroups.com, pynac...@googlegroups.com
Hi,

Here is a short example to replicate the first error mentioned below:

b = [var('b_%s'%i) for i in range(4)]

precomp = (2^b_2 + 2)*(2^b_1 + 2^(-b_1) + 2^b_1*2^b_0 - 2^b_1*2^(-b_0)
- 2^(-b_1)*2^b_0 - 2^(-b_1)*2^(-b_0) + 2^b_0 + 2^(-b_0) - 9) + (2^b_1 +
2^(-b_1) + 2^b_1*2^b_0 - 2^b_1*2^(-b_0) - 2^(-b_1)*2^b_0 -
2^(-b_1)*2^(-b_0) + 2^b_0 + 2^(-b_0) - 9)/2^b_2

repl_dict = {b_0: b_0, b_3: b_1, b_2: b_3, b_1: b_2}
P = precomp.substitute(repl_dict)
P.expand()


I will take a break now. I'd appreciate any help tracking down the
problem (in pynac) if anybody has the time.


Cheers,
Burcin


On Wed, 1 Sep 2010 08:25:07 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

kcrisman

unread,
Sep 7, 2010, 8:52:48 PM9/7/10
to sage-support


On Sep 7, 5:26 pm, Burcin Erocal <bur...@erocal.org> wrote:
> Hi,
>
> Here is a short example to replicate the first error mentioned below:
>
> b = [var('b_%s'%i) for i in range(4)]
>
> precomp = (2^b_2 + 2)*(2^b_1 + 2^(-b_1) + 2^b_1*2^b_0 - 2^b_1*2^(-b_0)
> - 2^(-b_1)*2^b_0 - 2^(-b_1)*2^(-b_0) + 2^b_0 + 2^(-b_0) - 9) + (2^b_1 +
> 2^(-b_1) + 2^b_1*2^b_0 - 2^b_1*2^(-b_0) - 2^(-b_1)*2^b_0 -
> 2^(-b_1)*2^(-b_0) + 2^b_0 + 2^(-b_0) - 9)/2^b_2
>
> repl_dict = {b_0: b_0, b_3: b_1, b_2: b_3, b_1: b_2}
> P = precomp.substitute(repl_dict)
> P.expand()
>
> I will take a break now. I'd appreciate any help tracking down the
> problem (in pynac) if anybody has the time.
>

Did you open a ticket for this?

- kcrisman

Jean-Pierre Flori

unread,
Sep 8, 2010, 4:04:48 AM9/8/10
to sage-support
Hi,

Thanks a lot for producing a small piece of code reproducing the bug !
I'll try to track something down from here with gdb.

By the way, the above code doesn't *always* produce a bug (in fact it
didn't the first time I tried it...).
You can verify it by putting it in a "bug.sage" file and running the
following bash script:

[jp@napoleon sage-current]\$ cat bug.sh
#!/bin/bash
for (( i=1; i<$1; i++ ))
do
./sage bug.sage &> /dev/null
if [ $? = 0 ]
then
echo $i
fi
done

It gave the following output:

[jp@napoleon sage-current]\$ ./bug.sh 100
48
57
75
79
83

then:

[jp@napoleon sage-current]\$ ./bug.sh 100
4
6
7
26
28
68
69
99


So there must be something strange going on.

>
> Did you open a ticket for this?
>

I did not.
I took a look on the recent tickets in Trac, and also used the search
function for Pynac, but could not find anything relevant.
So I guess Burcin did not either. I'll do it today.

> - kcrisman

Jean-Pierre Flori

unread,
Sep 9, 2010, 5:13:27 AM9/9/10
to sage-support
I created a Ticket on Sage's Trac:
http://trac.sagemath.org/sage_trac/ticket/9880

It could be a problem with the comparison function
"expair_rest_is_less" used by "std::sort" function like in the
following link:
http://gcc.gnu.org/ml/gcc-bugs/2010-08/msg01163.html

The backtrace of our bug look similar:
#0 GiNaC::power::compare (this=0x449be20, other=...) at power.cpp:899
#1 0x00007fffd7a9cc18 in void
std::__unguarded_linear_insert<__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair,
GiNaC::expair_rest_is_less>(__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair, GiNaC::expair_rest_is_less) ()
from /home/jp/boulot/thèse/sage/sage-4.5.2/local/lib/
libpynac-0.2.so.0
#2 0x00007fffd7a9cefa in void
std::__final_insertion_sort<__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair_rest_is_less>(__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair_rest_is_less) () from /home/jp/boulot/thèse/sage/
sage-4.5.2/local/lib/libpynac-0.2.so.0
#3 0x00007fffd7a95657 in
sort<__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair_rest_is_less> (this=<value optimized out>) at /usr/
include/c++/4.4/bits/stl_algo.h:5260
#4 GiNaC::expairseq::canonicalize (this=<value optimized out>) at
expairseq.cpp:1177
#5 0x00007fffd7a994a4 in GiNaC::expairseq::construct_from_epvector
(this=0x44a5070, v=<value optimized out>,
do_index_renaming=<value optimized out>) at expairseq.cpp:1055
#6 0x00007fffd7a58685 in GiNaC::add::add (this=0x44a5070, vp=...,
oc=<value optimized out>) at add.cpp:100
#7 0x00007fffd7a5871f in GiNaC::add::expand (this=0x449eb80,
options=<value optimized out>) at add.cpp:669
#8 0x00007fffd7a9231d in GiNaC::ex::expand (this=<value optimized
out>, options=3620337048) at ex.cpp:78
#9 0x00007fffd773e386 in
__pyx_pf_4sage_8symbolic_10expression_10Expression_expand
(__pyx_v_self=<value optimized out>,
__pyx_args=0x601160, __pyx_kwds=0x0) at sage/symbolic/
expression.cpp:13604
#10 0x00007ffff7b107ca in call_function (f=0x44a4580, throwflag=<value
optimized out>) at Python/ceval.c:3706
#11 PyEval_EvalFrameEx (f=0x44a4580, throwflag=<value optimized out>)
at Python/ceval.c:2389
#12 0x00007ffff7b124ad in PyEval_EvalCodeEx (co=0x79c4c0,
globals=<value optimized out>, locals=<value optimized out>,
args=0x0,
argcount=<value optimized out>, kws=0x0, kwcount=0, defs=0x0,
defcount=0, closure=0x0) at Python/ceval.c:2968
#13 0x00007ffff7b12582 in PyEval_EvalCode (co=0x449be20,
globals=0x1135200007879, locals=0x7fffd7c9f598) at Python/ceval.c:522
#14 0x00007ffff7b3014c in run_mod (
.
.
.

And valgrind output also:
==27910== Invalid read of size 8
==27910== at 0x282C4BF1: void
std::__unguarded_linear_insert<__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair,
GiNaC::expair_rest_is_less>(__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair, GiNaC::expair_rest_is_less) (ptr.h:99)
==27910== by 0x282C4EF9: void
std::__final_insertion_sort<__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair_rest_is_less>(__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
__gnu_cxx::__normal_iterator<GiNaC::expair*,
std::vector<GiNaC::expair, std::allocator<GiNaC::expair> > >,
GiNaC::expair_rest_is_less) (stl_algo.h:2161)
==27910== by 0x282BD656: GiNaC::expairseq::canonicalize()
(stl_algo.h:5260)
==27910== by 0x282C14A3:
GiNaC::expairseq::construct_from_epvector(std::vector<GiNaC::expair,
std::allocator<GiNaC::expair> > const&, bool) (expairseq.cpp:1055)
==27910== by 0x28280684:
GiNaC::add::add(std::auto_ptr<std::vector<GiNaC::expair,
std::allocator<GiNaC::expair> > >, GiNaC::ex const&) (add.cpp:100)
==27910== by 0x2828071E: GiNaC::add::expand(unsigned int) const
(add.cpp:669)
==27910== by 0x282BA31C: GiNaC::ex::expand(unsigned int) const
(ex.cpp:78)
==27910== by 0x28811385:
__pyx_pf_4sage_8symbolic_10expression_10Expression_expand(_object*,
_object*, _object*) (expression.cpp:13604)
==27910== by 0x4F137C9: PyEval_EvalFrameEx (ceval.c:3706)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F15581: PyEval_EvalCode (ceval.c:522)
==27910== by 0x4F3314B: PyRun_StringFlags (pythonrun.c:1335)
==27910== by 0x4F122DD: PyEval_EvalFrameEx (ceval.c:4437)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4E9BDDE: function_call (funcobject.c:524)
==27910== by 0x4E6FAA2: PyObject_Call (abstract.c:2492)
==27910== by 0xD7E0981:
__pyx_pf_4sage_9structure_11sage_object_load (sage_object.c:7304)
==27910== by 0x4F137C9: PyEval_EvalFrameEx (ceval.c:3706)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F15581: PyEval_EvalCode (ceval.c:522)
==27910== by 0x4F14853: PyEval_EvalFrameEx (ceval.c:4401)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F13844: PyEval_EvalFrameEx (ceval.c:3802)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F13844: PyEval_EvalFrameEx (ceval.c:3802)
==27910== Address 0x97e4280 is 16 bytes before a block of size 416
alloc'd
==27910== at 0x4C24CCC: operator new(unsigned long)
(vg_replace_malloc.c:220)
==27910== by 0x28286727: std::vector<GiNaC::expair,
std::allocator<GiNaC::expair> >::reserve(unsigned long)
(new_allocator.h:89)
==27910== by 0x282C0D4A:
GiNaC::expairseq::make_flat(std::vector<GiNaC::expair,
std::allocator<GiNaC::expair> > const&, bool) (expairseq.cpp:1138)
==27910== by 0x282C149B:
GiNaC::expairseq::construct_from_epvector(std::vector<GiNaC::expair,
std::allocator<GiNaC::expair> > const&, bool) (expairseq.cpp:1051)
==27910== by 0x28280684:
GiNaC::add::add(std::auto_ptr<std::vector<GiNaC::expair,
std::allocator<GiNaC::expair> > >, GiNaC::ex const&) (add.cpp:100)
==27910== by 0x2828071E: GiNaC::add::expand(unsigned int) const
(add.cpp:669)
==27910== by 0x282BA31C: GiNaC::ex::expand(unsigned int) const
(ex.cpp:78)
==27910== by 0x28811385:
__pyx_pf_4sage_8symbolic_10expression_10Expression_expand(_object*,
_object*, _object*) (expression.cpp:13604)
==27910== by 0x4F137C9: PyEval_EvalFrameEx (ceval.c:3706)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F15581: PyEval_EvalCode (ceval.c:522)
==27910== by 0x4F3314B: PyRun_StringFlags (pythonrun.c:1335)
==27910== by 0x4F122DD: PyEval_EvalFrameEx (ceval.c:4437)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4E9BDDE: function_call (funcobject.c:524)
==27910== by 0x4E6FAA2: PyObject_Call (abstract.c:2492)
==27910== by 0xD7E0981:
__pyx_pf_4sage_9structure_11sage_object_load (sage_object.c:7304)
==27910== by 0x4F137C9: PyEval_EvalFrameEx (ceval.c:3706)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F15581: PyEval_EvalCode (ceval.c:522)
==27910== by 0x4F14853: PyEval_EvalFrameEx (ceval.c:4401)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F13844: PyEval_EvalFrameEx (ceval.c:3802)
==27910== by 0x4F154AC: PyEval_EvalCodeEx (ceval.c:2968)
==27910== by 0x4F13844: PyEval_EvalFrameEx (ceval.c:3802)

GiNaC documentation also mentions that this comparison function is not
a "strict weak ordering" as required by "std::sort":
"Note that this does not define a strict weak ordering since for any
symbol x we have neither 3*x<2*x or 2*x<3*x. Handle with care!"
but that function is made to be used that way so...
Moreover replacing it by "expair_is_less" do not change anything.

However, something fishy is going on here.

Burcin Erocal

unread,
Sep 12, 2010, 9:40:29 AM9/12/10
to sage-s...@googlegroups.com, pynac...@googlegroups.com
Hi Jean-Pierre,

On Thu, 9 Sep 2010 02:13:27 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

> I created a Ticket on Sage's Trac:
> http://trac.sagemath.org/sage_trac/ticket/9880
>
> It could be a problem with the comparison function
> "expair_rest_is_less" used by "std::sort" function like in the
> following link:
> http://gcc.gnu.org/ml/gcc-bugs/2010-08/msg01163.html

This is indeed the problem. Thanks for tracking it down.

In the process of adapting GiNaC for use in Sage, I changed the
comparison functions to give a consistent order between different runs.
Apparently, the new comparison functions I wrote are not consistent.

For a long time, I wanted to rip out these new functions, and use them
only for printing. The following patch is a first step towards this:

http://sage.math.washington.edu/home/burcin/pynac/pynac_order.patch

After applying this patch to pynac (see here [1] for instructions on
how to set up a pynac development environment), the segfault you
reported, and the bug at #9046 goes away.

Note that this is very much a work in progress. I didn't take care to
make sure these comparison functions are consistent either. In the
process of copying things over, I even flipped a few signs in the
results.

This needs much work to be finished, we need to make sure the order
is consistent and correct, and it matches the one currently used by
Sage. The issue of accessing the operands of a symbolic expression also
needs to be sorted out, but that's merely a user interface issue.

I'd appreciate help with any of this, since I don't have much time to
spend on Sage/pynac these days. I'd be glad to answer any questions
though.


Cheers,
Burcin

Jean-Pierre Flori

unread,
Sep 12, 2010, 12:21:12 PM9/12/10
to sage-support
Hi,

Thanks for the patch, I'll have a look at the comparison functions
when I have some time, surely next week.
Where can I find the one defined by Sage ?

I have another question: are the source codes of GiNaC and pynac still
"synchronized" in some way ?
I saw that some recent changes are similar in both git repositories,
but as I was going through the source code, I saw that
numeric::do_print_csrc is defined in GiNaC but not in pynac.

Regards,

On 12 sep, 15:40, Burcin Erocal <bur...@erocal.org> wrote:
> Hi Jean-Pierre,
>
> On Thu, 9 Sep 2010 02:13:27 -0700 (PDT)
>

Burcin Erocal

unread,
Sep 12, 2010, 12:36:26 PM9/12/10
to sage-s...@googlegroups.com, pynac...@googlegroups.com
Hi,

On Sun, 12 Sep 2010 09:21:12 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

> Thanks for the patch, I'll have a look at the comparison functions
> when I have some time, surely next week.
> Where can I find the one defined by Sage ?

They are spread out through the classes in pynac. Look for compare()
and compare_same_type() functions in mul, power, function, etc.

> I have another question: are the source codes of GiNaC and pynac still
> "synchronized" in some way ?

Not really. I pull some patches from GiNaC when I think they are
relevant.

> I saw that some recent changes are similar in both git repositories,
> but as I was going through the source code, I saw that
> numeric::do_print_csrc is defined in GiNaC but not in pynac.

The numeric class in pynac is completely rewritten to use python
objects instead of CLN. It's normal that those methods are missing.

Cheers,
Burcin

Jean-Pierre Flori

unread,
Sep 13, 2010, 7:58:00 AM9/13/10
to sage-support
Hi Burcin,

Thanks for your quick answer!

I've got some other questions:
- Sage and pynac do not realize that 2^(-b_0) and (2^b_0)^(-1) are
equal.
I guess there is no canonical way of expanding 2^(a*b) into something
else, but it could be a good idea doing it when a or b is an integer
and the other one a variable ?
- Should symbolic sums be implemented into pynac (at some point...) to
avoid calling Maxima ?
By the way calling symbolic_sum(1,x,0,b_0) gives me an error whereas
calling sum(1,x,0,b_0) does not.

Cheers,
JP

Burcin Erocal

unread,
Sep 13, 2010, 8:38:26 AM9/13/10
to sage-s...@googlegroups.com
Hi Jean-Pierre,

On Mon, 13 Sep 2010 04:58:00 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

> I've got some other questions:

These should be on a separate thread, on sage-devel or pynac-devel.

> - Sage and pynac do not realize that 2^(-b_0) and (2^b_0)^(-1) are
> equal.
> I guess there is no canonical way of expanding 2^(a*b) into something
> else, but it could be a good idea doing it when a or b is an integer
> and the other one a variable ?

I don't understand your suggestion. Are you saying we should
* expand 2^(-b) to (2^b)^(-1) or
* simplify (2^b)^(-1) to 2^(-b)?

Apart from a few simple modifications, e.g., exp(a)^b -> exp(a*b) when
possible, pynac keeps the automatic evaluation rules from GiNaC.

These operations effect the performance of the library quite a bit. We
should be careful before adding more of them, and perhaps consult the
GiNaC developers in the process.

> - Should symbolic sums be implemented into pynac (at some point...) to
> avoid calling Maxima ?

They should be implemented in Sage. I don't think this requires any
modification at the C++ level in the pynac library.

Actually, I am working on this right now.

> By the way calling symbolic_sum(1,x,0,b_0) gives me an error whereas
> calling sum(1,x,0,b_0) does not.

This is the difference between sage.misc.functional.symbolic_sum and
sage.calculus.calculus.symbolic_sum. The latter tries to convert the
output received from maxima to the parent of the given expression. In
your example, this is ZZ. The conversion fails of course. The wrapper
sage.misc.functional.symbolic_sum first converts the first argument to
a symbolic expression, which prevents the error.

When we move to new code in sage/symbolic/summation/* (which is not in
the library yet), things will be much cleaner.


Cheers,
Burcin

kcrisman

unread,
Sep 13, 2010, 10:35:38 AM9/13/10
to sage-support

>
> When we move to new code in sage/symbolic/summation/* (which is not in
> the library yet), things will be much cleaner.

Can you give a ticket # for this? I would be interested in looking at
this.

- kcrisman

Jean-Pierre Flori

unread,
Sep 29, 2010, 10:00:39 AM9/29/10
to sage-support
Ok, I have finally looked at the comparison functions and exchanging :
cmpval = seq[0].coeff.compare(other.exponent);
by
cmpval = -seq[0].coeff.compare(other.exponent);
in mul::compare_pow (mul.cpp:1265) seems to prevent the above bug from
happening.
It seems to fit better with the change made by William Stein in
power::compare_same_type (power.cpp:951).
However it doesn't mean the problem is completely solved...
I'll try to take a deeper look at the comparison functions at some
point.

By the way I opened another thread on sage-devel:
http://groups.google.com/group/sage-devel/browse_thread/thread/41e3f276cdf8d08d#

Regards.

Burcin Erocal

unread,
Sep 29, 2010, 11:33:54 AM9/29/10
to sage-s...@googlegroups.com
Hi,

On Wed, 29 Sep 2010 07:00:39 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

> Ok, I have finally looked at the comparison functions and exchanging :
> cmpval = seq[0].coeff.compare(other.exponent);
> by
> cmpval = -seq[0].coeff.compare(other.exponent);
> in mul::compare_pow (mul.cpp:1265) seems to prevent the above bug from
> happening.
> It seems to fit better with the change made by William Stein in
> power::compare_same_type (power.cpp:951).
> However it doesn't mean the problem is completely solved...
> I'll try to take a deeper look at the comparison functions at some
> point.

Thank you for following up on this and taking the time to go through
the comparisons.

I was planning to build on the patch I posted before [1] to fix this
problem once and for all (perhaps rather optimistically).

[1] http://sage.math.washington.edu/home/burcin/pynac/pynac_order.patch

Note that your suggestions fixes your bug (#9880), but doesn't effect
#9046.


Cheers,
Burcin

Jean-Pierre Flori

unread,
Sep 29, 2010, 12:48:25 PM9/29/10
to sage-support
Hi,

Maybe it is a good thing to keep the same order as ginac internally
and your more usual ordering for printing.
However if you'd better not duplicate code, I can look at the "-
x^2+x^2" part of bug #9046.
Now I may understand a big enough part of pynac code to do that.
But if you'd better use the above patch, that won't be necessary.

Cheers,

On 29 sep, 17:33, Burcin Erocal <bur...@erocal.org> wrote:
> Hi,
>
> On Wed, 29 Sep 2010 07:00:39 -0700 (PDT)
>

Jean-Pierre Flori

unread,
Sep 29, 2010, 1:43:12 PM9/29/10
to sage-support
My first guess for #9046 is that the problem is located within the
last line of collect : return x + (*this-x).expand();
In that case (*this-x).expand() should be zero but looks like:
-x^2 - 2*y^3 - y^2*z + x^2 + 2*y^3 + y^2*z
I think that whithin the call to expand(),
expairseq::combine_same_terms_sorted_seq() is called after
canonicalize(); but with the new order x^2 and -x^2, etc. are not
adjacent and do not get simplified... whence the above expression for
(*this-x).expand() and finally for p.collect(x).
I'll look at the ordering problem causing that tomorrow.

Burcin Erocal

unread,
Sep 29, 2010, 2:28:14 PM9/29/10
to sage-s...@googlegroups.com, pynac...@googlegroups.com
On Wed, 29 Sep 2010 09:48:25 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

> Maybe it is a good thing to keep the same order as ginac internally
> and your more usual ordering for printing.

It is good to keep the ginac ordering internally. The user friendly
ordering is more expensive so it slows down computations, even when
you don't print intermediate expressions.

> However if you'd better not duplicate code, I can look at the "-
> x^2+x^2" part of bug #9046.
> Now I may understand a big enough part of pynac code to do that.

That's right. You're one of the few people who spent so much time
staring at pynac code. :)

> But if you'd better use the above patch, that won't be necessary.

I would like to replace the orders completely. I wanted to attack
this problem for a long time. Now that we started, we might as
well finish. :)

Can you help put the patch in usable shape? At the moment the ordering
it defines is completely different from the one we use in Sage. This
must be because I changed some signs while I was copying code.

If you can make the ordering consistent, I can fix the other relevant
places (op(), etc.).


Thanks.

Burcin

Jean-Pierre Flori

unread,
Sep 29, 2010, 3:57:54 PM9/29/10
to sage-support


On 29 sep, 20:28, Burcin Erocal <bur...@erocal.org> wrote:
> On Wed, 29 Sep 2010 09:48:25 -0700 (PDT)
>
> Jean-Pierre Flori <jpfl...@gmail.com> wrote:
> > Maybe it is a good thing to keep the same order as ginac internally
> > and your more usual ordering for printing.
>
> It is good to keep the ginac ordering internally. The user friendly
> ordering is more expensive so it slows down computations, even when
> you don't print intermediate expressions.
My bad, I didn't read your previous post correctly, I thought that
patch reverted the change you made in the ordering, but in fact it put
these changes in separate functions and use them only for printing...
That is why I tried to fix the bug in the code without your patch...
So I guess I should work WITH the patch !
>
> > However if you'd better not duplicate code, I can look at the "-
> > x^2+x^2" part of bug #9046.
> > Now I may understand a big enough part of pynac code to do that.
>
> That's right. You're one of the few people who spent so much time
> staring at pynac code. :)
>
> > But if you'd better use the above patch, that won't be necessary.
>
> I would like to replace the orders completely. I wanted to attack
> this problem for a long time. Now that we started, we might as
> well finish. :)
So once more I'm confused..
What do you mean by replace the orders completely ?
Having separate functions for internal use and printing as with what
your patch does ?
Or continue from the code before the above patch and make everything
back consistent so just one set of functions is present, even if it
can finally produce slower code ?
I guess that is because I was already confused before and what I said
may have seemed obscure..

>
> Can you help put the patch in usable shape? At the moment the ordering
> it defines is completely different from the one we use in Sage. This
> must be because I changed some signs while I was copying code.
I'll help with pleasure. Just want to be sure to take the right
direction !

>
> If you can make the ordering consistent, I can fix the other relevant
> places (op(), etc.).
>
> Thanks.
>
> Burcin

Cheers,

Burcin Erocal

unread,
Sep 29, 2010, 5:24:39 PM9/29/10
to sage-s...@googlegroups.com
On Wed, 29 Sep 2010 12:57:54 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

> On 29 sep, 20:28, Burcin Erocal <bur...@erocal.org> wrote:
> > On Wed, 29 Sep 2010 09:48:25 -0700 (PDT)
> >
> > Jean-Pierre Flori <jpfl...@gmail.com> wrote:
> > > Maybe it is a good thing to keep the same order as ginac
> > > internally and your more usual ordering for printing.
> >
> > It is good to keep the ginac ordering internally. The user friendly
> > ordering is more expensive so it slows down computations, even when
> > you don't print intermediate expressions.
> My bad, I didn't read your previous post correctly, I thought that
> patch reverted the change you made in the ordering, but in fact it put
> these changes in separate functions and use them only for printing...
> That is why I tried to fix the bug in the code without your patch...
> So I guess I should work WITH the patch !

Yes, please work with the patch. :)

<snip>


> > > But if you'd better use the above patch, that won't be necessary.
> >
> > I would like to replace the orders completely. I wanted to attack
> > this problem for a long time. Now that we started, we might as
> > well finish. :)
> So once more I'm confused..
> What do you mean by replace the orders completely ?
> Having separate functions for internal use and printing as with what
> your patch does ?

<snip>

I want to continue with the patch, to go back to using the original
comparison functions from ginac and use the new ones only for printing.


> >
> > Can you help put the patch in usable shape? At the moment the
> > ordering it defines is completely different from the one we use in
> > Sage. This must be because I changed some signs while I was copying
> > code.
> I'll help with pleasure. Just want to be sure to take the right
> direction !

Thank you. Please don't hesitate to ask if you need any help.

Cheers,
Burcin

Jean-Pierre Flori

unread,
Oct 1, 2010, 9:52:11 AM10/1/10
to sage-support
Hi,

Some remarks and questions...

I guess the order you want is 'degrevlex' as Sage's default order for
multivariate polynomial ring and as the name of your functions
suggest.
Am I wrong ?

At present, I did not modify a lot of things.
A few signs, replaced seq[0] by the smallest item of the sequence
sorted with the printing order.

I guess we want symbolic variables to be bigger than numerics.
So I changed the sign for different typeid's.

What is your opinion about symbolic exponents ?
(If I am not wrong, in mul object, they are encapsulated in power
object so it counts for 1.
So total_degree(a^b*x) == 2.)
In the compare_mul_power, you decided to compare the exponent of the
power object to the one of the smallest item (which could in fact be
the same power object) of the mul object if it is symbolic.
So 3*x^y will be smaller than x^y (considering that y > 1).
That is not a "real" problem because 3*x^y+x^y will be simplified into
4*x^y.
But we also get that a^b > a^b*x (if b > 1) with the current strategy
and it prints a^b + a^b*x...
Shouldn't we consider the degree of the power is 1, so at least the
above case is solved.
Or do something more complicated comparing exponents...

I do not understand the comment of William Stein about symbols
ordering.
When a polynomial ring is created in Sage, the variables are ordered
according to the order they were "created", i.e.
PR.<x,y> = PolynomialRing(QQ) will result in x > y
and PR.<y,x> = PolynomialRing(QQ) in y > x

Burcin Erocal

unread,
Oct 1, 2010, 3:02:41 PM10/1/10
to sage-s...@googlegroups.com
On Fri, 1 Oct 2010 06:52:11 -0700 (PDT)
Jean-Pierre Flori <jpf...@gmail.com> wrote:

> Some remarks and questions...
<snip>

I replied to this on pynac-devel:

http://groups.google.com/group/pynac-devel/msg/f3032856d1ed6953


Cheers,
Burcin

Reply all
Reply to author
Forward
0 new messages