Re: Help with algorithm to find rational roots of a
- To: mathgroup at smc.vnet.net
- Subject: [mg105207] Re: [mg105080] Help with algorithm to find rational roots of a
- From: Tito Piezas <tpiezas at gmail.com>
- Date: Tue, 24 Nov 2009 05:49:50 -0500 (EST)
- References: <200911201139.GAA03433@smc.vnet.net> <4B0B1AE9.1000004@wolfram.com>
Hello Daniel, Thanks for the reply. It had funny side comments, btw. :-) Anyway, this problem arose in the context of trying to solve the multi-grade eqn, x1^k + x2^k + ... + x6^k = y1^k + y2^k + ... + y6^k, for k = 2,4,6,8 One way to solve this is to use the highly symmetric form, (pa+qb+c)^k + (pa+qb-c)^k + (qa-pb+d)^k + (qa-pb-d)^k + (ra+sb)^k + (sa-rb)^k = (pa-qb+c)^k + (pa-qb-c)^k + (qa+qb+d)^k +(qa+pb-d)^k + (ra-sb)^k + (sa+rb)^k (Eq.1) where {c^2, d^2} = {ta^2+ub^2, tb^2+ua^2}. Note how "b" has just been negated in the RHS. For some constants {p,q,r,s,t,u}, this can have an _infinite_ number of non-trivial solns as the quadratic conditions on {c,d} imply an elliptic curve. This has only three known families: 1. {p,q,r,s,t,u} = {1, 3, 2, 8, 45, -11} (Piezas, 2009, for k = 2,4,6,8,10} 2. {p,q,r,s,t,u} = {1, 10, 1, 11, -27/5, 248/5} (Wroblewski, 2009) 3. {p,q,r,s,t,u} = {2, 5, 4, 6, -11/5, 64/5} (Wroblewski, 2009) One can assume p < q without loss of generality. I found the first one using one approach, while Wroblewski found the next two using a numerical search (something I should have done!). Solving for {c,d} by taking square roots, Eq. 1 is already true for k = 2. By expanding for k = 4,6, one can linearly express {t,u}, in terms of {p,q,r,s}, so those 4 are the only true unknowns. For k = 8, one gets a palindromic eqn of the form, P1a^4 + P2a^2b^2 + P1b^4 = 0 To make this true, one should find {p,q,r,s} such that P1 = P2 = 0. These are 15-deg eqns so to resolve them is horrendous. But we can use a trick to simplify matters. Note how r = np for n = {2,1,2} in the identities above. And since the system is homogeneous, it can be set s = 1 without loss of generality. So let, {r,s} = {pn, 1} The trick of making r = pn makes p have only even powers, so let p = z^(1/2) to reduce the degree even further. We end up with two 5th degs in z which is manageable. Let, Factor[Resultant[P1, P2, z]] and, after just a short while, Mathematica spits out an irreducible 22-deg eqn solely in n and q, call this E22. (Disregard the trivial linear factors in n,q.) Let n = 1, and E22 has one non-trivial linear factor, q = 10/11. Let n = 2, and E22 has two non-trivial linear factors, q = {3/8, 5/6} So far, I haven't been able to find another rational n that does not involve mere transposition of the known {p,q,r,s}. (For example, let n = 1/10, and it gives q = -1/11). Question: So what's the easier route: a) find n such that E22 has a linear factor?, or b) just brute-force search for non-trivial {p,q,r,s} that makes Eq.1 true for k = {2,4,6,8}? Tip 1: To speed up things, one can assume p < q. Tip 2: If Eq.1 = 0 for k = 12, then those {p,q,r,s} are trivial. Anyone can find a 4th infinite family? Any help is appreciated. - Titus On Mon, Nov 23, 2009 at 5:29 PM, Daniel Lichtblau <danl at wolfram.com> wrote: > TPiezas wrote: > >> Hello all, >> >> Does anyone know an efficient algorithm using Mathematica that can >> find _rational_ roots of a non-homogenous eqn in two variables with >> deg > 4? For ex, you want to find rational {x,y} such that, >> >> F(x,y) = x^n + (P_1)x^(n-1) + (P_2)x^(n-2) + .... = 0 >> >> where the P_i are polynomials in y. >> >> I _always_ come across this situation in the course of experimental >> mathematics, and it would be great if Mathematica had a built-in >> feature that solves _bivariate_ eqns in the rationals. Right now, I >> have a 22-deg eqn in x with coefficients in y. I know three non- >> trivial rational pairs {x,y} such that F(x,y) = 0, but I want to know >> if there are others. If there are, a certain family of identities >> would have more members. >> >> Any help will be appreciated. >> >> - Titus >> > > THis might give you cause to question that last remark. > > I can outline an approach that might get you rational solutions up to some > specified size (in numerator and denominator), but it is not terribly > efficient. > > The idea is to use homomorphic imaging, that is solve modulo a prime, and > lift to prime powers. Say size is the larger in absolute value of numerator > and denominator we will consider. In outline form, it goes like this. > > (1) First pick a modulus p of modest size, say 11. > > (2) For all 0<=x<=p-1, find solutions mod p for y. > > (3) Use rational recovery to go from solutions mod p to candidate solutions > over the rationals. > > (4) Select from the candidates all actual solutions. > > (5) Replace p by p^2, and use Use Hensel lifting to take solutions mod the > old p into solutions mod the new one. At this step an interesting thing > happens: individual solutions can split into multiple solutions. This is bad > in that it gives rise to a combinatorial explosion. But it is good because > it allows us to find more rational solution candidates. > > (6) If the new p is larger than twice size^2, we are done. Else go to step > (3). > > There are numerous caveats to this approach. First is that it is likely to > be slow in actual applications of interest. > > Second is that there are ways in which a prime can be "unlucky". Obviously > if a solution denominator is divisible by that prime, we'll not recover that > rational. But I believe other things can go wrong, just as they can when > using homomorphic images to factor a multivariate polynomial. > > A third is that I am not entirely certain the termination condition > guarantees we'll have all solutions up to the specified size (assuming none > were lost to the second caveat). I suspect it can be proven, in some manner > similar to how one shows rational recovery gives guaranteed results for > linear algebra problems after a certain amount of lifting (and assuming the > prime is not one of finitely many that are unlucky...). > > A fourth, perhaps part of the first, is that I've no idea how to cull out > lifted "solutions" that will not eventually give rise to actual solutions. I > could remove the ones that gave solutions at an earlier lifting step, but in > something like the example below this is but a small percentage. > > I'll illustrate with the usual Pythagorian triples problem, solving > x^2+y^2-1 == 0 over the rationals. Note that I am using Mathematica's > development version of Solve, so some amount of work would be needed to > translate the current version (though this can be done). > > rationalRecover[x_,pk_]:= > ((#[[2,2]]/#[[1,2,2]])&)[Internal`HGCD[pk,x]] > > bpoly = x^2+y^2-1; > > (* First we solve mod 11. *) > > firstsols = Flatten[Table[Solve[{bpoly,y-j}==0, {x,y}, > Modulus->11], {j,2,10}] /. {}:>Sequence[], 1]; > firstpairs = {x,y}/.firstsols; > > (* Now get candidates. *) > > firstrats = Map[rationalRecover[#,11^2]&, firstpairs, {2}]; > > (* Now lift and iterate the process, We'll harvest triples later. *) > > nextgbs = Table[GroebnerBasis[{11^2, bpoly, > (x-firstpairs[[j,1]])^2, (y-firstpairs[[j,2]])^2}, {x,y}, > CoefficientDomain->Integers], {j,Length[firstpairs]}]; > nextsols = Flatten[Map[Solve[Rest[#]==0, {x,y}, Modulus->First[#]]&, > nextgbs], 1]; > nextpairs = {x,y}/.nextsols; > nextrats = Map[rationalRecover[#,11^2]&, nextpairs, {2}]; > > (* Lift and iterate again. It starts to get slow at this point. I think > that's the sound of a combinatorial explosion. *) > > nextgbs2 = Table[GroebnerBasis[{11^4, bpoly, > (x-nextpairs[[j,1]])^2, (y-nextpairs[[j,2]])^2}, {x,y}, > CoefficientDomain->Integers], {j,Length[nextpairs]}]; > nextsols2 = Flatten[Map[Solve[Rest[#]==0, {x,y}, Modulus->First[#]]&, > nextgbs2], 1]; > nextpairs2 = {x,y}/.nextsols2; > nextrats2 = Map[rationalRecover[#,11^4]&, nextpairs2, {2}]; > > (* Let's see what solutions we obtained. *) > > In[39]:= InputForm[firsttriples = > Select[firstrats, (bpoly/.Thread[{x,y}->#])==0&]] > Out[39]//InputForm= {} > > (* Okay nothing there. *) > > In[40]:= InputForm[nexttriples = > Select[nextrats, (bpoly/.Thread[{x,y}->#])==0&]] > Out[40]//InputForm= > {{3/5, 4/5}, {-3/5, 4/5}, {4/5, 3/5}, {-4/5, 3/5}, {4/5, -3/5}, > {-4/5, -3/5}, {3/5, -4/5}, {-3/5, -4/5}, {0, -1}} > > (* Something, though not terribly exciting. *) > > In[41]:= InputForm[nexttriples2 = > Select[nextrats2, (bpoly/.Thread[{x,y}->#])==0&]] > > Out[41]//InputForm= > {{-20/29, 21/29}, {-12/13, -5/13}, {3/5, 4/5}, {3/5, 4/5}, {3/5, 4/5}, > {3/5, 4/5}, {3/5, 4/5}, {3/5, 4/5}, {3/5, 4/5}, {3/5, 4/5}, {3/5, 4/5}, > {3/5, 4/5}, {3/5, 4/5}, {40/41, -9/41}, {-35/37, 12/37}, {-7/25, -24/25}, > {8/17, -15/17}, {-16/65, 63/65}, {45/53, -28/53}, {-45/53, -28/53}, > {16/65, 63/65}, {-8/17, -15/17}, {7/25, -24/25}, {35/37, 12/37}, > {-40/41, -9/41}, {-3/5, 4/5}, {-3/5, 4/5}, {-3/5, 4/5}, {-3/5, 4/5}, > {-3/5, 4/5}, {-3/5, 4/5}, {-3/5, 4/5}, {-3/5, 4/5}, {-3/5, 4/5}, > {-3/5, 4/5}, {-3/5, 4/5}, {12/13, -5/13}, {20/29, 21/29}, {63/65, -16/65}, > {4/5, 3/5}, {4/5, 3/5}, {4/5, 3/5}, {4/5, 3/5}, {4/5, 3/5}, {4/5, 3/5}, > {4/5, 3/5}, {4/5, 3/5}, {4/5, 3/5}, {4/5, 3/5}, {4/5, 3/5}, {-28/53, > 45/53}, > {-9/41, 40/41}, {12/37, -35/37}, {21/29, -20/29}, {-24/25, -7/25}, > {-5/13, -12/13}, {-15/17, 8/17}, {15/17, 8/17}, {5/13, -12/13}, > {24/25, -7/25}, {-21/29, -20/29}, {-12/37, -35/37}, {9/41, 40/41}, > {28/53, 45/53}, {-4/5, 3/5}, {-4/5, 3/5}, {-4/5, 3/5}, {-4/5, 3/5}, > {-4/5, 3/5}, {-4/5, 3/5}, {-4/5, 3/5}, {-4/5, 3/5}, {-4/5, 3/5}, > {-4/5, 3/5}, {-4/5, 3/5}, {-63/65, -16/65}, {63/65, 16/65}, {4/5, -3/5}, > {4/5, -3/5}, {4/5, -3/5}, {4/5, -3/5}, {4/5, -3/5}, {4/5, -3/5}, > {4/5, -3/5}, {4/5, -3/5}, {4/5, -3/5}, {4/5, -3/5}, {4/5, -3/5}, > {-28/53, -45/53}, {-9/41, -40/41}, {12/37, 35/37}, {21/29, 20/29}, > {-24/25, 7/25}, {-5/13, 12/13}, {-15/17, -8/17}, {15/17, -8/17}, > {5/13, 12/13}, {24/25, 7/25}, {-21/29, 20/29}, {-12/37, 35/37}, > {9/41, -40/41}, {28/53, -45/53}, {-4/5, -3/5}, {-4/5, -3/5}, {-4/5, -3/5}, > {-4/5, -3/5}, {-4/5, -3/5}, {-4/5, -3/5}, {-4/5, -3/5}, {-4/5, -3/5}, > {-4/5, -3/5}, {-4/5, -3/5}, {-4/5, -3/5}, {-63/65, 16/65}, {-20/29, > -21/29}, > {-12/13, 5/13}, {3/5, -4/5}, {3/5, -4/5}, {3/5, -4/5}, {3/5, -4/5}, > {3/5, -4/5}, {3/5, -4/5}, {3/5, -4/5}, {3/5, -4/5}, {3/5, -4/5}, > {3/5, -4/5}, {3/5, -4/5}, {40/41, 9/41}, {-35/37, -12/37}, {-7/25, 24/25}, > {8/17, 15/17}, {-16/65, -63/65}, {45/53, 28/53}, {-45/53, 28/53}, > {16/65, -63/65}, {-8/17, 15/17}, {7/25, 24/25}, {35/37, -12/37}, > {-40/41, 9/41}, {-3/5, -4/5}, {-3/5, -4/5}, {-3/5, -4/5}, {-3/5, -4/5}, > {-3/5, -4/5}, {-3/5, -4/5}, {-3/5, -4/5}, {-3/5, -4/5}, {-3/5, -4/5}, > {-3/5, -4/5}, {-3/5, -4/5}, {12/13, 5/13}, {20/29, -21/29}, {0, -1}, > {0, -1}, {0, -1}, {0, -1}, {0, -1}, {0, -1}, {0, -1}, {0, -1}, {0, -1}, > {0, -1}, {0, -1}, {-99/101, 20/101}, {11/61, 60/61}, {-33/65, 56/65}, > {-55/73, 48/73}, {77/85, 36/85}, {-77/85, 36/85}, {55/73, 48/73}, > {33/65, 56/65}, {-11/61, 60/61}, {99/101, 20/101}} > > So that was a more substantial harvest. > > Clearly I skimped on a many details of why this might work at all. In part > this is because I don't know them all. But here are a few places one can get > background information on some of the tactics indicated above. > > For details on the upcoming Solve functionality: > http://library.wolfram.com/infocenter/Conferences/7552/ > > For some idea of how the lifting works, at least in the univariate case: > http://library.wolfram.com/infocenter/Conferences/7511/ > http://library.wolfram.com/infocenter/MathSource/7522/ > These show how one can use Groebner bases over the integers to lift p-adic > gcds in the univariate case. I do not claim that the same proof gives lifted > solutions in the bivariate case. But we do indeed get them by that method. > > For explanation of Rational Recovery, see the ARRA bill signed into law > last February. No, wait, that was a different recovery (though perhaps > rational all the same). Instead see: > http://library.wolfram.com/infocenter/Conferences/7534/ > > Daniel Lichtblau > Wolfram Research >
- Follow-Ups:
- Re: Re: Help with algorithm to find rational
- From: Daniel Lichtblau <danl@wolfram.com>
- Re: Re: Help with algorithm to find rational
- References:
- Help with algorithm to find rational roots of a bivariate equation?
- From: TPiezas <tpiezas@gmail.com>
- Help with algorithm to find rational roots of a bivariate equation?