Re: Diophantic equations
- To: mathgroup at smc.vnet.net
- Subject: [mg32777] Re: [mg32746] Diophantic equations
- From: Daniel Lichtblau <danl at wolfram.com>
- Date: Sat, 9 Feb 2002 05:11:52 -0500 (EST)
- References: <200202080849.DAA11925@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
Max Ulbrich wrote:
>
> Hi,
>
> I try to solve equations where only integers are allowed, like
>
> a + b + 10 c = 100
> 5 a + 2 b + c = 100
>
> a, b, c should be integers and positve.
>
> How can I do this?
>
> Max
There is a method for handling these by Groebner bases, due to
Conti and Traverso. Nice expositions may be found in
W. Adams and P. Loustaunau (1994). An Introduction to Gröbner Bases.
Graduate Studies in Mathematics 3, American Mathematical Society.
D. Cox, J. Little, and D. O'Shea (1998). Using Algebraic Geometry.
Graduate Texts in Mathematics 185, Springer-Verlag
The idea is as follows. We form two new sets of variables, one for
rows and one for columns. We show rows first. We will use x[j] for
row j and raise x[j] to a power given by each side of the jth
equation. For the example given we have:
x[1]^(a+b+10*c) == x[1]^100
x[2]^(5*a+2*b+c) == x[2]^100
Now multiply these left hand sides and right hand sides together
to obtain
x[1]^(a+b+10*c)*x[2]^(5*a+2*b+c) == x[1]^100*x[2]^100
We regroup so that our original variables are exponents of
certain bases.
(x[1]*x[2]^5)^a * (x[1]*x[2]^2)^b * (x[1]^10*x[2])^c ==
x[1]^100*x[2]^100
We now bring in the other set of variables. Specifically, we equate
them to the "bases" that appear on the left hand side above.
y[1] == (x[1]*x[2]^5)
y[2] == (x[1]*x[2]^2)
y[3] == (x[1]^10*x[2])
Here is where Groebner bases enter. The idea is that we form a basis
for the polynomials
polys = {y[1] - (x[1]*x[2]^5), y[2] - (x[1]*x[2]^2), y[3] -
(x[1]^10*x[2])}
with respect to a term ordering that places the x[j]'s higher than
the y[k]'s (we will use lexicographic).
gb = GroebnerBasis[polys, {x[1],x[2],y[1],y[2],y[3]}];
If there is a nonnegative integer combination that works for
{a,b,c} in order to satisfy the original system, then it will show
up when we reduce the right hand side, x[1]^100*x[2]^100, in terms
of that basis and the same term ordering. Specifically, if an
answer exists, then we will get a reductum entirely in terms of
the y[k]'s. If you require strict positivity, replace your original
variables with aprime==a-1, etc., and rewrite the equations in terms
of these.
InputForm[Last[
PolynomialReduce[x[1]^100*x[2]^100, gb, {x[1],x[2],y[1],y[2],y[3]}]]]
Out[191]//InputForm= y[1]^11*y[2]^19*y[3]^7
Our exponents {11,19,7} comprise a solution for the original
variables {a,b,c}. Note that such solutions need not be unique,
and others might be obtained by changing the term order used in
the computations. For example, if we make the original right hand
sides both 1000 then we can obtain:
In[217]:= gb = GroebnerBasis[polys, {x[1],x[2],y[1],y[2],y[3]}];
In[218]:= InputForm[Last[
PolynomialReduce[x[1]^1000*x[2]^1000, gb,
{x[1],x[2],y[1],y[2],y[3]}]]]
Out[218]//InputForm= y[1]^15*y[2]^435*y[3]^55
In[219]:= gb = GroebnerBasis[polys, {x[1],x[2],y[2],y[1],y[3]}];
In[220]:= InputForm[Last[
PolynomialReduce[x[1]^1000*x[2]^1000, gb,
{x[1],x[2],y[2],y[1],y[3]}]]]
Out[220]//InputForm= y[1]^167*y[2]^43*y[3]^79
The latter was performed alot faster, I might add (this sort of
thing is quite sensitive to term ordering).
To handle situations where we have negative coefficients one
may adapt by using extra variables and relations that in essence
force variables to be reciprocals.
The Groebner basis approach is more generally used in linear
integer programming. Details are in the references. One uses
the objective function cost to impose an appropriate term order.
A package to do this has been written by Devendra Kapadia here at
Wolfram Research (dkapadia at wolfram.com). I might add that it
demostrated some bottlenecks. Based on this I made modifications
to GroebnerBasis in our development version that are specific to
this sort of polynomial input (these are known as toric ideals).
I will also point out that the state of the art for this
general approach has advanced to ways that eradicate the
second set of variables and are in other respects specific to
toric ideals.
There are other methods, generally much faster, for attacking both
diophantine systems. There are also quite different methods for
handling linear integer programming. The former may be addressed
in a future release of Mathematica. The latter is also a topic of
some interest. One might see preliminary work in this direction at
http://library.wolfram.com/conferences/devconf2001/#numerics
The first bullet item discusses a heuristic approach to discrete
optimization.
Daniel Lichtblau
Wolfram Research
- References:
- Diophantic equations
- From: Max Ulbrich <ulbrich@biochem.mpg.de>
- Diophantic equations