source: git/Singular/LIB/toric.lib

spielwiese
Last change on this file was 3f7e01a, checked in by Hans Schoenemann <hannes@…>, 9 months ago
remnemae: sat -> sat, sat_with_exp
  • Property mode set to 100644
File size: 27.6 KB
Line 
1///////////////////////////////////////////////////////////////////////////
2version="version toric.lib 4.3.1.3 Feb_2023 "; // $Id$
3category="Commutative Algebra";
4info="
5LIBRARY: toric.lib   Standard Basis of Toric Ideals
6AUTHOR:  Christine Theis, email: ctheis@math.uni-sb.de
7
8PROCEDURES:
9toric_ideal(A,..);   computes the toric ideal of A
10toric_std(ideal I);  standard basis of I by a specialized Buchberger algorithm
11";
12
13///////////////////////////////////////////////////////////////////////////////
14
15static proc toric_ideal_1(intmat A, string alg)
16{
17  ideal I;
18  // to be returned
19
20  // check suitability of actual basering
21  if(nvars(basering)<ncols(A))
22  {
23    ERROR("The number of matrix columns is smaller than the number of ring variables.");
24  }
25
26  // check suitability of actual term ordering
27  // the "module ordering" c or C is ignored
28  string singular_ord=ordstr(basering);
29  string test_ord;
30  string external_ord="";
31  int i,j;
32  intvec weightvec;
33
34  if(find(singular_ord,"lp")==1)
35  {
36    external_ord="W_LEX";
37    for(i=1;i<=nvars(basering);i++)
38    {
39      weightvec[i]=0;
40    }
41    test_ord="lp("+string(nvars(basering))+"),";
42    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
43    {
44      "Warning: Block orderings are not supported; lp used for computation.";
45    }
46  }
47  if(external_ord=="" && find(singular_ord,"lp")==3)
48  {
49    external_ord="W_LEX";
50    for(i=1;i<=nvars(basering);i++)
51    {
52      weightvec[i]=0;
53    }
54    test_ord="lp("+string(nvars(basering))+"),";
55    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
56    {
57      "Warning: Block orderings are not supported; lp used for computation.";
58    }
59  }
60
61  if(external_ord=="" && find(singular_ord,"dp")==1)
62  {
63    external_ord="W_REV_LEX";
64    for(i=1;i<=nvars(basering);i++)
65    {
66      weightvec[i]=1;
67    }
68    test_ord="dp("+string(nvars(basering))+"),";
69    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
70    {
71      "Warning: Block orderings are not supported; dp used for computation.";
72    }
73  }
74  if(external_ord=="" && find(singular_ord,"dp")==3)
75  {
76    external_ord="W_REV_LEX";
77    for(i=1;i<=nvars(basering);i++)
78    {
79      weightvec[i]=1;
80    }
81    test_ord="dp("+string(nvars(basering))+"),";
82    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
83    {
84      "Warning: Block orderings are not supported; dp used for computation.";
85    }
86  }
87
88  if(external_ord=="" && find(singular_ord,"Dp")==1)
89  {
90    external_ord="W_LEX";
91    for(i=1;i<=nvars(basering);i++)
92    {
93      weightvec[i]=1;
94    }
95    test_ord="Dp("+string(nvars(basering))+"),";
96    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
97    {
98      "Warning: Block orderings are not supported; Dp used for computation.";
99    }
100  }
101  if(external_ord=="" && find(singular_ord,"Dp")==3)
102  {
103    external_ord="W_LEX";
104    for(i=1;i<=nvars(basering);i++)
105    {
106      weightvec[i]=1;
107    }
108    test_ord="Dp("+string(nvars(basering))+"),";
109    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
110    {
111      "Warning: Block orderings are not supported; Dp used for computation.";
112    }
113  }
114
115  int pos;
116  string number_string;
117
118  if(external_ord=="" && find(singular_ord,"wp")==1)
119  {
120    external_ord="W_REV_LEX";
121    pos=3;
122    for(i=1;i<=nvars(basering);i++)
123    {
124      pos++;
125      number_string="";
126      while(singular_ord[pos]!=",")
127      {
128        number_string=number_string+singular_ord[pos];
129        pos++;
130      }
131      weightvec[i]=`number_string`;
132    }
133    test_ord="wp("+string(weightvec)+"),";
134    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
135    {
136      "Warning: Block orderings are not supported; wp("+string(weightvec)+") used for computation.";
137    }
138  }
139  if(external_ord=="" && find(singular_ord,"wp")==3)
140  {
141    external_ord="W_REV_LEX";
142    pos=5;
143    for(i=1;i<=nvars(basering);i++)
144    {
145      pos++;
146      number_string="";
147      while(singular_ord[pos]!=",")
148      {
149        number_string=number_string+singular_ord[pos];
150        pos++;
151      }
152      weightvec[i]=`number_string`;
153    }
154    test_ord="wp("+string(weightvec)+"),";
155    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
156    {
157      "Warning: Block orderings are not supported; wp("+string(weightvec)+") used for computation.";
158    }
159  }
160
161  if(external_ord=="" && find(singular_ord,"Wp")==1)
162  {
163    external_ord="W_LEX";
164    pos=3;
165    for(i=1;i<=nvars(basering);i++)
166    {
167      pos++;
168      number_string="";
169      while(singular_ord[pos]!=",")
170      {
171        number_string=number_string+singular_ord[pos];
172        pos++;
173      }
174      weightvec[i]=`number_string`;
175    }
176    test_ord="Wp("+string(weightvec)+"),";
177    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
178    {
179      "Warning: Block orderings are not supported; Wp("+string(weightvec)+") used for computation.";
180    }
181  }
182  if(external_ord=="" && find(singular_ord,"Wp")==3)
183  {
184    external_ord="W_LEX";
185    pos=5;
186    for(i=1;i<=nvars(basering);i++)
187    {
188      pos++;
189      number_string="";
190      while(singular_ord[pos]!=",")
191      {
192        number_string=number_string+singular_ord[pos];
193        pos++;
194      }
195      weightvec[i]=`number_string`;
196    }
197    test_ord="Wp("+string(weightvec)+"),";
198    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
199    {
200      "Warning: Block orderings are not supported; Wp("+string(weightvec)+") used for computation.";
201    }
202  }
203
204  if(external_ord=="")
205  {
206    ERROR("The term ordering of the actual basering is not supported.");
207  }
208
209  // check algorithm
210  if(alg=="ct" || alg=="pct")
211    // these algorithms will not cause an error in the external program;
212    // however, they do not compute the toric ideal of A, but of an
213    // extended matrix
214  {
215    ERROR("The chosen algorithm is not suitable.");
216  }
217
218  // create temporary file with which the external program is called
219
220  int dummy;
221  int process=system("pid");
222  string matrixfile="temp_MATRIX"+string(process);
223  link MATRIX=":w "+matrixfile;
224  open(MATRIX);
225
226  write(MATRIX,"MATRIX","columns:",ncols(A),"cost vector:");
227  for(j=1;j<=ncols(A);j++)
228  {
229    write(MATRIX,weightvec[j]);
230  }
231  write(MATRIX,"rows:",nrows(A),"matrix:");
232  for(i=1;i<=nrows(A);i++)
233  {
234    for(j=1;j<=ncols(A);j++)
235    {
236      write(MATRIX,A[i,j]);
237    }
238  }
239
240  // search for positive row space vector, if required by the
241  // algorithm
242  int found=0;
243  if((alg=="blr") || (alg=="hs"))
244  {
245    for(i=1;i<=nrows(A);i++)
246    {
247      found=i;
248      for(j=1;j<=ncols(A);j++)
249      {
250        if(A[i,j]<=0)
251        {
252          found=0;
253        }
254      }
255      if(found>0)
256      {
257        break;
258      }
259    }
260    if(found==0)
261    {
262      close(MATRIX);
263      dummy=system("sh","rm -f "+matrixfile);
264      ERROR("The chosen algorithm needs a positive vector in the row space of the matrix.");
265    }
266    write(MATRIX,"positive row space vector:");
267    for(j=1;j<=ncols(A);j++)
268    {
269      write(MATRIX,A[found,j]);
270    }
271  }
272  close(MATRIX);
273
274
275  // call external program
276  if (status(system("SingularBin")+"toric_ideal","exists")=="yes")
277  {
278    dummy=system("sh",system("SingularBin")+"toric_ideal -alg "+alg+" "+matrixfile);
279  }
280  else
281  {
282    dummy=system("sh","toric_ideal -alg "+alg+" "+matrixfile);
283  }
284  if (dummy!=0) { ERROR("toric_ideal failed with error code "+string(dummy)); }
285
286  // read toric ideal from created file
287  link TORIC_IDEAL=":r "+matrixfile+".GB."+alg;
288  string toric_id=read(TORIC_IDEAL);
289
290  int generators;
291  pos=find(toric_id,"size");
292  pos=find(toric_id,":",pos);
293  pos++;
294
295  while(toric_id[pos]==" " || toric_id[pos]==newline)
296  {
297    pos++;
298  }
299  number_string="";
300  while(toric_id[pos]!=" " && toric_id[pos]!=newline)
301  {
302    number_string=number_string+toric_id[pos];
303    pos++;
304  }
305  generators=`number_string`;
306
307  intvec v;
308  poly head;
309  poly tail;
310
311  pos=find(toric_id,"basis");
312  pos=find(toric_id,":",pos);
313  pos++;
314
315  for(i=1;i<=generators;i++)
316  {
317    head=1;
318    tail=1;
319
320    for(j=1;j<=ncols(A);j++)
321    {
322      while(toric_id[pos]==" " || toric_id[pos]==newline)
323      {
324        pos++;
325      }
326      number_string="";
327      while(toric_id[pos]!=" " && toric_id[pos]!=newline)
328      {
329        number_string=number_string+toric_id[pos];
330        pos++;
331      }
332      v[j]=`number_string`;
333      if(v[j]<0)
334      {
335        tail=tail*var(j)^(-v[j]);
336      }
337      if(v[j]>0)
338      {
339        head=head*var(j)^v[j];
340      }
341    }
342    I[i]=head-tail;
343  }
344
345  // delete all created files
346  dummy=system("sh","rm -f "+matrixfile);
347  dummy=system("sh","rm -f "+matrixfile+".GB."+alg);
348
349  return(I);
350}
351///////////////////////////////////////////////////////////////////////////////
352
353static proc toric_ideal_2(intmat A, string alg, intvec prsv)
354{
355  ideal I;
356  // to be returned
357
358  // check arguments
359  if(size(prsv)<ncols(A))
360  {
361    ERROR("The number of matrix columns does not equal the size of the positive row space vector.");
362  }
363
364  // check suitability of actual basering
365  if(nvars(basering)!=ncols(A))
366  {
367    ERROR("The number of matrix columns is smaller than the number of ring variables.");
368  }
369
370  // check suitability of actual term ordering
371  // the "module ordering" c or C is ignored
372  string singular_ord=ordstr(basering);
373  string test_ord;
374  string external_ord="";
375  int i,j;
376  intvec weightvec;
377
378  if(find(singular_ord,"lp")==1)
379  {
380    external_ord="W_LEX";
381    for(i=1;i<=nvars(basering);i++)
382    {
383      weightvec[i]=0;
384    }
385    test_ord="lp("+string(nvars(basering))+"),";
386    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
387    {
388      "Warning: Block orderings are not supported; lp used for computation.";
389    }
390  }
391  if(external_ord=="" && find(singular_ord,"lp")==3)
392  {
393    external_ord="W_LEX";
394    for(i=1;i<=nvars(basering);i++)
395    {
396      weightvec[i]=0;
397    }
398    test_ord="lp("+string(nvars(basering))+"),";
399    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
400    {
401      "Warning: Block orderings are not supported; lp used for computation.";
402    }
403  }
404
405  if(external_ord=="" && find(singular_ord,"dp")==1)
406  {
407    external_ord="W_REV_LEX";
408    for(i=1;i<=nvars(basering);i++)
409    {
410      weightvec[i]=1;
411    }
412    test_ord="dp("+string(nvars(basering))+"),";
413    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
414    {
415      "Warning: Block orderings are not supported; dp used for computation.";
416    }
417  }
418  if(external_ord=="" && find(singular_ord,"dp")==3)
419  {
420    external_ord="W_REV_LEX";
421    for(i=1;i<=nvars(basering);i++)
422    {
423      weightvec[i]=1;
424    }
425    test_ord="dp("+string(nvars(basering))+"),";
426    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
427    {
428      "Warning: Block orderings are not supported; dp used for computation.";
429    }
430  }
431
432  if(external_ord=="" && find(singular_ord,"Dp")==1)
433  {
434    external_ord="W_LEX";
435    for(i=1;i<=nvars(basering);i++)
436    {
437      weightvec[i]=1;
438    }
439    test_ord="Dp("+string(nvars(basering))+"),";
440    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
441    {
442      "Warning: Block orderings are not supported; Dp used for computation.";
443    }
444  }
445  if(external_ord=="" && find(singular_ord,"Dp")==3)
446  {
447    external_ord="W_LEX";
448    for(i=1;i<=nvars(basering);i++)
449    {
450      weightvec[i]=1;
451    }
452    test_ord="Dp("+string(nvars(basering))+"),";
453    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
454    {
455      "Warning: Block orderings are not supported; Dp used for computation.";
456    }
457  }
458
459  int pos;
460  string number_string;
461
462  if(external_ord=="" && find(singular_ord,"wp")==1)
463  {
464    external_ord="W_REV_LEX";
465    pos=3;
466    for(i=1;i<=nvars(basering);i++)
467    {
468      pos++;
469      number_string="";
470      while(singular_ord[pos]!=",")
471      {
472        number_string=number_string+singular_ord[pos];
473        pos++;
474      }
475      weightvec[i]=`number_string`;
476    }
477    test_ord="wp("+string(weightvec)+"),";
478    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
479    {
480      "Warning: Block orderings are not supported; wp("+string(weightvec)+") used for computation.";
481    }
482  }
483  if(external_ord=="" && find(singular_ord,"wp")==3)
484  {
485    external_ord="W_REV_LEX";
486    pos=5;
487    for(i=1;i<=nvars(basering);i++)
488    {
489      pos++;
490      number_string="";
491      while(singular_ord[pos]!=",")
492      {
493        number_string=number_string+singular_ord[pos];
494        pos++;
495      }
496      weightvec[i]=`number_string`;
497    }
498    test_ord="wp("+string(weightvec)+"),";
499    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
500    {
501      "Warning: Block orderings are not supported; wp("+string(weightvec)+") used for computation.";
502    }
503  }
504
505  if(external_ord=="" && find(singular_ord,"Wp")==1)
506  {
507    external_ord="W_LEX";
508    pos=3;
509    for(i=1;i<=nvars(basering);i++)
510    {
511      pos++;
512      number_string="";
513      while(singular_ord[pos]!=",")
514      {
515        number_string=number_string+singular_ord[pos];
516        pos++;
517      }
518      weightvec[i]=`number_string`;
519    }
520    test_ord="Wp("+string(weightvec)+"),";
521    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
522    {
523      "Warning: Block orderings are not supported; Wp("+string(weightvec)+") used for computation.";
524    }
525  }
526  if(external_ord=="" && find(singular_ord,"Wp")==3)
527  {
528    external_ord="W_LEX";
529    pos=5;
530    for(i=1;i<=nvars(basering);i++)
531    {
532      pos++;
533      number_string="";
534      while(singular_ord[pos]!=",")
535      {
536        number_string=number_string+singular_ord[pos];
537        pos++;
538      }
539      weightvec[i]=`number_string`;
540    }
541    test_ord="Wp("+string(weightvec)+"),";
542    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
543    {
544      "Warning: Block orderings are not supported; Wp("+string(weightvec)+") used for computation.";
545    }
546  }
547
548  if(external_ord=="")
549  {
550    ERROR("The term ordering of the actual basering is not supported.");
551  }
552
553  // check algorithm
554  if(alg=="ct" || alg=="pct")
555    // these algorithms will not cause an error in the external program;
556    // however, they do not compute the toric ideal of A, but of an
557    // extended matrix
558  {
559    ERROR("The chosen algorithm is not suitable.");
560  }
561
562  // create temporary file with that the external program is called
563
564  int dummy;
565  int process=system("pid");
566  string matrixfile="temp_MATRIX"+string(process);
567  link MATRIX=":w "+matrixfile;
568  open(MATRIX);
569
570  write(MATRIX,"MATRIX","columns:",ncols(A),"cost vector:");
571  for(j=1;j<=ncols(A);j++)
572  {
573    write(MATRIX,weightvec[j]);
574  }
575  write(MATRIX,"rows:",nrows(A),"matrix:");
576  for(i=1;i<=nrows(A);i++)
577  {
578    for(j=1;j<=ncols(A);j++)
579    {
580      write(MATRIX,A[i,j]);
581    }
582  }
583
584  // enter positive row space vector, if required by the algorithm
585  if((alg=="blr") || (alg=="hs"))
586  {
587    write(MATRIX,"positive row space vector:");
588    for(j=1;j<=ncols(A);j++)
589    {
590      write(MATRIX,prsv[j]);
591    }
592  }
593  close(MATRIX);
594
595  // call external program
596  if (status(system("SingularBin")+"toric_ideal","exists")=="yes")
597  {
598    dummy=system("sh",system("SingularBin")+"toric_ideal -alg "+alg+" "+matrixfile);
599  }
600  else
601  {
602    dummy=system("sh","toric_ideal -alg "+alg+" "+matrixfile);
603  }
604  if (dummy!=0) { ERROR("toric_ideal failed with error code "+string(dummy)); }
605
606  // read toric ideal from created file
607  link TORIC_IDEAL=":r "+matrixfile+".GB."+alg;
608  string toric_id=read(TORIC_IDEAL);
609
610  int generators;
611  pos=find(toric_id,"size");
612  pos=find(toric_id,":",pos);
613  pos++;
614
615  while(toric_id[pos]==" " || toric_id[pos]==newline)
616  {
617    pos++;
618  }
619  number_string="";
620  while(toric_id[pos]!=" " && toric_id[pos]!=newline)
621  {
622    number_string=number_string+toric_id[pos];
623    pos++;
624  }
625  generators=`number_string`;
626
627  intvec v;
628  poly head;
629  poly tail;
630
631  pos=find(toric_id,"basis");
632  pos=find(toric_id,":",pos);
633  pos++;
634
635  for(i=1;i<=generators;i++)
636  {
637    head=1;
638    tail=1;
639
640    for(j=1;j<=ncols(A);j++)
641    {
642      while(toric_id[pos]==" " || toric_id[pos]==newline)
643      {
644        pos++;
645      }
646      number_string="";
647      while(toric_id[pos]!=" " && toric_id[pos]!=newline)
648      {
649        number_string=number_string+toric_id[pos];
650        pos++;
651      }
652      v[j]=`number_string`;
653      if(v[j]<0)
654      {
655        tail=tail*var(j)^(-v[j]);
656      }
657      if(v[j]>0)
658      {
659        head=head*var(j)^v[j];
660      }
661    }
662    I[i]=head-tail;
663  }
664
665  // delete all created files
666  dummy=system("sh","rm -f "+matrixfile);
667  dummy=system("sh","rm -f "+matrixfile+".GB."+alg);
668
669  return(I);
670}
671///////////////////////////////////////////////////////////////////////////////
672
673proc toric_ideal
674"USAGE: toric_ideal(A,alg);      A intmat, alg string
675        toric_ideal(A,alg,prsv); A intmat, alg string, prsv intvec
676RETURN: ideal: standard basis of the toric ideal of A
677NOTE:   These procedures return the standard basis of the toric ideal of A
678        with respect to the term ordering in the current basering. Not all
679        term orderings are supported: The usual global term orderings may be
680        used, but no block orderings combining them.
681        One may call the procedure with several different algorithms: @*
682        - the algorithm of Conti/Traverso using elimination (ect), @*
683        - the algorithm of Pottier (pt),
684        - an algorithm of Bigatti/La Scala/Robbiano (blr),
685        - the algorithm of Hosten/Sturmfels (hs),
686        - the algorithm of DiBiase/Urbanke (du).
687       The argument `alg' should be the abbreviation for an algorithm as
688       above: ect, pt, blr, hs or du.
689
690      If `alg' is chosen to be `blr' or `hs', the algorithm needs a vector
691      with positive coefficients in the row space of A.
692      If no row of A contains only positive entries, one has to use the
693      second version of toric_ideal which takes such a vector as its third
694      argument.
695      For the mathematical background, see
696  @texinfo
697  @ref{Toric ideals and integer programming}.
698  @end texinfo
699EXAMPLE:  example toric_ideal; shows an example
700SEE ALSO: toric_std, intprog_lib, Toric ideals
701"
702{
703  if(size(#)==2)
704  {
705    return(toric_ideal_1(#[1],#[2]));
706  }
707  else
708  {
709    return(toric_ideal_2(#[1],#[2],#[3]));
710  }
711}
712
713
714
715example
716{
717"EXAMPLE"; echo=2;
718
719ring r=0,(x,y,z),dp;
720
721// call with two arguments
722intmat A[2][3]=1,1,0,0,1,1;
723A;
724
725ideal I=toric_ideal(A,"du");
726I;
727
728I=toric_ideal(A,"blr");
729I;
730
731// call with three arguments
732intvec prsv=1,2,1;
733I=toric_ideal(A,"blr",prsv);
734I;
735
736}
737///////////////////////////////////////////////////////////////////////////////
738
739proc toric_std(ideal I)
740"USAGE:   toric_std(I);      I ideal
741RETURN:   ideal: standard basis of I
742NOTE:     This procedure computes the standard basis of I using a specialized
743          Buchberger algorithm. The generating system by which I is given has
744          to consist of binomials of the form x^u-x^v. There is no real check
745          if I is toric. If I is generated by binomials of the above form,
746          but not toric, toric_std computes an ideal `between' I and its
747          saturation with respect to all variables.
748          For the mathematical background, see
749   @texinfo
750   @ref{Toric ideals and integer programming}.
751   @end texinfo
752EXAMPLE:  example toric_std; shows an example
753SEE ALSO: toric_ideal, toric_lib, intprog_lib, Toric ideals
754"
755{
756  ideal J;
757  // to be returned
758
759  // check suitability of actual term ordering
760  // the "module ordering" c or C is ignored
761  string singular_ord=ordstr(basering);
762  string test_ord;
763  string external_ord="";
764  int i,j;
765  intvec weightvec;
766
767  if(find(singular_ord,"lp")==1)
768  {
769    external_ord="W_LEX";
770    for(i=1;i<=nvars(basering);i++)
771    {
772      weightvec[i]=0;
773    }
774    test_ord="lp("+string(nvars(basering))+"),";
775    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
776    {
777      "Warning: Block orderings are not supported; lp used for computation.";
778    }
779  }
780  if(external_ord=="" && find(singular_ord,"lp")==3)
781  {
782    external_ord="W_LEX";
783    for(i=1;i<=nvars(basering);i++)
784    {
785      weightvec[i]=0;
786    }
787    test_ord="lp("+string(nvars(basering))+"),";
788    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
789    {
790      "Warning: Block orderings are not supported; lp used for computation.";
791    }
792  }
793
794  if(external_ord=="" && find(singular_ord,"dp")==1)
795  {
796    external_ord="W_REV_LEX";
797    for(i=1;i<=nvars(basering);i++)
798    {
799      weightvec[i]=1;
800    }
801    test_ord="dp("+string(nvars(basering))+"),";
802    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
803    {
804      "Warning: Block orderings are not supported; dp used for computation.";
805    }
806  }
807  if(external_ord=="" && find(singular_ord,"dp")==3)
808  {
809    external_ord="W_REV_LEX";
810    for(i=1;i<=nvars(basering);i++)
811    {
812      weightvec[i]=1;
813    }
814    test_ord="dp("+string(nvars(basering))+"),";
815    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
816    {
817      "Warning: Block orderings are not supported; dp used for computation.";
818    }
819  }
820
821  if(external_ord=="" && find(singular_ord,"Dp")==1)
822  {
823    external_ord="W_LEX";
824    for(i=1;i<=nvars(basering);i++)
825    {
826      weightvec[i]=1;
827    }
828    test_ord="Dp("+string(nvars(basering))+"),";
829    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
830    {
831      "Warning: Block orderings are not supported; Dp used for computation.";
832    }
833  }
834  if(external_ord=="" && find(singular_ord,"Dp")==3)
835  {
836    external_ord="W_LEX";
837    for(i=1;i<=nvars(basering);i++)
838    {
839      weightvec[i]=1;
840    }
841    test_ord="Dp("+string(nvars(basering))+"),";
842    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
843    {
844      "Warning: Block orderings are not supported; Dp used for computation.";
845    }
846  }
847
848  int pos;
849  string number_string;
850
851  if(external_ord=="" && find(singular_ord,"wp")==1)
852  {
853    external_ord="W_REV_LEX";
854    pos=3;
855    for(i=1;i<=nvars(basering);i++)
856    {
857      pos++;
858      number_string="";
859      while(singular_ord[pos]!="," && singular_ord[pos]!=")")
860      {
861        number_string=number_string+singular_ord[pos];
862        pos++;
863      }
864      weightvec[i]=`number_string`;
865    }
866    test_ord="wp("+string(weightvec)+"),";
867    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
868    {
869      "Warning: Block orderings are not supported; wp("+string(weightvec)+") used for computation.";
870    }
871  }
872
873  if(external_ord=="" && find(singular_ord,"wp")==3)
874  {
875    external_ord="W_REV_LEX";
876    pos=5;
877    for(i=1;i<=nvars(basering);i++)
878    {
879      pos++;
880      number_string="";
881      while(singular_ord[pos]!=",")
882      {
883        number_string=number_string+singular_ord[pos];
884        pos++;
885      }
886      weightvec[i]=`number_string`;
887    }
888    test_ord="wp("+string(weightvec)+"),";
889    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
890    {
891      "Warning: Block orderings are not supported; wp("+string(weightvec)+") used for computation.";
892    }
893  }
894
895  if(external_ord=="" && find(singular_ord,"Wp")==1)
896  {
897    external_ord="W_LEX";
898    pos=3;
899    for(i=1;i<=nvars(basering);i++)
900    {
901      pos++;
902      number_string="";
903      while(singular_ord[pos]!=",")
904      {
905        number_string=number_string+singular_ord[pos];
906        pos++;
907      }
908      weightvec[i]=`number_string`;
909    }
910    test_ord="Wp("+string(weightvec)+"),";
911    if(singular_ord!=(test_ord+"C") && singular_ord!=(test_ord+"c"))
912    {
913      "Warning: Block orderings are not supported; Wp("+string(weightvec)+") used for computation.";
914    }
915  }
916  if(external_ord=="" && find(singular_ord,"Wp")==3)
917  {
918    external_ord="W_LEX";
919    pos=5;
920    for(i=1;i<=nvars(basering);i++)
921    {
922      pos++;
923      number_string="";
924      while(singular_ord[pos]!=",")
925      {
926        number_string=number_string+singular_ord[pos];
927        pos++;
928      }
929      weightvec[i]=`number_string`;
930    }
931    test_ord="Wp("+string(weightvec)+"),";
932    if(singular_ord!=("C"+test_ord) && singular_ord!=("c"+test_ord))
933    {
934      "Warning: Block orderings are not supported; Wp("+string(weightvec)+") used for computation.";
935    }
936  }
937
938  if(external_ord=="")
939  {
940    ERROR("The term ordering of the actual basering is not supported.");
941  }
942
943  // create first temporary file with which the external program is called
944
945  int dummy;
946  int process=system("pid");
947  string groebnerfile="temp_GROEBNER"+string(process);
948  link GROEBNER=":w "+groebnerfile;
949  open(GROEBNER);
950
951  write(GROEBNER,"GROEBNER","computed with algorithm:","pt","term ordering:","elimination block",0,"weighted block",nvars(basering),external_ord);
952  // algorithm is totally unimportant, only required by the external program
953
954  for(i=1;i<=nvars(basering);i++)
955  {
956    write(GROEBNER,weightvec[i]);
957  }
958
959  write(GROEBNER,"size:",size(I),"Groebner basis:");
960  poly head;
961  poly tail;
962  poly rest;
963  intvec v;
964
965  for(j=1;j<=size(I);j++)
966  {
967    // test suitability of generator j
968    rest=I[j];
969    head=lead(rest);
970    rest=rest-head;
971    tail=lead(rest);
972    rest=rest-tail;
973
974    if(head==0 && tail==0 && rest!=0)
975    {
976      close(GROEBNER);
977      dummy=system("sh","rm -f "+groebnerfile);
978      ERROR("Generator "+string(j)+" of the input ideal is no binomial.");
979    }
980
981    if(leadcoef(tail)!=-leadcoef(head))
982      // generator is no difference of monomials (or a constant multiple)
983    {
984      close(GROEBNER);
985      dummy=system("sh","rm -f "+groebnerfile);
986      ERROR("Generator "+string(j)+" of the input ideal is no difference of monomials.");
987    }
988
989    if(gcd(head,tail)!=1)
990    {
991      "Warning: The monomials of generator "+string(j)+" of the input ideal are not relatively prime.";
992    }
993
994    // write vector representation of generator j into the file
995    v=leadexp(head)-leadexp(tail);
996    for(i=1;i<=nvars(basering);i++)
997    {
998      write(GROEBNER,v[i]);
999    }
1000  }
1001  close(GROEBNER);
1002
1003  // create second temporary file
1004
1005  string newcostfile="temp_NEW_COST"+string(process);
1006  link NEW_COST=":w "+newcostfile;
1007  open(NEW_COST);
1008
1009  write(NEW_COST,"NEW_COST","variables:",nvars(basering),"cost vector:");
1010  for(i=1;i<=nvars(basering);i++)
1011  {
1012    write(NEW_COST,weightvec[i]);
1013  }
1014
1015  // call external program
1016  if (status(system("SingularBin")+"change_cost","exists")=="yes")
1017  {
1018    dummy=system("sh",system("SingularBin")+"change_cost "+groebnerfile+" "+newcostfile);
1019  }
1020  else
1021  {
1022    dummy=system("sh","change_cost "+groebnerfile+" "+newcostfile);
1023  }
1024  if (dummy!=0) { ERROR("change_cost failed with error code "+string(dummy)); }
1025
1026  // read toric standard basis from created file
1027  link TORIC_IDEAL=":r "+newcostfile+".GB.pt";
1028  string toric_id=read(TORIC_IDEAL);
1029
1030  int generators;
1031  pos=find(toric_id,"size");
1032  pos=find(toric_id,":",pos);
1033  pos++;
1034
1035  while(toric_id[pos]==" " || toric_id[pos]==newline)
1036  {
1037    pos++;
1038  }
1039  number_string="";
1040  while(toric_id[pos]!=" " && toric_id[pos]!=newline)
1041  {
1042    number_string=number_string+toric_id[pos];
1043    pos++;
1044  }
1045  generators=`number_string`;
1046
1047  pos=find(toric_id,"basis");
1048  pos=find(toric_id,":",pos);
1049  pos++;
1050
1051  for(j=1;j<=generators;j++)
1052  {
1053    head=1;
1054    tail=1;
1055
1056    for(i=1;i<=nvars(basering);i++)
1057    {
1058      while(toric_id[pos]==" " || toric_id[pos]==newline)
1059      {
1060        pos++;
1061      }
1062      number_string="";
1063      while(toric_id[pos]!=" " && toric_id[pos]!=newline)
1064      {
1065        number_string=number_string+toric_id[pos];
1066        pos++;
1067      }
1068      v[i]=`number_string`;
1069      if(v[i]<0)
1070      {
1071        tail=tail*var(i)^(-v[i]);
1072      }
1073      if(v[i]>0)
1074      {
1075        head=head*var(i)^v[i];
1076      }
1077    }
1078    J[j]=head-tail;
1079  }
1080
1081  // delete all created files
1082  dummy=system("sh","rm -f "+groebnerfile);
1083  dummy=system("sh","rm -f "+groebnerfile+".GB.pt");
1084  dummy=system("sh","rm -f "+newcostfile);
1085
1086  return(J);
1087}
1088
1089example
1090{
1091"EXAMPLE"; echo=2;
1092
1093ring r=0,(x,y,z),wp(3,2,1);
1094
1095// call with toric ideal (of the matrix A=(1,1,1))
1096ideal I=x-y,x-z;
1097ideal J=toric_std(I);
1098J;
1099
1100// call with the same ideal, but badly chosen generators:
1101// 1) not only binomials
1102I=x-y,2x-y-z;
1103J=toric_std(I);
1104// 2) binomials whose monomials are not relatively prime
1105I=x-y,xy-yz,y-z;
1106J=toric_std(I);
1107J;
1108
1109// call with a non-toric ideal that seems to be toric
1110I=x-yz,xy-z;
1111J=toric_std(I);
1112J;
1113// comparison with real standard basis and saturation
1114ideal H=std(I);
1115H;
1116LIB "elim.lib";
1117sat_with_exp(H,xyz);
1118}
1119///////////////////////////////////////////////////////////////////////////////
Note: See TracBrowser for help on using the repository browser.