Services & Resources / Wolfram Forums / MathGroup Archive
-----

MathGroup Archive 2009

[Date Index] [Thread Index] [Author Index]

Search the Archive

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
>



  • Prev by Date: Re: Help with algorithm to find rational roots of a bivariate
  • Next by Date: Re: newbie q-n about FinancialData
  • Previous by thread: Re: Help with algorithm to find rational roots of a bivariate
  • Next by thread: Re: Re: Help with algorithm to find rational