source: git/Singular/LIB/groups.lib @ 6d37e8

spielwiese
Last change on this file since 6d37e8 was 5b3302c, checked in by Gerhard Pfister <pfister@…>, 23 years ago
lib to solve a problem in finite group theory git-svn-id: file:///usr/local/Singular/svn/trunk@5666 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 28.8 KB
Line 
1// $Id: groups.lib,v 1.1 2001-11-05 16:04:11 pfister Exp $
2//(GP, last modified 26.10.01)
3///////////////////////////////////////////////////////////////////////////////
4version="$Id: groups.lib,v 1.1 2001-11-05 16:04:11 pfister Exp $";
5category="Applications";
6info="
7LIBRARY:  groups.lib      Finite Group Theory
8AUTHORS:  Gert-Martin Greuel, email: greuel@mathematik.uni-kl.de @*
9          Gerhard Pfister, email: pfister@mathematik.uni-kl.de
10
11PROCEDURES:
12
13   noSolution(I)         I an ideal in a polynomial ring over Z[x1..xn]
14                         returns a list l of primes <=32003 such that 1
15                         is in IZ/p[x1..xn] for p not in l or an ERROR
16                         if this is wrong or not decided.
17"
18/*
19// ===================== A Problem in Finite Group Theory ===================
20// Posed by Boris Kunyavskii, Bar-Ilan University, Tel Aviv
21//
22// For any word w in X,Y,X^(-1),Y^(-1) consider the sequence U_n
23// of words (depending on w) inductively
24//           U_1 = w
25//         U_n+1 = [X*U_n*X^(-1),Y*U_n*Y^(-1)]
26// with
27//         [X,Y] = X*Y*X^(-1)*Y^(-1)  (the commutator)
28
29
30// Conjecture (1):
31// (by B. Plotkin, slightly modified):
32// A finite group G is solvable <==>
33// there is an n >= 1 such that U_n(x,y) = 1 for all x,y in G
34// and for any of the 4 following words
35//         w1 = X^(-1)*Y*X*Y^(-1)*X
36//         w2 = X^(-2)*Y^(-1)*X
37//         w3 = Y^(-2)*X^(-1)*X
38//         w4 = X*Y^(-2)*X^(-1)*Y*X^(-1)
39//
40// (These words remained as possibly good words by a computer search
41// through about 10 000 candidates, by Fritz Grunewald)
42
43
44//  ==>  is clear
45// The minimal finite non-solvable groups (i.e. every proper subgroup is
46// solvable) have been classified by Thompson in 1968, they are:
47//
48// 1. PSL(2,p)    (p=5 or p=+-2 (mod 5), p!=3)
49// 2. PSL(2,2^p)
50// 3. PSL(2,3^p)  (p odd)
51// 4. Sz(2^p)     (p odd)      (Suzuki group)
52// 5. PSL(3,3)
53//
54// In view of this result Conjecture (1) is equivalent to
55//
56// Conjecture (2):
57// Let G be one of the groups above, then, for at least one of the 4 words w
58// above there are x,y in G such that 1 != U_n(x,y) = U_n+1(x,y) for some n.
59// (then U_n(x,y) != 1 for all n by definitionof U_n)
60
61
62// We give a computer aided proof, using SINGULAR, of Conjecture (2) for
63// the groups 1. - 3. (up to checking for a small, explicit number of primes)
64// Hence only the Suzuki groups have to be checked.
65//
66// We show:  1 != U_1(x,y) = U_2(x,y) for word w1 and some x,y in G
67// We need:
68// Theory
69//   - simple facts from algebraic geometry, singularity theory, finite fields
70//   - theorem of Lang-Weil, estimating the number of rational points on an
71//     absolutely irreducible projective curve C defined over Z (g=genus):
72//     if q >= N(g) ==> #C(F_q) != 0
73//     (N(g) = 4g^2 -2, g arithmetic genus)
74// SINGULAR
75//   - Groebner basis (elimination), multivariate factorization, resolution of
76//     plane curve singularities (Hamburger Noether developement), primary
77//     decomposition
78
79*/
80
81LIB "standard.lib";
82LIB "general.lib";
83LIB "matrix.lib";
84///////////////////////////////////////////////////////////////////////////////
85
86static proc splitS1(ideal I,int s,list #)
87"USAGE: splitS1(I,s[,l]) I ideal, s integer, l list of ideals
88COMPUTE: Factorizes the generators of I. If for one generator f=f1*f2 and
89         fi=ni*gi^ri, ni an integer and gi a polynomial, then
90         compute a standardbasis of I1=I,g1 and I2=I,g2, if this can be
91         done in < s seconds. Then apply splitS to I1 and I2 and continue
92         in the same way. The procedure stops if no generator can be
93         factorized.
94RETURN:  A list L of ideals and prime numbers
95         L[1]: A list of ideals such that the radical of the intersection
96               of these ideals coincides with the radical of I
97               If the optional list l of the input is not empty then
98               the ideals of L[1] which contain an ideal of l are canceled.
99         L[2]: A list of prime numbers appearing as factors of the ni.
100NOTE:   The computation avoids division by integers (by using
101        option(contentSB)) hence the result is correct modulo any prime
102        number which does not appear in the list L[2].
103EXAMPLE: example splitS1; shows an example
104"
105{
106   option(redSB);
107   option(contentSB);
108   int j,k,e;
109   int i=1;
110   int l=attrib(I,"isSB");
111   ideal J;
112   list re,fac,te,pr,qr,w;
113   number n;
114   poly p;
115   re=#;
116     
117   if(deg(I[1])==0){return(list(re+list(std(I)),qr));}
118 
119   fac=factorize(I[1]);
120 
121   while((size(fac[1])==2)&&(i<size(I)))
122   {
123      I[i]=fac[1][2]*fac[1][1];   //not in splitS
124      i++;
125      fac=factorize(I[i]);
126   }
127   if(size(fac[1])==2){I[size(I)]=fac[1][2]*fac[1][1];}  //not in splitS
128   if(size(fac[1])>2)
129   {
130      w=squarefreeP(number(fac[1][1]));
131      n=w[1];
132      qr=insResult(qr,w[2]);
133      for(j=2;j<=size(fac[1]);j++)
134      {
135         I[i]=fac[1][j];
136         attrib(I,"isSB",1);
137         e=1;
138         k=0;
139         while(k<size(re))
140         {
141            k++;
142            if(size(reduce(re[k],I))==0){e=0;break;}
143            attrib(re[k],"isSB",1);
144            if(size(reduce(I,re[k]))==0){re=delete(re,k);k--;}
145         }
146         if(j==2){I[i]=I[i]*n;}
147         if(e)
148         {
149            if(l)
150            {
151               J=I;
152               p=I[i];
153               J[i]=0;
154               J=simplify(J,2);
155               attrib(J,"isSB",1);
156               pr=splitS(std(J,p),s,re);
157               re=pr[1];
158               qr=insResult(qr,pr[2]);
159            }
160            else
161            {
162               J=interred(I);
163               pr=splitS(timeStd(J,s),s,re);
164               re=pr[1];
165               qr=insResult(qr,pr[2]);
166            }
167         }
168      }
169      return(list(re,qr));
170   }
171   J=timeStd(I,s);       //J=std(I) in splitS
172   attrib(I,"isSB",1);
173   if(size(reduce(J,I))==0){return(list(re+list(I),qr));}
174   pr=splitS(J,s,re);
175   return(list(re+pr[1],pr[2]));
176}
177example
178{ "EXAMPLE:"; echo = 2;
179   ring r=0,(b,s,t,u,v,w,x,y,z),dp;
180   ideal i=
181   bv+su,
182   bw+tu,
183   sw+tv,
184   by+sx,
185   bz+tx,
186   sz+ty,
187   uy+vx,
188   uz+wx,
189   vz+wy,
190   bvz;
191   splitS1(i,5);
192}
193
194///////////////////////////////////////////////////////////////////////////////
195static proc splitS(ideal I,int s,list #)
196"USAGE: splitS(I,s[,l]) I ideal, s integer, l list of ideals
197COMPUTE: Factorizes the generators of I. If for one generator f=f1*f2 and
198         fi=ni*gi^ri, ni an integer and gi a polynomial, then
199         compute a standardbasis of I1=I,g1 and I2=I,g2, if this can be
200         done in < s seconds. Then apply splitS to I1 and I2 and continue
201         in the same way. The procedure stops if no generator can be
202         factorized.
203RETURN:  A list L of ideals and prime numbers
204         L[1]: A list of ideals such that the radical of the intersection
205               of these ideals coincides with the radical of I
206               If the optional list l of the input is not empty then
207               the ideals of L[1] which contain an ideal of l are canceled.
208         L[2]: A list of prime numbers appearing as factors of the ni.
209NOTE:   The computation avoids division by integers (by using
210        option(contentSB)) hence the result is correct modulo any prime
211        number which does not appear in the list L[2].
212EXAMPLE: example splitS; shows an example
213"
214{
215   option(redSB);
216   option(contentSB);
217   int j,k,e;
218   int i=1;
219   int l=attrib(I,"isSB");
220   ideal J;
221   list re,fac,te,pr,qr,w;
222   number n;
223   poly p;
224   re=#;
225   
226   if(deg(I[1])==0){return(list(re+list(std(I)),qr));}
227 
228   fac=factorize(I[1]);
229 
230   while((size(fac[1])==2)&&(i<size(I)))
231   {
232      i++;
233      fac=factorize(I[i]);
234   }
235   if(size(fac[1])>2)
236   {
237      w=squarefreeP(number(fac[1][1]));
238      n=w[1];
239      qr=insResult(qr,w[2]);
240      for(j=2;j<=size(fac[1]);j++)
241      {
242         I[i]=fac[1][j];
243         attrib(I,"isSB",1);
244         e=1;
245         k=0;
246         while(k<size(re))
247         {
248            k++;
249            if(size(reduce(re[k],I))==0){e=0;break;}
250            attrib(re[k],"isSB",1);
251            if(size(reduce(I,re[k]))==0){re=delete(re,k);k--;}
252         }
253         if(j==2){I[i]=I[i]*n;}
254         if(e)
255         {
256            if(l)
257            {
258               J=I;
259               p=I[i];
260               J[i]=0;
261               J=simplify(J,2);
262               attrib(J,"isSB",1);
263               pr=splitS(std(J,p),s,re);
264               re=pr[1];
265               qr=insResult(qr,pr[2]);
266            }
267            else
268            {
269               J=interred(I);
270               pr=splitS(timeStd(J,s),s,re);
271               re=pr[1];
272               qr=insResult(qr,pr[2]);
273            }
274         }
275      }
276      return(list(re,qr));
277   }
278   J=std(I);
279   attrib(I,"isSB",1);
280   if(size(reduce(J,I))==0){return(list(re+list(I),qr));}
281   pr=splitS(J,s,re);
282   return(list(re+pr[1],pr[2]));
283}
284example
285{ "EXAMPLE:"; echo = 2;
286   ring r=0,(b,s,t,u,v,w,x,y,z),dp;
287   ideal i=
288   bv+su,
289   bw+tu,
290   sw+tv,
291   by+sx,
292   bz+tx,
293   sz+ty,
294   uy+vx,
295   uz+wx,
296   vz+wy,
297   bvz;
298   splitS(i,5);
299}
300
301///////////////////////////////////////////////////////////////////////////////
302static proc finalSplit(list I,list qr)
303"USAGE: finalSplit(I,qr) I list of ideals, qr list of primes
304RETURN: list l of primes <=32003 such that 1 is in I[j]Z/p[x1..xn]
305        for p not in l and all j or an ERROR if this is wrong or
306        not decided.
307EXAMPLE: example finalSplit; shows an example
308"
309{
310   option(redThrough);
311   option(contentSB);
312   ideal J,K;
313   int i,j,k,n;
314   int count=1;
315   list q,pr,l;
316
317   "trivial splitting";
318   pr=trivialSplit(I,4);
319   l=pr[1];
320   qr=insResult(qr,pr[2]);
321   while(count)
322   {
323      k++;
324      "loop";k;size(l);"======";
325      count=0;list p;
326      for(i=1;i<=size(l);i++)
327      {
328         i;
329         K=changeOrdTest(l[i]);
330         if(deg(K[1])!=0)
331         {
332            "split";
333            pr=splitS1(shortid(K,3),2);
334            q=pr[1];
335            qr=insResult(qr,pr[2]);
336            size(q);"out";
337            for(j=1;j<=size(q);j++)
338            {
339               if(deg(q[j][1])==0)
340               {
341                  pr=contentS(q[j]);
342                  q[j]=pr[1];
343                  qr=insResult(qr,pr[2]);
344               }
345               else
346               {
347                  pr=simpliFy(q[j]+K);
348                  q[j]=pr[1];
349                  qr=insResult(qr,pr[2]);
350                  pr=contentS(q[j]);
351                  q[j]=pr[1];
352                  qr=insResult(qr,pr[2]);
353               }
354            }
355            if((size(q)==1)&&(deg(q[1][1])>0))
356            {
357               "split again";
358               K=sort(q[1])[1];
359               J=K[1..20+4*k];
360               pr=splitS(J,5);
361               q=pr[1];
362               qr=insResult(qr,pr[2]); 
363               size(q);"out";
364               for(j=1;j<=size(q);j++)
365               {
366                  if(deg(q[j][1])==0)
367                  {
368                     pr=contentS(q[j]);
369                     q[j]=pr[1];
370                     qr=insResult(qr,pr[2]);
371                  }
372                  else
373                  {
374                     pr=simpliFy(q[j]+K);
375                     q[j]=pr[1];
376                     qr=insResult(qr,pr[2]);
377                     pr=contentS(q[j]);
378                     q[j]=pr[1];
379                     qr=insResult(qr,pr[2]);
380                  }
381               }
382 
383            }
384            if((size(q)==1)&&(deg(q[1][1])>0))
385            {
386               n++;
387               "split again";
388               K=q[1];
389               J=shortid(K,3+n);
390               pr=splitS1(J,5);
391               if(size(pr[1])>1)
392               {
393                  q=pr[1];             
394                  qr=insResult(qr,pr[2]);
395                  size(q);"out";
396                  for(j=1;j<=size(q);j++)
397                  {
398                     if(deg(q[j][1])==0)
399                     {
400                       pr=contentS(q[j]);
401                       q[j]=pr[1];
402                       qr=insResult(qr,pr[2]);
403                     }
404                     else
405                     {
406                       pr=simpliFy(q[j]+K);
407                       q[j]=pr[1];
408                       qr=insResult(qr,pr[2]);
409                       pr=contentS(q[j]);
410                       q[j]=pr[1];
411                       qr=insResult(qr,pr[2]);
412                     }
413                  }
414               }
415            }
416            if(size(q)>1)
417            {
418                count++;
419            }
420            else
421            {
422               if(deg(q[1][1])>0)
423               {
424                  q[1]=changeOrdTest(q[1]);
425               }
426            }
427            p=p+q;
428         }
429         else
430         {
431           pr=primefactors(number(K[1]));
432           qr=insResult(qr,pr[1]);
433         }
434      }
435      l=p;
436      kill p;
437   }
438l;
439   for(i=1;i<=size(l);i++)
440   {
441     
442      if(deg(l[i][1])>0){l;ERROR("not ready");}
443      pr=primefactors(number(l[i][1]));
444      if(pr[3]!=1){pr;ERROR("not ready");}
445      qr=insResult(qr,pr[1]);
446   }
447   return(qr);
448}
449example
450{ "EXAMPLE:"; echo = 2;
451   ring r=0,x,dp;
452   list qr=2,5,7;
453   ideal i=181x-181,11x+11;
454   list pr=i;
455   finalSplit(pr,qr);
456}
457
458///////////////////////////////////////////////////////////////////////////////
459proc noSolution(ideal I)
460"USAGE: noSolution(I) I ideal
461RETURN: list l of primes <=32003 such that 1 is in IZ/p[x1..xn]
462        for p not in l or an ERROR if this is wrong or not decided.
463EXAMPLE: example noSolution; shows an example
464"
465{
466   int s=1;
467   int t=30;
468   option(redThrough);
469   option(contentSB);
470   ideal J=shortid(I,3);
471   ideal K;
472   number n;
473
474   "first splitting";
475   int i,j,k;
476   list l,p,q,re,pr,qr;
477   pr=splitS(J,s);
478   l=pr[1];
479   qr=pr[2];
480   size(l);
481
482   for(i=1;i<=size(l);i++)
483   {
484      if(deg(l[i][1])==0)
485      {
486         K=l[i][1];
487      }
488      else
489      {
490         pr=simpliFy(I+l[i]);
491         J=pr[1];
492         qr=insResult(qr,pr[2]);
493         pr=contentS(J);
494         J=pr[1];
495         qr=insResult(qr,pr[2]);
496         K=timeStd(J,1);
497      }
498      if(deg(K[1])==0)
499      {
500         n=number(K[1]);
501         pr=primefactors(n);
502         if(pr[3]!=1)
503         {
504             J=changeOrdTest(J);
505             p=p+list(J);
506         }
507         else
508         {
509             qr=insResult(qr,pr[1]);
510         }
511      }
512      else
513      {
514         pr=contentS(K);
515         p=p+list(pr[1]);
516         qr=insResult(qr,pr[2]);
517      }
518   }
519   "trivial splitting 1";
520   pr=trivialSplit(p,3);
521   p=pr[1];
522   qr=insResult(qr,pr[2]);
523   size(p);
524   "second splitting";
525
526   for(i=1;i<=size(p);i++)
527   {
528      i;
529      if(size(p[i])<t)
530      {
531         K=timeStd(p[i],10);
532         if(deg(K[1])==0)
533         {
534            n=number(K[1]);
535            pr=primefactors(n);
536            if(pr[3]!=1)
537            {
538                p[i]=changeOrdTest(p[i]);
539                re=re+list(p[i]);
540            }
541            else
542            {
543               qr=insResult(qr,pr[1]);
544            }
545         }
546         else
547         {
548            pr=contentS(K);
549            re=re+list(pr[1]);
550            qr=insResult(qr,pr[2]);
551         }
552      }
553      else
554      {
555         J=p[i];
556         J=J[1..t];
557         "in splitting";
558         pr=splitS(J,s);
559         l=pr[1];
560         qr=insResult(qr,pr[2]);
561         size(l);
562         "out";
563         for(k=1;k<=size(l);k++)
564         {
565            pr=contentS(l[k]);
566            l[k]= pr[1];
567            qr=insResult(qr,pr[2]);
568            pr=simpliFy(l[k]+p[i]);
569            J=pr[1];
570            qr=insResult(qr,pr[2]);
571            K=timeStd(J,10);
572            if(deg(K[1])==0)
573            {
574               n=number(K[1]);
575               pr=primefactors(n);
576               if(pr[3]!=1)
577               {
578                  pr=contentS(J);
579                  qr=insResult(qr,pr[2]);
580                  J=pr[1];
581                  J=changeOrdTest(J);
582                  re=re+list(J);
583               }
584               else
585               {
586                  qr=insResult(qr,pr[1]);
587               }
588            }
589            else
590            {
591               pr=contentS(K);
592               qr=insResult(qr,pr[2]);
593               re=re+list(pr[1]);
594            }
595         }
596      }
597   }
598   "trivial splitting 2";
599   size(re);
600   pr=trivialSplit(re,2);
601   re=pr[1];
602   qr=insResult(qr,pr[2]);
603   "third splitting";
604   size(re);
605   for(i=1;i<=size(re);i++)
606   {
607      if(deg(re[i][1])>0)
608      {
609         i;
610         pr=simpliFy(re[i]);
611         J=pr[1];
612         qr=insResult(qr,pr[2]);
613         "in splitting";
614         J=shortid(J,3);
615         pr=splitS(J,s);
616         l=pr[1];
617         if(size(l)==1)
618         {
619            "split again";
620            pr=simpliFy(re[i]);
621            J=pr[1];
622            qr=insResult(qr,pr[2]);
623            J=J[1..t];
624            pr=splitS(J,s);
625            l=pr[1];
626            qr=insResult(qr,pr[2]);
627         }
628         else
629         {
630            qr=insResult(qr,pr[2]);
631         }
632         size(l);
633         "out";
634         for(j=1;j<=size(l);j++)
635         {
636            pr=contentS(l[j]);
637            l[j]=pr[1];
638            qr=insResult(qr,pr[2]);
639            pr=simpliFy(re[i]+l[j]);
640            J=pr[1];
641            qr=insResult(qr,pr[2]);
642            K=timeStd(J,10);
643            if(deg(K[1])==0)
644            {
645               n=number(K[1]);
646               pr=primefactors(n);
647               if(pr[3]!=1)
648               {
649                  pr=contentS(J);
650                  J=pr[1];
651                  qr=insResult(qr,pr[2]);
652                  J=changeOrdTest(J);
653                  q=q+list(J);
654               }
655               else
656               {
657                 qr=insResult(qr,pr[1]);
658               }
659            }
660            else
661            {
662               pr=contentS(K);
663               qr=insResult(qr,pr[2]);
664               q=q+list(pr[1]);
665            }
666         }
667      }
668      else
669      {
670         i;
671         pr=primefactors(number(re[i][1]));
672         qr=insResult(qr,pr[1]);
673      }
674   }
675"jetzt geht es los";
676   return(finalSplit(q,qr));
677}
678example
679{ "EXAMPLE:"; echo = 2;
680   ring r=0,x,dp;
681   ideal i=181x-181,11x+11;
682   noSolution(i);
683}
684
685///////////////////////////////////////////////////////////////////////////////
686
687static proc changeOrdTest(ideal I)
688{
689   def R=basering;
690   if(deg(I[1])==0){return(I);}
691   changeord("RR","dp");
692   ideal I=imap(R,I);
693   ideal K=timeStd(I,5);
694   if(deg(K[1])==0)
695   {
696      number n=number(K[1]);
697      if(primefactors(n)[3]==1)
698      {
699         I=K;
700      }
701   }
702   setring R;
703   I=imap(RR,I);
704   kill RR;
705   return(I);
706}
707///////////////////////////////////////////////////////////////////////////////
708
709static proc trivialSplit(list p,int depth, list #)
710"USAGE: trivialSplit(p,s) p list of ideals, s integer the number of iterations
711COMPUTE: Factorizes the monomials among the generators of I.
712         If  one monomial contains the variables x1,..xr, then
713         I1=I(x1=0),...,Ir=I(xr=0) is considered. Then apply
714         trivialSplit to I1 ...Ir and continue in the same way s times.
715RETURN:  A list L of ideals and prime numbers
716         L[1]: A list of ideals such that the radical of the intersection
717               of these ideals at x1=...=xr=0 coincides with the radical of
718               I at x1=...=xr=0.
719         L[2]: A list of prime numbers appearing as factors of the monomials.
720NOTE:   The computation avoids division by integers  hence the result
721        is correct modulo any prime
722        number which does not appear in the list L[2].
723EXAMPLE: example trivialSplit; shows an example
724"
725{
726   list re,l,pr,qr;
727   int i,k;
728   ideal J,K,T,Ke;
729   number n;
730   
731   if(size(p)==0){return(list(p,qr));}
732   if(depth<=0){return(list(p,qr));}
733   if(size(#)>0)
734   { 
735      T=#[1];
736   }
737   else
738   {
739      T=1;
740   }
741   for(k=1;k<=size(p);k++)
742   {
743      pr=simpliFy(p[k]);
744      p[k]=pr[1];
745      qr=insResult(qr,pr[2]);
746      J=shortid(p[k],1);
747      if((size(J)>0)&&(deg(J[1])>=1))
748      {
749         pr=splitS1(J,10);
750         l=pr[1];
751         qr=insResult(qr,pr[2]);
752         for(i=1;i<=size(l);i++)
753         {
754            Ke=l[i];
755            l[i]=trivialSimplify(p[k],l[i]);
756            pr=simpliFy(l[i]);
757            l[i]=pr[1];
758            qr=insResult(qr,pr[2]);
759            K=timeStd(l[i],1);
760            attrib(K,"isSB",1);
761            if(deg(K[1])==0)
762            {
763               n=number(K[1]);
764               if(primefactors(n)[3]!=1)
765               {
766                  l[i]=changeOrdTest(l[i]);
767                  pr=trivialSplit(l[i],depth-1,T);
768                  re=re+pr[1];
769                  qr=insResult(qr,pr[2]);
770               }
771               else
772               {
773                  l[i]=K; //neu
774               }
775            }
776            else
777            {
778               if(size(reduce(trivialSimplify(T,Ke),K))!=0)
779               {
780                  pr=trivialSplit(K,depth-1,trivialSimplify(T,Ke));
781                  re=re+pr[1];
782                  qr=insResult(qr,pr[2]);
783               }
784            }
785            T=intersect(T,Ke);
786         }
787      }
788      else
789      {
790         J=timeStd(p[k],5);
791         if(deg(J[1])==0)
792         {
793            n=number(J[1]);
794            if(primefactors(n)[3]!=1)
795            {
796               p[k]=changeOrdTest(p[k]);
797               re=re+list(p[k]);
798            }
799            else
800            {
801                re=re+list(J);
802            }
803         }
804         else
805         {
806            re=re+list(J);
807         }
808      }
809   }
810   return(list(re,qr));
811}
812example
813{ "EXAMPLE:"; echo = 2;
814   ring r=0,(b,s,t,u,v,w,x,y,z),dp;
815   ideal i=
816   bv+su,
817   bw+tu,
818   sw+tv,
819   by+sx,
820   bz+tx,
821   sz+ty,
822   uy+vx,
823   uz+wx,
824   vz+wy,
825   bvz;
826   trivialSplit(i,2);
827}
828
829///////////////////////////////////////////////////////////////////////////////
830
831static proc trivialSimplify(ideal I, ideal J)
832"USAGE:  trivialSimplify(I,J);  I,J ideals, J generated by variables
833RETURN:  ideal K = eliminate(I,m) m the product of variables of J
834EXAMPLE: example trivialSimplify; shows an example
835"
836{
837   int i;
838   for(i=1;i<=size(J);i++){I=subst(I,J[i],0);}
839   return(simplify(I,2));
840}
841example
842{ "EXAMPLE:"; echo = 2;
843   ring r=0,(x,y,z,w,t,u),dp;
844   ideal i=
845   t,u,
846   x3+y2+2z,
847   x2+3y,
848   x2+y2+z2,
849   w2+z+u;
850   ideal j=t,u;
851   trivialSimplify(i,j);
852}
853
854///////////////////////////////////////////////////////////////////////////////
855
856static proc simpliFy(ideal J,list #)
857"USAGE:  simpliFy(id);  id ideal
858RETURN:  ideal I = eliminate(Id,m) m is a product of variables
859         which are in linear equations involved in the generators of id
860EXAMPLE: example simpliFy; shows an example
861"
862{
863   ideal I=J;
864   if(size(#)!=0){I=#[1];}
865   def R=basering;
866   matrix M=jacob(I);
867   ideal ma=maxideal(1);
868   int i,j,k;
869   map phi;
870   list pr,re;
871
872   for(i=1;i<=nrows(M);i++)
873   {
874      for(j=1;j<=ncols(M);j++)
875      {
876         if((deg(M[i,j])==0)&&(deg(I[i])==1))       
877         {
878            ma[j]=(-1/M[i,j])*(I[i]-M[i,j]*var(j));
879            pr=primefactors(number(M[i,j]));
880            if(pr[3]>1){pr[3];ERROR("number too big in simpliFy");}
881            re=insResult(re,pr[1]);
882            phi=R,ma;
883            I=phi(I);
884            J=phi(J);
885            for(k=1;k<=ncols(I);k++)
886            {
887               I[k]=cleardenom(I[k]);
888            }
889            M=jacob(I);
890         }
891      }
892   }
893   J=simplify(J,2);
894   for(i=1;i<=size(J);i++)
895   {
896      if(deg(J[i])==0){J=std(J);break;}
897      J[i]=cleardenom(J[i]);
898   }
899   return(list(J,re));
900}
901example
902{ "EXAMPLE:"; echo = 2;
903   ring r=0,(x,y,z,w),dp;
904   ideal i=
905   x3+y2+2z,
906   x2+3y,
907   x2+y2+z2,
908   w2+z;
909   simpliFy(i);
910}
911///////////////////////////////////////////////////////////////////////////////
912
913static proc squarefreeP(number n)
914"USAGE:  squarefreeP(n);  n number
915RETURN:  list l of numbers. l[2] contains all prime factors of n less
916         then 32003, l[1] is the part of m which has prime factors greater
917         then 32003
918EXAMPLE: example squarefreeP; shows an example
919"
920{
921   list re;
922   if(n<0){n=-n;}
923   if(n==1){return(list(n,re));}
924   list pr=primefactors(n);
925   int i;
926   number m=number(pr[3][1]);
927   re=insResult(re,pr[1]);
928   return(list(m,re));
929}
930example
931{ "EXAMPLE:"; echo = 2;
932   ring r = 0,x,dp;
933   squarefreeP(123456789100);
934
935
936///////////////////////////////////////////////////////////////////////////////
937
938static proc contentS(ideal I)
939"USAGE:  contentS(I);  I ideal
940RETURN:  list l
941         l[1] the ideal I. Te generators of I are the generators of the
942         input ideal devided by the part of the content which have
943         prime factors less then 32003.
944         l[2] contains the prime numbers which occured in the division
945EXAMPLE: example contentS; shows an example
946"
947
948{
949   option(contentSB);
950   int i,k;
951   number n,m;
952   list pr,re;
953   for(i=1;i<=size(I);i++)
954   {
955      n=leadcoef(I[i]);
956      if(n<0){n=-n;}
957      if(n>1)
958      {
959         pr=primefactors(n);
960         if(n<=32003)
961         {
962            m=n;
963         }
964         else
965         {
966            m=number(pr[1][1])^pr[2][1];
967            for(k=2;k<=size(pr[1]);k++)
968            {
969               m=m*number(pr[1][k])^pr[2][k];
970            }
971         }
972         I[i]=cleardenom(I[i]/m);
973         re=insResult(re,pr[1]); 
974      }
975   }
976   return(list(I,re)); 
977}
978example
979{ "EXAMPLE:"; echo = 2;
980   ring r = 0,(x,y),dp;
981   ideal I=2x+2,3y+3x,1234567891x+1234567891;
982   contentS(I);
983
984
985///////////////////////////////////////////////////////////////////////////////
986
987static proc insResult(list r,#)
988{
989   if(size(#)==0){return(r);}
990   if(size(r)==0){r[1]=1;}
991   int i,j;
992   for(i=1;i<=size(#);i++)
993   {
994      j=1;
995      while((#[i]>r[j])&&(j<size(r))){j++;}
996      if(#[i]>r[j]){r=insert(r,#[i],j);}
997      if(#[i]<r[j]){r=insert(r,#[i],j-1);}     
998   }
999   return(r);
1000}
1001example
1002{ "EXAMPLE:"; echo = 2;
1003   list r=list(3,7,13);
1004   intvec v=2,3,5,11,17;
1005   insResult(r,v);
1006}
1007///////////////////////////////////////////////////////////////////////////////
1008
1009static proc shortid (ideal id,int n)
1010{
1011  ideal Lc;
1012  intvec v;
1013  int ii;
1014  for(ii= 1; ii<=ncols(id); ii++)
1015  {
1016   if (size(id[ii])<=n and id[ii]!=0)
1017   {
1018       Lc = Lc,id[ii];
1019       v=v,ii;
1020   }
1021  }
1022  Lc = simplify(Lc,2);
1023  list L = Lc,v;
1024  return(Lc);
1025}
1026
1027
1028/*
1029//-------------- the solution of the problem -----------------
1030ring r = 0,(b,c,t),dp;   //global (affine) ring
1031matrix X[2][2] = 0, -1,
1032                 1, t;
1033matrix Y[2][2] = 1, b,
1034                 c, 1+bc;
1035;
1036
1037//------------ create word W1 and ideal I of U1-U2 ------------
1038matrix iX = inverse(X);
1039matrix iY = inverse(Y);
1040matrix U1 = iX*Y*X*iY*X;          //the word w1
1041//matrix U1=iX*iX*iY*X;           //das neue Wort
1042matrix N  = X*U1*iX;
1043matrix M  = Y*U1*iY;
1044matrix iN = inverse(N);
1045matrix iM = inverse(M);
1046matrix U2 = N*M*iN*iM;
1047ideal I = ideal(U1-U2);     
1048
1049//list qr=primdecGTZ(I);
1050//I=qr[1][2];
1051
1052ring rh = 0,(b,c,t,h),dp;
1053ideal I = imap(r,I);
1054ideal sI = groebner(I);
1055ideal hI = homog(sI,h);
1056ideal shI =std(hI);   
1057ideal J = eliminate(shI,c);       //eliminate c
1058poly f = J[1];
1059
1060ring r1 = 0,(b,t,h),dp;
1061poly hf=b6t6-2b5t7+b4t8+b8t3h-4b7t4h+7b6t5h-6b5t6h+2b4t7h-7b6t4h2+12b5t5h2+b4t6h2-6b3t7h2-3b8th3+12b7t2h3-16b6t3h3+19b4t5h3-12b3t6h3-2b8h4+8b7th4-3b6t2h4+2b5t3h4-45b4t4h4+32b3t5h4+12b2t6h4-12b6th5+50b5t2h5-64b4t3h5+4b3t4h5+26b2t5h5-12b6h6+24b5th6+22b4t2h6+12b3t3h6-73b2t4h6-10bt5h6-8b5h7-6b4th7+88b3t2h7-68b2t3h7-26bt4h7-29b4h8+16b3th8+42b2t2h8+54bt3h8+3t4h8-28b3h9-12b2th9+88bt2h9+10t3h9-38b2h10-8bth10-11t2h10-28bh11-34th11-17h12;                                               
1062//poly hf=b3t4-b2t5+b4t2h-2b3t3h+2b2t4h-bt5h-2b3t2h2+4bt4h2+b2t2h3-bt3h3+t4h3
1063//+2b2th4-6bt2h4-4t3h4+b2h5-2bth5+2t2h5+4th6+h7;
1064
1065int n,m,sA,sB;
1066 n=6;
1067//n=3;
1068//m=7-n;
1069
1070 m = 12-n;
1071 ideal A = maxideal(n);       ideal B = maxideal(m);     
1072 sA =size(A);
1073 sB = size(B);                                                               
1074             
1075ring R = 0,(b(1..sB),a(1..sA),b,t,h),dp;
1076poly hf = imap(r1,hf);
1077ideal A = imap(r1,A);
1078ideal B = imap(r1,B);
1079matrix aa[sA][1]= a(1..sA);
1080matrix bb[sB][1] = b(1..sB);
1081poly f1= (matrix(A)*aa)[1,1];   //Ansatz for degree n
1082poly f2= (matrix(B)*bb)[1,1];   //Ansatz for degree m=12-n
1083poly g = f1*f2-hf;              //assume hf factors
1084matrix C = coef(g,bth);
1085ideal co = submat(C,2,1..ncols(C));//condition for decomposition, size 91
1086 
1087ring R1 = 0,(b(1..sB),a(1..sA)),lp;
1088ideal co = imap(R,co);                       
1089co=subst(co,a(1),1);          //OE a1=1
1090co = subst(co,b(1),-17);      //co1[91]=b(1)+17
1091//co = subst(co,b(1),1);
1092
1093int aa=timer;
1094list pr=noSolution(co);
1095timer-aa;
1096
1097
1098//n=1 0 sec
1099// keine Primzahl
1100//n=2   628 sec
1101// 2,3,5,7,11,13,17,19,23,29,31,37,43,47,61,89,293,347,367,487,491,3463,7498
1102//n=3  1604 sec
1103// 2,3,5,7,11,13,17,19,23,29,31,59,71,79,101,163,197,211,269,281,541,647,863,
1104// 7129,9041,10343,18413,20857
1105//n=4  8296 sec
1106// 2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,61,67,71,73,83,89,101,109,127,
1107   131,149,151,157,163,167,173,179,193,199,211,223,229,251,263,283,307,313,331,
1108   347,373,401,601,607,631,643,701,719,727,829,857,887,971,1279,1361,1453,1721,
1109   1847,2003,2069,2213,2671,2693,3299,3373,3391,3517,3593,3701,3779,4111,4409,
1110   4423,4561,4657,4793,5273,5399,5659,5987,6949,7487,8011,8243,8887,9769,10159,
1111   10177,12007,26347
1112//n=5  4715 sec
1113// 2,3,5,7,11,13,17,19,23,31,41,43,61,283,421,631,1609
1114//n=6 18317 sec
1115// 2,3,5,7,11,13,17,19,29,31,37,41,47,61,71,73,79,89,97,127,131,173,181,223,
1116// 269,953,1039,6151,6343,7823
1117
1118
1119//fuer das andere Wort
1120//n=1  keine Primzahl
1121//n=2  2,3,7,109
1122//n=3  2,3,5,7,11,13,19,29,31,41,43,47,89,139,149,167,173,991,1381,1637,27367
1123
1124*/
Note: See TracBrowser for help on using the repository browser.