Help with algorithm to find rational roots of a bivariate equation?
TPiezas
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
Andrzej Kozlowski
points is that due to Poincare, who showed hoe to find all rational
points on a curve of genus 1. Since the genus is (n-1)(n-2)/2 - d, where
n is the degree and d the number of double points, your curve must have
2 double points. In fact, if I remember correctly, even the Poincare
method is not really an algorithm but a collection of special techniques
(but I may be wrong).
Anyway, I am sure that not only there is no such Mathematica program but
that there is no known algorithm for solving this type of problem at
all.
Andrzej Kozlowski
Tito Piezas
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 <da...@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
>
Daniel Lichtblau
Daniel Lichtblau
About the only pointer I can give is that the method I outlined is not
likely to be of significant help, unless you get lucky and for some
prime most substitutions in one variable give no solution in the other.
The more I look at it, the more that method appears to be brute force,
albeit well disguised. I should have caught that earlier (like, before I
posted...).
Daniel Lichtblau
Wolfram Research