# source:git/Singular/LIB/primitiv.lib@229948

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