1
/* -*- Mode: Maxima -*- */
5
** $Id: grobner.demo,v 1.2 2003/05/03 11:40:00 starseeker Exp $
6
** Copyright (C) 1999, 2002 Marek Rychlik <rychlik@u.arizona.edu>
8
** This program is free software; you can redistribute it and/or modify
9
** it under the terms of the GNU General Public License as published by
10
** the Free Software Foundation; either version 2 of the License, or
11
** (at your option) any later version.
13
** This program is distributed in the hope that it will be useful,
14
** but WITHOUT ANY WARRANTY; without even the implied warranty of
15
** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16
** GNU General Public License for more details.
18
** You should have received a copy of the GNU General Public License
19
** along with this program; if not, write to the Free Software
20
** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
25
/* POLY_MONOMIAL_ORDER switch represents the monomial order that will globally be in effect
26
for the succeeding operations. */
28
poly_monomial_order:'lex;
30
/* POLY_EXPAND parses polynomials to internal form and back. It can be used to test
31
whether an expression correctly parses to the internal representation.
32
The following examples illustrate that indexed and transcendental function variables
36
poly_expand(x+y,[x,y]);
37
poly_expand(x-y,[x,y]);
38
poly_expand((x-y)*(x+y),[x,y]);
39
poly_expand((x+y)^2,[x,y]);
40
poly_expand((x+y)^5,[x,y]);
41
poly_expand(x/y-1,[x]);
42
poly_expand(x^2/sqrt(y)-x*exp(y)-1,[x]);
43
poly_expand(sin(x)-sin(x)^2-1,[sin(x)]);
44
poly_expand((x[2]/sin(y[3])-1)^5,[x[2]]),poly_return_term_list:true;
46
/* POLY_ADD, POLY_SUBTRACT, POLY_MULTIPLY and POLY_EXPT are the arithmetical operations on polynomials.
47
These are performed using the internal representation, but the results are converted back to the
48
Maxima general form */
50
poly_add(x^2*y+z,x-z,[x,y,z]);
51
poly_subtract(x^2*y+z,x-z,[x,y,z]);
52
poly_multiply(x^2*y+z,x-z,[x,y,z]) - (x^2*y+z)*(x-z), expand;
53
poly_expt(x-y, 3, [x,y]) - (x-y)^3, expand;
55
/* POLY_CONTENT extracts the GCD of its coefficients */
56
poly_content(21*x+35*y,[x,y]);
58
/* POLY_PRIMITIVE_PART divides a polynomial by the GCD of its coefficients */
59
poly_primitive_part(21*x+35*y,[x,y]);
61
/* POLY_S_POLYNOMIAL computest the syzygy polynomial (S-polynomial) of two polynomials */
62
poly_s_polynomial(x+y,x-y,[x,y]);
65
/* POLY_NORMAL_FORM finds the normal form of a polynomial with respect to a set of polynomials. */
66
poly_normal_form(x^2+y^2,[x-y,x+y],[x,y]);
67
poly_pseudo_divide(2*x^2+3*y^2,[7*x-y^2,11*x+y],[x,y]);
68
poly_exact_divide((x+y)^2,x+y,[x,y]);
70
/* POLY_BUCHBERGER performs the Buchberger algorithm on a list of polynomials and returns
71
the resulting Grobner basis */
72
poly_buchberger([x^2-y*x,x^2+y+x*y^2],[x,y]);
74
/* POLY_REDUCTION reduces a set of polynomials, so that
75
each polynomial is fully reduced with respect to the other polynomials */
77
poly_reduction([x^2-x*y,x*y^2+y+x^2,x*y^2+x*y+y,x*y-y^2,y^3+y^2+y],[x,y]);
79
/* POLY_MINIMIZATION selects a subset of a set of polynomials, so that no leading monomial is divisible by
80
another leading monomial */
82
poly_minimization([x^2-x*y,x*y^2+y+x^2,x*y^2+x*y+y,x*y-y^2,y^3+y^2+y],[x,y]);
84
/* POLY_REDUCED_GROBNER returns a reduced Grobner basis */
85
poly_reduced_grobner([x^2-y*x,x^2+y+x*y^2],[x,y]);
87
/* POLY_NORMALIZE divides a polynomial by its leading coefficient */
88
poly_normalize(2*x+y,[x,y]);
90
/* POLY_NORMALIZE_LIST applies POLY_NORMALIZE to each polynomial in the list */
92
poly_normalize_list([2*x+y,3*x^2+7],[x,y]);
94
/* POLY_DEPENDS_P tests whether a polynomial depends on a variable */
96
poly_depends_p(x^2+y,x,[x,y,z]);
97
poly_depends_p(x^2+y,z,[x,y,z]);
100
/* POLY_ELIMINATION_IDEAL returns the grobner basis of the K-th elimination ideal of an
101
ideal specified as a list of generating polynomials (not necessarily Grobner basis */
103
poly_elimination_ideal([x+y,x-y],0,[x,y]);
104
poly_elimination_ideal([x+y,x-y],1,[x,y]);
105
poly_elimination_ideal([x+y,x-y],2,[x,y]);
107
/* POLY_IDEAL_INTERSECTION returns the intersection of two ideals */
108
poly_ideal_intersection([x^2+y,x^2-y],[x,y^2],[x,y]);
110
/* POLY_LCM and POLY_GCD are the Grobner versions of LCM and GCD */
112
poly_lcm(x*y^2-x,x^2*y+x,[x,y]);
113
poly_gcd(x*y^2-x,x^2*y+x,[x,y]);
115
/* POLY_GROBNER_MEMBER tests whether a polynomial belongs to an ideal generated by a list of polynomials,
116
which is assumed to be a Grobner basis. Equivalent to NORMAL_FORM being 0. */
118
poly_grobner_member(x+y,[x,y],[x,y]);
120
/* POLY_GROBNER_EQUAL tests whether two Grobner bases generate the same ideal.
121
This is equivalent to checking that every polynomial of the first basis reduces to 0
122
modulo the second basis and vice versa. Note that in the example below the
123
first list is not a Grobner basis, and thus the result is FALSE. */
125
poly_grobner_equal([x+y,x-y],[x,y],[x,y]);
127
/* POLY_GROBNER_SUBSETP tests whether an ideal generated by the first list of polynomials
128
is contained in the ideal generated by the second list. For this test to always succeed,
129
the second list must be a Grobner basis */
131
poly_grobner_subsetp([x+y,x-y],[x,y],[x,y]);
133
/* POLY_POLYSATURATION_EXTENSION implements the famous Rabinowitz trick. */
134
poly_polysaturation_extension([x,y],[x^2,y^3],[x,y],[u,v]);
136
poly_saturation_extension([x,y],[x^2,y^3],[x,y],[u,v]);
138
/* POLY_IDEAL_POLYSATURATION1 for a given ideal I and polynomials f, g, ..., finds
139
the colon ideal I : f^inf : g^inf : ... */
140
poly_ideal_polysaturation1([x,y],[x^2,y^3],[x,y]);
142
/* POLY_IDEAL_SATURATION for given ideals I and J computes the ideal I : J^inf. */
143
poly_ideal_saturation([x,y],[x^2,y^3],[x,y]);
145
/* POLY_IDEAL_POLYSATURATION for a given ideal I and a sequence of ideals J1, J2, J3, ...,
146
finds the ideal I : J1^inf : J2^inf : J3^inf : ... */
147
poly_ideal_polysaturation([x,y],[[x^2],[y^3]],[x,y]);
148
poly_ideal_polysaturation([x^4-y^4], [[x-y],[x^2+y^2, x+y]],[x,y]);
150
/* POLY_COLON_IDEAL finds the reduced Grobner basis of the colon ideal I:J, i.e. the set of polynomials h
151
such that there is a polynomial F in J for which H*F is in I */
153
poly_colon_ideal([x^2*y],[y],[x,y]);
155
/* POLY_BUCHBERGER_CRITERION verifies whether a given set of polynomials is a Grobner basis with respect
156
to the current term order */
157
poly_buchberger_criterion([x,y],[x,y]);
158
poly_buchberger_criterion([x-y,x+y],[x,y]);
160
/* Grobner basis associated with Enneper minimal surface */
161
poly_grobner([x-3*u-3*u*v^2+u^3,y-3*v-3*u^2*v+v^3,z-3*u^2+3*v^2],[u,v,x,y,z]);
162
poly_reduced_grobner([x-3*u-3*u*v^2+u^3,y-3*v-3*u^2*v+v^3,z-3*u^2+3*v^2],[u,v,x,y,z]);
164
/* Cyclic roots of degree 5 */
165
poly_reduced_grobner([x+y+z+u+v,x*y+y*z+z*u+u*v+v*x,x*y*z+y*z*u+z*u*v+u*v*x+v*x*y,x*y*z*u+y*z*u*v+z*u*v*x+u*v*x*y+v*x*y*z,x*y*z*u*v-1],[u,v,x,y,z]);
167
/* The next example demonstrates the use of the switch
168
POLY_RETURN_TERM_LIST, which, if set to TRUE, makes the results to
169
appear as lists of terms listed in the current monomial order rather
170
than a general form expression */
172
block([orders:[lex,grlex,grevlex,invlex]],
173
for i:1 thru length(orders) do (
174
print(ev([orders[i], poly_expand((x^2+x+y)^3,[x,y])], poly_monomial_order=orders[i]))
176
), poly_return_term_list=true;
178
/* Grobner bases can be computed over coefficient ring of maxima general expressions */
179
poly_grobner([x*y-1,x+y],[x]);
181
/* A tough example learned from Cox */
182
poly_grobner([x^5+y^4+z^3-1,x^3+y^3+z^2-1], [x,y,z]);
184
/* An even tougher example of Cox */
185
poly_grobner([x^5+y^4+z^3-1,x^3+y^3+z^2-1], [x,y,z]);
187
/* We can also perform Grobner basis calculations modulo prime */
188
poly_grobner([x^5+y^4+z^3-1,x^3+y^3+z^2-1], [x,y,z]), modulus=3;
190
/* We can also explicitly ask for the Grobner basis to be calculated using only
191
integer coefficients. An error will result if this assertion is not satisfied. */
192
poly_grobner([x^5+y^4+z^3-1,x^3+y^3+z^2-1], [x,y,z]), poly_coefficient_ring='ring_of_integers;
194
/* The following several tests demonstrate the use of jet variables useful in processing differential equations */
197
/* Clear some variables */
198
kill(ode,t,x,y,u,v,r);
200
/* Set up dependencies */
201
depends([x,y,u,v,r],t);
203
/* These are equations representing mathematical pendulum */
204
ode:[x^2+y^2-c,'diff(x,t)-u,'diff(y,t)-v,'diff(u,t)+r*x,'diff(v,t)+r*y+1];
206
jet_vars(k):=apply(append,reverse(makelist(['diff(x,t,j),'diff(y,t,j),'diff(u,t,j),'diff(v,t,j),'diff(r,t,j)],j,0,k+1)));
208
/* Define k-fold prolongation */
209
prolongate(k):=apply(append,makelist(diff(ode,t,j),j,0,k));
211
/* Define Grobner basis of k-fold prolongation */
212
gb(k):=poly_reduced_grobner(prolongate(k),jet_vars(k));
214
/* Define the l-th projection of the k-th prolongation */
215
projection(l, k):=poly_elimination_ideal(prolongate(k),5*l,jet_vars(k));
217
/* Compute some projections */
added
removed
share/contrib/Grobner/grobner.demo
