1 | /////////////////////////////////////////////////////////////////////////////// |
---|
2 | version="version primitiv.lib 4.0.0.0 Jun_2013 "; // $Id$ |
---|
3 | category="Commutative Algebra"; |
---|
4 | info=" |
---|
5 | LIBRARY: primitiv.lib Computing a Primitive Element |
---|
6 | AUTHOR: Martin Lamm, email: lamm@mathematik.uni-kl.de |
---|
7 | |
---|
8 | PROCEDURES: |
---|
9 | primitive(ideal i); find minimal polynomial for a primitive element |
---|
10 | primitive_extra(i); find primitive element for two generators |
---|
11 | splitring(f,R[,L]); define ring extension with name R and switch to it |
---|
12 | "; |
---|
13 | |
---|
14 | LIB "random.lib"; |
---|
15 | /////////////////////////////////////////////////////////////////////////////// |
---|
16 | |
---|
17 | proc primitive(ideal i) |
---|
18 | "USAGE: primitive(i); i ideal |
---|
19 | ASSUME: i is given by generators m[1],...,m[n] such that for j=1,...,n @* |
---|
20 | - m[j] is a polynomial in k[x(1),...,x(j)] @* |
---|
21 | - m[j](a[1],...,a[j-1],x(j)) is the minimal polynomial for a[j] over |
---|
22 | k(a[1],...,a[j-1]) @* |
---|
23 | (k the ground field of the current basering and x(1),...,x(n) |
---|
24 | the ring variables). |
---|
25 | RETURN: ideal j in k[x(n)] with |
---|
26 | - j[1] a minimal polynomial for a primitive element b of |
---|
27 | k(a[1],...,a[n]) over k, |
---|
28 | - j[2],...,j[n+1] polynomials in k[x(n)] such that j[i+1](b)=a[i] |
---|
29 | for i=1,...,n. |
---|
30 | NOTE: the number of variables in the basering has to be exactly n, |
---|
31 | the number of given generators (i.e., minimal polynomials).@* |
---|
32 | If the ground field k has only a few elements it may happen that no |
---|
33 | linear combination of a[1],...,a[n] is a primitive element. In this |
---|
34 | case @code{primitive(i)} returns the zero ideal, and one should use |
---|
35 | @code{primitive_extra(i)} instead. |
---|
36 | SEE ALSO: primitive_extra |
---|
37 | KEYWORDS: primitive element |
---|
38 | EXAMPLE: example primitive; shows an example |
---|
39 | " |
---|
40 | { |
---|
41 | def altring=basering; |
---|
42 | execute("ring deglexring=("+charstr(altring)+"),("+varstr(altring)+"),dp;"); |
---|
43 | ideal j; |
---|
44 | execute("ring lexring=("+charstr(altring)+"),("+varstr(altring)+"),lp;"); |
---|
45 | ideal i=fetch(altring,i); |
---|
46 | |
---|
47 | int k,schlecht,Fehlversuche,maxtry; |
---|
48 | int nva = nvars(basering); |
---|
49 | int p=char(basering); |
---|
50 | if (p==0) { |
---|
51 | p=100000; |
---|
52 | if (nva<3) { maxtry= 100000000; } |
---|
53 | else { maxtry=2147483647; } |
---|
54 | } |
---|
55 | else { |
---|
56 | if ((nva<4) || (p<60)) { |
---|
57 | maxtry=p^(nva-1); } |
---|
58 | else { |
---|
59 | maxtry=2147483647; // int overflow(^) vermeiden |
---|
60 | } |
---|
61 | } |
---|
62 | ideal jmap,j; |
---|
63 | map phi; |
---|
64 | option(redSB); |
---|
65 | |
---|
66 | //-------- Mache so lange Random-Koord.wechsel, bis letztes Polynom ------------- |
---|
67 | //--------------- das Minpoly eines primitiven Elements ist : ---------------- |
---|
68 | for (Fehlversuche=0; Fehlversuche<maxtry; Fehlversuche++) { |
---|
69 | schlecht=0; |
---|
70 | if ((p<60) && (nva==2)) { // systematische Suche statt random |
---|
71 | jmap=ideal(var(1),var(2)+Fehlversuche*var(1)); |
---|
72 | } |
---|
73 | else { |
---|
74 | if (Fehlversuche==0) { jmap=maxideal(1);} |
---|
75 | else { |
---|
76 | if (Fehlversuche<5) { jmap=randomLast(10);} |
---|
77 | else { |
---|
78 | if (Fehlversuche<20) { jmap=randomLast(100);} |
---|
79 | else { jmap=randomLast(100000000);} |
---|
80 | }} // groessere Werte machen keinen Sinn |
---|
81 | } |
---|
82 | phi=lexring,jmap; |
---|
83 | j=phi(i); |
---|
84 | setring deglexring; |
---|
85 | //--------------- Berechne reduzierte Standardbasis mit fglm: ---------------- |
---|
86 | j=std(fetch(lexring,j)); |
---|
87 | setring lexring; |
---|
88 | j=fglm(deglexring,j); |
---|
89 | //-- teste, ob SB n Elemente enthaelt (falls ja, ob lead(Fi)=xi i=1... n-1): - |
---|
90 | if (size(j)==nva) { |
---|
91 | for (k=1; k<nva; k++) { |
---|
92 | j[k+1]=j[k+1]/leadcoef(j[k+1]); // normiere die Erzeuger |
---|
93 | if (lead(j[k+1]) != var(nva-k)) { schlecht=1;} |
---|
94 | } |
---|
95 | if (schlecht==0) { |
---|
96 | //--- Random-Koord.wechsel war gut: Berechne das zurueckzugebende Ideal: ----- |
---|
97 | ideal erg; |
---|
98 | for (k=1; k<nva; k++) { erg[k]=var(k)-j[nva-k+1]; } |
---|
99 | // =g_k(x_n) mit a_k=g_k(a_n) |
---|
100 | erg[nva]=var(nva); |
---|
101 | map chi=lexring,erg; |
---|
102 | ideal extra=maxideal(1);extra=phi(extra); |
---|
103 | // sonst: "argument of a map must have a name" |
---|
104 | erg=j[1],chi(extra); // j[1] = Minimalpolynom |
---|
105 | setring altring; |
---|
106 | return(fetch(lexring,erg)); |
---|
107 | } |
---|
108 | } |
---|
109 | dbprint("The random coordinate change was bad!"); |
---|
110 | } |
---|
111 | if (voice==2) { |
---|
112 | "// ** Warning: No primitive element could be found."; |
---|
113 | "// If the given ideal really describes the minimal polynomials of"; |
---|
114 | "// a series of algebraic elements (cf. `help primitive;') then"; |
---|
115 | "// try `primitive_extra'."; |
---|
116 | } |
---|
117 | setring altring; |
---|
118 | return(ideal(0)); |
---|
119 | } |
---|
120 | example |
---|
121 | { "EXAMPLE:"; echo = 2; |
---|
122 | ring exring=0,(x,y),dp; |
---|
123 | ideal i=x2+1,y2-x; // compute Q(i,i^(1/2))=:L |
---|
124 | ideal j=primitive(i); |
---|
125 | j[1]; // L=Q(a) with a=(-1)^(1/4) |
---|
126 | j[2]; // i=a^2 |
---|
127 | j[3]; // i^(1/2)=a |
---|
128 | // the 2nd element was already primitive! |
---|
129 | j=primitive(ideal(x2-2,y2-3)); // compute Q(sqrt(2),sqrt(3)) |
---|
130 | j[1]; |
---|
131 | j[2]; |
---|
132 | j[3]; |
---|
133 | // no element was primitive -- the calculation of primitive elements |
---|
134 | // is based on a random choice. |
---|
135 | } |
---|
136 | /////////////////////////////////////////////////////////////////////////////// |
---|
137 | |
---|
138 | proc primitive_extra(ideal i) |
---|
139 | "USAGE: primitive_extra(i); i ideal |
---|
140 | ASSUME: The ground field of the basering is k=Q or k=Z/pZ and the ideal |
---|
141 | i is given by 2 generators f,g with the following properties: |
---|
142 | @format |
---|
143 | f is the minimal polynomial of a in k[x], |
---|
144 | g is a polynomial in k[x,y] s.th. g(a,y) is the minpoly of b in k(a)[y]. |
---|
145 | @end format |
---|
146 | Here, x is the name of the first ring variable, y the name of the |
---|
147 | second. |
---|
148 | RETURN: ideal j in k[y] such that |
---|
149 | @format |
---|
150 | j[1] is the minimal polynomial for a primitive element c of k(a,b) over k, |
---|
151 | j[2] is a polynomial s.th. j[2](c)=a. |
---|
152 | @end format |
---|
153 | NOTE: While @code{primitive(i)} may fail for finite fields, |
---|
154 | @code{primitive_extra(i)} tries all elements of k(a,b) and, hence, |
---|
155 | always finds a primitive element. @* |
---|
156 | In order to do this (try all elements), field extensions like Z/pZ(a) |
---|
157 | are not allowed for the ground field k. @* |
---|
158 | @code{primitive_extra(i)} assumes that the second generator, g, is |
---|
159 | monic as polynomial in (k[x])[y]. |
---|
160 | EXAMPLE: example primitive_extra; shows an example |
---|
161 | " |
---|
162 | { |
---|
163 | def altring=basering; |
---|
164 | int grad1=deg(i[1]); |
---|
165 | int grad2=deg(jet(i[2],0,intvec(1,0))); |
---|
166 | if (grad2==0) { ERROR("i[2] is not monic"); } |
---|
167 | int countx,countz; |
---|
168 | if (size(variables(i[1]))!=1) { ERROR("i[1] must be poly in x"); } |
---|
169 | if (size(variables(i[2]))>2) { ERROR("i[2] must be poly in x,a"); } |
---|
170 | //if (variables(i[2])[2]!=a) { ERROR("i[2] must be poly in x,a"); } |
---|
171 | ring deglexring=char(altring),(x,y,z),dp; |
---|
172 | map transfer=altring,x,z; |
---|
173 | ideal i=transfer(i); |
---|
174 | if (size(i)!=2) |
---|
175 | { |
---|
176 | ERROR("either wrong number of given minimal polynomials"+newline+ |
---|
177 | "or wrong choice of ring variables (must use the first two)"); |
---|
178 | } |
---|
179 | matrix mat; |
---|
180 | ring lexring=char(altring),(x,y),lp; |
---|
181 | ideal j; |
---|
182 | ring deglex2ring=char(altring),(x,y),dp; |
---|
183 | ideal j; |
---|
184 | setring deglexring; |
---|
185 | ideal j; |
---|
186 | option(redSB); |
---|
187 | poly g=z; |
---|
188 | int found=0; |
---|
189 | |
---|
190 | //---------------- Schleife zum Finden des primitiven Elements --------------- |
---|
191 | //--- Schleife ist so angordnet, dass g in Charakteristik 0 linear bleibt ---- |
---|
192 | while (found==0) |
---|
193 | { |
---|
194 | j=eliminate(i+ideal(g-y),z); |
---|
195 | setring deglex2ring; |
---|
196 | j=std(imap(deglexring,j)); |
---|
197 | setring lexring; |
---|
198 | j=fglm(deglex2ring,j); |
---|
199 | if (size(j)==2) |
---|
200 | { |
---|
201 | if (deg(j[1])==grad1*grad2) |
---|
202 | { |
---|
203 | j[2]=j[2]/leadcoef(j[2]); // Normierung |
---|
204 | if (lead(j[2])==x) |
---|
205 | { // Alles ok |
---|
206 | found=1; |
---|
207 | } |
---|
208 | } |
---|
209 | } |
---|
210 | setring deglexring; |
---|
211 | if (found==0) |
---|
212 | { |
---|
213 | //------------------ waehle ein neues Polynom g ------------------------------ |
---|
214 | dbprint("Still searching for primitive element..."); |
---|
215 | countx=0; |
---|
216 | countz=0; |
---|
217 | while (found==0) |
---|
218 | { |
---|
219 | countx++; |
---|
220 | if (countx>=grad1) |
---|
221 | { |
---|
222 | countx=0; |
---|
223 | countz++; |
---|
224 | if (countz>=grad2) |
---|
225 | { ERROR("No primitive element found!! This should NEVER happen!"); } |
---|
226 | } |
---|
227 | g = g +x^countx *z^countz; |
---|
228 | mat=coeffs(g,z); |
---|
229 | if (size(mat)>countz) |
---|
230 | { |
---|
231 | mat=coeffs(mat[countz+1,1],x); |
---|
232 | if (size(mat)>countx) |
---|
233 | { |
---|
234 | if (mat[countx+1,1] != 0) |
---|
235 | { |
---|
236 | found=1; // d.h. hier: neues g gefunden |
---|
237 | }}} |
---|
238 | } |
---|
239 | found=0; |
---|
240 | } |
---|
241 | } |
---|
242 | //------------------- primitives Element gefunden; Rueckgabe ----------------- |
---|
243 | setring lexring; |
---|
244 | j[2]=x-j[2]; |
---|
245 | setring altring; |
---|
246 | map transfer=lexring,var(1),var(2); |
---|
247 | return(transfer(j)); |
---|
248 | } |
---|
249 | example |
---|
250 | { "EXAMPLE:"; echo = 2; |
---|
251 | ring exring=3,(x,y),dp; |
---|
252 | ideal i=x2+1,y3+y2-1; |
---|
253 | primitive_extra(i); |
---|
254 | ring extension=(3,y),x,dp; |
---|
255 | minpoly=y6-y5+y4-y3-y-1; |
---|
256 | number a=y5+y4+y2+y+1; |
---|
257 | a^2; |
---|
258 | factorize(x2+1); |
---|
259 | factorize(x3+x2-1); |
---|
260 | } |
---|
261 | /////////////////////////////////////////////////////////////////////////////// |
---|
262 | |
---|
263 | proc splitring(poly f,list #) |
---|
264 | "USAGE: splitring(f[,L]); f poly, L list of polys and/or ideals |
---|
265 | (optional) |
---|
266 | ASSUME: f is univariate and irreducible over the active ring. @* |
---|
267 | The active ring must allow an algebraic extension (e.g., it cannot |
---|
268 | be a transcendent ring extension of Q or Z/p). |
---|
269 | RETURN: ring; @* |
---|
270 | if called with a nonempty second parameter L, then in the output |
---|
271 | ring there is defined a list erg ( =L mapped to the new ring); |
---|
272 | if the minpoly of the active ring is non-zero, then the image of |
---|
273 | the primitive root of f in the output ring is appended as last |
---|
274 | entry of the list erg. |
---|
275 | NOTE: If the old ring has no parameter, the name @code{a} is chosen for the |
---|
276 | parameter of R (if @code{a} is no ring variable; if it is, @code{b} is |
---|
277 | chosen, etc.; if @code{a,b,c,o} are ring variables, |
---|
278 | @code{splitring(f[,L])} produces an error message), otherwise the |
---|
279 | name of the parameter is kept and only the minimal polynomial is |
---|
280 | changed. @* |
---|
281 | The names of the ring variables and the orderings are not affected. @* |
---|
282 | KEYWORDS: algebraic field extension; extension of rings |
---|
283 | EXAMPLE: example splitring; shows an example |
---|
284 | " |
---|
285 | { |
---|
286 | //----------------- split ist bereits eine proc in 'inout.lib' ! ------------- |
---|
287 | if (size(#)>=1) { |
---|
288 | list L=#; |
---|
289 | int L_groesse=size(L); |
---|
290 | } |
---|
291 | else { int L_groesse=-1; } |
---|
292 | //-------------- ermittle das Minimalpolynom des aktuellen Rings: ------------ |
---|
293 | string minp=string(minpoly); |
---|
294 | |
---|
295 | def altring=basering; |
---|
296 | string charakt=string(char(altring)); |
---|
297 | string varnames=varstr(altring); |
---|
298 | string algname; |
---|
299 | int i; |
---|
300 | int anzvar=size(maxideal(1)); |
---|
301 | //--------------- Fall 1: Bisheriger Ring hatte kein Minimalpolynom ---------- |
---|
302 | if (minp=="0") { // only possible without parameters (by assumption) |
---|
303 | if (find(varnames,"a")==0) { algname="a";} |
---|
304 | else { if (find(varnames,"b")==0) { algname="b";} |
---|
305 | else { if (find(varnames,"c")==0) |
---|
306 | { algname="c";} |
---|
307 | else { if (find(varnames,"o")==0) |
---|
308 | { algname="o";} |
---|
309 | else { |
---|
310 | "** Sorry -- could not find a free name for the primitive element."; |
---|
311 | "** Try e.g. a ring without 'a' or 'b' as variable."; |
---|
312 | return(); |
---|
313 | }} |
---|
314 | } |
---|
315 | } |
---|
316 | //-- erzeuge einen String, der das Minimalpolynom des neuen Rings enthaelt: - |
---|
317 | execute("ring splt1="+charakt+","+algname+",dp;"); |
---|
318 | ideal abbnach=var(1); |
---|
319 | for (i=1; i<anzvar; i++) { abbnach=abbnach,var(1); } |
---|
320 | map nach_splt1=altring,abbnach; |
---|
321 | execute("poly mipol="+string(nach_splt1(f))+";"); |
---|
322 | string Rminp=string(mipol); |
---|
323 | |
---|
324 | //--------------------- definiere den neuen Ring: --------------------------- |
---|
325 | execute("ring neuring = ("+charakt+","+algname+"),("+varnames+"),(" |
---|
326 | +ordstr(altring)+");"); |
---|
327 | execute("minpoly="+Rminp+";"); |
---|
328 | |
---|
329 | //---------------------- Berechne die zurueckzugebende Liste: --------------- |
---|
330 | if (L_groesse>0) { |
---|
331 | list erg; |
---|
332 | map take=altring,maxideal(1); |
---|
333 | erg=take(L); |
---|
334 | } |
---|
335 | } |
---|
336 | else { |
---|
337 | |
---|
338 | //------------- Fall 2: Bisheriger Ring hatte ein Minimalpolynom: ----------- |
---|
339 | algname=parstr(altring); // Name des algebraischen Elements |
---|
340 | if (npars(altring)>1) {"only one Parameter is allowed!!"; return(altring);} |
---|
341 | |
---|
342 | //---------------- Minimalpolynom in ein Polynom umwandeln: ----------------- |
---|
343 | execute("ring splt2="+charakt+","+algname+",dp;"); |
---|
344 | execute("poly mipol="+minp+";"); |
---|
345 | // f ist Polynom in algname und einer weiteren Variablen -> mache f bivariat: |
---|
346 | execute("ring splt3="+charakt+",("+algname+","+varnames+"),dp;"); |
---|
347 | poly f=imap(altring,f); |
---|
348 | |
---|
349 | //-------------- Vorbereitung des Aufrufes von primitive: ------------------- |
---|
350 | execute("ring splt1="+charakt+",(x,y),dp;"); |
---|
351 | ideal abbnach=x; |
---|
352 | for (i=1; i<=anzvar; i++) { abbnach=abbnach,y; } |
---|
353 | map nach_splt1_3=splt3,abbnach; |
---|
354 | map nach_splt1_2=splt2,x; |
---|
355 | ideal maxid=nach_splt1_2(mipol),nach_splt1_3(f); |
---|
356 | ideal primit=primitive(maxid); |
---|
357 | if (size(primit)==0) { // Suche mit 1. Proc erfolglos |
---|
358 | primit=primitive_extra(maxid); |
---|
359 | } |
---|
360 | //-- erzeuge einen String, der das Minimalpolynom des neuen Rings enthaelt: - |
---|
361 | setring splt2; |
---|
362 | map nach_splt2=splt1,0,var(1); // x->0, y->a |
---|
363 | minp=string(nach_splt2(primit)[1]); |
---|
364 | if (printlevel > -1) { "// new minimal polynomial:",minp; } |
---|
365 | //--------------------- definiere den neuen Ring: --------------------------- |
---|
366 | execute("ring neuring = ("+charakt+","+algname+"),("+varnames+"),(" |
---|
367 | +ordstr(altring)+");"); |
---|
368 | execute("minpoly="+minp+";"); |
---|
369 | |
---|
370 | if (L_groesse>0) { |
---|
371 | //---------------------- Berechne die zurueckzugebende Liste: ------------- |
---|
372 | list erg; |
---|
373 | setring splt3; |
---|
374 | list zwi=imap(altring,L); |
---|
375 | map nach_splt3_1=splt1,0,var(1); // x->0, y->a |
---|
376 | //----- rechne das primitive Element von altring in das von neuring um: --- |
---|
377 | ideal convid=maxideal(1); |
---|
378 | convid[1]=nach_splt3_1(primit)[2]; |
---|
379 | poly new_b=nach_splt3_1(primit)[3]; |
---|
380 | map convert=splt3,convid; |
---|
381 | zwi=convert(zwi); |
---|
382 | setring neuring; |
---|
383 | erg=imap(splt3,zwi); |
---|
384 | erg[size(erg)+1]=imap(splt3,new_b); |
---|
385 | } |
---|
386 | } |
---|
387 | if (defined(erg)){export erg;} |
---|
388 | return(neuring); |
---|
389 | } |
---|
390 | example |
---|
391 | { "EXAMPLE:"; echo = 2; |
---|
392 | ring r=0,(x,y),dp; |
---|
393 | def r1=splitring(x2-2); |
---|
394 | setring r1; basering; // change to Q(sqrt(2)) |
---|
395 | // change to Q(sqrt(2),sqrt(sqrt(2)))=Q(a) and return the transformed |
---|
396 | // old parameter: |
---|
397 | def r2=splitring(x2-a,a); |
---|
398 | setring r2; basering; erg; |
---|
399 | // the result is (a)^2 = (sqrt(sqrt(2)))^2 |
---|
400 | kill r1; kill r2; |
---|
401 | } |
---|
402 | /////////////////////////////////////////////////////////////////////////////// |
---|