factor() hangs on a simple expression, randomly
Hi, folks. I have this problem: $ ginsh ginsh - GiNaC Interactive Shell (GiNaC V1.7.2) __, _______ Copyright (C) 1999-2017 Johannes Gutenberg University Mainz, (__) * | Germany. This is free software with ABSOLUTELY NO WARRANTY. ._) i N a C | You are welcome to redistribute it under certain conditions. <-------------' For details type `warranty;'. Type ?? for a list of help topics. > factor((2*x-l)^2*(-1+6*x-l)*(-1+3*x-l)*(-2+4*x-l)*(1+4*x-l)*(4*x-l)*(2+4*x-l)); At this point ginsh either hangs (80% of cases), or gives the correct result (20% of cases). As you can see the expression is already factored, so a hang is pretty unexpected here. Note that if you'll drop two of the factors from that expression, GiNaC won't hang anymore. Here are the cumulative timings you'll typically see via 'perf top' during this hang: 91.54% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::factor_multivariate 81.32% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::hensel_multivar 79.67% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::multivar_diophant 79.67% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::univar_diophant 79.67% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::multiterm_eea_lift 72.77% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::factor_sqrfree 66.18% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::eea_lift 64.28% libginac.so.6.0.2 [.] GiNaC::factor 26.99% libginac.so.6.0.2 [.] GiNaC::(anonymous namespace)::operator* 23.67% ginsh [.] yyparse 23.67% ginsh [.] f_factor 16.91% libginac.so.6.0.2 [.] cln::cl_MI::~cl_MI 15.26% libginac.so.6.0.2 [.] cln::cl_I::~cl_I 13.66% libginac.so.6.0.2 [.] cln::cl_MI::cl_MI 13.55% libginac.so.6.0.2 [.] cln::cl_RA::~cl_RA So, any ideas on the nature of this problem, and how to fix it?
Hi Vitaly, On 11/30/2017 05:16 PM, Vitaly Magerya wrote:
$ ginsh ginsh - GiNaC Interactive Shell (GiNaC V1.7.2) __, _______ Copyright (C) 1999-2017 Johannes Gutenberg University Mainz, (__) * | Germany. This is free software with ABSOLUTELY NO WARRANTY. ._) i N a C | You are welcome to redistribute it under certain conditions. <-------------' For details type `warranty;'.
Type ?? for a list of help topics. > factor((2*x-l)^2*(-1+6*x-l)*(-1+3*x-l)*(-2+4*x-l)*(1+4*x-l)*(4*x-l)*(2+4*x-l));
At this point ginsh either hangs (80% of cases), or gives the correct result (20% of cases).
As you can see the expression is already factored, so a hang is pretty unexpected here. Note that if you'll drop two of the factors from that expression, GiNaC won't hang anymore.
I'm not an expert in factorization, but apparently multivariate Hensel lifting keeps failing (in file factor.cpp). You're right that this is unexpected. GiNaC first does a square-free factorization and finds the terms (2*x-l)^2 and the expanded form of all the rest, i.e. of (-1+6*x-l)*(-1+3*x-l)*(-2+4*x-l)*(1+4*x-l)*(4*x-l)*(2+4*x-l). The factorization of the rest then fails, depending on the order it finds the symbols x and l. Try the easier test case with the (2*x-l)^2 removed. No idea why it doesn't factor the product term by term, first. I suppose this is an easy and risk-free optimization! Can you try to send a patch? I'ld be happy to review it and help, but have not the time to do it myself, I'm afraid. All my best, -richy. -- Richard B. Kreckel <https://in.terlu.de/~kreckel/>
On 12/08/2017 11:22 PM, Richard B. Kreckel wrote:
You're right that this is unexpected. GiNaC first does a square-free factorization and finds the terms (2*x-l)^2 and the expanded form of all the rest, i.e. of (-1+6*x-l)*(-1+3*x-l)*(-2+4*x-l)*(1+4*x-l)*(4*x-l)*(2+4*x-l). The factorization of the rest then fails, depending on the order it finds the symbols x and l.
Yes, this expand()ing seems unfortunate.
No idea why it doesn't factor the product term by term, first. I suppose this is an easy and risk-free optimization!
One potential risk I see is if you factor() each of the pre-existing factors separately, then identical factors of each part may not get re-combined properly. Not that current factor() ever provided such a guarantee. I also can't seem to find a good counter-example, so maybe it's a non-issue.
Can you try to send a patch? I'ld be happy to review it and help, but have not the time to do it myself, I'm afraid.
I'm attaching a brute force approach, see if you like it. (Does this list accept attachments?).
Hi Vitaly, On 12/13/2017 05:05 PM, Vitaly Magerya wrote:
One potential risk I see is if you factor() each of the pre-existing factors separately, then identical factors of each part may not get re-combined properly. Not that current factor() ever provided such a guarantee. I also can't seem to find a good counter-example, so maybe it's a non-issue.
The automatic .eval() thing should take care of this situation and re-combine identical factors.
Can you try to send a patch? I'ld be happy to review it and help, but have not the time to do it myself, I'm afraid.
I'm attaching a brute force approach, see if you like it. (Does this list accept attachments?).
Thanks a lot for the patch! I see that the comment of the factor() function is now outdated. Would you like to fix the comment before I push the patch? All my best, -richy. -- Richard B. Kreckel <https://in.terlu.de/~kreckel/>
Hi, On 12/18/2017 07:32 PM, Richard B. Kreckel wrote:
I've pushed your patch. Thanks a lot! Summary: Vitaly's original problem is solved now.
But still, factor(expand((-1+6*x-l)*(-1+3*x-l)*(-2+4*x-l)*(-1+4*x-l)*(2+4*x-l))); does not finish in half of the cases. Hensel lifting just keeps failing. Anyone with better grasp of this is more than welcome to suggest a fix! Bye, -richy. -- Richard B. Kreckel <https://in.terlu.de/~kreckel/>
Hi Vitaly, On 11/30/2017 05:16 PM, Vitaly Magerya wrote:
> factor((2*x-l)^2*(-1+6*x-l)*(-1+3*x-l)*(-2+4*x-l)*(1+4*x-l)*(4*x-l)*(2+4*x-l));
At this point ginsh either hangs (80% of cases), or gives the correct result (20% of cases).
[...]
So, any ideas on the nature of this problem, and how to fix it?
How did you find this polynomial which makes Hensel lifting keep failing? How many of them are there? Are they always multivariate? Do they always fail sometimes, work at other times? If the answer to the last two questions is yes, then we might add a stop condition and when that has been reached, try factorizing with an alternative permutation of symbols. just an idea -richy. -- Richard B. Kreckel <https://in.terlu.de/~kreckel/>
On 12/22/2017 10:09 AM, Richard B. Kreckel wrote:
On 11/30/2017 05:16 PM, Vitaly Magerya wrote:
> factor((2*x-l)^2*(-1+6*x-l)*(-1+3*x-l)*(-2+4*x-l)*(1+4*x-l)*(4*x-l)*(2+4*x-l));
At this point ginsh either hangs (80% of cases), or gives the correct result (20% of cases).
How did you find this polynomial which makes Hensel lifting keep failing?
I'm currently writing a particular tool related to high energy physics with GiNaC, and that polynomial is a charpoly of one of the matrices I'm working with.
How many of them are there?
As many as you want, but not as many as you'd think: most of similar polynomials are factor()ized easily, and yet you can construct enough of hanging examples too. E.g.: (7+4*x+l)*(-5+2*x-l)*(-6+6*x+l)*l*(10+2*x+l) (8+6*x-l)*(-5+4*x-l)*(-5+6*x+l)*(-2+2*x-l)*(2+4*x+l) (-9+7*x+l)*(-5+6*x+l)*(4+2*x+l)*(5+2*x-l)*(7+9*x-l)*(8+6*x-l) I haven't looked into the factorization algorithm to know what makes these particular polys special.
Are they always multivariate?
Seems that way. I only really care for the bivariate case though, so who knows.
Do they always fail sometimes, work at other times?
Seems so, yes.
If the answer to the last two questions is yes, then we might add a stop condition and when that has been reached, try factorizing with an alternative permutation of symbols.
This may be a reasonable workaround (if it doesn't affect the non-hanging cases), but I'd personally be wary of treating a symptom without understanding the underlying cause.
On 12/22/2017 07:31 PM, Vitaly Magerya wrote:
On 12/22/2017 10:09 AM, Richard B. Kreckel wrote:
If the answer to the last two questions is yes, then we might add a stop condition and when that has been reached, try factorizing with an alternative permutation of symbols.
This may be a reasonable workaround (if it doesn't affect the non-hanging cases), but I'd personally be wary of treating a symptom without understanding the underlying cause.
Well, I'm unconvinced that this failure is some kind of bug in the implementation. If it's a problem with the algorithm, it could be worth trying alternative permutations of symbols. (And, of course, it shouldn't cause regressions.) -richy. -- Richard B. Kreckel <https://in.terlu.de/~kreckel/>
participants (2)
-
Richard B. Kreckel
-
Vitaly Magerya