source: git/Singular/ipshell.cc @ e6f1e6

spielwiese
Last change on this file since e6f1e6 was e6f1e6, checked in by Hans Schoenemann <hannes@…>, 8 years ago
code cleanup: HAVE_RINGS is mostly not needed
  • Property mode set to 100644
File size: 159.6 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT:
6*/
7
8#include <kernel/mod2.h>
9
10#include <omalloc/omalloc.h>
11
12#include <factory/factory.h>
13
14#include <misc/auxiliary.h>
15#include <misc/options.h>
16#include <misc/mylimits.h>
17#include <misc/intvec.h>
18#include <misc/prime.h>
19
20#include <coeffs/numbers.h>
21#include <coeffs/coeffs.h>
22
23#include <coeffs/rmodulon.h>
24#include <coeffs/longrat.h>
25
26#include <polys/monomials/ring.h>
27#include <polys/monomials/maps.h>
28
29#include <polys/prCopy.h>
30#include <polys/matpol.h>
31
32#include <polys/weight.h>
33#include <polys/clapsing.h>
34
35
36#include <polys/ext_fields/algext.h>
37#include <polys/ext_fields/transext.h>
38
39#include <kernel/polys.h>
40#include <kernel/ideals.h>
41
42#include <kernel/numeric/mpr_base.h>
43#include <kernel/numeric/mpr_numeric.h>
44
45#include <kernel/GBEngine/syz.h>
46#include <kernel/GBEngine/kstd1.h>
47#include <kernel/GBEngine/kutil.h> // denominator_list
48
49#include <kernel/combinatorics/stairc.h>
50#include <kernel/combinatorics/hutil.h>
51
52#include <kernel/spectrum/semic.h>
53#include <kernel/spectrum/splist.h>
54#include <kernel/spectrum/spectrum.h>
55
56#include <kernel/oswrapper/feread.h>
57
58#include <Singular/lists.h>
59#include <Singular/attrib.h>
60#include <Singular/ipconv.h>
61#include <Singular/links/silink.h>
62#include <Singular/ipshell.h>
63#include <Singular/maps_ip.h>
64#include <Singular/tok.h>
65#include <Singular/ipid.h>
66#include <Singular/subexpr.h>
67#include <Singular/fevoices.h>
68#include <Singular/sdb.h>
69
70#include <math.h>
71#include <ctype.h>
72
73#include <kernel/maps/gen_maps.h>
74
75#ifdef SINGULAR_4_1
76#include <Singular/number2.h>
77#include <coeffs/bigintmat.h>
78#endif
79leftv iiCurrArgs=NULL;
80idhdl iiCurrProc=NULL;
81const char *lastreserved=NULL;
82
83static BOOLEAN iiNoKeepRing=TRUE;
84
85/*0 implementation*/
86
87const char * iiTwoOps(int t)
88{
89  if (t<127)
90  {
91    static char ch[2];
92    switch (t)
93    {
94      case '&':
95        return "and";
96      case '|':
97        return "or";
98      default:
99        ch[0]=t;
100        ch[1]='\0';
101        return ch;
102    }
103  }
104  switch (t)
105  {
106    case COLONCOLON:  return "::";
107    case DOTDOT:      return "..";
108    //case PLUSEQUAL:   return "+=";
109    //case MINUSEQUAL:  return "-=";
110    case MINUSMINUS:  return "--";
111    case PLUSPLUS:    return "++";
112    case EQUAL_EQUAL: return "==";
113    case LE:          return "<=";
114    case GE:          return ">=";
115    case NOTEQUAL:    return "<>";
116    default:          return Tok2Cmdname(t);
117  }
118}
119
120int iiOpsTwoChar(const char *s)
121{
122/* not handling: &&, ||, ** */
123  if (s[1]=='\0') return s[0];
124  else if (s[2]!='\0') return 0;
125  switch(s[0])
126  {
127    case '.': if (s[1]=='.') return DOTDOT;
128              else           return 0;
129    case ':': if (s[1]==':') return COLONCOLON;
130              else           return 0;
131    case '-': if (s[1]=='-') return MINUSMINUS;
132              else           return 0;
133    case '+': if (s[1]=='+') return PLUSPLUS;
134              else           return 0;
135    case '=': if (s[1]=='=') return EQUAL_EQUAL;
136              else           return 0;
137    case '<': if (s[1]=='=') return LE;
138              else if (s[1]=='>') return NOTEQUAL;
139              else           return 0;
140    case '>': if (s[1]=='=') return GE;
141              else           return 0;
142    case '!': if (s[1]=='=') return NOTEQUAL;
143              else           return 0;
144  }
145  return 0;
146}
147
148static void list1(const char* s, idhdl h,BOOLEAN c, BOOLEAN fullname)
149{
150  char buffer[22];
151  int l;
152  char buf2[128];
153
154  if(fullname) sprintf(buf2, "%s::%s", "", IDID(h));
155  else sprintf(buf2, "%s", IDID(h));
156
157  Print("%s%-30.30s [%d]  ",s,buf2,IDLEV(h));
158  if (h == currRingHdl) PrintS("*");
159  PrintS(Tok2Cmdname((int)IDTYP(h)));
160
161  ipListFlag(h);
162  switch(IDTYP(h))
163  {
164    case ALIAS_CMD: Print(" for %s",IDID((idhdl)IDDATA(h))); break;
165    case INT_CMD:   Print(" %d",IDINT(h)); break;
166    case INTVEC_CMD:Print(" (%d)",IDINTVEC(h)->length()); break;
167    case INTMAT_CMD:Print(" %d x %d",IDINTVEC(h)->rows(),IDINTVEC(h)->cols());
168                    break;
169    case POLY_CMD:
170    case VECTOR_CMD:if (c)
171                    {
172                      PrintS(" ");wrp(IDPOLY(h));
173                      if(IDPOLY(h) != NULL)
174                      {
175                        Print(", %d monomial(s)",pLength(IDPOLY(h)));
176                      }
177                    }
178                    break;
179    case MODUL_CMD: Print(", rk %d", (int)(IDIDEAL(h)->rank));
180    case IDEAL_CMD: Print(", %u generator(s)",
181                    IDELEMS(IDIDEAL(h))); break;
182    case MAP_CMD:
183                    Print(" from %s",IDMAP(h)->preimage); break;
184    case MATRIX_CMD:Print(" %u x %u"
185                      ,MATROWS(IDMATRIX(h))
186                      ,MATCOLS(IDMATRIX(h))
187                    );
188                    break;
189    case PACKAGE_CMD:
190                    paPrint(IDID(h),IDPACKAGE(h));
191                    break;
192    case PROC_CMD: if((IDPROC(h)->libname!=NULL)
193                   && (strlen(IDPROC(h)->libname)>0))
194                     Print(" from %s",IDPROC(h)->libname);
195                   if(IDPROC(h)->is_static)
196                     PrintS(" (static)");
197                   break;
198    case STRING_CMD:
199                   {
200                     char *s;
201                     l=strlen(IDSTRING(h));
202                     memset(buffer,0,22);
203                     strncpy(buffer,IDSTRING(h),si_min(l,20));
204                     if ((s=strchr(buffer,'\n'))!=NULL)
205                     {
206                       *s='\0';
207                     }
208                     PrintS(" ");
209                     PrintS(buffer);
210                     if((s!=NULL) ||(l>20))
211                     {
212                       Print("..., %d char(s)",l);
213                     }
214                     break;
215                   }
216    case LIST_CMD: Print(", size: %d",IDLIST(h)->nr+1);
217                   break;
218    case QRING_CMD:
219    case RING_CMD:
220                   if ((IDRING(h)==currRing) && (currRingHdl!=h))
221                     PrintS("(*)"); /* this is an alias to currRing */
222#ifdef RDEBUG
223                   if (traceit &TRACE_SHOW_RINGS)
224                     Print(" <%lx>",(long)(IDRING(h)));
225#endif
226                   break;
227#ifdef SINGULAR_4_1
228    case CNUMBER_CMD:
229                   {  number2 n=(number2)IDDATA(h);
230                      Print(" (%s)",nCoeffName(n->cf));
231                      break;
232                   }
233    case CMATRIX_CMD:
234                   {  bigintmat *b=(bigintmat*)IDDATA(h);
235                      Print(" %d x %d (%s)",
236                        b->rows(),b->cols(),
237                        nCoeffName(b->basecoeffs()));
238                      break;
239                   }
240#endif
241    /*default:     break;*/
242  }
243  PrintLn();
244}
245
246void type_cmd(leftv v)
247{
248  BOOLEAN oldShortOut = FALSE;
249
250  if (currRing != NULL)
251  {
252    oldShortOut = currRing->ShortOut;
253    currRing->ShortOut = 1;
254  }
255  int t=v->Typ();
256  Print("// %s %s ",v->Name(),Tok2Cmdname(t));
257  switch (t)
258  {
259    case MAP_CMD:Print(" from %s\n",((map)(v->Data()))->preimage); break;
260    case INTMAT_CMD: Print(" %d x %d\n",((intvec*)(v->Data()))->rows(),
261                                      ((intvec*)(v->Data()))->cols()); break;
262    case MATRIX_CMD:Print(" %u x %u\n" ,
263       MATROWS((matrix)(v->Data())),
264       MATCOLS((matrix)(v->Data())));break;
265    case MODUL_CMD: Print(", rk %d\n", (int)(((ideal)(v->Data()))->rank));break;
266    case LIST_CMD: Print(", size %d\n",((lists)(v->Data()))->nr+1); break;
267
268    case PROC_CMD:
269    case RING_CMD:
270    case IDEAL_CMD:
271    case QRING_CMD: PrintLn(); break;
272
273    //case INT_CMD:
274    //case STRING_CMD:
275    //case INTVEC_CMD:
276    //case POLY_CMD:
277    //case VECTOR_CMD:
278    //case PACKAGE_CMD:
279
280    default:
281      break;
282  }
283  v->Print();
284  if (currRing != NULL)
285    currRing->ShortOut = oldShortOut;
286}
287
288static void killlocals0(int v, idhdl * localhdl, const ring r)
289{
290  idhdl h = *localhdl;
291  while (h!=NULL)
292  {
293    int vv;
294    //Print("consider %s, lev: %d:",IDID(h),IDLEV(h));
295    if ((vv=IDLEV(h))>0)
296    {
297      if (vv < v)
298      {
299        if (iiNoKeepRing)
300        {
301          //PrintS(" break\n");
302          return;
303        }
304        h = IDNEXT(h);
305        //PrintLn();
306      }
307      else //if (vv >= v)
308      {
309        idhdl nexth = IDNEXT(h);
310        killhdl2(h,localhdl,r);
311        h = nexth;
312        //PrintS("kill\n");
313      }
314    }
315    else
316    {
317      h = IDNEXT(h);
318      //PrintLn();
319    }
320  }
321}
322
323void killlocals_rec(idhdl *root,int v, ring r)
324{
325  idhdl h=*root;
326  while (h!=NULL)
327  {
328    if (IDLEV(h)>=v)
329    {
330//      Print("kill %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
331      idhdl n=IDNEXT(h);
332      killhdl2(h,root,r);
333      h=n;
334    }
335    else if (IDTYP(h)==PACKAGE_CMD)
336    {
337 //     Print("into pack %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
338      if (IDPACKAGE(h)!=basePack)
339        killlocals_rec(&(IDRING(h)->idroot),v,r);
340      h=IDNEXT(h);
341    }
342    else if ((IDTYP(h)==RING_CMD)
343    ||(IDTYP(h)==QRING_CMD))
344    {
345      if ((IDRING(h)!=NULL) && (IDRING(h)->idroot!=NULL))
346      // we have to test IDRING(h)!=NULL: qring Q=groebner(...): killlocals
347      {
348  //    Print("into ring %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
349        killlocals_rec(&(IDRING(h)->idroot),v,IDRING(h));
350      }
351      h=IDNEXT(h);
352    }
353    else
354    {
355//      Print("skip %s lev %d for lev %d\n",IDID(h),IDLEV(h),v);
356      h=IDNEXT(h);
357    }
358  }
359}
360BOOLEAN killlocals_list(int v, lists L)
361{
362  if (L==NULL) return FALSE;
363  BOOLEAN changed=FALSE;
364  int n=L->nr;
365  for(;n>=0;n--)
366  {
367    leftv h=&(L->m[n]);
368    void *d=h->data;
369    if (((h->rtyp==RING_CMD) || (h->rtyp==QRING_CMD))
370    && (((ring)d)->idroot!=NULL))
371    {
372      if (d!=currRing) {changed=TRUE;rChangeCurrRing((ring)d);}
373      killlocals0(v,&(((ring)h->data)->idroot),(ring)h->data);
374    }
375    else if (h->rtyp==LIST_CMD)
376      changed|=killlocals_list(v,(lists)d);
377  }
378  return changed;
379}
380void killlocals(int v)
381{
382  BOOLEAN changed=FALSE;
383  idhdl sh=currRingHdl;
384  ring cr=currRing;
385  if (sh!=NULL) changed=((IDLEV(sh)<v) || (IDRING(sh)->ref>0));
386  //if (changed) Print("currRing=%s(%x), lev=%d,ref=%d\n",IDID(sh),IDRING(sh),IDLEV(sh),IDRING(sh)->ref);
387
388  killlocals_rec(&(basePack->idroot),v,currRing);
389
390  if (iiRETURNEXPR_len > myynest)
391  {
392    int t=iiRETURNEXPR.Typ();
393    if ((/*iiRETURNEXPR.Typ()*/ t==RING_CMD)
394    || (/*iiRETURNEXPR.Typ()*/ t==QRING_CMD))
395    {
396      leftv h=&iiRETURNEXPR;
397      if (((ring)h->data)->idroot!=NULL)
398        killlocals0(v,&(((ring)h->data)->idroot),(ring)h->data);
399    }
400    else if (/*iiRETURNEXPR.Typ()*/ t==LIST_CMD)
401    {
402      leftv h=&iiRETURNEXPR;
403      changed |=killlocals_list(v,(lists)h->data);
404    }
405  }
406  if (changed)
407  {
408    currRingHdl=rFindHdl(cr,NULL);
409    if (currRingHdl==NULL)
410      currRing=NULL;
411    else if(cr!=currRing)
412      rChangeCurrRing(cr);
413  }
414
415  if (myynest<=1) iiNoKeepRing=TRUE;
416  //Print("end killlocals  >= %d\n",v);
417  //listall();
418}
419
420void list_cmd(int typ, const char* what, const char *prefix,BOOLEAN iterate, BOOLEAN fullname)
421{
422  package savePack=currPack;
423  idhdl h,start;
424  BOOLEAN all = typ<0;
425  BOOLEAN really_all=FALSE;
426
427  if ( typ==0 )
428  {
429    if (strcmp(what,"all")==0)
430    {
431      if (currPack!=basePack)
432        list_cmd(-1,NULL,prefix,iterate,fullname); // list current package
433      really_all=TRUE;
434      h=basePack->idroot;
435    }
436    else
437    {
438      h = ggetid(what);
439      if (h!=NULL)
440      {
441        if (iterate) list1(prefix,h,TRUE,fullname);
442        if (IDTYP(h)==ALIAS_CMD) PrintS("A");
443        if ((IDTYP(h)==RING_CMD)
444            || (IDTYP(h)==QRING_CMD)
445            //|| (IDTYP(h)==PACKE_CMD)
446        )
447        {
448          h=IDRING(h)->idroot;
449        }
450        else if(IDTYP(h)==PACKAGE_CMD)
451        {
452          currPack=IDPACKAGE(h);
453          //Print("list_cmd:package\n");
454          all=TRUE;typ=PROC_CMD;fullname=TRUE;really_all=TRUE;
455          h=IDPACKAGE(h)->idroot;
456        }
457        else
458        {
459          currPack=savePack;
460          return;
461        }
462      }
463      else
464      {
465        Werror("%s is undefined",what);
466        currPack=savePack;
467        return;
468      }
469    }
470    all=TRUE;
471  }
472  else if (RingDependend(typ))
473  {
474    h = currRing->idroot;
475  }
476  else
477    h = IDROOT;
478  start=h;
479  while (h!=NULL)
480  {
481    if ((all
482      && (IDTYP(h)!=PROC_CMD)
483      &&(IDTYP(h)!=PACKAGE_CMD)
484      #ifdef SINGULAR_4_1
485      &&(IDTYP(h)!=CRING_CMD)
486      #endif
487      )
488    || (typ == IDTYP(h))
489    #ifdef SINGULAR_4_1
490    || ((IDTYP(h)==CRING_CMD) && (typ==RING_CMD))
491    #else
492    || ((IDTYP(h)==QRING_CMD) && (typ==RING_CMD))
493    #endif
494    )
495    {
496      list1(prefix,h,start==currRingHdl, fullname);
497      if (((IDTYP(h)==RING_CMD)||(IDTYP(h)==QRING_CMD))
498        && (really_all || (all && (h==currRingHdl)))
499        && ((IDLEV(h)==0)||(IDLEV(h)==myynest)))
500      {
501        list_cmd(0,IDID(h),"//      ",FALSE);
502      }
503      if (IDTYP(h)==PACKAGE_CMD && really_all)
504      {
505        package save_p=currPack;
506        currPack=IDPACKAGE(h);
507        list_cmd(0,IDID(h),"//      ",FALSE);
508        currPack=save_p;
509      }
510    }
511    h = IDNEXT(h);
512  }
513  currPack=savePack;
514}
515
516void test_cmd(int i)
517{
518  int ii;
519
520  if (i<0)
521  {
522    ii= -i;
523    if (ii < 32)
524    {
525      si_opt_1 &= ~Sy_bit(ii);
526    }
527    else if (ii < 64)
528    {
529      si_opt_2 &= ~Sy_bit(ii-32);
530    }
531    else
532      WerrorS("out of bounds\n");
533  }
534  else if (i<32)
535  {
536    ii=i;
537    if (Sy_bit(ii) & kOptions)
538    {
539      Warn("Gerhard, use the option command");
540      si_opt_1 |= Sy_bit(ii);
541    }
542    else if (Sy_bit(ii) & validOpts)
543      si_opt_1 |= Sy_bit(ii);
544  }
545  else if (i<64)
546  {
547    ii=i-32;
548    si_opt_2 |= Sy_bit(ii);
549  }
550  else
551    WerrorS("out of bounds\n");
552}
553
554int exprlist_length(leftv v)
555{
556  int rc = 0;
557  while (v!=NULL)
558  {
559    switch (v->Typ())
560    {
561      case INT_CMD:
562      case POLY_CMD:
563      case VECTOR_CMD:
564      case NUMBER_CMD:
565        rc++;
566        break;
567      case INTVEC_CMD:
568      case INTMAT_CMD:
569        rc += ((intvec *)(v->Data()))->length();
570        break;
571      case MATRIX_CMD:
572      case IDEAL_CMD:
573      case MODUL_CMD:
574        {
575          matrix mm = (matrix)(v->Data());
576          rc += mm->rows() * mm->cols();
577        }
578        break;
579      case LIST_CMD:
580        rc+=((lists)v->Data())->nr+1;
581        break;
582      default:
583        rc++;
584    }
585    v = v->next;
586  }
587  return rc;
588}
589
590BOOLEAN iiWRITE(leftv,leftv v)
591{
592  sleftv vf;
593  if (iiConvert(v->Typ(),LINK_CMD,iiTestConvert(v->Typ(),LINK_CMD),v,&vf))
594  {
595    WerrorS("link expected");
596    return TRUE;
597  }
598  si_link l=(si_link)vf.Data();
599  if (vf.next == NULL)
600  {
601    WerrorS("write: need at least two arguments");
602    return TRUE;
603  }
604
605  BOOLEAN b=slWrite(l,vf.next); /* iiConvert preserves next */
606  if (b)
607  {
608    const char *s;
609    if ((l!=NULL)&&(l->name!=NULL)) s=l->name;
610    else                            s=sNoName;
611    Werror("cannot write to %s",s);
612  }
613  vf.CleanUp();
614  return b;
615}
616
617leftv iiMap(map theMap, const char * what)
618{
619  idhdl w,r;
620  leftv v;
621  int i;
622  nMapFunc nMap;
623
624  r=IDROOT->get(theMap->preimage,myynest);
625  if ((currPack!=basePack)
626  &&((r==NULL) || ((r->typ != RING_CMD) && (r->typ != QRING_CMD))))
627    r=basePack->idroot->get(theMap->preimage,myynest);
628  if ((r==NULL) && (currRingHdl!=NULL)
629  && (strcmp(theMap->preimage,IDID(currRingHdl))==0))
630  {
631    r=currRingHdl;
632  }
633  if ((r!=NULL) && ((r->typ == RING_CMD) || (r->typ== QRING_CMD)))
634  {
635    ring src_ring=IDRING(r);
636    if ((nMap=n_SetMap(src_ring->cf, currRing->cf))==NULL)
637    {
638      Werror("can not map from ground field of %s to current ground field",
639          theMap->preimage);
640      return NULL;
641    }
642    if (IDELEMS(theMap)<src_ring->N)
643    {
644      theMap->m=(polyset)omReallocSize((ADDRESS)theMap->m,
645                                 IDELEMS(theMap)*sizeof(poly),
646                                 (src_ring->N)*sizeof(poly));
647      for(i=IDELEMS(theMap);i<src_ring->N;i++)
648        theMap->m[i]=NULL;
649      IDELEMS(theMap)=src_ring->N;
650    }
651    if (what==NULL)
652    {
653      WerrorS("argument of a map must have a name");
654    }
655    else if ((w=src_ring->idroot->get(what,myynest))!=NULL)
656    {
657      char *save_r=NULL;
658      v=(leftv)omAlloc0Bin(sleftv_bin);
659      sleftv tmpW;
660      memset(&tmpW,0,sizeof(sleftv));
661      tmpW.rtyp=IDTYP(w);
662      if (tmpW.rtyp==MAP_CMD)
663      {
664        tmpW.rtyp=IDEAL_CMD;
665        save_r=IDMAP(w)->preimage;
666        IDMAP(w)->preimage=0;
667      }
668      tmpW.data=IDDATA(w);
669      // check overflow
670      BOOLEAN overflow=FALSE;
671      if ((tmpW.rtyp==IDEAL_CMD)
672        || (tmpW.rtyp==MODUL_CMD)
673        || (tmpW.rtyp==MAP_CMD))
674      {
675        ideal id=(ideal)tmpW.data;
676        long *degs=(long*)omAlloc(IDELEMS(id)*sizeof(long));
677        for(int i=IDELEMS(id)-1;i>=0;i--)
678        {
679          poly p=id->m[i];
680          if (p!=NULL) degs[i]=p_Totaldegree(p,src_ring);
681          else         degs[i]=0;
682        }
683        for(int j=IDELEMS(theMap)-1;j>=0 && !overflow;j--)
684        {
685          if (theMap->m[j]!=NULL)
686          {
687            long deg_monexp=pTotaldegree(theMap->m[j]);
688
689            for(int i=IDELEMS(id)-1;i>=0;i--)
690            {
691              poly p=id->m[i];
692              if ((p!=NULL) && (degs[i]!=0) &&
693              ((unsigned long)deg_monexp > (currRing->bitmask / ((unsigned long)degs[i])/2)))
694              {
695                overflow=TRUE;
696                break;
697              }
698            }
699          }
700        }
701        omFreeSize(degs,IDELEMS(id)*sizeof(long));
702      }
703      else if (tmpW.rtyp==POLY_CMD)
704      {
705        for(int j=IDELEMS(theMap)-1;j>=0 && !overflow;j--)
706        {
707          if (theMap->m[j]!=NULL)
708          {
709            long deg_monexp=pTotaldegree(theMap->m[j]);
710            poly p=(poly)tmpW.data;
711            long deg=0;
712            if ((p!=NULL) && ((deg=p_Totaldegree(p,src_ring))!=0) &&
713            ((unsigned long)deg_monexp > (currRing->bitmask / ((unsigned long)deg)/2)))
714            {
715              overflow=TRUE;
716              break;
717            }
718          }
719        }
720      }
721      if (overflow)
722        Warn("possible OVERFLOW in map, max exponent is %ld",currRing->bitmask/2);
723#if 0
724      if (((tmpW.rtyp==IDEAL_CMD)||(tmpW.rtyp==MODUL_CMD)) && idIs0(IDIDEAL(w)))
725      {
726        v->rtyp=tmpW.rtyp;
727        v->data=idInit(IDELEMS(IDIDEAL(w)),IDIDEAL(w)->rank);
728      }
729      else
730#endif
731      {
732        if ((tmpW.rtyp==IDEAL_CMD)
733        ||(tmpW.rtyp==MODUL_CMD)
734        ||(tmpW.rtyp==MATRIX_CMD)
735        ||(tmpW.rtyp==MAP_CMD))
736        {
737          v->rtyp=tmpW.rtyp;
738          char *tmp = theMap->preimage;
739          theMap->preimage=(char*)1L;
740          // map gets 1 as its rank (as an ideal)
741          v->data=maMapIdeal(IDIDEAL(w), src_ring, (ideal)theMap, currRing,nMap);
742          theMap->preimage=tmp; // map gets its preimage back
743        }
744        if (v->data==NULL) /*i.e. not IDEAL_CMD/MODUL_CMD/MATRIX_CMD/MAP */
745        {
746          if (maApplyFetch(MAP_CMD,theMap,v,&tmpW,src_ring,NULL,NULL,0,nMap))
747          {
748            Werror("cannot map %s(%d)",Tok2Cmdname(w->typ),w->typ);
749            omFreeBin((ADDRESS)v, sleftv_bin);
750            if (save_r!=NULL) IDMAP(w)->preimage=save_r;
751            return NULL;
752          }
753        }
754      }
755      if (save_r!=NULL)
756      {
757        IDMAP(w)->preimage=save_r;
758        IDMAP((idhdl)v)->preimage=omStrDup(save_r);
759        v->rtyp=MAP_CMD;
760      }
761      return v;
762    }
763    else
764    {
765      Werror("%s undefined in %s",what,theMap->preimage);
766    }
767  }
768  else
769  {
770    Werror("cannot find preimage %s",theMap->preimage);
771  }
772  return NULL;
773}
774
775#ifdef OLD_RES
776void  iiMakeResolv(resolvente r, int length, int rlen, char * name, int typ0,
777                   intvec ** weights)
778{
779  lists L=liMakeResolv(r,length,rlen,typ0,weights);
780  int i=0;
781  idhdl h;
782  char * s=(char *)omAlloc(strlen(name)+5);
783
784  while (i<=L->nr)
785  {
786    sprintf(s,"%s(%d)",name,i+1);
787    if (i==0)
788      h=enterid(s,myynest,typ0,&(currRing->idroot), FALSE);
789    else
790      h=enterid(s,myynest,MODUL_CMD,&(currRing->idroot), FALSE);
791    if (h!=NULL)
792    {
793      h->data.uideal=(ideal)L->m[i].data;
794      h->attribute=L->m[i].attribute;
795      if (BVERBOSE(V_DEF_RES))
796        Print("//defining: %s as %d-th syzygy module\n",s,i+1);
797    }
798    else
799    {
800      idDelete((ideal *)&(L->m[i].data));
801      Warn("cannot define %s",s);
802    }
803    //L->m[i].data=NULL;
804    //L->m[i].rtyp=0;
805    //L->m[i].attribute=NULL;
806    i++;
807  }
808  omFreeSize((ADDRESS)L->m,(L->nr+1)*sizeof(sleftv));
809  omFreeBin((ADDRESS)L, slists_bin);
810  omFreeSize((ADDRESS)s,strlen(name)+5);
811}
812#endif
813
814//resolvente iiFindRes(char * name, int * len, int *typ0)
815//{
816//  char *s=(char *)omAlloc(strlen(name)+5);
817//  int i=-1;
818//  resolvente r;
819//  idhdl h;
820//
821//  do
822//  {
823//    i++;
824//    sprintf(s,"%s(%d)",name,i+1);
825//    h=currRing->idroot->get(s,myynest);
826//  } while (h!=NULL);
827//  *len=i-1;
828//  if (*len<=0)
829//  {
830//    Werror("no objects %s(1),.. found",name);
831//    omFreeSize((ADDRESS)s,strlen(name)+5);
832//    return NULL;
833//  }
834//  r=(ideal *)omAlloc(/*(len+1)*/ i*sizeof(ideal));
835//  memset(r,0,(*len)*sizeof(ideal));
836//  i=-1;
837//  *typ0=MODUL_CMD;
838//  while (i<(*len))
839//  {
840//    i++;
841//    sprintf(s,"%s(%d)",name,i+1);
842//    h=currRing->idroot->get(s,myynest);
843//    if (h->typ != MODUL_CMD)
844//    {
845//      if ((i!=0) || (h->typ!=IDEAL_CMD))
846//      {
847//        Werror("%s is not of type module",s);
848//        omFreeSize((ADDRESS)r,(*len)*sizeof(ideal));
849//        omFreeSize((ADDRESS)s,strlen(name)+5);
850//        return NULL;
851//      }
852//      *typ0=IDEAL_CMD;
853//    }
854//    if ((i>0) && (idIs0(r[i-1])))
855//    {
856//      *len=i-1;
857//      break;
858//    }
859//    r[i]=IDIDEAL(h);
860//  }
861//  omFreeSize((ADDRESS)s,strlen(name)+5);
862//  return r;
863//}
864
865static resolvente iiCopyRes(resolvente r, int l)
866{
867  int i;
868  resolvente res=(ideal *)omAlloc0((l+1)*sizeof(ideal));
869
870  for (i=0; i<l; i++)
871    if (r[i]!=NULL) res[i]=idCopy(r[i]);
872  return res;
873}
874
875BOOLEAN jjMINRES(leftv res, leftv v)
876{
877  int len=0;
878  int typ0;
879  lists L=(lists)v->Data();
880  intvec *weights=(intvec*)atGet(v,"isHomog",INTVEC_CMD);
881  int add_row_shift = 0;
882  if (weights==NULL)
883    weights=(intvec*)atGet(&(L->m[0]),"isHomog",INTVEC_CMD);
884  if (weights!=NULL)  add_row_shift=weights->min_in();
885  resolvente rr=liFindRes(L,&len,&typ0);
886  if (rr==NULL) return TRUE;
887  resolvente r=iiCopyRes(rr,len);
888
889  syMinimizeResolvente(r,len,0);
890  omFreeSize((ADDRESS)rr,len*sizeof(ideal));
891  len++;
892  res->data=(char *)liMakeResolv(r,len,-1,typ0,NULL,add_row_shift);
893  return FALSE;
894}
895
896BOOLEAN jjBETTI(leftv res, leftv u)
897{
898  sleftv tmp;
899  memset(&tmp,0,sizeof(tmp));
900  tmp.rtyp=INT_CMD;
901  tmp.data=(void *)1;
902  if ((u->Typ()==IDEAL_CMD)
903  || (u->Typ()==MODUL_CMD))
904    return jjBETTI2_ID(res,u,&tmp);
905  else
906    return jjBETTI2(res,u,&tmp);
907}
908
909BOOLEAN jjBETTI2_ID(leftv res, leftv u, leftv v)
910{
911  lists l=(lists) omAllocBin(slists_bin);
912  l->Init(1);
913  l->m[0].rtyp=u->Typ();
914  l->m[0].data=u->Data();
915  attr *a=u->Attribute();
916  if (a!=NULL)
917  l->m[0].attribute=*a;
918  sleftv tmp2;
919  memset(&tmp2,0,sizeof(tmp2));
920  tmp2.rtyp=LIST_CMD;
921  tmp2.data=(void *)l;
922  BOOLEAN r=jjBETTI2(res,&tmp2,v);
923  l->m[0].data=NULL;
924  l->m[0].attribute=NULL;
925  l->m[0].rtyp=DEF_CMD;
926  l->Clean();
927  return r;
928}
929
930BOOLEAN jjBETTI2(leftv res, leftv u, leftv v)
931{
932  resolvente r;
933  int len;
934  int reg,typ0;
935  lists l=(lists)u->Data();
936
937  intvec *weights=NULL;
938  int add_row_shift=0;
939  intvec *ww=NULL;
940  if (l->nr>=0) ww=(intvec *)atGet(&(l->m[0]),"isHomog",INTVEC_CMD);
941  if (ww!=NULL)
942  {
943     weights=ivCopy(ww);
944     add_row_shift = ww->min_in();
945     (*weights) -= add_row_shift;
946  }
947  //Print("attr:%x\n",weights);
948
949  r=liFindRes(l,&len,&typ0);
950  if (r==NULL) return TRUE;
951  intvec* res_im=syBetti(r,len,&reg,weights,(int)(long)v->Data());
952  res->data=(void*)res_im;
953  omFreeSize((ADDRESS)r,(len)*sizeof(ideal));
954  //Print("rowShift: %d ",add_row_shift);
955  for(int i=1;i<=res_im->rows();i++)
956  {
957    if (IMATELEM(*res_im,1,i)==0) { add_row_shift--; }
958    else break;
959  }
960  //Print(" %d\n",add_row_shift);
961  atSet(res,omStrDup("rowShift"),(void*)(long)add_row_shift,INT_CMD);
962  if (weights!=NULL) delete weights;
963  return FALSE;
964}
965
966int iiRegularity(lists L)
967{
968  int len,reg,typ0;
969
970  resolvente r=liFindRes(L,&len,&typ0);
971
972  if (r==NULL)
973    return -2;
974  intvec *weights=NULL;
975  int add_row_shift=0;
976  intvec *ww=(intvec *)atGet(&(L->m[0]),"isHomog",INTVEC_CMD);
977  if (ww!=NULL)
978  {
979     weights=ivCopy(ww);
980     add_row_shift = ww->min_in();
981     (*weights) -= add_row_shift;
982  }
983  //Print("attr:%x\n",weights);
984
985  intvec *dummy=syBetti(r,len,&reg,weights);
986  if (weights!=NULL) delete weights;
987  delete dummy;
988  omFreeSize((ADDRESS)r,len*sizeof(ideal));
989  return reg+1+add_row_shift;
990}
991
992BOOLEAN iiDebugMarker=TRUE;
993#define BREAK_LINE_LENGTH 80
994void iiDebug()
995{
996#ifdef HAVE_SDB
997  sdb_flags=1;
998#endif
999  Print("\n-- break point in %s --\n",VoiceName());
1000  if (iiDebugMarker) VoiceBackTrack();
1001  char * s;
1002  iiDebugMarker=FALSE;
1003  s = (char *)omAlloc(BREAK_LINE_LENGTH+4);
1004  loop
1005  {
1006    memset(s,0,80);
1007    fe_fgets_stdin("",s,BREAK_LINE_LENGTH);
1008    if (s[BREAK_LINE_LENGTH-1]!='\0')
1009    {
1010      Print("line too long, max is %d chars\n",BREAK_LINE_LENGTH);
1011    }
1012    else
1013      break;
1014  }
1015  if (*s=='\n')
1016  {
1017    iiDebugMarker=TRUE;
1018  }
1019#if MDEBUG
1020  else if(strncmp(s,"cont;",5)==0)
1021  {
1022    iiDebugMarker=TRUE;
1023  }
1024#endif /* MDEBUG */
1025  else
1026  {
1027    strcat( s, "\n;~\n");
1028    newBuffer(s,BT_execute);
1029  }
1030}
1031
1032lists scIndIndset(ideal S, BOOLEAN all, ideal Q)
1033{
1034  int i;
1035  indset save;
1036  lists res=(lists)omAlloc0Bin(slists_bin);
1037
1038  hexist = hInit(S, Q, &hNexist, currRing);
1039  if (hNexist == 0)
1040  {
1041    intvec *iv=new intvec(rVar(currRing));
1042    for(i=0; i<rVar(currRing); i++) (*iv)[i]=1;
1043    res->Init(1);
1044    res->m[0].rtyp=INTVEC_CMD;
1045    res->m[0].data=(intvec*)iv;
1046    return res;
1047  }
1048  else if (hisModule!=0)
1049  {
1050    res->Init(0);
1051    return res;
1052  }
1053  save = ISet = (indset)omAlloc0Bin(indlist_bin);
1054  hMu = 0;
1055  hwork = (scfmon)omAlloc(hNexist * sizeof(scmon));
1056  hvar = (varset)omAlloc((rVar(currRing) + 1) * sizeof(int));
1057  hpure = (scmon)omAlloc((1 + (rVar(currRing) * rVar(currRing))) * sizeof(long));
1058  hrad = hexist;
1059  hNrad = hNexist;
1060  radmem = hCreate(rVar(currRing) - 1);
1061  hCo = rVar(currRing) + 1;
1062  hNvar = rVar(currRing);
1063  hRadical(hrad, &hNrad, hNvar);
1064  hSupp(hrad, hNrad, hvar, &hNvar);
1065  if (hNvar)
1066  {
1067    hCo = hNvar;
1068    memset(hpure, 0, (rVar(currRing) + 1) * sizeof(long));
1069    hPure(hrad, 0, &hNrad, hvar, hNvar, hpure, &hNpure);
1070    hLexR(hrad, hNrad, hvar, hNvar);
1071    hDimSolve(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1072  }
1073  if (hCo && (hCo < rVar(currRing)))
1074  {
1075    hIndMult(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1076  }
1077  if (hMu!=0)
1078  {
1079    ISet = save;
1080    hMu2 = 0;
1081    if (all && (hCo+1 < rVar(currRing)))
1082    {
1083      JSet = (indset)omAlloc0Bin(indlist_bin);
1084      hIndAllMult(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1085      i=hMu+hMu2;
1086      res->Init(i);
1087      if (hMu2 == 0)
1088      {
1089        omFreeBin((ADDRESS)JSet, indlist_bin);
1090      }
1091    }
1092    else
1093    {
1094      res->Init(hMu);
1095    }
1096    for (i=0;i<hMu;i++)
1097    {
1098      res->m[i].data = (void *)save->set;
1099      res->m[i].rtyp = INTVEC_CMD;
1100      ISet = save;
1101      save = save->nx;
1102      omFreeBin((ADDRESS)ISet, indlist_bin);
1103    }
1104    omFreeBin((ADDRESS)save, indlist_bin);
1105    if (hMu2 != 0)
1106    {
1107      save = JSet;
1108      for (i=hMu;i<hMu+hMu2;i++)
1109      {
1110        res->m[i].data = (void *)save->set;
1111        res->m[i].rtyp = INTVEC_CMD;
1112        JSet = save;
1113        save = save->nx;
1114        omFreeBin((ADDRESS)JSet, indlist_bin);
1115      }
1116      omFreeBin((ADDRESS)save, indlist_bin);
1117    }
1118  }
1119  else
1120  {
1121    res->Init(0);
1122    omFreeBin((ADDRESS)ISet,  indlist_bin);
1123  }
1124  hKill(radmem, rVar(currRing) - 1);
1125  omFreeSize((ADDRESS)hpure, (1 + (rVar(currRing) * rVar(currRing))) * sizeof(long));
1126  omFreeSize((ADDRESS)hvar, (rVar(currRing) + 1) * sizeof(int));
1127  omFreeSize((ADDRESS)hwork, hNexist * sizeof(scmon));
1128  hDelete(hexist, hNexist);
1129  return res;
1130}
1131
1132int iiDeclCommand(leftv sy, leftv name, int lev,int t, idhdl* root,BOOLEAN isring, BOOLEAN init_b)
1133{
1134  BOOLEAN res=FALSE;
1135  const char *id = name->name;
1136
1137  memset(sy,0,sizeof(sleftv));
1138  if ((name->name==NULL)||(isdigit(name->name[0])))
1139  {
1140    WerrorS("object to declare is not a name");
1141    res=TRUE;
1142  }
1143  else
1144  {
1145    if (TEST_V_ALLWARN
1146    && (name->rtyp!=0)
1147    && (name->rtyp!=IDHDL)
1148    && (currRingHdl!=NULL) && (IDLEV(currRingHdl)==myynest))
1149    {
1150      Warn("`%s` is %s in %s:%d:%s",name->name,Tok2Cmdname(name->rtyp),
1151      currentVoice->filename,yylineno,my_yylinebuf);
1152    }
1153    {
1154      sy->data = (char *)enterid(id,lev,t,root,init_b);
1155    }
1156    if (sy->data!=NULL)
1157    {
1158      sy->rtyp=IDHDL;
1159      currid=sy->name=IDID((idhdl)sy->data);
1160      // name->name=NULL; /* used in enterid */
1161      //sy->e = NULL;
1162      if (name->next!=NULL)
1163      {
1164        sy->next=(leftv)omAllocBin(sleftv_bin);
1165        res=iiDeclCommand(sy->next,name->next,lev,t,root, isring);
1166      }
1167    }
1168    else res=TRUE;
1169  }
1170  name->CleanUp();
1171  return res;
1172}
1173
1174BOOLEAN iiDefaultParameter(leftv p)
1175{
1176  attr at=NULL;
1177  if (iiCurrProc!=NULL)
1178     at=iiCurrProc->attribute->get("default_arg");
1179  if (at==NULL)
1180    return FALSE;
1181  sleftv tmp;
1182  memset(&tmp,0,sizeof(sleftv));
1183  tmp.rtyp=at->atyp;
1184  tmp.data=at->CopyA();
1185  return iiAssign(p,&tmp);
1186}
1187BOOLEAN iiBranchTo(leftv, leftv args)
1188{
1189  // <string1...stringN>,<proc>
1190  // known: args!=NULL, l>=1
1191  int l=args->listLength();
1192  int ll=0;
1193  if (iiCurrArgs!=NULL) ll=iiCurrArgs->listLength();
1194  if (ll!=(l-1)) return FALSE;
1195  leftv h=args;
1196  short *t=(short*)omAlloc(l*sizeof(short));
1197  t[0]=l-1;
1198  int b;
1199  int i;
1200  for(i=1;i<l;i++,h=h->next)
1201  {
1202    if (h->Typ()!=STRING_CMD)
1203    {
1204      omFree(t);
1205      Werror("arg %d is not a string",i);
1206      return TRUE;
1207    }
1208    int tt;
1209    b=IsCmd((char *)h->Data(),tt);
1210    if(b) t[i]=tt;
1211    else
1212    {
1213      omFree(t);
1214      Werror("arg %d is not a type name",i);
1215      return TRUE;
1216    }
1217  }
1218  if (h->Typ()!=PROC_CMD)
1219  {
1220    omFree(t);
1221    Werror("last arg (%d) is not a proc",i);
1222    return TRUE;
1223  }
1224  b=iiCheckTypes(iiCurrArgs,t,0);
1225  omFree(t);
1226  if (b && (h->rtyp==IDHDL) && (h->e==NULL))
1227  {
1228    BOOLEAN err;
1229    //Print("branchTo: %s\n",h->Name());
1230    iiCurrProc=(idhdl)h->data;
1231    procinfo * pi=IDPROC(iiCurrProc);
1232    if( pi->data.s.body==NULL )
1233    {
1234      iiGetLibProcBuffer(pi);
1235      if (pi->data.s.body==NULL) return TRUE;
1236    }
1237    if ((pi->pack!=NULL)&&(currPack!=pi->pack))
1238    {
1239      currPack=pi->pack;
1240      iiCheckPack(currPack);
1241      currPackHdl=packFindHdl(currPack);
1242      //Print("set pack=%s\n",IDID(currPackHdl));
1243    }
1244    err=iiAllStart(pi,pi->data.s.body,BT_proc,pi->data.s.body_lineno-(iiCurrArgs==NULL));
1245    exitBuffer(BT_proc);
1246    if (iiCurrArgs!=NULL)
1247    {
1248      if (!err) Warn("too many arguments for %s",IDID(iiCurrProc));
1249      iiCurrArgs->CleanUp();
1250      omFreeBin((ADDRESS)iiCurrArgs, sleftv_bin);
1251      iiCurrArgs=NULL;
1252    }
1253    return 2-err;
1254  }
1255  return FALSE;
1256}
1257BOOLEAN iiParameter(leftv p)
1258{
1259  if (iiCurrArgs==NULL)
1260  {
1261    if (strcmp(p->name,"#")==0)
1262      return iiDefaultParameter(p);
1263    Werror("not enough arguments for proc %s",VoiceName());
1264    p->CleanUp();
1265    return TRUE;
1266  }
1267  leftv h=iiCurrArgs;
1268  leftv rest=h->next; /*iiCurrArgs is not NULL here*/
1269  BOOLEAN is_default_list=FALSE;
1270  if (strcmp(p->name,"#")==0)
1271  {
1272    is_default_list=TRUE;
1273    rest=NULL;
1274  }
1275  else
1276  {
1277    h->next=NULL;
1278  }
1279  BOOLEAN res=iiAssign(p,h);
1280  if (is_default_list)
1281  {
1282    iiCurrArgs=NULL;
1283  }
1284  else
1285  {
1286    iiCurrArgs=rest;
1287  }
1288  h->CleanUp();
1289  omFreeBin((ADDRESS)h, sleftv_bin);
1290  return res;
1291}
1292
1293static BOOLEAN iiInternalExport (leftv v, int toLev)
1294{
1295  idhdl h=(idhdl)v->data;
1296  //Print("iiInternalExport('%s',%d)%s\n", v->name, toLev,"");
1297  if (IDLEV(h)==0)
1298  {
1299    if (BVERBOSE(V_REDEFINE)) Warn("`%s` is already global",IDID(h));
1300  }
1301  else
1302  {
1303    h=IDROOT->get(v->name,toLev);
1304    idhdl *root=&IDROOT;
1305    if ((h==NULL)&&(currRing!=NULL))
1306    {
1307      h=currRing->idroot->get(v->name,toLev);
1308      root=&currRing->idroot;
1309    }
1310    BOOLEAN keepring=FALSE;
1311    if ((h!=NULL)&&(IDLEV(h)==toLev))
1312    {
1313      if (IDTYP(h)==v->Typ())
1314      {
1315        if (((IDTYP(h)==RING_CMD)||(IDTYP(h)==QRING_CMD))
1316        && (v->Data()==IDDATA(h)))
1317        {
1318          IDRING(h)->ref++;
1319          keepring=TRUE;
1320          IDLEV(h)=toLev;
1321          //WarnS("keepring");
1322          return FALSE;
1323        }
1324        if (BVERBOSE(V_REDEFINE))
1325        {
1326          Warn("redefining %s (%s)",IDID(h),my_yylinebuf);
1327        }
1328#ifdef USE_IILOCALRING
1329        if (iiLocalRing[0]==IDRING(h) && (!keepring)) iiLocalRing[0]=NULL;
1330#else
1331        proclevel *p=procstack;
1332        while (p->next!=NULL) p=p->next;
1333        if ((p->cRing==IDRING(h)) && (!keepring))
1334        {
1335          p->cRing=NULL;
1336          p->cRingHdl=NULL;
1337        }
1338#endif
1339        killhdl2(h,root,currRing);
1340      }
1341      else
1342      {
1343        return TRUE;
1344      }
1345    }
1346    h=(idhdl)v->data;
1347    IDLEV(h)=toLev;
1348    if (keepring) IDRING(h)->ref--;
1349    iiNoKeepRing=FALSE;
1350    //Print("export %s\n",IDID(h));
1351  }
1352  return FALSE;
1353}
1354
1355BOOLEAN iiInternalExport (leftv v, int toLev, package rootpack)
1356{
1357  idhdl h=(idhdl)v->data;
1358  if(h==NULL)
1359  {
1360    Warn("'%s': no such identifier\n", v->name);
1361    return FALSE;
1362  }
1363  package frompack=v->req_packhdl;
1364  if (frompack==NULL) frompack=currPack;
1365  if ((RingDependend(IDTYP(h)))
1366  || ((IDTYP(h)==LIST_CMD)
1367     && (lRingDependend(IDLIST(h)))
1368     )
1369  )
1370  {
1371    //Print("// ==> Ringdependent set nesting to 0\n");
1372    return (iiInternalExport(v, toLev));
1373  }
1374  else
1375  {
1376    IDLEV(h)=toLev;
1377    v->req_packhdl=rootpack;
1378    if (h==frompack->idroot)
1379    {
1380      frompack->idroot=h->next;
1381    }
1382    else
1383    {
1384      idhdl hh=frompack->idroot;
1385      while ((hh!=NULL) && (hh->next!=h))
1386        hh=hh->next;
1387      if ((hh!=NULL) && (hh->next==h))
1388        hh->next=h->next;
1389      else
1390      {
1391        Werror("`%s` not found",v->Name());
1392        return TRUE;
1393      }
1394    }
1395    h->next=rootpack->idroot;
1396    rootpack->idroot=h;
1397  }
1398  return FALSE;
1399}
1400
1401BOOLEAN iiExport (leftv v, int toLev)
1402{
1403  BOOLEAN nok=FALSE;
1404  leftv r=v;
1405  while (v!=NULL)
1406  {
1407    if ((v->name==NULL)||(v->rtyp==0)||(v->e!=NULL))
1408    {
1409      Werror("cannot export:%s of internal type %d",v->name,v->rtyp);
1410      nok=TRUE;
1411    }
1412    else
1413    {
1414      if(iiInternalExport(v, toLev))
1415      {
1416        r->CleanUp();
1417        return TRUE;
1418      }
1419    }
1420    v=v->next;
1421  }
1422  r->CleanUp();
1423  return nok;
1424}
1425
1426/*assume root!=idroot*/
1427BOOLEAN iiExport (leftv v, int toLev, package pack)
1428{
1429#ifdef SINGULAR_4_1
1430  if ((pack==basePack)&&(pack!=currPack))
1431  { Warn("'exportto' to Top is depreciated in >>%s<<",my_yylinebuf);}
1432#endif
1433  BOOLEAN nok=FALSE;
1434  leftv rv=v;
1435  while (v!=NULL)
1436  {
1437    if ((v->name==NULL)||(v->rtyp==0)||(v->e!=NULL)
1438    )
1439    {
1440      Werror("cannot export:%s of internal type %d",v->name,v->rtyp);
1441      nok=TRUE;
1442    }
1443    else
1444    {
1445      idhdl old=pack->idroot->get( v->name,toLev);
1446      if (old!=NULL)
1447      {
1448        if ((pack==currPack) && (old==(idhdl)v->data))
1449        {
1450          if (BVERBOSE(V_REDEFINE)) Warn("`%s` is already global",IDID(old));
1451          break;
1452        }
1453        else if (IDTYP(old)==v->Typ())
1454        {
1455          if (BVERBOSE(V_REDEFINE))
1456          {
1457            Warn("redefining %s (%s)",IDID(old),my_yylinebuf);
1458          }
1459          v->name=omStrDup(v->name);
1460          killhdl2(old,&(pack->idroot),currRing);
1461        }
1462        else
1463        {
1464          rv->CleanUp();
1465          return TRUE;
1466        }
1467      }
1468      //Print("iiExport: pack=%s\n",IDID(root));
1469      if(iiInternalExport(v, toLev, pack))
1470      {
1471        rv->CleanUp();
1472        return TRUE;
1473      }
1474    }
1475    v=v->next;
1476  }
1477  rv->CleanUp();
1478  return nok;
1479}
1480
1481BOOLEAN iiCheckRing(int i)
1482{
1483  if (currRing==NULL)
1484  {
1485    #ifdef SIQ
1486    if (siq<=0)
1487    {
1488    #endif
1489      if (RingDependend(i))
1490      {
1491        WerrorS("no ring active");
1492        return TRUE;
1493      }
1494    #ifdef SIQ
1495    }
1496    #endif
1497  }
1498  return FALSE;
1499}
1500
1501poly    iiHighCorner(ideal I, int ak)
1502{
1503  int i;
1504  if(!idIsZeroDim(I)) return NULL; // not zero-dim.
1505  poly po=NULL;
1506  if (rHasLocalOrMixedOrdering_currRing())
1507  {
1508    scComputeHC(I,currRing->qideal,ak,po);
1509    if (po!=NULL)
1510    {
1511      pGetCoeff(po)=nInit(1);
1512      for (i=rVar(currRing); i>0; i--)
1513      {
1514        if (pGetExp(po, i) > 0) pDecrExp(po,i);
1515      }
1516      pSetComp(po,ak);
1517      pSetm(po);
1518    }
1519  }
1520  else
1521    po=pOne();
1522  return po;
1523}
1524
1525void iiCheckPack(package &p)
1526{
1527  if (p!=basePack)
1528  {
1529    idhdl t=basePack->idroot;
1530    while ((t!=NULL) && (IDTYP(t)!=PACKAGE_CMD) && (IDPACKAGE(t)!=p)) t=t->next;
1531    if (t==NULL)
1532    {
1533      WarnS("package not found\n");
1534      p=basePack;
1535    }
1536  }
1537}
1538
1539idhdl rDefault(const char *s)
1540{
1541  idhdl tmp=NULL;
1542
1543  if (s!=NULL) tmp = enterid(s, myynest, RING_CMD, &IDROOT);
1544  if (tmp==NULL) return NULL;
1545
1546// if ((currRing->ppNoether)!=NULL) pDelete(&(currRing->ppNoether));
1547  if (sLastPrinted.RingDependend())
1548  {
1549    sLastPrinted.CleanUp();
1550    memset(&sLastPrinted,0,sizeof(sleftv));
1551  }
1552
1553  ring r = IDRING(tmp);
1554
1555  r->cf = nInitChar(n_Zp, (void*)32003); //   r->cf->ch = 32003;
1556  r->N      = 3;
1557  /*r->P     = 0; Alloc0 in idhdl::set, ipid.cc*/
1558  /*names*/
1559  r->names = (char **) omAlloc0(3 * sizeof(char_ptr));
1560  r->names[0]  = omStrDup("x");
1561  r->names[1]  = omStrDup("y");
1562  r->names[2]  = omStrDup("z");
1563  /*weights: entries for 3 blocks: NULL*/
1564  r->wvhdl = (int **)omAlloc0(3 * sizeof(int_ptr));
1565  /*order: dp,C,0*/
1566  r->order = (int *) omAlloc(3 * sizeof(int *));
1567  r->block0 = (int *)omAlloc0(3 * sizeof(int *));
1568  r->block1 = (int *)omAlloc0(3 * sizeof(int *));
1569  /* ringorder dp for the first block: var 1..3 */
1570  r->order[0]  = ringorder_dp;
1571  r->block0[0] = 1;
1572  r->block1[0] = 3;
1573  /* ringorder C for the second block: no vars */
1574  r->order[1]  = ringorder_C;
1575  /* the last block: everything is 0 */
1576  r->order[2]  = 0;
1577
1578  /* complete ring intializations */
1579  rComplete(r);
1580  rSetHdl(tmp);
1581  return currRingHdl;
1582}
1583
1584idhdl rFindHdl(ring r, idhdl n)
1585{
1586  idhdl h=rSimpleFindHdl(r,IDROOT,n);
1587  if (h!=NULL)  return h;
1588  if (IDROOT!=basePack->idroot) h=rSimpleFindHdl(r,basePack->idroot,n);
1589  if (h!=NULL)  return h;
1590  proclevel *p=procstack;
1591  while(p!=NULL)
1592  {
1593    if ((p->cPack!=basePack)
1594    && (p->cPack!=currPack))
1595      h=rSimpleFindHdl(r,p->cPack->idroot,n);
1596    if (h!=NULL)  return h;
1597    p=p->next;
1598  }
1599  idhdl tmp=basePack->idroot;
1600  while (tmp!=NULL)
1601  {
1602    if (IDTYP(tmp)==PACKAGE_CMD)
1603      h=rSimpleFindHdl(r,IDPACKAGE(tmp)->idroot,n);
1604    if (h!=NULL)  return h;
1605    tmp=IDNEXT(tmp);
1606  }
1607  return NULL;
1608}
1609
1610void rDecomposeCF(leftv h,const ring r,const ring R)
1611{
1612  lists L=(lists)omAlloc0Bin(slists_bin);
1613  L->Init(4);
1614  h->rtyp=LIST_CMD;
1615  h->data=(void *)L;
1616  // 0: char/ cf - ring
1617  // 1: list (var)
1618  // 2: list (ord)
1619  // 3: qideal
1620  // ----------------------------------------
1621  // 0: char/ cf - ring
1622  L->m[0].rtyp=INT_CMD;
1623  L->m[0].data=(void *)(long)r->cf->ch;
1624  // ----------------------------------------
1625  // 1: list (var)
1626  lists LL=(lists)omAlloc0Bin(slists_bin);
1627  LL->Init(r->N);
1628  int i;
1629  for(i=0; i<r->N; i++)
1630  {
1631    LL->m[i].rtyp=STRING_CMD;
1632    LL->m[i].data=(void *)omStrDup(r->names[i]);
1633  }
1634  L->m[1].rtyp=LIST_CMD;
1635  L->m[1].data=(void *)LL;
1636  // ----------------------------------------
1637  // 2: list (ord)
1638  LL=(lists)omAlloc0Bin(slists_bin);
1639  i=rBlocks(r)-1;
1640  LL->Init(i);
1641  i--;
1642  lists LLL;
1643  for(; i>=0; i--)
1644  {
1645    intvec *iv;
1646    int j;
1647    LL->m[i].rtyp=LIST_CMD;
1648    LLL=(lists)omAlloc0Bin(slists_bin);
1649    LLL->Init(2);
1650    LLL->m[0].rtyp=STRING_CMD;
1651    LLL->m[0].data=(void *)omStrDup(rSimpleOrdStr(r->order[i]));
1652    if (r->block1[i]-r->block0[i] >=0 )
1653    {
1654      j=r->block1[i]-r->block0[i];
1655      if(r->order[i]==ringorder_M) j=(j+1)*(j+1)-1;
1656      iv=new intvec(j+1);
1657      if ((r->wvhdl!=NULL) && (r->wvhdl[i]!=NULL))
1658      {
1659        for(;j>=0; j--) (*iv)[j]=r->wvhdl[i][j];
1660      }
1661      else switch (r->order[i])
1662      {
1663        case ringorder_dp:
1664        case ringorder_Dp:
1665        case ringorder_ds:
1666        case ringorder_Ds:
1667        case ringorder_lp:
1668          for(;j>=0; j--) (*iv)[j]=1;
1669          break;
1670        default: /* do nothing */;
1671      }
1672    }
1673    else
1674    {
1675      iv=new intvec(1);
1676    }
1677    LLL->m[1].rtyp=INTVEC_CMD;
1678    LLL->m[1].data=(void *)iv;
1679    LL->m[i].data=(void *)LLL;
1680  }
1681  L->m[2].rtyp=LIST_CMD;
1682  L->m[2].data=(void *)LL;
1683  // ----------------------------------------
1684  // 3: qideal
1685  L->m[3].rtyp=IDEAL_CMD;
1686  if (nCoeff_is_transExt(R->cf))
1687    L->m[3].data=(void *)idInit(1,1);
1688  else
1689  {
1690    ideal q=idInit(IDELEMS(r->qideal));
1691    q->m[0]=p_Init(R);
1692    pSetCoeff0(q->m[0],(number)(r->qideal->m[0]));
1693    L->m[3].data=(void *)q;
1694//    I->m[0] = pNSet(R->minpoly);
1695  }
1696  // ----------------------------------------
1697}
1698#ifdef SINGULAR_4_1
1699static void rDecomposeC_41(leftv h,const coeffs C)
1700/* field is R or C */
1701{
1702  lists L=(lists)omAlloc0Bin(slists_bin);
1703  if (nCoeff_is_long_C(C)) L->Init(3);
1704  else                     L->Init(2);
1705  h->rtyp=LIST_CMD;
1706  h->data=(void *)L;
1707  // 0: char/ cf - ring
1708  // 1: list (var)
1709  // 2: list (ord)
1710  // ----------------------------------------
1711  // 0: char/ cf - ring
1712  L->m[0].rtyp=INT_CMD;
1713  L->m[0].data=(void *)0;
1714  // ----------------------------------------
1715  // 1:
1716  lists LL=(lists)omAlloc0Bin(slists_bin);
1717  LL->Init(2);
1718    LL->m[0].rtyp=INT_CMD;
1719    LL->m[0].data=(void *)(long)si_max(C->float_len,SHORT_REAL_LENGTH/2);
1720    LL->m[1].rtyp=INT_CMD;
1721    LL->m[1].data=(void *)(long)si_max(C->float_len2,SHORT_REAL_LENGTH);
1722  L->m[1].rtyp=LIST_CMD;
1723  L->m[1].data=(void *)LL;
1724  // ----------------------------------------
1725  // 2: list (par)
1726  if (nCoeff_is_long_C(C))
1727  {
1728    L->m[2].rtyp=STRING_CMD;
1729    L->m[2].data=(void *)omStrDup(*n_ParameterNames(C));
1730  }
1731  // ----------------------------------------
1732}
1733#endif
1734static void rDecomposeC(leftv h,const ring R)
1735/* field is R or C */
1736{
1737  lists L=(lists)omAlloc0Bin(slists_bin);
1738  if (rField_is_long_C(R)) L->Init(3);
1739  else                     L->Init(2);
1740  h->rtyp=LIST_CMD;
1741  h->data=(void *)L;
1742  // 0: char/ cf - ring
1743  // 1: list (var)
1744  // 2: list (ord)
1745  // ----------------------------------------
1746  // 0: char/ cf - ring
1747  L->m[0].rtyp=INT_CMD;
1748  L->m[0].data=(void *)0;
1749  // ----------------------------------------
1750  // 1:
1751  lists LL=(lists)omAlloc0Bin(slists_bin);
1752  LL->Init(2);
1753    LL->m[0].rtyp=INT_CMD;
1754    LL->m[0].data=(void *)(long)si_max(R->cf->float_len,SHORT_REAL_LENGTH/2);
1755    LL->m[1].rtyp=INT_CMD;
1756    LL->m[1].data=(void *)(long)si_max(R->cf->float_len2,SHORT_REAL_LENGTH);
1757  L->m[1].rtyp=LIST_CMD;
1758  L->m[1].data=(void *)LL;
1759  // ----------------------------------------
1760  // 2: list (par)
1761  if (rField_is_long_C(R))
1762  {
1763    L->m[2].rtyp=STRING_CMD;
1764    L->m[2].data=(void *)omStrDup(*rParameter(R));
1765  }
1766  // ----------------------------------------
1767}
1768
1769#ifdef SINGULAR_4_1
1770#ifdef HAVE_RINGS
1771void rDecomposeRing_41(leftv h,const coeffs C)
1772/* field is R or C */
1773{
1774  lists L=(lists)omAlloc0Bin(slists_bin);
1775  if (nCoeff_is_Ring(C)) L->Init(1);
1776  else                   L->Init(2);
1777  h->rtyp=LIST_CMD;
1778  h->data=(void *)L;
1779  // 0: char/ cf - ring
1780  // 1: list (module)
1781  // ----------------------------------------
1782  // 0: char/ cf - ring
1783  L->m[0].rtyp=STRING_CMD;
1784  L->m[0].data=(void *)omStrDup("integer");
1785  // ----------------------------------------
1786  // 1: modulo
1787  if (nCoeff_is_Ring_Z(C)) return;
1788  lists LL=(lists)omAlloc0Bin(slists_bin);
1789  LL->Init(2);
1790  LL->m[0].rtyp=BIGINT_CMD;
1791  LL->m[0].data=nlMapGMP((number) C->modBase, C, coeffs_BIGINT);
1792  LL->m[1].rtyp=INT_CMD;
1793  LL->m[1].data=(void *) C->modExponent;
1794  L->m[1].rtyp=LIST_CMD;
1795  L->m[1].data=(void *)LL;
1796}
1797#endif
1798#endif
1799
1800void rDecomposeRing(leftv h,const ring R)
1801/* field is R or C */
1802{
1803#ifdef HAVE_RINGS
1804  lists L=(lists)omAlloc0Bin(slists_bin);
1805  if (rField_is_Ring_Z(R)) L->Init(1);
1806  else                     L->Init(2);
1807  h->rtyp=LIST_CMD;
1808  h->data=(void *)L;
1809  // 0: char/ cf - ring
1810  // 1: list (module)
1811  // ----------------------------------------
1812  // 0: char/ cf - ring
1813  L->m[0].rtyp=STRING_CMD;
1814  L->m[0].data=(void *)omStrDup("integer");
1815  // ----------------------------------------
1816  // 1: module
1817  if (rField_is_Ring_Z(R)) return;
1818  lists LL=(lists)omAlloc0Bin(slists_bin);
1819  LL->Init(2);
1820  LL->m[0].rtyp=BIGINT_CMD;
1821  LL->m[0].data=nlMapGMP((number) R->cf->modBase, R->cf, R->cf); // TODO: what is this?? // extern number nlMapGMP(number from, const coeffs src, const coeffs dst); // FIXME: replace with n_InitMPZ(R->cf->modBase, coeffs_BIGINT); ?
1822  LL->m[1].rtyp=INT_CMD;
1823  LL->m[1].data=(void *) R->cf->modExponent;
1824  L->m[1].rtyp=LIST_CMD;
1825  L->m[1].data=(void *)LL;
1826#else
1827  WerrorS("rDecomposeRing");
1828#endif
1829}
1830
1831
1832#ifdef SINGULAR_4_1
1833BOOLEAN rDecompose_CF(leftv res,const coeffs C)
1834{
1835  assume( C != NULL );
1836
1837  // sanity check: require currRing==r for rings with polynomial data
1838  if ( nCoeff_is_algExt(C) && (C != currRing->cf))
1839  {
1840    WerrorS("ring with polynomial data must be the base ring or compatible");
1841    return TRUE;
1842  }
1843  if (nCoeff_is_numeric(C))
1844  {
1845    rDecomposeC_41(res,C);
1846  }
1847#ifdef HAVE_RINGS
1848  else if (nCoeff_is_Ring(C))
1849  {
1850    rDecomposeRing_41(res,C);
1851  }
1852#endif
1853  else if ( C->extRing!=NULL )// nCoeff_is_algExt(r->cf))
1854  {
1855    rDecomposeCF(res, C->extRing, currRing);
1856  }
1857  else if(nCoeff_is_GF(C))
1858  {
1859    lists Lc=(lists)omAlloc0Bin(slists_bin);
1860    Lc->Init(4);
1861    // char:
1862    Lc->m[0].rtyp=INT_CMD;
1863    Lc->m[0].data=(void*)(long)C->m_nfCharQ;
1864    // var:
1865    lists Lv=(lists)omAlloc0Bin(slists_bin);
1866    Lv->Init(1);
1867    Lv->m[0].rtyp=STRING_CMD;
1868    Lv->m[0].data=(void *)omStrDup(*n_ParameterNames(C));
1869    Lc->m[1].rtyp=LIST_CMD;
1870    Lc->m[1].data=(void*)Lv;
1871    // ord:
1872    lists Lo=(lists)omAlloc0Bin(slists_bin);
1873    Lo->Init(1);
1874    lists Loo=(lists)omAlloc0Bin(slists_bin);
1875    Loo->Init(2);
1876    Loo->m[0].rtyp=STRING_CMD;
1877    Loo->m[0].data=(void *)omStrDup(rSimpleOrdStr(ringorder_lp));
1878
1879    intvec *iv=new intvec(1); (*iv)[0]=1;
1880    Loo->m[1].rtyp=INTVEC_CMD;
1881    Loo->m[1].data=(void *)iv;
1882
1883    Lo->m[0].rtyp=LIST_CMD;
1884    Lo->m[0].data=(void*)Loo;
1885
1886    Lc->m[2].rtyp=LIST_CMD;
1887    Lc->m[2].data=(void*)Lo;
1888    // q-ideal:
1889    Lc->m[3].rtyp=IDEAL_CMD;
1890    Lc->m[3].data=(void *)idInit(1,1);
1891    // ----------------------
1892    res->rtyp=LIST_CMD;
1893    res->data=(void*)Lc;
1894  }
1895  else
1896  {
1897    res->rtyp=INT_CMD;
1898    res->data=(void *)(long)C->ch;
1899  }
1900  // ----------------------------------------
1901  return FALSE;
1902}
1903#endif
1904
1905#ifdef SINGULAR_4_1
1906lists rDecompose_list_cf(const ring r)
1907{
1908  assume( r != NULL );
1909  const coeffs C = r->cf;
1910  assume( C != NULL );
1911
1912  // sanity check: require currRing==r for rings with polynomial data
1913  if ( (r!=currRing) && (
1914           (nCoeff_is_algExt(C) && (C != currRing->cf))
1915        || (r->qideal != NULL)
1916#ifdef HAVE_PLURAL
1917        || (rIsPluralRing(r))
1918#endif
1919                        )
1920     )
1921  {
1922    WerrorS("ring with polynomial data must be the base ring or compatible");
1923    return NULL;
1924  }
1925  // 0: char/ cf - ring
1926  // 1: list (var)
1927  // 2: list (ord)
1928  // 3: qideal
1929  // possibly:
1930  // 4: C
1931  // 5: D
1932  lists L=(lists)omAlloc0Bin(slists_bin);
1933  if (rIsPluralRing(r))
1934    L->Init(6);
1935  else
1936    L->Init(4);
1937  // ----------------------------------------
1938  // 0: char/ cf - ring
1939  L->m[0].rtyp=CRING_CMD;
1940  L->m[0].data=(char*)r->cf; r->cf->ref++;
1941  // ----------------------------------------
1942  // 1: list (var)
1943  lists LL=(lists)omAlloc0Bin(slists_bin);
1944  LL->Init(r->N);
1945  int i;
1946  for(i=0; i<r->N; i++)
1947  {
1948    LL->m[i].rtyp=STRING_CMD;
1949    LL->m[i].data=(void *)omStrDup(r->names[i]);
1950  }
1951  L->m[1].rtyp=LIST_CMD;
1952  L->m[1].data=(void *)LL;
1953  // ----------------------------------------
1954  // 2: list (ord)
1955  LL=(lists)omAlloc0Bin(slists_bin);
1956  i=rBlocks(r)-1;
1957  LL->Init(i);
1958  i--;
1959  lists LLL;
1960  for(; i>=0; i--)
1961  {
1962    intvec *iv;
1963    int j;
1964    LL->m[i].rtyp=LIST_CMD;
1965    LLL=(lists)omAlloc0Bin(slists_bin);
1966    LLL->Init(2);
1967    LLL->m[0].rtyp=STRING_CMD;
1968    LLL->m[0].data=(void *)omStrDup(rSimpleOrdStr(r->order[i]));
1969
1970    if(r->order[i] == ringorder_IS) //  || r->order[i] == ringorder_s || r->order[i] == ringorder_S)
1971    {
1972      assume( r->block0[i] == r->block1[i] );
1973      const int s = r->block0[i];
1974      assume( -2 < s && s < 2);
1975
1976      iv=new intvec(1);
1977      (*iv)[0] = s;
1978    }
1979    else if (r->block1[i]-r->block0[i] >=0 )
1980    {
1981      int bl=j=r->block1[i]-r->block0[i];
1982      if (r->order[i]==ringorder_M)
1983      {
1984        j=(j+1)*(j+1)-1;
1985        bl=j+1;
1986      }
1987      else if (r->order[i]==ringorder_am)
1988      {
1989        j+=r->wvhdl[i][bl+1];
1990      }
1991      iv=new intvec(j+1);
1992      if ((r->wvhdl!=NULL) && (r->wvhdl[i]!=NULL))
1993      {
1994        for(;j>=0; j--) (*iv)[j]=r->wvhdl[i][j+(j>bl)];
1995      }
1996      else switch (r->order[i])
1997      {
1998        case ringorder_dp:
1999        case ringorder_Dp:
2000        case ringorder_ds:
2001        case ringorder_Ds:
2002        case ringorder_lp:
2003          for(;j>=0; j--) (*iv)[j]=1;
2004          break;
2005        default: /* do nothing */;
2006      }
2007    }
2008    else
2009    {
2010      iv=new intvec(1);
2011    }
2012    LLL->m[1].rtyp=INTVEC_CMD;
2013    LLL->m[1].data=(void *)iv;
2014    LL->m[i].data=(void *)LLL;
2015  }
2016  L->m[2].rtyp=LIST_CMD;
2017  L->m[2].data=(void *)LL;
2018  // ----------------------------------------
2019  // 3: qideal
2020  L->m[3].rtyp=IDEAL_CMD;
2021  if (r->qideal==NULL)
2022    L->m[3].data=(void *)idInit(1,1);
2023  else
2024    L->m[3].data=(void *)idCopy(r->qideal);
2025  // ----------------------------------------
2026#ifdef HAVE_PLURAL // NC! in rDecompose
2027  if (rIsPluralRing(r))
2028  {
2029    L->m[4].rtyp=MATRIX_CMD;
2030    L->m[4].data=(void *)mp_Copy(r->GetNC()->C, r, r);
2031    L->m[5].rtyp=MATRIX_CMD;
2032    L->m[5].data=(void *)mp_Copy(r->GetNC()->D, r, r);
2033  }
2034#endif
2035  return L;
2036}
2037#endif
2038
2039lists rDecompose(const ring r)
2040{
2041  assume( r != NULL );
2042  const coeffs C = r->cf;
2043  assume( C != NULL );
2044
2045  // sanity check: require currRing==r for rings with polynomial data
2046  if ( (r!=currRing) && (
2047           (nCoeff_is_algExt(C) && (C != currRing->cf))
2048        || (r->qideal != NULL)
2049#ifdef HAVE_PLURAL
2050        || (rIsPluralRing(r))
2051#endif
2052                        )
2053     )
2054  {
2055    WerrorS("ring with polynomial data must be the base ring or compatible");
2056    return NULL;
2057  }
2058  // 0: char/ cf - ring
2059  // 1: list (var)
2060  // 2: list (ord)
2061  // 3: qideal
2062  // possibly:
2063  // 4: C
2064  // 5: D
2065  lists L=(lists)omAlloc0Bin(slists_bin);
2066  if (rIsPluralRing(r))
2067    L->Init(6);
2068  else
2069    L->Init(4);
2070  // ----------------------------------------
2071  // 0: char/ cf - ring
2072  if (rField_is_numeric(r))
2073  {
2074    rDecomposeC(&(L->m[0]),r);
2075  }
2076  else if (rField_is_Ring(r))
2077  {
2078    rDecomposeRing(&(L->m[0]),r);
2079  }
2080  else if ( r->cf->extRing!=NULL )// nCoeff_is_algExt(r->cf))
2081  {
2082    rDecomposeCF(&(L->m[0]), r->cf->extRing, r);
2083  }
2084  else if(rField_is_GF(r))
2085  {
2086    lists Lc=(lists)omAlloc0Bin(slists_bin);
2087    Lc->Init(4);
2088    // char:
2089    Lc->m[0].rtyp=INT_CMD;
2090    Lc->m[0].data=(void*)(long)r->cf->m_nfCharQ;
2091    // var:
2092    lists Lv=(lists)omAlloc0Bin(slists_bin);
2093    Lv->Init(1);
2094    Lv->m[0].rtyp=STRING_CMD;
2095    Lv->m[0].data=(void *)omStrDup(*rParameter(r));
2096    Lc->m[1].rtyp=LIST_CMD;
2097    Lc->m[1].data=(void*)Lv;
2098    // ord:
2099    lists Lo=(lists)omAlloc0Bin(slists_bin);
2100    Lo->Init(1);
2101    lists Loo=(lists)omAlloc0Bin(slists_bin);
2102    Loo->Init(2);
2103    Loo->m[0].rtyp=STRING_CMD;
2104    Loo->m[0].data=(void *)omStrDup(rSimpleOrdStr(ringorder_lp));
2105
2106    intvec *iv=new intvec(1); (*iv)[0]=1;
2107    Loo->m[1].rtyp=INTVEC_CMD;
2108    Loo->m[1].data=(void *)iv;
2109
2110    Lo->m[0].rtyp=LIST_CMD;
2111    Lo->m[0].data=(void*)Loo;
2112
2113    Lc->m[2].rtyp=LIST_CMD;
2114    Lc->m[2].data=(void*)Lo;
2115    // q-ideal:
2116    Lc->m[3].rtyp=IDEAL_CMD;
2117    Lc->m[3].data=(void *)idInit(1,1);
2118    // ----------------------
2119    L->m[0].rtyp=LIST_CMD;
2120    L->m[0].data=(void*)Lc;
2121  }
2122  else
2123  {
2124    L->m[0].rtyp=INT_CMD;
2125    L->m[0].data=(void *)(long)r->cf->ch;
2126  }
2127  // ----------------------------------------
2128  // 1: list (var)
2129  lists LL=(lists)omAlloc0Bin(slists_bin);
2130  LL->Init(r->N);
2131  int i;
2132  for(i=0; i<r->N; i++)
2133  {
2134    LL->m[i].rtyp=STRING_CMD;
2135    LL->m[i].data=(void *)omStrDup(r->names[i]);
2136  }
2137  L->m[1].rtyp=LIST_CMD;
2138  L->m[1].data=(void *)LL;
2139  // ----------------------------------------
2140  // 2: list (ord)
2141  LL=(lists)omAlloc0Bin(slists_bin);
2142  i=rBlocks(r)-1;
2143  LL->Init(i);
2144  i--;
2145  lists LLL;
2146  for(; i>=0; i--)
2147  {
2148    intvec *iv;
2149    int j;
2150    LL->m[i].rtyp=LIST_CMD;
2151    LLL=(lists)omAlloc0Bin(slists_bin);
2152    LLL->Init(2);
2153    LLL->m[0].rtyp=STRING_CMD;
2154    LLL->m[0].data=(void *)omStrDup(rSimpleOrdStr(r->order[i]));
2155
2156    if(r->order[i] == ringorder_IS) //  || r->order[i] == ringorder_s || r->order[i] == ringorder_S)
2157    {
2158      assume( r->block0[i] == r->block1[i] );
2159      const int s = r->block0[i];
2160      assume( -2 < s && s < 2);
2161
2162      iv=new intvec(1);
2163      (*iv)[0] = s;
2164    }
2165    else if (r->block1[i]-r->block0[i] >=0 )
2166    {
2167      int bl=j=r->block1[i]-r->block0[i];
2168      if (r->order[i]==ringorder_M)
2169      {
2170        j=(j+1)*(j+1)-1;
2171        bl=j+1;
2172      }
2173      else if (r->order[i]==ringorder_am)
2174      {
2175        j+=r->wvhdl[i][bl+1];
2176      }
2177      iv=new intvec(j+1);
2178      if ((r->wvhdl!=NULL) && (r->wvhdl[i]!=NULL))
2179      {
2180        for(;j>=0; j--) (*iv)[j]=r->wvhdl[i][j+(j>bl)];
2181      }
2182      else switch (r->order[i])
2183      {
2184        case ringorder_dp:
2185        case ringorder_Dp:
2186        case ringorder_ds:
2187        case ringorder_Ds:
2188        case ringorder_lp:
2189          for(;j>=0; j--) (*iv)[j]=1;
2190          break;
2191        default: /* do nothing */;
2192      }
2193    }
2194    else
2195    {
2196      iv=new intvec(1);
2197    }
2198    LLL->m[1].rtyp=INTVEC_CMD;
2199    LLL->m[1].data=(void *)iv;
2200    LL->m[i].data=(void *)LLL;
2201  }
2202  L->m[2].rtyp=LIST_CMD;
2203  L->m[2].data=(void *)LL;
2204  // ----------------------------------------
2205  // 3: qideal
2206  L->m[3].rtyp=IDEAL_CMD;
2207  if (r->qideal==NULL)
2208    L->m[3].data=(void *)idInit(1,1);
2209  else
2210    L->m[3].data=(void *)idCopy(r->qideal);
2211  // ----------------------------------------
2212#ifdef HAVE_PLURAL // NC! in rDecompose
2213  if (rIsPluralRing(r))
2214  {
2215    L->m[4].rtyp=MATRIX_CMD;
2216    L->m[4].data=(void *)mp_Copy(r->GetNC()->C, r, r);
2217    L->m[5].rtyp=MATRIX_CMD;
2218    L->m[5].data=(void *)mp_Copy(r->GetNC()->D, r, r);
2219  }
2220#endif
2221  return L;
2222}
2223
2224void rComposeC(lists L, ring R)
2225/* field is R or C */
2226{
2227  // ----------------------------------------
2228  // 0: char/ cf - ring
2229  if ((L->m[0].rtyp!=INT_CMD) || (L->m[0].data!=(char *)0))
2230  {
2231    Werror("invald coeff. field description, expecting 0");
2232    return;
2233  }
2234//  R->cf->ch=0;
2235  // ----------------------------------------
2236  // 1:
2237  if (L->m[1].rtyp!=LIST_CMD)
2238  {
2239    Werror("invald coeff. field description, expecting precision list");
2240    return;
2241  }
2242  lists LL=(lists)L->m[1].data;
2243  if (((LL->nr!=2)
2244    || (LL->m[0].rtyp!=INT_CMD)
2245    || (LL->m[1].rtyp!=INT_CMD))
2246  && ((LL->nr!=1)
2247    || (LL->m[0].rtyp!=INT_CMD)))
2248  {
2249    Werror("invald coeff. field description list");
2250    return;
2251  }
2252  int r1=(int)(long)LL->m[0].data;
2253  int r2=(int)(long)LL->m[1].data;
2254  if (L->nr==2) // complex
2255    R->cf = nInitChar(n_long_C, NULL);
2256  else if ((r1<=SHORT_REAL_LENGTH)
2257  && (r2=SHORT_REAL_LENGTH))
2258    R->cf = nInitChar(n_R, NULL);
2259  else
2260  {
2261    LongComplexInfo* p = (LongComplexInfo *)omAlloc0(sizeof(LongComplexInfo));
2262    p->float_len=r1;
2263    p->float_len2=r2;
2264    R->cf = nInitChar(n_long_R, NULL);
2265  }
2266
2267  if ((r1<=SHORT_REAL_LENGTH)   // should go into nInitChar
2268  && (r2=SHORT_REAL_LENGTH))
2269  {
2270    R->cf->float_len=SHORT_REAL_LENGTH/2;
2271    R->cf->float_len2=SHORT_REAL_LENGTH;
2272  }
2273  else
2274  {
2275    R->cf->float_len=si_min(r1,32767);
2276    R->cf->float_len2=si_min(r2,32767);
2277  }
2278  // ----------------------------------------
2279  // 2: list (par)
2280  if (L->nr==2)
2281  {
2282    //R->cf->extRing->N=1;
2283    if (L->m[2].rtyp!=STRING_CMD)
2284    {
2285      Werror("invald coeff. field description, expecting parameter name");
2286      return;
2287    }
2288    //(rParameter(R))=(char**)omAlloc0(rPar(R)*sizeof(char_ptr));
2289    rParameter(R)[0]=omStrDup((char *)L->m[2].data);
2290  }
2291  // ----------------------------------------
2292}
2293
2294#ifdef HAVE_RINGS
2295void rComposeRing(lists L, ring R)
2296/* field is R or C */
2297{
2298  // ----------------------------------------
2299  // 0: string: integer
2300  // no further entries --> Z
2301  mpz_ptr modBase = NULL;
2302  unsigned int modExponent = 1;
2303
2304  modBase = (mpz_ptr) omAlloc(sizeof(mpz_t));
2305  if (L->nr == 0)
2306  {
2307    mpz_init_set_ui(modBase,0);
2308    modExponent = 1;
2309  }
2310  // ----------------------------------------
2311  // 1:
2312  else
2313  {
2314    if (L->m[1].rtyp!=LIST_CMD) Werror("invald data, expecting list of numbers");
2315    lists LL=(lists)L->m[1].data;
2316    if ((LL->nr >= 0) && LL->m[0].rtyp == BIGINT_CMD)
2317    {
2318      number tmp= (number) LL->m[0].data; // never use CopyD() on list elements
2319                                    // assume that tmp is integer, not rational
2320      n_MPZ (modBase, tmp, coeffs_BIGINT);
2321    }
2322    else if (LL->nr >= 0 && LL->m[0].rtyp == INT_CMD)
2323    {
2324      mpz_init_set_ui(modBase,(unsigned long) LL->m[0].data);
2325    }
2326    else
2327    {
2328      mpz_init_set_ui(modBase,0);
2329    }
2330    if (LL->nr >= 1)
2331    {
2332      modExponent = (unsigned long) LL->m[1].data;
2333    }
2334    else
2335    {
2336      modExponent = 1;
2337    }
2338  }
2339  // ----------------------------------------
2340  if ((mpz_cmp_ui(modBase, 1) == 0) && (mpz_cmp_ui(modBase, 0) < 0))
2341  {
2342    Werror("Wrong ground ring specification (module is 1)");
2343    return;
2344  }
2345  if (modExponent < 1)
2346  {
2347    Werror("Wrong ground ring specification (exponent smaller than 1");
2348    return;
2349  }
2350  // module is 0 ---> integers
2351  if (mpz_cmp_ui(modBase, 0) == 0)
2352  {
2353    R->cf=nInitChar(n_Z,NULL);
2354  }
2355  // we have an exponent
2356  else if (modExponent > 1)
2357  {
2358    //R->cf->ch = R->cf->modExponent;
2359    if ((mpz_cmp_ui(modBase, 2) == 0) && (modExponent <= 8*sizeof(unsigned long)))
2360    {
2361      /* this branch should be active for modExponent = 2..32 resp. 2..64,
2362           depending on the size of a long on the respective platform */
2363      R->cf=nInitChar(n_Z2m,(void*)(long)modExponent);       // Use Z/2^ch
2364      omFreeSize (modBase, sizeof(mpz_t));
2365    }
2366    else
2367    {
2368      //ringtype 3
2369      ZnmInfo info;
2370      info.base= modBase;
2371      info.exp= modExponent;
2372      R->cf=nInitChar(n_Znm,(void*) &info);
2373    }
2374  }
2375  // just a module m > 1
2376  else
2377  {
2378    //ringtype = 2;
2379    //const int ch = mpz_get_ui(modBase);
2380    ZnmInfo info;
2381    info.base= modBase;
2382    info.exp= modExponent;
2383    R->cf=nInitChar(n_Zn,(void*) &info);
2384  }
2385}
2386#endif
2387
2388static void rRenameVars(ring R)
2389{
2390  int i,j;
2391  BOOLEAN ch;
2392  do
2393  {
2394    ch=0;
2395    for(i=0;i<R->N-1;i++)
2396    {
2397      for(j=i+1;j<R->N;j++)
2398      {
2399        if (strcmp(R->names[i],R->names[j])==0)
2400        {
2401          ch=TRUE;
2402          Warn("name conflict var(%d) and var(%d): `%s`, rename to `@%s`",i+1,j+1,R->names[i],R->names[i]);
2403          omFree(R->names[j]);
2404          R->names[j]=(char *)omAlloc(2+strlen(R->names[i]));
2405          sprintf(R->names[j],"@%s",R->names[i]);
2406        }
2407      }
2408    }
2409  }
2410  while (ch);
2411  for(i=0;i<rPar(R); i++)
2412  {
2413    for(j=0;j<R->N;j++)
2414    {
2415      if (strcmp(rParameter(R)[i],R->names[j])==0)
2416      {
2417        Warn("name conflict par(%d) and var(%d): `%s`, renaming the VARIABLE to `@@(%d)`",i+1,j+1,R->names[j],i+1);
2418//        omFree(rParameter(R)[i]);
2419//        rParameter(R)[i]=(char *)omAlloc(10);
2420//        sprintf(rParameter(R)[i],"@@(%d)",i+1);
2421        omFree(R->names[j]);
2422        R->names[j]=(char *)omAlloc(10);
2423        sprintf(R->names[j],"@@(%d)",i+1);
2424      }
2425    }
2426  }
2427}
2428
2429static inline BOOLEAN rComposeVar(const lists  L, ring R)
2430{
2431  assume(R!=NULL);
2432  if (L->m[1].Typ()==LIST_CMD)
2433  {
2434    lists v=(lists)L->m[1].Data();
2435    R->N = v->nr+1;
2436    if (R->N<=0)
2437    {
2438      WerrorS("no ring variables");
2439      return TRUE;
2440    }
2441    R->names   = (char **)omAlloc0(R->N * sizeof(char_ptr));
2442    int i;
2443    for(i=0;i<R->N;i++)
2444    {
2445      if (v->m[i].Typ()==STRING_CMD)
2446        R->names[i]=omStrDup((char *)v->m[i].Data());
2447      else if (v->m[i].Typ()==POLY_CMD)
2448      {
2449        poly p=(poly)v->m[i].Data();
2450        int nr=pIsPurePower(p);
2451        if (nr>0)
2452          R->names[i]=omStrDup(currRing->names[nr-1]);
2453        else
2454        {
2455          Werror("var name %d must be a string or a ring variable",i+1);
2456          return TRUE;
2457        }
2458      }
2459      else
2460      {
2461        Werror("var name %d must be `string`",i+1);
2462        return TRUE;
2463      }
2464    }
2465  }
2466  else
2467  {
2468    WerrorS("variable must be given as `list`");
2469    return TRUE;
2470  }
2471  return FALSE;
2472}
2473
2474static inline BOOLEAN rComposeOrder(const lists  L, const BOOLEAN check_comp, ring R)
2475{
2476  assume(R!=NULL);
2477  long bitmask=0L;
2478  if (L->m[2].Typ()==LIST_CMD)
2479  {
2480    lists v=(lists)L->m[2].Data();
2481    int n= v->nr+2;
2482    int j_in_R,j_in_L;
2483    // do we have an entry "L",... ?: set bitmask
2484    for (int j=0; j < n-1; j++)
2485    {
2486      if (v->m[j].Typ()==LIST_CMD)
2487      {
2488        lists vv=(lists)v->m[j].Data();
2489        if ((vv->nr==1)
2490        &&(vv->m[0].Typ()==STRING_CMD)
2491        &&(strcmp((char*)vv->m[0].Data(),"L")==0))
2492        {
2493          number nn=(number)vv->m[1].Data();
2494          if (vv->m[1].Typ()==BIGINT_CMD)
2495            bitmask=n_Int(nn,coeffs_BIGINT);
2496          else if (vv->m[1].Typ()==INT_CMD)
2497            bitmask=(long)nn;
2498          else
2499          {
2500            Werror("illegal argument for pseudo ordering L: %d",vv->m[1].Typ());
2501            return TRUE;
2502          }
2503          break;
2504        }
2505      }
2506    }
2507    if (bitmask!=0) n--;
2508
2509    // initialize fields of R
2510    R->order=(int *)omAlloc0(n*sizeof(int));
2511    R->block0=(int *)omAlloc0(n*sizeof(int));
2512    R->block1=(int *)omAlloc0(n*sizeof(int));
2513    R->wvhdl=(int**)omAlloc0(n*sizeof(int_ptr));
2514    // init order, so that rBlocks works correctly
2515    for (j_in_R= n-2; j_in_R>=0; j_in_R--)
2516      R->order[j_in_R] = (int) ringorder_unspec;
2517    // orderings
2518    for(j_in_R=0,j_in_L=0;j_in_R<n-1;j_in_R++,j_in_L++)
2519    {
2520    // todo: a(..), M
2521      if (v->m[j_in_L].Typ()!=LIST_CMD)
2522      {
2523        WerrorS("ordering must be list of lists");
2524        return TRUE;
2525      }
2526      lists vv=(lists)v->m[j_in_L].Data();
2527      if ((vv->nr==1)
2528      && (vv->m[0].Typ()==STRING_CMD))
2529      {
2530        if (strcmp((char*)vv->m[0].Data(),"L")==0)
2531        {
2532          j_in_R--;
2533          continue;
2534        }
2535        if ((vv->m[1].Typ()!=INTVEC_CMD) && (vv->m[1].Typ()!=INT_CMD))
2536        {
2537          PrintS(lString(vv));
2538          WerrorS("ordering name must be a (string,intvec)(1)");
2539          return TRUE;
2540        }
2541        R->order[j_in_R]=rOrderName(omStrDup((char*)vv->m[0].Data())); // assume STRING
2542
2543        if (j_in_R==0) R->block0[0]=1;
2544        else
2545        {
2546           int jj=j_in_R-1;
2547           while((jj>=0)
2548           && ((R->order[jj]== ringorder_a)
2549              || (R->order[jj]== ringorder_aa)
2550              || (R->order[jj]== ringorder_am)
2551              || (R->order[jj]== ringorder_c)
2552              || (R->order[jj]== ringorder_C)
2553              || (R->order[jj]== ringorder_s)
2554              || (R->order[jj]== ringorder_S)
2555           ))
2556           {
2557             //Print("jj=%, skip %s\n",rSimpleOrdStr(R->order[jj]));
2558             jj--;
2559           }
2560           if (jj<0) R->block0[j_in_R]=1;
2561           else       R->block0[j_in_R]=R->block1[jj]+1;
2562        }
2563        intvec *iv;
2564        if (vv->m[1].Typ()==INT_CMD)
2565          iv=new intvec((int)(long)vv->m[1].Data(),(int)(long)vv->m[1].Data());
2566        else
2567          iv=ivCopy((intvec*)vv->m[1].Data()); //assume INTVEC
2568        int iv_len=iv->length();
2569        R->block1[j_in_R]=si_max(R->block0[j_in_R],R->block0[j_in_R]+iv_len-1);
2570        if (R->block1[j_in_R]>R->N)
2571        {
2572          R->block1[j_in_R]=R->N;
2573          iv_len=R->block1[j_in_R]-R->block0[j_in_R]+1;
2574        }
2575        //Print("block %d from %d to %d\n",j,R->block0[j], R->block1[j]);
2576        int i;
2577        switch (R->order[j_in_R])
2578        {
2579           case ringorder_ws:
2580           case ringorder_Ws:
2581              R->OrdSgn=-1;
2582           case ringorder_aa:
2583           case ringorder_a:
2584           case ringorder_wp:
2585           case ringorder_Wp:
2586             R->wvhdl[j_in_R] =( int *)omAlloc(iv_len*sizeof(int));
2587             for (i=0; i<iv_len;i++)
2588             {
2589               R->wvhdl[j_in_R][i]=(*iv)[i];
2590             }
2591             break;
2592           case ringorder_am:
2593             R->wvhdl[j_in_R] =( int *)omAlloc((iv->length()+1)*sizeof(int));
2594             for (i=0; i<iv_len;i++)
2595             {
2596               R->wvhdl[j_in_R][i]=(*iv)[i];
2597             }
2598             R->wvhdl[j_in_R][i]=iv->length() - iv_len;
2599             //printf("ivlen:%d,iv->len:%d,mod:%d\n",iv_len,iv->length(),R->wvhdl[j][i]);
2600             for (; i<iv->length(); i++)
2601             {
2602                R->wvhdl[j_in_R][i+1]=(*iv)[i];
2603             }
2604             break;
2605           case ringorder_M:
2606             R->wvhdl[j_in_R] =( int *)omAlloc((iv->length())*sizeof(int));
2607             for (i=0; i<iv->length();i++) R->wvhdl[j_in_R][i]=(*iv)[i];
2608             R->block1[j_in_R]=si_max(R->block0[j_in_R],R->block0[j_in_R]+(int)sqrt((double)(iv->length()-1)));
2609             if (R->block1[j_in_R]>R->N)
2610             {
2611               WerrorS("ordering matrix too big");
2612               return TRUE;
2613             }
2614             break;
2615           case ringorder_ls:
2616           case ringorder_ds:
2617           case ringorder_Ds:
2618           case ringorder_rs:
2619             R->OrdSgn=-1;
2620           case ringorder_lp:
2621           case ringorder_dp:
2622           case ringorder_Dp:
2623           case ringorder_rp:
2624             break;
2625           case ringorder_S:
2626             break;
2627           case ringorder_c:
2628           case ringorder_C:
2629             R->block1[j_in_R]=R->block0[j_in_R]=0;
2630             break;
2631
2632           case ringorder_s:
2633             break;
2634
2635           case ringorder_IS:
2636           {
2637             R->block1[j_in_R] = R->block0[j_in_R] = 0;
2638             if( iv->length() > 0 )
2639             {
2640               const int s = (*iv)[0];
2641               assume( -2 < s && s < 2 );
2642               R->block1[j_in_R] = R->block0[j_in_R] = s;
2643             }
2644             break;
2645           }
2646           case 0:
2647           case ringorder_unspec:
2648             break;
2649        }
2650        delete iv;
2651      }
2652      else
2653      {
2654        PrintS(lString(vv));
2655        WerrorS("ordering name must be a (string,intvec)");
2656        return TRUE;
2657      }
2658    }
2659    // sanity check
2660    j_in_R=n-2;
2661    if ((R->order[j_in_R]==ringorder_c)
2662    || (R->order[j_in_R]==ringorder_C)
2663    || (R->order[j_in_R]==ringorder_unspec)) j_in_R--;
2664    if (R->block1[j_in_R] != R->N)
2665    {
2666      if (((R->order[j_in_R]==ringorder_dp) ||
2667           (R->order[j_in_R]==ringorder_ds) ||
2668           (R->order[j_in_R]==ringorder_Dp) ||
2669           (R->order[j_in_R]==ringorder_Ds) ||
2670           (R->order[j_in_R]==ringorder_rp) ||
2671           (R->order[j_in_R]==ringorder_rs) ||
2672           (R->order[j_in_R]==ringorder_lp) ||
2673           (R->order[j_in_R]==ringorder_ls))
2674          &&
2675            R->block0[j_in_R] <= R->N)
2676      {
2677        R->block1[j_in_R] = R->N;
2678      }
2679      else
2680      {
2681        Werror("ordering incomplete: size (%d) should be %d",R->block1[j_in_R],R->N);
2682        return TRUE;
2683      }
2684    }
2685    if (R->block0[j_in_R]>R->N)
2686    {
2687      Werror("not enough variables (%d) for ordering block %d, scanned so far:",R->N,j_in_R+1);
2688      for(int ii=0;ii<=j_in_R;ii++)
2689        Werror("ord[%d]: %s from v%d to v%d",ii+1,rSimpleOrdStr(R->order[ii]),R->block0[ii],R->block1[ii]);
2690      return TRUE;
2691    }
2692    if (check_comp)
2693    {
2694      BOOLEAN comp_order=FALSE;
2695      int jj;
2696      for(jj=0;jj<n;jj++)
2697      {
2698        if ((R->order[jj]==ringorder_c) ||
2699            (R->order[jj]==ringorder_C)) { comp_order=TRUE; break; }
2700      }
2701      if (!comp_order)
2702      {
2703        R->order=(int*)omRealloc0Size(R->order,n*sizeof(int),(n+1)*sizeof(int));
2704        R->block0=(int*)omRealloc0Size(R->block0,n*sizeof(int),(n+1)*sizeof(int));
2705        R->block1=(int*)omRealloc0Size(R->block1,n*sizeof(int),(n+1)*sizeof(int));
2706        R->wvhdl=(int**)omRealloc0Size(R->wvhdl,n*sizeof(int_ptr),(n+1)*sizeof(int_ptr));
2707        R->order[n-1]=ringorder_C;
2708        R->block0[n-1]=0;
2709        R->block1[n-1]=0;
2710        R->wvhdl[n-1]=NULL;
2711        n++;
2712      }
2713    }
2714  }
2715  else
2716  {
2717    WerrorS("ordering must be given as `list`");
2718    return TRUE;
2719  }
2720  if (bitmask!=0) R->bitmask=bitmask*2;
2721  return FALSE;
2722}
2723
2724ring rCompose(const lists  L, const BOOLEAN check_comp)
2725{
2726  if ((L->nr!=3)
2727#ifdef HAVE_PLURAL
2728  &&(L->nr!=5)
2729#endif
2730  )
2731    return NULL;
2732  int is_gf_char=0;
2733  // 0: char/ cf - ring
2734  // 1: list (var)
2735  // 2: list (ord)
2736  // 3: qideal
2737  // possibly:
2738  // 4: C
2739  // 5: D
2740
2741  ring R = (ring) omAlloc0Bin(sip_sring_bin);
2742
2743  // ------------------------------------------------------------------
2744  // 0: char:
2745#ifdef SINGULAR_4_1
2746  if (L->m[0].Typ()==CRING_CMD)
2747  {
2748    R->cf=(coeffs)L->m[0].Data();
2749    R->cf->ref++;
2750  }
2751  else
2752#endif
2753  if (L->m[0].Typ()==INT_CMD)
2754  {
2755    int ch = (int)(long)L->m[0].Data();
2756    assume( ch >= 0 );
2757
2758    if (ch == 0) // Q?
2759      R->cf = nInitChar(n_Q, NULL);
2760    else
2761    {
2762      int l = IsPrime(ch); // Zp?
2763      if( l != ch )
2764      {
2765        Warn("%d is invalid characteristic of ground field. %d is used.", ch, l);
2766        ch = l;
2767      }
2768      R->cf = nInitChar(n_Zp, (void*)(long)ch);
2769    }
2770  }
2771  else if (L->m[0].Typ()==LIST_CMD) // something complicated...
2772  {
2773    lists LL=(lists)L->m[0].Data();
2774
2775#ifdef HAVE_RINGS
2776    if (LL->m[0].Typ() == STRING_CMD) // 1st comes a string?
2777    {
2778      rComposeRing(LL, R); // Ring!?
2779    }
2780    else
2781#endif
2782    if (LL->nr < 3)
2783      rComposeC(LL,R); // R, long_R, long_C
2784    else
2785    {
2786      if (LL->m[0].Typ()==INT_CMD)
2787      {
2788        int ch = (int)(long)LL->m[0].Data();
2789        while ((ch!=fftable[is_gf_char]) && (fftable[is_gf_char])) is_gf_char++;
2790        if (fftable[is_gf_char]==0) is_gf_char=-1;
2791
2792        if(is_gf_char!= -1)
2793        {
2794          GFInfo param;
2795
2796          param.GFChar = ch;
2797          param.GFDegree = 1;
2798          param.GFPar_name = (const char*)(((lists)(LL->m[1].Data()))->m[0].Data());
2799
2800          // nfInitChar should be able to handle the case when ch is in fftables!
2801          R->cf = nInitChar(n_GF, (void*)&param);
2802        }
2803      }
2804
2805      if( R->cf == NULL )
2806      {
2807        ring extRing = rCompose((lists)L->m[0].Data(),FALSE);
2808
2809        if (extRing==NULL)
2810        {
2811          WerrorS("could not create the specified coefficient field");
2812          goto rCompose_err;
2813        }
2814
2815        if( extRing->qideal != NULL ) // Algebraic extension
2816        {
2817          AlgExtInfo extParam;
2818
2819          extParam.r = extRing;
2820
2821          R->cf = nInitChar(n_algExt, (void*)&extParam);
2822        }
2823        else // Transcendental extension
2824        {
2825          TransExtInfo extParam;
2826          extParam.r = extRing;
2827          assume( extRing->qideal == NULL );
2828
2829          R->cf = nInitChar(n_transExt, &extParam);
2830        }
2831      }
2832    }
2833  }
2834  else
2835  {
2836    WerrorS("coefficient field must be described by `int` or `list`");
2837    goto rCompose_err;
2838  }
2839
2840  if( R->cf == NULL )
2841  {
2842    WerrorS("could not create coefficient field described by the input!");
2843    goto rCompose_err;
2844  }
2845
2846  // ------------------------- VARS ---------------------------
2847  if (rComposeVar(L,R)) goto rCompose_err;
2848  // ------------------------ ORDER ------------------------------
2849  if (rComposeOrder(L,check_comp,R)) goto rCompose_err;
2850
2851  // ------------------------ ??????? --------------------
2852
2853  rRenameVars(R);
2854  rComplete(R);
2855
2856  // ------------------------ Q-IDEAL ------------------------
2857
2858  if (L->m[3].Typ()==IDEAL_CMD)
2859  {
2860    ideal q=(ideal)L->m[3].Data();
2861    if (q->m[0]!=NULL)
2862    {
2863      if (R->cf != currRing->cf) //->cf->ch!=currRing->cf->ch)
2864      {
2865      #if 0
2866            WerrorS("coefficient fields must be equal if q-ideal !=0");
2867            goto rCompose_err;
2868      #else
2869        ring orig_ring=currRing;
2870        rChangeCurrRing(R);
2871        int *perm=NULL;
2872        int *par_perm=NULL;
2873        int par_perm_size=0;
2874        nMapFunc nMap;
2875
2876        if ((nMap=nSetMap(orig_ring->cf))==NULL)
2877        {
2878          if (rEqual(orig_ring,currRing))
2879          {
2880            nMap=n_SetMap(currRing->cf, currRing->cf);
2881          }
2882          else
2883          // Allow imap/fetch to be make an exception only for:
2884          if ( (rField_is_Q_a(orig_ring) &&  // Q(a..) -> Q(a..) || Q || Zp || Zp(a)
2885            (rField_is_Q(currRing) || rField_is_Q_a(currRing) ||
2886             (rField_is_Zp(currRing) || rField_is_Zp_a(currRing))))
2887           ||
2888           (rField_is_Zp_a(orig_ring) &&  // Zp(a..) -> Zp(a..) || Zp
2889            (rField_is_Zp(currRing, rInternalChar(orig_ring)) ||
2890             rField_is_Zp_a(currRing, rInternalChar(orig_ring)))) )
2891          {
2892            par_perm_size=rPar(orig_ring);
2893
2894//            if ((orig_ring->minpoly != NULL) || (orig_ring->qideal != NULL))
2895//              naSetChar(rInternalChar(orig_ring),orig_ring);
2896//            else ntSetChar(rInternalChar(orig_ring),orig_ring);
2897
2898            nSetChar(currRing->cf);
2899          }
2900          else
2901          {
2902            WerrorS("coefficient fields must be equal if q-ideal !=0");
2903            goto rCompose_err;
2904          }
2905        }
2906        perm=(int *)omAlloc0((orig_ring->N+1)*sizeof(int));
2907        if (par_perm_size!=0)
2908          par_perm=(int *)omAlloc0(par_perm_size*sizeof(int));
2909        int i;
2910        #if 0
2911        // use imap:
2912        maFindPerm(orig_ring->names,orig_ring->N,orig_ring->parameter,orig_ring->P,
2913          currRing->names,currRing->N,currRing->parameter, currRing->P,
2914          perm,par_perm, currRing->ch);
2915        #else
2916        // use fetch
2917        if ((rPar(orig_ring)>0) && (rPar(currRing)==0))
2918        {
2919          for(i=si_min(rPar(orig_ring),rVar(currRing))-1;i>=0;i--) par_perm[i]=i+1;
2920        }
2921        else if (par_perm_size!=0)
2922          for(i=si_min(rPar(orig_ring),rPar(currRing))-1;i>=0;i--) par_perm[i]=-(i+1);
2923        for(i=si_min(orig_ring->N,rVar(currRing));i>0;i--) perm[i]=i;
2924        #endif
2925        ideal dest_id=idInit(IDELEMS(q),1);
2926        for(i=IDELEMS(q)-1; i>=0; i--)
2927        {
2928          dest_id->m[i]=p_PermPoly(q->m[i],perm,orig_ring, currRing,nMap,
2929                                  par_perm,par_perm_size);
2930          //  PrintS("map:");pWrite(dest_id->m[i]);PrintLn();
2931          pTest(dest_id->m[i]);
2932        }
2933        R->qideal=dest_id;
2934        if (perm!=NULL)
2935          omFreeSize((ADDRESS)perm,(orig_ring->N+1)*sizeof(int));
2936        if (par_perm!=NULL)
2937          omFreeSize((ADDRESS)par_perm,par_perm_size*sizeof(int));
2938        rChangeCurrRing(orig_ring);
2939      #endif
2940      }
2941      else
2942        R->qideal=idrCopyR(q,currRing,R);
2943    }
2944  }
2945  else
2946  {
2947    WerrorS("q-ideal must be given as `ideal`");
2948    goto rCompose_err;
2949  }
2950
2951
2952  // ---------------------------------------------------------------
2953  #ifdef HAVE_PLURAL
2954  if (L->nr==5)
2955  {
2956    if (nc_CallPlural((matrix)L->m[4].Data(),
2957                      (matrix)L->m[5].Data(),
2958                      NULL,NULL,
2959                      R,
2960                      true, // !!!
2961                      true, false,
2962                      currRing, FALSE)) goto rCompose_err;
2963    // takes care about non-comm. quotient! i.e. calls "nc_SetupQuotient" due to last true
2964  }
2965  #endif
2966  return R;
2967
2968rCompose_err:
2969  if (R->N>0)
2970  {
2971    int i;
2972    if (R->names!=NULL)
2973    {
2974      i=R->N-1;
2975      while (i>=0) { if (R->names[i]!=NULL) omFree(R->names[i]); i--; }
2976      omFree(R->names);
2977    }
2978  }
2979  if (R->order!=NULL) omFree(R->order);
2980  if (R->block0!=NULL) omFree(R->block0);
2981  if (R->block1!=NULL) omFree(R->block1);
2982  if (R->wvhdl!=NULL) omFree(R->wvhdl);
2983  omFree(R);
2984  return NULL;
2985}
2986
2987// from matpol.cc
2988
2989/*2
2990* compute the jacobi matrix of an ideal
2991*/
2992BOOLEAN mpJacobi(leftv res,leftv a)
2993{
2994  int     i,j;
2995  matrix result;
2996  ideal id=(ideal)a->Data();
2997
2998  result =mpNew(IDELEMS(id),rVar(currRing));
2999  for (i=1; i<=IDELEMS(id); i++)
3000  {
3001    for (j=1; j<=rVar(currRing); j++)
3002    {
3003      MATELEM(result,i,j) = pDiff(id->m[i-1],j);
3004    }
3005  }
3006  res->data=(char *)result;
3007  return FALSE;
3008}
3009
3010/*2
3011* returns the Koszul-matrix of degree d of a vectorspace with dimension n
3012* uses the first n entrees of id, if id <> NULL
3013*/
3014BOOLEAN mpKoszul(leftv res,leftv c/*ip*/, leftv b/*in*/, leftv id)
3015{
3016  int n=(int)(long)b->Data();
3017  int d=(int)(long)c->Data();
3018  int     k,l,sign,row,col;
3019  matrix  result;
3020  ideal temp;
3021  BOOLEAN bo;
3022  poly    p;
3023
3024  if ((d>n) || (d<1) || (n<1))
3025  {
3026    res->data=(char *)mpNew(1,1);
3027    return FALSE;
3028  }
3029  int *choise = (int*)omAlloc(d*sizeof(int));
3030  if (id==NULL)
3031    temp=idMaxIdeal(1);
3032  else
3033    temp=(ideal)id->Data();
3034
3035  k = binom(n,d);
3036  l = k*d;
3037  l /= n-d+1;
3038  result =mpNew(l,k);
3039  col = 1;
3040  idInitChoise(d,1,n,&bo,choise);
3041  while (!bo)
3042  {
3043    sign = 1;
3044    for (l=1;l<=d;l++)
3045    {
3046      if (choise[l-1]<=IDELEMS(temp))
3047      {
3048        p = pCopy(temp->m[choise[l-1]-1]);
3049        if (sign == -1) p = pNeg(p);
3050        sign *= -1;
3051        row = idGetNumberOfChoise(l-1,d,1,n,choise);
3052        MATELEM(result,row,col) = p;
3053      }
3054    }
3055    col++;
3056    idGetNextChoise(d,n,&bo,choise);
3057  }
3058  if (id==NULL) idDelete(&temp);
3059
3060  res->data=(char *)result;
3061  return FALSE;
3062}
3063
3064// from syz1.cc
3065/*2
3066* read out the Betti numbers from resolution
3067* (interpreter interface)
3068*/
3069BOOLEAN syBetti2(leftv res, leftv u, leftv w)
3070{
3071  syStrategy syzstr=(syStrategy)u->Data();
3072
3073  BOOLEAN minim=(int)(long)w->Data();
3074  int row_shift=0;
3075  int add_row_shift=0;
3076  intvec *weights=NULL;
3077  intvec *ww=(intvec *)atGet(u,"isHomog",INTVEC_CMD);
3078  if (ww!=NULL)
3079  {
3080     weights=ivCopy(ww);
3081     add_row_shift = ww->min_in();
3082     (*weights) -= add_row_shift;
3083  }
3084
3085  res->data=(void *)syBettiOfComputation(syzstr,minim,&row_shift,weights);
3086  //row_shift += add_row_shift;
3087  //Print("row_shift=%d, add_row_shift=%d\n",row_shift,add_row_shift);
3088  atSet(res,omStrDup("rowShift"),(void*)(long)add_row_shift,INT_CMD);
3089
3090  return FALSE;
3091}
3092BOOLEAN syBetti1(leftv res, leftv u)
3093{
3094  sleftv tmp;
3095  memset(&tmp,0,sizeof(tmp));
3096  tmp.rtyp=INT_CMD;
3097  tmp.data=(void *)1;
3098  return syBetti2(res,u,&tmp);
3099}
3100
3101/*3
3102* converts a resolution into a list of modules
3103*/
3104lists syConvRes(syStrategy syzstr,BOOLEAN toDel,int add_row_shift)
3105{
3106  resolvente fullres = syzstr->fullres;
3107  resolvente minres = syzstr->minres;
3108
3109  const int length = syzstr->length;
3110
3111  if ((fullres==NULL) && (minres==NULL))
3112  {
3113    if (syzstr->hilb_coeffs==NULL)
3114    { // La Scala
3115      fullres = syReorder(syzstr->res, length, syzstr);
3116    }
3117    else
3118    { // HRES
3119      minres = syReorder(syzstr->orderedRes, length, syzstr);
3120      syKillEmptyEntres(minres, length);
3121    }
3122  }
3123
3124  resolvente tr;
3125  int typ0=IDEAL_CMD;
3126
3127  if (minres!=NULL)
3128    tr = minres;
3129  else
3130    tr = fullres;
3131
3132  resolvente trueres=NULL; intvec ** w=NULL;
3133
3134  if (length>0)
3135  {
3136    trueres = (resolvente)omAlloc0((length)*sizeof(ideal));
3137    for (int i=(length)-1;i>=0;i--)
3138    {
3139      if (tr[i]!=NULL)
3140      {
3141        trueres[i] = idCopy(tr[i]);
3142      }
3143    }
3144    if ( id_RankFreeModule(trueres[0], currRing) > 0)
3145      typ0 = MODUL_CMD;
3146    if (syzstr->weights!=NULL)
3147    {
3148      w = (intvec**)omAlloc0(length*sizeof(intvec*));
3149      for (int i=length-1;i>=0;i--)
3150      {
3151        if (syzstr->weights[i]!=NULL) w[i] = ivCopy(syzstr->weights[i]);
3152      }
3153    }
3154  }
3155
3156  lists li = liMakeResolv(trueres, length, syzstr->list_length,typ0,
3157                          w, add_row_shift);
3158
3159  if (w != NULL) omFreeSize(w, length*sizeof(intvec*));
3160
3161  if (toDel)
3162    syKillComputation(syzstr);
3163  else
3164  {
3165    if( fullres != NULL && syzstr->fullres == NULL )
3166      syzstr->fullres = fullres;
3167
3168    if( minres != NULL && syzstr->minres == NULL )
3169      syzstr->minres = minres;
3170  }
3171  return li;
3172}
3173
3174/*3
3175* converts a list of modules into a resolution
3176*/
3177syStrategy syConvList(lists li)
3178{
3179  int typ0;
3180  syStrategy result=(syStrategy)omAlloc0(sizeof(ssyStrategy));
3181
3182  resolvente fr = liFindRes(li,&(result->length),&typ0,&(result->weights));
3183  if (fr != NULL)
3184  {
3185
3186    result->fullres = (resolvente)omAlloc0((result->length+1)*sizeof(ideal));
3187    for (int i=result->length-1;i>=0;i--)
3188    {
3189      if (fr[i]!=NULL)
3190        result->fullres[i] = idCopy(fr[i]);
3191    }
3192    result->list_length=result->length;
3193    omFreeSize((ADDRESS)fr,(result->length)*sizeof(ideal));
3194  }
3195  else
3196  {
3197    omFreeSize(result, sizeof(ssyStrategy));
3198    result = NULL;
3199  }
3200  return result;
3201}
3202
3203/*3
3204* converts a list of modules into a minimal resolution
3205*/
3206syStrategy syForceMin(lists li)
3207{
3208  int typ0;
3209  syStrategy result=(syStrategy)omAlloc0(sizeof(ssyStrategy));
3210
3211  resolvente fr = liFindRes(li,&(result->length),&typ0);
3212  result->minres = (resolvente)omAlloc0((result->length+1)*sizeof(ideal));
3213  for (int i=result->length-1;i>=0;i--)
3214  {
3215    if (fr[i]!=NULL)
3216      result->minres[i] = idCopy(fr[i]);
3217  }
3218  omFreeSize((ADDRESS)fr,(result->length)*sizeof(ideal));
3219  return result;
3220}
3221// from weight.cc
3222BOOLEAN kWeight(leftv res,leftv id)
3223{
3224  ideal F=(ideal)id->Data();
3225  intvec * iv = new intvec(rVar(currRing));
3226  polyset s;
3227  int  sl, n, i;
3228  int  *x;
3229
3230  res->data=(char *)iv;
3231  s = F->m;
3232  sl = IDELEMS(F) - 1;
3233  n = rVar(currRing);
3234  double wNsqr = (double)2.0 / (double)n;
3235  wFunctional = wFunctionalBuch;
3236  x = (int * )omAlloc(2 * (n + 1) * sizeof(int));
3237  wCall(s, sl, x, wNsqr, currRing);
3238  for (i = n; i!=0; i--)
3239    (*iv)[i-1] = x[i + n + 1];
3240  omFreeSize((ADDRESS)x, 2 * (n + 1) * sizeof(int));
3241  return FALSE;
3242}
3243
3244BOOLEAN kQHWeight(leftv res,leftv v)
3245{
3246  res->data=(char *)id_QHomWeight((ideal)v->Data(), currRing);
3247  if (res->data==NULL)
3248    res->data=(char *)new intvec(rVar(currRing));
3249  return FALSE;
3250}
3251/*==============================================================*/
3252// from clapsing.cc
3253#if 0
3254BOOLEAN jjIS_SQR_FREE(leftv res, leftv u)
3255{
3256  BOOLEAN b=singclap_factorize((poly)(u->CopyD()), &v, 0);
3257  res->data=(void *)b;
3258}
3259#endif
3260
3261BOOLEAN jjRESULTANT(leftv res, leftv u, leftv v, leftv w)
3262{
3263  res->data=singclap_resultant((poly)u->CopyD(),(poly)v->CopyD(),
3264                  (poly)w->CopyD(), currRing);
3265  return errorreported;
3266}
3267
3268BOOLEAN jjCHARSERIES(leftv res, leftv u)
3269{
3270  res->data=singclap_irrCharSeries((ideal)u->Data(), currRing);
3271  return (res->data==NULL);
3272}
3273
3274// from semic.cc
3275#ifdef HAVE_SPECTRUM
3276
3277// ----------------------------------------------------------------------------
3278//  Initialize a  spectrum  deep from a  singular  lists
3279// ----------------------------------------------------------------------------
3280
3281void copy_deep( spectrum& spec, lists l )
3282{
3283    spec.mu = (int)(long)(l->m[0].Data( ));
3284    spec.pg = (int)(long)(l->m[1].Data( ));
3285    spec.= (int)(long)(l->m[2].Data( ));
3286
3287    spec.copy_new( spec.n );
3288
3289    intvec  *num = (intvec*)l->m[3].Data( );
3290    intvec  *den = (intvec*)l->m[4].Data( );
3291    intvec  *mul = (intvec*)l->m[5].Data( );
3292
3293    for( int i=0; i<spec.n; i++ )
3294    {
3295        spec.s[i] = (Rational)((*num)[i])/(Rational)((*den)[i]);
3296        spec.w[i] = (*mul)[i];
3297    }
3298}
3299
3300// ----------------------------------------------------------------------------
3301//  singular lists  constructor for  spectrum
3302// ----------------------------------------------------------------------------
3303
3304spectrum /*former spectrum::spectrum ( lists l )*/
3305spectrumFromList( lists l )
3306{
3307    spectrum result;
3308    copy_deep( result, l );
3309    return result;
3310}
3311
3312// ----------------------------------------------------------------------------
3313//  generate a Singular  lists  from a spectrum
3314// ----------------------------------------------------------------------------
3315
3316/* former spectrum::thelist ( void )*/
3317lists   getList( spectrum& spec )
3318{
3319    lists   L  = (lists)omAllocBin( slists_bin);
3320
3321    L->Init( 6 );
3322
3323    intvec            *num  = new intvec( spec.n );
3324    intvec            *den  = new intvec( spec.n );
3325    intvec            *mult = new intvec( spec.n );
3326
3327    for( int i=0; i<spec.n; i++ )
3328    {
3329        (*num) [i] = spec.s[i].get_num_si( );
3330        (*den) [i] = spec.s[i].get_den_si( );
3331        (*mult)[i] = spec.w[i];
3332    }
3333
3334    L->m[0].rtyp = INT_CMD;    //  milnor number
3335    L->m[1].rtyp = INT_CMD;    //  geometrical genus
3336    L->m[2].rtyp = INT_CMD;    //  # of spectrum numbers
3337    L->m[3].rtyp = INTVEC_CMD; //  numerators
3338    L->m[4].rtyp = INTVEC_CMD; //  denomiantors
3339    L->m[5].rtyp = INTVEC_CMD; //  multiplicities
3340
3341    L->m[0].data = (void*)(long)spec.mu;
3342    L->m[1].data = (void*)(long)spec.pg;
3343    L->m[2].data = (void*)(long)spec.n;
3344    L->m[3].data = (void*)num;
3345    L->m[4].data = (void*)den;
3346    L->m[5].data = (void*)mult;
3347
3348    return  L;
3349}
3350// from spectrum.cc
3351// ----------------------------------------------------------------------------
3352//  print out an error message for a spectrum list
3353// ----------------------------------------------------------------------------
3354
3355typedef enum
3356{
3357    semicOK,
3358    semicMulNegative,
3359
3360    semicListTooShort,
3361    semicListTooLong,
3362
3363    semicListFirstElementWrongType,
3364    semicListSecondElementWrongType,
3365    semicListThirdElementWrongType,
3366    semicListFourthElementWrongType,
3367    semicListFifthElementWrongType,
3368    semicListSixthElementWrongType,
3369
3370    semicListNNegative,
3371    semicListWrongNumberOfNumerators,
3372    semicListWrongNumberOfDenominators,
3373    semicListWrongNumberOfMultiplicities,
3374
3375    semicListMuNegative,
3376    semicListPgNegative,
3377    semicListNumNegative,
3378    semicListDenNegative,
3379    semicListMulNegative,
3380
3381    semicListNotSymmetric,
3382    semicListNotMonotonous,
3383
3384    semicListMilnorWrong,
3385    semicListPGWrong
3386
3387} semicState;
3388
3389void    list_error( semicState state )
3390{
3391    switch( state )
3392    {
3393        case semicListTooShort:
3394            WerrorS( "the list is too short" );
3395            break;
3396        case semicListTooLong:
3397            WerrorS( "the list is too long" );
3398            break;
3399
3400        case semicListFirstElementWrongType:
3401            WerrorS( "first element of the list should be int" );
3402            break;
3403        case semicListSecondElementWrongType:
3404            WerrorS( "second element of the list should be int" );
3405            break;
3406        case semicListThirdElementWrongType:
3407            WerrorS( "third element of the list should be int" );
3408            break;
3409        case semicListFourthElementWrongType:
3410            WerrorS( "fourth element of the list should be intvec" );
3411            break;
3412        case semicListFifthElementWrongType:
3413            WerrorS( "fifth element of the list should be intvec" );
3414            break;
3415        case semicListSixthElementWrongType:
3416            WerrorS( "sixth element of the list should be intvec" );
3417            break;
3418
3419        case semicListNNegative:
3420            WerrorS( "first element of the list should be positive" );
3421            break;
3422        case semicListWrongNumberOfNumerators:
3423            WerrorS( "wrong number of numerators" );
3424            break;
3425        case semicListWrongNumberOfDenominators:
3426            WerrorS( "wrong number of denominators" );
3427            break;
3428        case semicListWrongNumberOfMultiplicities:
3429            WerrorS( "wrong number of multiplicities" );
3430            break;
3431
3432        case semicListMuNegative:
3433            WerrorS( "the Milnor number should be positive" );
3434            break;
3435        case semicListPgNegative:
3436            WerrorS( "the geometrical genus should be nonnegative" );
3437            break;
3438        case semicListNumNegative:
3439            WerrorS( "all numerators should be positive" );
3440            break;
3441        case semicListDenNegative:
3442            WerrorS( "all denominators should be positive" );
3443            break;
3444        case semicListMulNegative:
3445            WerrorS( "all multiplicities should be positive" );
3446            break;
3447
3448        case semicListNotSymmetric:
3449            WerrorS( "it is not symmetric" );
3450            break;
3451        case semicListNotMonotonous:
3452            WerrorS( "it is not monotonous" );
3453            break;
3454
3455        case semicListMilnorWrong:
3456            WerrorS( "the Milnor number is wrong" );
3457            break;
3458        case semicListPGWrong:
3459            WerrorS( "the geometrical genus is wrong" );
3460            break;
3461
3462        default:
3463            WerrorS( "unspecific error" );
3464            break;
3465    }
3466}
3467// ----------------------------------------------------------------------------
3468//  this is the main spectrum computation function
3469// ----------------------------------------------------------------------------
3470
3471enum    spectrumState
3472{
3473    spectrumOK,
3474    spectrumZero,
3475    spectrumBadPoly,
3476    spectrumNoSingularity,
3477    spectrumNotIsolated,
3478    spectrumDegenerate,
3479    spectrumWrongRing,
3480    spectrumNoHC,
3481    spectrumUnspecErr
3482};
3483
3484// from splist.cc
3485// ----------------------------------------------------------------------------
3486//  Compute the spectrum of a  spectrumPolyList
3487// ----------------------------------------------------------------------------
3488
3489/* former spectrumPolyList::spectrum ( lists*, int) */
3490spectrumState   spectrumStateFromList( spectrumPolyList& speclist, lists *L,int fast )
3491{
3492  spectrumPolyNode  **node = &speclist.root;
3493  spectrumPolyNode  *search;
3494
3495  poly              f,tmp;
3496  int               found,cmp;
3497
3498  Rational smax( ( fast==0 ? 0 : rVar(currRing) ),
3499                 ( fast==2 ? 2 : 1 ) );
3500
3501  Rational weight_prev( 0,1 );
3502
3503  int     mu = 0;          // the milnor number
3504  int     pg = 0;          // the geometrical genus
3505  int     n  = 0;          // number of different spectral numbers
3506  int     z  = 0;          // number of spectral number equal to smax
3507
3508  while( (*node)!=(spectrumPolyNode*)NULL &&
3509         ( fast==0 || (*node)->weight<=smax ) )
3510  {
3511        // ---------------------------------------
3512        //  determine the first normal form which
3513        //  contains the monomial  node->mon
3514        // ---------------------------------------
3515
3516    found  = FALSE;
3517    search = *node;
3518
3519    while( search!=(spectrumPolyNode*)NULL && found==FALSE )
3520    {
3521      if( search->nf!=(poly)NULL )
3522      {
3523        f = search->nf;
3524
3525        do
3526        {
3527                    // --------------------------------
3528                    //  look for  (*node)->mon  in   f
3529                    // --------------------------------
3530
3531          cmp = pCmp( (*node)->mon,f );
3532
3533          if( cmp<0 )
3534          {
3535            f = pNext( f );
3536          }
3537          else if( cmp==0 )
3538          {
3539                        // -----------------------------
3540                        //  we have found a normal form
3541                        // -----------------------------
3542
3543            found = TRUE;
3544
3545                        //  normalize coefficient
3546
3547            number inv = nInvers( pGetCoeff( f ) );
3548            pMult_nn( search->nf,inv );
3549            nDelete( &inv );
3550
3551                        //  exchange  normal forms
3552
3553            tmp         = (*node)->nf;
3554            (*node)->nf = search->nf;
3555            search->nf  = tmp;
3556          }
3557        }
3558        while( cmp<0 && f!=(poly)NULL );
3559      }
3560      search = search->next;
3561    }
3562
3563    if( found==FALSE )
3564    {
3565            // ------------------------------------------------
3566            //  the weight of  node->mon  is a spectrum number
3567            // ------------------------------------------------
3568
3569      mu++;
3570
3571      if( (*node)->weight<=(Rational)1 )              pg++;
3572      if( (*node)->weight==smax )           z++;
3573      if( (*node)->weight>weight_prev )     n++;
3574
3575      weight_prev = (*node)->weight;
3576      node = &((*node)->next);
3577    }
3578    else
3579    {
3580            // -----------------------------------------------
3581            //  determine all other normal form which contain
3582            //  the monomial  node->mon
3583            //  replace for  node->mon  its normal form
3584            // -----------------------------------------------
3585
3586      while( search!=(spectrumPolyNode*)NULL )
3587      {
3588        if( search->nf!=(poly)NULL )
3589        {
3590          f = search->nf;
3591
3592          do
3593          {
3594                        // --------------------------------
3595                        //  look for  (*node)->mon  in   f
3596                        // --------------------------------
3597
3598            cmp = pCmp( (*node)->mon,f );
3599
3600            if( cmp<0 )
3601            {
3602              f = pNext( f );
3603            }
3604            else if( cmp==0 )
3605            {
3606              search->nf = pSub( search->nf,
3607                                 ppMult_nn( (*node)->nf,pGetCoeff( f ) ) );
3608              pNorm( search->nf );
3609            }
3610          }
3611          while( cmp<0 && f!=(poly)NULL );
3612        }
3613        search = search->next;
3614      }
3615      speclist.delete_node( node );
3616    }
3617
3618  }
3619
3620    // --------------------------------------------------------
3621    //  fast computation exploits the symmetry of the spectrum
3622    // --------------------------------------------------------
3623
3624  if( fast==2 )
3625  {
3626    mu = 2*mu - z;
3627    n  = ( z > 0 ? 2*n - 1 : 2*n );
3628  }
3629
3630    // --------------------------------------------------------
3631    //  compute the spectrum numbers with their multiplicities
3632    // --------------------------------------------------------
3633
3634  intvec            *nom  = new intvec( n );
3635  intvec            *den  = new intvec( n );
3636  intvec            *mult = new intvec( n );
3637
3638  int count         = 0;
3639  int multiplicity  = 1;
3640
3641  for( search=speclist.root; search!=(spectrumPolyNode*)NULL &&
3642              ( fast==0 || search->weight<=smax );
3643       search=search->next )
3644  {
3645    if( search->next==(spectrumPolyNode*)NULL ||
3646        search->weight<search->next->weight )
3647    {
3648      (*nom) [count] = search->weight.get_num_si( );
3649      (*den) [count] = search->weight.get_den_si( );
3650      (*mult)[count] = multiplicity;
3651
3652      multiplicity=1;
3653      count++;
3654    }
3655    else
3656    {
3657      multiplicity++;
3658    }
3659  }
3660
3661    // --------------------------------------------------------
3662    //  fast computation exploits the symmetry of the spectrum
3663    // --------------------------------------------------------
3664
3665  if( fast==2 )
3666  {
3667    int n1,n2;
3668    for( n1=0, n2=n-1; n1<n2; n1++, n2-- )
3669    {
3670      (*nom) [n2] = rVar(currRing)*(*den)[n1]-(*nom)[n1];
3671      (*den) [n2] = (*den)[n1];
3672      (*mult)[n2] = (*mult)[n1];
3673    }
3674  }
3675
3676    // -----------------------------------
3677    //  test if the spectrum is symmetric
3678    // -----------------------------------
3679
3680  if( fast==0 || fast==1 )
3681  {
3682    int symmetric=TRUE;
3683
3684    for( int n1=0, n2=n-1 ; n1<n2 && symmetric==TRUE; n1++, n2-- )
3685    {
3686      if( (*mult)[n1]!=(*mult)[n2] ||
3687          (*den) [n1]!= (*den)[n2] ||
3688          (*nom)[n1]+(*nom)[n2]!=rVar(currRing)*(*den) [n1] )
3689      {
3690        symmetric = FALSE;
3691      }
3692    }
3693
3694    if( symmetric==FALSE )
3695    {
3696            // ---------------------------------------------
3697            //  the spectrum is not symmetric => degenerate
3698            //  principal part
3699            // ---------------------------------------------
3700
3701      *L = (lists)omAllocBin( slists_bin);
3702      (*L)->Init( 1 );
3703      (*L)->m[0].rtyp = INT_CMD;    //  milnor number
3704      (*L)->m[0].data = (void*)(long)mu;
3705
3706      return spectrumDegenerate;
3707    }
3708  }
3709
3710  *L = (lists)omAllocBin( slists_bin);
3711
3712  (*L)->Init( 6 );
3713
3714  (*L)->m[0].rtyp = INT_CMD;    //  milnor number
3715  (*L)->m[1].rtyp = INT_CMD;    //  geometrical genus
3716  (*L)->m[2].rtyp = INT_CMD;    //  number of spectrum values
3717  (*L)->m[3].rtyp = INTVEC_CMD; //  nominators
3718  (*L)->m[4].rtyp = INTVEC_CMD; //  denomiantors
3719  (*L)->m[5].rtyp = INTVEC_CMD; //  multiplicities
3720
3721  (*L)->m[0].data = (void*)(long)mu;
3722  (*L)->m[1].data = (void*)(long)pg;
3723  (*L)->m[2].data = (void*)(long)n;
3724  (*L)->m[3].data = (void*)nom;
3725  (*L)->m[4].data = (void*)den;
3726  (*L)->m[5].data = (void*)mult;
3727
3728  return  spectrumOK;
3729}
3730
3731spectrumState   spectrumCompute( poly h,lists *L,int fast )
3732{
3733  int i;
3734
3735  #ifdef SPECTRUM_DEBUG
3736  #ifdef SPECTRUM_PRINT
3737  #ifdef SPECTRUM_IOSTREAM
3738    cout << "spectrumCompute\n";
3739    if( fast==0 ) cout << "    no optimization" << endl;
3740    if( fast==1 ) cout << "    weight optimization" << endl;
3741    if( fast==2 ) cout << "    symmetry optimization" << endl;
3742  #else
3743    fprintf( stdout,"spectrumCompute\n" );
3744    if( fast==0 ) fprintf( stdout,"    no optimization\n" );
3745    if( fast==1 ) fprintf( stdout,"    weight optimization\n" );
3746    if( fast==2 ) fprintf( stdout,"    symmetry optimization\n" );
3747  #endif
3748  #endif
3749  #endif
3750
3751  // ----------------------
3752  //  check if  h  is zero
3753  // ----------------------
3754
3755  if( h==(poly)NULL )
3756  {
3757    return  spectrumZero;
3758  }
3759
3760  // ----------------------------------
3761  //  check if  h  has a constant term
3762  // ----------------------------------
3763
3764  if( hasConstTerm( h, currRing ) )
3765  {
3766    return  spectrumBadPoly;
3767  }
3768
3769  // --------------------------------
3770  //  check if  h  has a linear term
3771  // --------------------------------
3772
3773  if( hasLinearTerm( h, currRing ) )
3774  {
3775    *L = (lists)omAllocBin( slists_bin);
3776    (*L)->Init( 1 );
3777    (*L)->m[0].rtyp = INT_CMD;    //  milnor number
3778    /* (*L)->m[0].data = (void*)0;a  -- done by Init */
3779
3780    return  spectrumNoSingularity;
3781  }
3782
3783  // ----------------------------------
3784  //  compute the jacobi ideal of  (h)
3785  // ----------------------------------
3786
3787  ideal J = NULL;
3788  J = idInit( rVar(currRing),1 );
3789
3790  #ifdef SPECTRUM_DEBUG
3791  #ifdef SPECTRUM_PRINT
3792  #ifdef SPECTRUM_IOSTREAM
3793    cout << "\n   computing the Jacobi ideal...\n";
3794  #else
3795    fprintf( stdout,"\n   computing the Jacobi ideal...\n" );
3796  #endif
3797  #endif
3798  #endif
3799
3800  for( i=0; i<rVar(currRing); i++ )
3801  {
3802    J->m[i] = pDiff( h,i+1); //j );
3803
3804    #ifdef SPECTRUM_DEBUG
3805    #ifdef SPECTRUM_PRINT
3806    #ifdef SPECTRUM_IOSTREAM
3807      cout << "        ";
3808    #else
3809      fprintf( stdout,"        " );
3810    #endif
3811      pWrite( J->m[i] );
3812    #endif
3813    #endif
3814  }
3815
3816  // --------------------------------------------
3817  //  compute a standard basis  stdJ  of  jac(h)
3818  // --------------------------------------------
3819
3820  #ifdef SPECTRUM_DEBUG
3821  #ifdef SPECTRUM_PRINT
3822  #ifdef SPECTRUM_IOSTREAM
3823    cout << endl;
3824    cout << "    computing a standard basis..." << endl;
3825  #else
3826    fprintf( stdout,"\n" );
3827    fprintf( stdout,"    computing a standard basis...\n" );
3828  #endif
3829  #endif
3830  #endif
3831
3832  ideal stdJ = kStd(J,currRing->qideal,isNotHomog,NULL);
3833  idSkipZeroes( stdJ );
3834
3835  #ifdef SPECTRUM_DEBUG
3836  #ifdef SPECTRUM_PRINT
3837    for( i=0; i<IDELEMS(stdJ); i++ )
3838    {
3839      #ifdef SPECTRUM_IOSTREAM
3840        cout << "        ";
3841      #else
3842        fprintf( stdout,"        " );
3843      #endif
3844
3845      pWrite( stdJ->m[i] );
3846    }
3847  #endif
3848  #endif
3849
3850  idDelete( &J );
3851
3852  // ------------------------------------------
3853  //  check if the  h  has a singularity
3854  // ------------------------------------------
3855
3856  if( hasOne( stdJ, currRing ) )
3857  {
3858    // -------------------------------
3859    //  h is smooth in the origin
3860    //  return only the Milnor number
3861    // -------------------------------
3862
3863    *L = (lists)omAllocBin( slists_bin);
3864    (*L)->Init( 1 );
3865    (*L)->m[0].rtyp = INT_CMD;    //  milnor number
3866    /* (*L)->m[0].data = (void*)0;a  -- done by Init */
3867
3868    return  spectrumNoSingularity;
3869  }
3870
3871  // ------------------------------------------
3872  //  check if the singularity  h  is isolated
3873  // ------------------------------------------
3874
3875  for( i=rVar(currRing); i>0; i-- )
3876  {
3877    if( hasAxis( stdJ,i, currRing )==FALSE )
3878    {
3879      return  spectrumNotIsolated;
3880    }
3881  }
3882
3883  // ------------------------------------------
3884  //  compute the highest corner  hc  of  stdJ
3885  // ------------------------------------------
3886
3887  #ifdef SPECTRUM_DEBUG
3888  #ifdef SPECTRUM_PRINT
3889  #ifdef SPECTRUM_IOSTREAM
3890    cout << "\n    computing the highest corner...\n";
3891  #else
3892    fprintf( stdout,"\n    computing the highest corner...\n" );
3893  #endif
3894  #endif
3895  #endif
3896
3897  poly hc = (poly)NULL;
3898
3899  scComputeHC( stdJ,currRing->qideal, 0,hc );
3900
3901  if( hc!=(poly)NULL )
3902  {
3903    pGetCoeff(hc) = nInit(1);
3904
3905    for( i=rVar(currRing); i>0; i-- )
3906    {
3907      if( pGetExp( hc,i )>0 ) pDecrExp( hc,i );
3908    }
3909    pSetm( hc );
3910  }
3911  else
3912  {
3913    return  spectrumNoHC;
3914  }
3915
3916  #ifdef SPECTRUM_DEBUG
3917  #ifdef SPECTRUM_PRINT
3918  #ifdef SPECTRUM_IOSTREAM
3919    cout << "       ";
3920  #else
3921    fprintf( stdout,"       " );
3922  #endif
3923    pWrite( hc );
3924  #endif
3925  #endif
3926
3927  // ----------------------------------------
3928  //  compute the Newton polygon  nph  of  h
3929  // ----------------------------------------
3930
3931  #ifdef SPECTRUM_DEBUG
3932  #ifdef SPECTRUM_PRINT
3933  #ifdef SPECTRUM_IOSTREAM
3934    cout << "\n    computing the newton polygon...\n";
3935  #else
3936    fprintf( stdout,"\n    computing the newton polygon...\n" );
3937  #endif
3938  #endif
3939  #endif
3940
3941  newtonPolygon nph( h, currRing );
3942
3943  #ifdef SPECTRUM_DEBUG
3944  #ifdef SPECTRUM_PRINT
3945    cout << nph;
3946  #endif
3947  #endif
3948
3949  // -----------------------------------------------
3950  //  compute the weight corner  wc  of  (stdj,nph)
3951  // -----------------------------------------------
3952
3953  #ifdef SPECTRUM_DEBUG
3954  #ifdef SPECTRUM_PRINT
3955  #ifdef SPECTRUM_IOSTREAM
3956    cout << "\n    computing the weight corner...\n";
3957  #else
3958    fprintf( stdout,"\n    computing the weight corner...\n" );
3959  #endif
3960  #endif
3961  #endif
3962
3963  poly    wc = ( fast==0 ? pCopy( hc ) :
3964               ( fast==1 ? computeWC( nph,(Rational)rVar(currRing), currRing ) :
3965              /* fast==2 */computeWC( nph,
3966                      ((Rational)rVar(currRing))/(Rational)2, currRing ) ) );
3967
3968  #ifdef SPECTRUM_DEBUG
3969  #ifdef SPECTRUM_PRINT
3970  #ifdef SPECTRUM_IOSTREAM
3971    cout << "        ";
3972  #else
3973    fprintf( stdout,"        " );
3974  #endif
3975    pWrite( wc );
3976  #endif
3977  #endif
3978
3979  // -------------
3980  //  compute  NF
3981  // -------------
3982
3983  #ifdef SPECTRUM_DEBUG
3984  #ifdef SPECTRUM_PRINT
3985  #ifdef SPECTRUM_IOSTREAM
3986    cout << "\n    computing NF...\n" << endl;
3987  #else
3988    fprintf( stdout,"\n    computing NF...\n" );
3989  #endif
3990  #endif
3991  #endif
3992
3993  spectrumPolyList NF( &nph );
3994
3995  computeNF( stdJ,hc,wc,&NF, currRing );
3996
3997  #ifdef SPECTRUM_DEBUG
3998  #ifdef SPECTRUM_PRINT
3999    cout << NF;
4000  #ifdef SPECTRUM_IOSTREAM
4001    cout << endl;
4002  #else
4003    fprintf( stdout,"\n" );
4004  #endif
4005  #endif
4006  #endif
4007
4008  // ----------------------------
4009  //  compute the spectrum of  h
4010  // ----------------------------
4011//  spectrumState spectrumStateFromList( spectrumPolyList& speclist, lists *L, int fast );
4012
4013  return spectrumStateFromList(NF, L, fast );
4014}
4015
4016// ----------------------------------------------------------------------------
4017//  this procedure is called from the interpreter
4018// ----------------------------------------------------------------------------
4019//  first  = polynomial
4020//  result = list of spectrum numbers
4021// ----------------------------------------------------------------------------
4022
4023void spectrumPrintError(spectrumState state)
4024{
4025  switch( state )
4026  {
4027    case spectrumZero:
4028      WerrorS( "polynomial is zero" );
4029      break;
4030    case spectrumBadPoly:
4031      WerrorS( "polynomial has constant term" );
4032      break;
4033    case spectrumNoSingularity:
4034      WerrorS( "not a singularity" );
4035      break;
4036    case spectrumNotIsolated:
4037      WerrorS( "the singularity is not isolated" );
4038      break;
4039    case spectrumNoHC:
4040      WerrorS( "highest corner cannot be computed" );
4041      break;
4042    case spectrumDegenerate:
4043      WerrorS( "principal part is degenerate" );
4044      break;
4045    case spectrumOK:
4046      break;
4047
4048    default:
4049      WerrorS( "unknown error occurred" );
4050      break;
4051  }
4052}
4053
4054BOOLEAN spectrumProc( leftv result,leftv first )
4055{
4056  spectrumState state = spectrumOK;
4057
4058  // -------------------
4059  //  check consistency
4060  // -------------------
4061
4062  //  check for a local ring
4063
4064  if( !ringIsLocal(currRing ) )
4065  {
4066    WerrorS( "only works for local orderings" );
4067    state = spectrumWrongRing;
4068  }
4069
4070  //  no quotient rings are allowed
4071
4072  else if( currRing->qideal != NULL )
4073  {
4074    WerrorS( "does not work in quotient rings" );
4075    state = spectrumWrongRing;
4076  }
4077  else
4078  {
4079    lists   L    = (lists)NULL;
4080    int     flag = 1; // weight corner optimization is safe
4081
4082    state = spectrumCompute( (poly)first->Data( ),&L,flag );
4083
4084    if( state==spectrumOK )
4085    {
4086      result->rtyp = LIST_CMD;
4087      result->data = (char*)L;
4088    }
4089    else
4090    {
4091      spectrumPrintError(state);
4092    }
4093  }
4094
4095  return  (state!=spectrumOK);
4096}
4097
4098// ----------------------------------------------------------------------------
4099//  this procedure is called from the interpreter
4100// ----------------------------------------------------------------------------
4101//  first  = polynomial
4102//  result = list of spectrum numbers
4103// ----------------------------------------------------------------------------
4104
4105BOOLEAN spectrumfProc( leftv result,leftv first )
4106{
4107  spectrumState state = spectrumOK;
4108
4109  // -------------------
4110  //  check consistency
4111  // -------------------
4112
4113  //  check for a local polynomial ring
4114
4115  if( currRing->OrdSgn != -1 )
4116  // ?? HS: the test above is also true for k[x][[y]], k[[x]][y]
4117  // or should we use:
4118  //if( !ringIsLocal( ) )
4119  {
4120    WerrorS( "only works for local orderings" );
4121    state = spectrumWrongRing;
4122  }
4123  else if( currRing->qideal != NULL )
4124  {
4125    WerrorS( "does not work in quotient rings" );
4126    state = spectrumWrongRing;
4127  }
4128  else
4129  {
4130    lists   L    = (lists)NULL;
4131    int     flag = 2; // symmetric optimization
4132
4133    state = spectrumCompute( (poly)first->Data( ),&L,flag );
4134
4135    if( state==spectrumOK )
4136    {
4137      result->rtyp = LIST_CMD;
4138      result->data = (char*)L;
4139    }
4140    else
4141    {
4142      spectrumPrintError(state);
4143    }
4144  }
4145
4146  return  (state!=spectrumOK);
4147}
4148
4149// ----------------------------------------------------------------------------
4150//  check if a list is a spectrum
4151//  check for:
4152//      list has 6 elements
4153//      1st element is int (mu=Milnor number)
4154//      2nd element is int (pg=geometrical genus)
4155//      3rd element is int (n =number of different spectrum numbers)
4156//      4th element is intvec (num=numerators)
4157//      5th element is intvec (den=denomiantors)
4158//      6th element is intvec (mul=multiplicities)
4159//      exactly n numerators
4160//      exactly n denominators
4161//      exactly n multiplicities
4162//      mu>0
4163//      pg>=0
4164//      n>0
4165//      num>0
4166//      den>0
4167//      mul>0
4168//      symmetriy with respect to numberofvariables/2
4169//      monotony
4170//      mu = sum of all multiplicities
4171//      pg = sum of all multiplicities where num/den<=1
4172// ----------------------------------------------------------------------------
4173
4174semicState  list_is_spectrum( lists l )
4175{
4176    // -------------------
4177    //  check list length
4178    // -------------------
4179
4180    if( l->nr < 5 )
4181    {
4182        return  semicListTooShort;
4183    }
4184    else if( l->nr > 5 )
4185    {
4186        return  semicListTooLong;
4187    }
4188
4189    // -------------
4190    //  check types
4191    // -------------
4192
4193    if( l->m[0].rtyp != INT_CMD )
4194    {
4195        return  semicListFirstElementWrongType;
4196    }
4197    else if( l->m[1].rtyp != INT_CMD )
4198    {
4199        return  semicListSecondElementWrongType;
4200    }
4201    else if( l->m[2].rtyp != INT_CMD )
4202    {
4203        return  semicListThirdElementWrongType;
4204    }
4205    else if( l->m[3].rtyp != INTVEC_CMD )
4206    {
4207        return  semicListFourthElementWrongType;
4208    }
4209    else if( l->m[4].rtyp != INTVEC_CMD )
4210    {
4211        return  semicListFifthElementWrongType;
4212    }
4213    else if( l->m[5].rtyp != INTVEC_CMD )
4214    {
4215        return  semicListSixthElementWrongType;
4216    }
4217
4218    // -------------------------
4219    //  check number of entries
4220    // -------------------------
4221
4222    int     mu = (int)(long)(l->m[0].Data( ));
4223    int     pg = (int)(long)(l->m[1].Data( ));
4224    int     n  = (int)(long)(l->m[2].Data( ));
4225
4226    if( n <= 0 )
4227    {
4228        return  semicListNNegative;
4229    }
4230
4231    intvec  *num = (intvec*)l->m[3].Data( );
4232    intvec  *den = (intvec*)l->m[4].Data( );
4233    intvec  *mul = (intvec*)l->m[5].Data( );
4234
4235    if( n != num->length( ) )
4236    {
4237        return  semicListWrongNumberOfNumerators;
4238    }
4239    else if( n != den->length( ) )
4240    {
4241        return  semicListWrongNumberOfDenominators;
4242    }
4243    else if( n != mul->length( ) )
4244    {
4245        return  semicListWrongNumberOfMultiplicities;
4246    }
4247
4248    // --------
4249    //  values
4250    // --------
4251
4252    if( mu <= 0 )
4253    {
4254        return  semicListMuNegative;
4255    }
4256    if( pg < 0 )
4257    {
4258        return  semicListPgNegative;
4259    }
4260
4261    int i;
4262
4263    for( i=0; i<n; i++ )
4264    {
4265        if( (*num)[i] <= 0 )
4266        {
4267            return  semicListNumNegative;
4268        }
4269        if( (*den)[i] <= 0 )
4270        {
4271            return  semicListDenNegative;
4272        }
4273        if( (*mul)[i] <= 0 )
4274        {
4275            return  semicListMulNegative;
4276        }
4277    }
4278
4279    // ----------------
4280    //  check symmetry
4281    // ----------------
4282
4283    int     j;
4284
4285    for( i=0, j=n-1; i<=j; i++,j-- )
4286    {
4287        if( (*num)[i] != rVar(currRing)*((*den)[i]) - (*num)[j] ||
4288            (*den)[i] != (*den)[j] ||
4289            (*mul)[i] != (*mul)[j] )
4290        {
4291            return  semicListNotSymmetric;
4292        }
4293    }
4294
4295    // ----------------
4296    //  check monotony
4297    // ----------------
4298
4299    for( i=0, j=1; i<n/2; i++,j++ )
4300    {
4301        if( (*num)[i]*(*den)[j] >= (*num)[j]*(*den)[i] )
4302        {
4303            return  semicListNotMonotonous;
4304        }
4305    }
4306
4307    // ---------------------
4308    //  check Milnor number
4309    // ---------------------
4310
4311    for( mu=0, i=0; i<n; i++ )
4312    {
4313        mu += (*mul)[i];
4314    }
4315
4316    if( mu != (int)(long)(l->m[0].Data( )) )
4317    {
4318        return  semicListMilnorWrong;
4319    }
4320
4321    // -------------------------
4322    //  check geometrical genus
4323    // -------------------------
4324
4325    for( pg=0, i=0; i<n; i++ )
4326    {
4327        if( (*num)[i]<=(*den)[i] )
4328        {
4329            pg += (*mul)[i];
4330        }
4331    }
4332
4333    if( pg != (int)(long)(l->m[1].Data( )) )
4334    {
4335        return  semicListPGWrong;
4336    }
4337
4338    return  semicOK;
4339}
4340
4341// ----------------------------------------------------------------------------
4342//  this procedure is called from the interpreter
4343// ----------------------------------------------------------------------------
4344//  first  = list of spectrum numbers
4345//  second = list of spectrum numbers
4346//  result = sum of the two lists
4347// ----------------------------------------------------------------------------
4348
4349BOOLEAN spaddProc( leftv result,leftv first,leftv second )
4350{
4351    semicState  state;
4352
4353    // -----------------
4354    //  check arguments
4355    // -----------------
4356
4357    lists l1 = (lists)first->Data( );
4358    lists l2 = (lists)second->Data( );
4359
4360    if( (state=list_is_spectrum( l1 )) != semicOK )
4361    {
4362        WerrorS( "first argument is not a spectrum:" );
4363        list_error( state );
4364    }
4365    else if( (state=list_is_spectrum( l2 )) != semicOK )
4366    {
4367        WerrorS( "second argument is not a spectrum:" );
4368        list_error( state );
4369    }
4370    else
4371    {
4372        spectrum s1= spectrumFromList ( l1 );
4373        spectrum s2= spectrumFromList ( l2 );
4374        spectrum sum( s1+s2 );
4375
4376        result->rtyp = LIST_CMD;
4377        result->data = (char*)(getList(sum));
4378    }
4379
4380    return  (state!=semicOK);
4381}
4382
4383// ----------------------------------------------------------------------------
4384//  this procedure is called from the interpreter
4385// ----------------------------------------------------------------------------
4386//  first  = list of spectrum numbers
4387//  second = integer
4388//  result = the multiple of the first list by the second factor
4389// ----------------------------------------------------------------------------
4390
4391BOOLEAN spmulProc( leftv result,leftv first,leftv second )
4392{
4393    semicState  state;
4394
4395    // -----------------
4396    //  check arguments
4397    // -----------------
4398
4399    lists   l = (lists)first->Data( );
4400    int     k = (int)(long)second->Data( );
4401
4402    if( (state=list_is_spectrum( l ))!=semicOK )
4403    {
4404        WerrorS( "first argument is not a spectrum" );
4405        list_error( state );
4406    }
4407    else if( k < 0 )
4408    {
4409        WerrorS( "second argument should be positive" );
4410        state = semicMulNegative;
4411    }
4412    else
4413    {
4414        spectrum s= spectrumFromList( l );
4415        spectrum product( k*s );
4416
4417        result->rtyp = LIST_CMD;
4418        result->data = (char*)getList(product);
4419    }
4420
4421    return  (state!=semicOK);
4422}
4423
4424// ----------------------------------------------------------------------------
4425//  this procedure is called from the interpreter
4426// ----------------------------------------------------------------------------
4427//  first  = list of spectrum numbers
4428//  second = list of spectrum numbers
4429//  result = semicontinuity index
4430// ----------------------------------------------------------------------------
4431
4432BOOLEAN    semicProc3   ( leftv res,leftv u,leftv v,leftv w )
4433{
4434  semicState  state;
4435  BOOLEAN qh=(((int)(long)w->Data())==1);
4436
4437  // -----------------
4438  //  check arguments
4439  // -----------------
4440
4441  lists l1 = (lists)u->Data( );
4442  lists l2 = (lists)v->Data( );
4443
4444  if( (state=list_is_spectrum( l1 ))!=semicOK )
4445  {
4446    WerrorS( "first argument is not a spectrum" );
4447    list_error( state );
4448  }
4449  else if( (state=list_is_spectrum( l2 ))!=semicOK )
4450  {
4451    WerrorS( "second argument is not a spectrum" );
4452    list_error( state );
4453  }
4454  else
4455  {
4456    spectrum s1= spectrumFromList( l1 );
4457    spectrum s2= spectrumFromList( l2 );
4458
4459    res->rtyp = INT_CMD;
4460    if (qh)
4461      res->data = (void*)(long)(s1.mult_spectrumh( s2 ));
4462    else
4463      res->data = (void*)(long)(s1.mult_spectrum( s2 ));
4464  }
4465
4466  // -----------------
4467  //  check status
4468  // -----------------
4469
4470  return  (state!=semicOK);
4471}
4472BOOLEAN    semicProc   ( leftv res,leftv u,leftv v )
4473{
4474  sleftv tmp;
4475  memset(&tmp,0,sizeof(tmp));
4476  tmp.rtyp=INT_CMD;
4477  /* tmp.data = (void *)0;  -- done by memset */
4478
4479  return  semicProc3(res,u,v,&tmp);
4480}
4481
4482#endif
4483
4484BOOLEAN loNewtonP( leftv res, leftv arg1 )
4485{
4486  res->data= (void*)loNewtonPolytope( (ideal)arg1->Data() );
4487  return FALSE;
4488}
4489
4490BOOLEAN loSimplex( leftv res, leftv args )
4491{
4492  if ( !(rField_is_long_R(currRing)) )
4493  {
4494    WerrorS("Ground field not implemented!");
4495    return TRUE;
4496  }
4497
4498  simplex * LP;
4499  matrix m;
4500
4501  leftv v= args;
4502  if ( v->Typ() != MATRIX_CMD ) // 1: matrix
4503    return TRUE;
4504  else
4505    m= (matrix)(v->CopyD());
4506
4507  LP = new simplex(MATROWS(m),MATCOLS(m));
4508  LP->mapFromMatrix(m);
4509
4510  v= v->next;
4511  if ( v->Typ() != INT_CMD )    // 2: m = number of constraints
4512    return TRUE;
4513  else
4514    LP->m= (int)(long)(v->Data());
4515
4516  v= v->next;
4517  if ( v->Typ() != INT_CMD )    // 3: n = number of variables
4518    return TRUE;
4519  else
4520    LP->n= (int)(long)(v->Data());
4521
4522  v= v->next;
4523  if ( v->Typ() != INT_CMD )    // 4: m1 = number of <= constraints
4524    return TRUE;
4525  else
4526    LP->m1= (int)(long)(v->Data());
4527
4528  v= v->next;
4529  if ( v->Typ() != INT_CMD )    // 5: m2 = number of >= constraints
4530    return TRUE;
4531  else
4532    LP->m2= (int)(long)(v->Data());
4533
4534  v= v->next;
4535  if ( v->Typ() != INT_CMD )    // 6: m3 = number of == constraints
4536    return TRUE;
4537  else
4538    LP->m3= (int)(long)(v->Data());
4539
4540#ifdef mprDEBUG_PROT
4541  Print("m (constraints) %d\n",LP->m);
4542  Print("n (columns) %d\n",LP->n);
4543  Print("m1 (<=) %d\n",LP->m1);
4544  Print("m2 (>=) %d\n",LP->m2);
4545  Print("m3 (==) %d\n",LP->m3);
4546#endif
4547
4548  LP->compute();
4549
4550  lists lres= (lists)omAlloc( sizeof(slists) );
4551  lres->Init( 6 );
4552
4553  lres->m[0].rtyp= MATRIX_CMD; // output matrix
4554  lres->m[0].data=(void*)LP->mapToMatrix(m);
4555
4556  lres->m[1].rtyp= INT_CMD;   // found a solution?
4557  lres->m[1].data=(void*)(long)LP->icase;
4558
4559  lres->m[2].rtyp= INTVEC_CMD;
4560  lres->m[2].data=(void*)LP->posvToIV();
4561
4562  lres->m[3].rtyp= INTVEC_CMD;
4563  lres->m[3].data=(void*)LP->zrovToIV();
4564
4565  lres->m[4].rtyp= INT_CMD;
4566  lres->m[4].data=(void*)(long)LP->m;
4567
4568  lres->m[5].rtyp= INT_CMD;
4569  lres->m[5].data=(void*)(long)LP->n;
4570
4571  res->data= (void*)lres;
4572
4573  return FALSE;
4574}
4575
4576BOOLEAN nuMPResMat( leftv res, leftv arg1, leftv arg2 )
4577{
4578  ideal gls = (ideal)(arg1->Data());
4579  int imtype= (int)(long)arg2->Data();
4580
4581  uResultant::resMatType mtype= determineMType( imtype );
4582
4583  // check input ideal ( = polynomial system )
4584  if ( mprIdealCheck( gls, arg1->Name(), mtype, true ) != mprOk )
4585  {
4586    return TRUE;
4587  }
4588
4589  uResultant *resMat= new uResultant( gls, mtype, false );
4590  if (resMat!=NULL)
4591  {
4592    res->rtyp = MODUL_CMD;
4593    res->data= (void*)resMat->accessResMat()->getMatrix();
4594    if (!errorreported) delete resMat;
4595  }
4596  return errorreported;
4597}
4598
4599BOOLEAN nuLagSolve( leftv res, leftv arg1, leftv arg2, leftv arg3 )
4600{
4601
4602  poly gls;
4603  gls= (poly)(arg1->Data());
4604  int howclean= (int)(long)arg3->Data();
4605
4606  if ( !(rField_is_R(currRing) ||
4607         rField_is_Q(currRing) ||
4608         rField_is_long_R(currRing) ||
4609         rField_is_long_C(currRing)) )
4610  {
4611    WerrorS("Ground field not implemented!");
4612    return TRUE;
4613  }
4614
4615  if ( !(rField_is_R(currRing) || rField_is_long_R(currRing) || \
4616                          rField_is_long_C(currRing)) )
4617  {
4618    unsigned long int ii = (unsigned long int)arg2->Data();
4619    setGMPFloatDigits( ii, ii );
4620  }
4621
4622  if ( gls == NULL || pIsConstant( gls ) )
4623  {
4624    WerrorS("Input polynomial is constant!");
4625    return TRUE;
4626  }
4627
4628  int ldummy;
4629  int deg= currRing->pLDeg( gls, &ldummy, currRing );
4630  //  int deg= pDeg( gls );
4631  //  int len= pLength( gls );
4632  int i,vpos=0;
4633  poly piter;
4634  lists elist;
4635  lists rlist;
4636
4637  elist= (lists)omAlloc( sizeof(slists) );
4638  elist->Init( 0 );
4639
4640  if ( rVar(currRing) > 1 )
4641  {
4642    piter= gls;
4643    for ( i= 1; i <= rVar(currRing); i++ )
4644      if ( pGetExp( piter, i ) )
4645      {
4646        vpos= i;
4647        break;
4648      }
4649    while ( piter )
4650    {
4651      for ( i= 1; i <= rVar(currRing); i++ )
4652        if ( (vpos != i) && (pGetExp( piter, i ) != 0) )
4653        {
4654          WerrorS("The input polynomial must be univariate!");
4655          return TRUE;
4656        }
4657      pIter( piter );
4658    }
4659  }
4660
4661  rootContainer * roots= new rootContainer();
4662  number * pcoeffs= (number *)omAlloc( (deg+1) * sizeof( number ) );
4663  piter= gls;
4664  for ( i= deg; i >= 0; i-- )
4665  {
4666    //if ( piter ) Print("deg %d, pDeg(piter) %d\n",i,pTotaldegree(piter));
4667    if ( piter && pTotaldegree(piter) == i )
4668    {
4669      pcoeffs[i]= nCopy( pGetCoeff( piter ) );
4670      //nPrint( pcoeffs[i] );PrintS("  ");
4671      pIter( piter );
4672    }
4673    else
4674    {
4675      pcoeffs[i]= nInit(0);
4676    }
4677  }
4678
4679#ifdef mprDEBUG_PROT
4680  for (i=deg; i >= 0; i--)
4681  {
4682    nPrint( pcoeffs[i] );PrintS("  ");
4683  }
4684  PrintLn();
4685#endif
4686
4687  roots->fillContainer( pcoeffs, NULL, 1, deg, rootContainer::onepoly, 1 );
4688  roots->solver( howclean );
4689
4690  int elem= roots->getAnzRoots();
4691  char *dummy;
4692  int j;
4693
4694  rlist= (lists)omAlloc( sizeof(slists) );
4695  rlist->Init( elem );
4696
4697  if (rField_is_long_C(currRing))
4698  {
4699    for ( j= 0; j < elem; j++ )
4700    {
4701      rlist->m[j].rtyp=NUMBER_CMD;
4702      rlist->m[j].data=(void *)nCopy((number)(roots->getRoot(j)));
4703      //rlist->m[j].data=(void *)(number)(roots->getRoot(j));
4704    }
4705  }
4706  else
4707  {
4708    for ( j= 0; j < elem; j++ )
4709    {
4710      dummy = complexToStr( (*roots)[j], gmp_output_digits, currRing->cf );
4711      rlist->m[j].rtyp=STRING_CMD;
4712      rlist->m[j].data=(void *)dummy;
4713    }
4714  }
4715
4716  elist->Clean();
4717  //omFreeSize( (ADDRESS) elist, sizeof(slists) );
4718
4719  // this is (via fillContainer) the same data as in root
4720  //for ( i= deg; i >= 0; i-- ) nDelete( &pcoeffs[i] );
4721  //omFreeSize( (ADDRESS) pcoeffs, (deg+1) * sizeof( number ) );
4722
4723  delete roots;
4724
4725  res->rtyp= LIST_CMD;
4726  res->data= (void*)rlist;
4727
4728  return FALSE;
4729}
4730
4731BOOLEAN nuVanderSys( leftv res, leftv arg1, leftv arg2, leftv arg3)
4732{
4733  int i;
4734  ideal p,w;
4735  p= (ideal)arg1->Data();
4736  w= (ideal)arg2->Data();
4737
4738  // w[0] = f(p^0)
4739  // w[1] = f(p^1)
4740  // ...
4741  // p can be a vector of numbers (multivariate polynom)
4742  //   or one number (univariate polynom)
4743  // tdg = deg(f)
4744
4745  int n= IDELEMS( p );
4746  int m= IDELEMS( w );
4747  int tdg= (int)(long)arg3->Data();
4748
4749  res->data= (void*)NULL;
4750
4751  // check the input
4752  if ( tdg < 1 )
4753  {
4754    WerrorS("Last input parameter must be > 0!");
4755    return TRUE;
4756  }
4757  if ( n != rVar(currRing) )
4758  {
4759    Werror("Size of first input ideal must be equal to %d!",rVar(currRing));
4760    return TRUE;
4761  }
4762  if ( m != (int)pow((double)tdg+1,(double)n) )
4763  {
4764    Werror("Size of second input ideal must be equal to %d!",
4765      (int)pow((double)tdg+1,(double)n));
4766    return TRUE;
4767  }
4768  if ( !(rField_is_Q(currRing) /* ||
4769         rField_is_R() || rField_is_long_R() ||
4770         rField_is_long_C()*/ ) )
4771         {
4772    WerrorS("Ground field not implemented!");
4773    return TRUE;
4774  }
4775
4776  number tmp;
4777  number *pevpoint= (number *)omAlloc( n * sizeof( number ) );
4778  for ( i= 0; i < n; i++ )
4779  {
4780    pevpoint[i]=nInit(0);
4781    if (  (p->m)[i] )
4782    {
4783      tmp = pGetCoeff( (p->m)[i] );
4784      if ( nIsZero(tmp) || nIsOne(tmp) || nIsMOne(tmp) )
4785      {
4786        omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4787        WerrorS("Elements of first input ideal must not be equal to -1, 0, 1!");
4788        return TRUE;
4789      }
4790    } else tmp= NULL;
4791    if ( !nIsZero(tmp) )
4792    {
4793      if ( !pIsConstant((p->m)[i]))
4794      {
4795        omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4796        WerrorS("Elements of first input ideal must be numbers!");
4797        return TRUE;
4798      }
4799      pevpoint[i]= nCopy( tmp );
4800    }
4801  }
4802
4803  number *wresults= (number *)omAlloc( m * sizeof( number ) );
4804  for ( i= 0; i < m; i++ )
4805  {
4806    wresults[i]= nInit(0);
4807    if ( (w->m)[i] && !nIsZero(pGetCoeff((w->m)[i])) )
4808    {
4809      if ( !pIsConstant((w->m)[i]))
4810      {
4811        omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4812        omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
4813        WerrorS("Elements of second input ideal must be numbers!");
4814        return TRUE;
4815      }
4816      wresults[i]= nCopy(pGetCoeff((w->m)[i]));
4817    }
4818  }
4819
4820  vandermonde vm( m, n, tdg, pevpoint, FALSE );
4821  number *ncpoly= vm.interpolateDense( wresults );
4822  // do not free ncpoly[]!!
4823  poly rpoly= vm.numvec2poly( ncpoly );
4824
4825  omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4826  omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
4827
4828  res->data= (void*)rpoly;
4829  return FALSE;
4830}
4831
4832BOOLEAN nuUResSolve( leftv res, leftv args )
4833{
4834  leftv v= args;
4835
4836  ideal gls;
4837  int imtype;
4838  int howclean;
4839
4840  // get ideal
4841  if ( v->Typ() != IDEAL_CMD )
4842    return TRUE;
4843  else gls= (ideal)(v->Data());
4844  v= v->next;
4845
4846  // get resultant matrix type to use (0,1)
4847  if ( v->Typ() != INT_CMD )
4848    return TRUE;
4849  else imtype= (int)(long)v->Data();
4850  v= v->next;
4851
4852  if (imtype==0)
4853  {
4854    ideal test_id=idInit(1,1);
4855    int j;
4856    for(j=IDELEMS(gls)-1;j>=0;j--)
4857    {
4858      if (gls->m[j]!=NULL)
4859      {
4860        test_id->m[0]=gls->m[j];
4861        intvec *dummy_w=id_QHomWeight(test_id, currRing);
4862        if (dummy_w!=NULL)
4863        {
4864          WerrorS("Newton polytope not of expected dimension");
4865          delete dummy_w;
4866          return TRUE;
4867        }
4868      }
4869    }
4870  }
4871
4872  // get and set precision in digits ( > 0 )
4873  if ( v->Typ() != INT_CMD )
4874    return TRUE;
4875  else if ( !(rField_is_R(currRing) || rField_is_long_R(currRing) || \
4876                          rField_is_long_C(currRing)) )
4877  {
4878    unsigned long int ii=(unsigned long int)v->Data();
4879    setGMPFloatDigits( ii, ii );
4880  }
4881  v= v->next;
4882
4883  // get interpolation steps (0,1,2)
4884  if ( v->Typ() != INT_CMD )
4885    return TRUE;
4886  else howclean= (int)(long)v->Data();
4887
4888  uResultant::resMatType mtype= determineMType( imtype );
4889  int i,count;
4890  lists listofroots= NULL;
4891  number smv= NULL;
4892  BOOLEAN interpolate_det= (mtype==uResultant::denseResMat)?TRUE:FALSE;
4893
4894  //emptylist= (lists)omAlloc( sizeof(slists) );
4895  //emptylist->Init( 0 );
4896
4897  //res->rtyp = LIST_CMD;
4898  //res->data= (void *)emptylist;
4899
4900  // check input ideal ( = polynomial system )
4901  if ( mprIdealCheck( gls, args->Name(), mtype ) != mprOk )
4902  {
4903    return TRUE;
4904  }
4905
4906  uResultant * ures;
4907  rootContainer ** iproots;
4908  rootContainer ** muiproots;
4909  rootArranger * arranger;
4910
4911  // main task 1: setup of resultant matrix
4912  ures= new uResultant( gls, mtype );
4913  if ( ures->accessResMat()->initState() != resMatrixBase::ready )
4914  {
4915    WerrorS("Error occurred during matrix setup!");
4916    return TRUE;
4917  }
4918
4919  // if dense resultant, check if minor nonsingular
4920  if ( mtype == uResultant::denseResMat )
4921  {
4922    smv= ures->accessResMat()->getSubDet();
4923#ifdef mprDEBUG_PROT
4924    PrintS("// Determinant of submatrix: ");nPrint(smv);PrintLn();
4925#endif
4926    if ( nIsZero(smv) )
4927    {
4928      WerrorS("Unsuitable input ideal: Minor of resultant matrix is singular!");
4929      return TRUE;
4930    }
4931  }
4932
4933  // main task 2: Interpolate specialized resultant polynomials
4934  if ( interpolate_det )
4935    iproots= ures->interpolateDenseSP( false, smv );
4936  else
4937    iproots= ures->specializeInU( false, smv );
4938
4939  // main task 3: Interpolate specialized resultant polynomials
4940  if ( interpolate_det )
4941    muiproots= ures->interpolateDenseSP( true, smv );
4942  else
4943    muiproots= ures->specializeInU( true, smv );
4944
4945#ifdef mprDEBUG_PROT
4946  int c= iproots[0]->getAnzElems();
4947  for (i=0; i < c; i++) pWrite(iproots[i]->getPoly());
4948  c= muiproots[0]->getAnzElems();
4949  for (i=0; i < c; i++) pWrite(muiproots[i]->getPoly());
4950#endif
4951
4952  // main task 4: Compute roots of specialized polys and match them up
4953  arranger= new rootArranger( iproots, muiproots, howclean );
4954  arranger->solve_all();
4955
4956  // get list of roots
4957  if ( arranger->success() )
4958  {
4959    arranger->arrange();
4960    listofroots= listOfRoots(arranger, gmp_output_digits );
4961  }
4962  else
4963  {
4964    WerrorS("Solver was unable to find any roots!");
4965    return TRUE;
4966  }
4967
4968  // free everything
4969  count= iproots[0]->getAnzElems();
4970  for (i=0; i < count; i++) delete iproots[i];
4971  omFreeSize( (ADDRESS) iproots, count * sizeof(rootContainer*) );
4972  count= muiproots[0]->getAnzElems();
4973  for (i=0; i < count; i++) delete muiproots[i];
4974  omFreeSize( (ADDRESS) muiproots, count * sizeof(rootContainer*) );
4975
4976  delete ures;
4977  delete arranger;
4978  nDelete( &smv );
4979
4980  res->data= (void *)listofroots;
4981
4982  //emptylist->Clean();
4983  //  omFreeSize( (ADDRESS) emptylist, sizeof(slists) );
4984
4985  return FALSE;
4986}
4987
4988// from mpr_numeric.cc
4989lists listOfRoots( rootArranger* self, const unsigned int oprec )
4990{
4991  int i,j;
4992  int count= self->roots[0]->getAnzRoots(); // number of roots
4993  int elem= self->roots[0]->getAnzElems();  // number of koordinates per root
4994
4995  lists listofroots= (lists)omAlloc( sizeof(slists) ); // must be done this way!
4996
4997  if ( self->found_roots )
4998  {
4999    listofroots->Init( count );
5000
5001    for (i=0; i < count; i++)
5002    {
5003      lists onepoint= (lists)omAlloc(sizeof(slists)); // must be done this way!
5004      onepoint->Init(elem);
5005      for ( j= 0; j < elem; j++ )
5006      {
5007        if ( !rField_is_long_C(currRing) )
5008        {
5009          onepoint->m[j].rtyp=STRING_CMD;
5010          onepoint->m[j].data=(void *)complexToStr((*self->roots[j])[i],oprec, currRing->cf);
5011        }
5012        else
5013        {
5014          onepoint->m[j].rtyp=NUMBER_CMD;
5015          onepoint->m[j].data=(void *)n_Copy((number)(self->roots[j]->getRoot(i)), currRing->cf);
5016        }
5017        onepoint->m[j].next= NULL;
5018        onepoint->m[j].name= NULL;
5019      }
5020      listofroots->m[i].rtyp=LIST_CMD;
5021      listofroots->m[i].data=(void *)onepoint;
5022      listofroots->m[j].next= NULL;
5023      listofroots->m[j].name= NULL;
5024    }
5025
5026  }
5027  else
5028  {
5029    listofroots->Init( 0 );
5030  }
5031
5032  return listofroots;
5033}
5034
5035// from ring.cc
5036void rSetHdl(idhdl h)
5037{
5038  ring rg = NULL;
5039  if (h!=NULL)
5040  {
5041//   Print(" new ring:%s (l:%d)\n",IDID(h),IDLEV(h));
5042    rg = IDRING(h);
5043    if (rg==NULL) return; //id <>NULL, ring==NULL
5044    omCheckAddrSize((ADDRESS)h,sizeof(idrec));
5045    if (IDID(h))  // OB: ????
5046      omCheckAddr((ADDRESS)IDID(h));
5047    rTest(rg);
5048  }
5049
5050  // clean up history
5051  if (sLastPrinted.RingDependend())
5052  {
5053    sLastPrinted.CleanUp();
5054    memset(&sLastPrinted,0,sizeof(sleftv));
5055  }
5056
5057  if ((rg!=currRing)&&(currRing!=NULL))
5058  {
5059    denominator_list dd=DENOMINATOR_LIST;
5060    if (DENOMINATOR_LIST!=NULL)
5061    {
5062      if (TEST_V_ALLWARN)
5063        Warn("deleting denom_list for ring change to %s",IDID(h));
5064      do
5065      {
5066        n_Delete(&(dd->n),currRing->cf);
5067        dd=dd->next;
5068        omFree(DENOMINATOR_LIST);
5069        DENOMINATOR_LIST=dd;
5070      } while(DENOMINATOR_LIST!=NULL);
5071    }
5072  }
5073
5074  // test for valid "currRing":
5075  if ((rg!=NULL) && (rg->idroot==NULL))
5076  {
5077    ring old=rg;
5078    rg=rAssure_HasComp(rg);
5079    if (old!=rg)
5080    {
5081      rKill(old);
5082      IDRING(h)=rg;
5083    }
5084  }
5085   /*------------ change the global ring -----------------------*/
5086  rChangeCurrRing(rg);
5087  currRingHdl = h;
5088}
5089
5090static leftv rOptimizeOrdAsSleftv(leftv ord)
5091{
5092  // change some bad orderings/combination into better ones
5093  leftv h=ord;
5094  while(h!=NULL)
5095  {
5096    BOOLEAN change=FALSE;
5097    intvec *iv = (intvec *)(h->data);
5098 // ws(-i) -> wp(i)
5099    if ((*iv)[1]==ringorder_ws)
5100    {
5101      BOOLEAN neg=TRUE;
5102      for(int i=2;i<iv->length();i++)
5103        if((*iv)[i]>=0) { neg=FALSE; break; }
5104      if (neg)
5105      {
5106        (*iv)[1]=ringorder_wp;
5107        for(int i=2;i<iv->length();i++)
5108          (*iv)[i]= - (*iv)[i];
5109        change=TRUE;
5110      }
5111    }
5112 // Ws(-i) -> Wp(i)
5113    if ((*iv)[1]==ringorder_Ws)
5114    {
5115      BOOLEAN neg=TRUE;
5116      for(int i=2;i<iv->length();i++)
5117        if((*iv)[i]>=0) { neg=FALSE; break; }
5118      if (neg)
5119      {
5120        (*iv)[1]=ringorder_Wp;
5121        for(int i=2;i<iv->length();i++)
5122          (*iv)[i]= -(*iv)[i];
5123        change=TRUE;
5124      }
5125    }
5126 // wp(1) -> dp
5127    if ((*iv)[1]==ringorder_wp)
5128    {
5129      BOOLEAN all_one=TRUE;
5130      for(int i=2;i<iv->length();i++)
5131        if((*iv)[i]!=1) { all_one=FALSE; break; }
5132      if (all_one)
5133      {
5134        intvec *iv2=new intvec(3);
5135        (*iv2)[0]=1;
5136        (*iv2)[1]=ringorder_dp;
5137        (*iv2)[2]=iv->length()-2;
5138        delete iv;
5139        iv=iv2;
5140        h->data=iv2;
5141        change=TRUE;
5142      }
5143    }
5144 // Wp(1) -> Dp
5145    if ((*iv)[1]==ringorder_Wp)
5146    {
5147      BOOLEAN all_one=TRUE;
5148      for(int i=2;i<iv->length();i++)
5149        if((*iv)[i]!=1) { all_one=FALSE; break; }
5150      if (all_one)
5151      {
5152        intvec *iv2=new intvec(3);
5153        (*iv2)[0]=1;
5154        (*iv2)[1]=ringorder_Dp;
5155        (*iv2)[2]=iv->length()-2;
5156        delete iv;
5157        iv=iv2;
5158        h->data=iv2;
5159        change=TRUE;
5160      }
5161    }
5162 // dp(1)/Dp(1)/rp(1) -> lp(1)
5163    if (((*iv)[1]==ringorder_dp)
5164    || ((*iv)[1]==ringorder_Dp)
5165    || ((*iv)[1]==ringorder_rp))
5166    {
5167      if (iv->length()==3)
5168      {
5169        if ((*iv)[2]==1)
5170        {
5171          (*iv)[1]=ringorder_lp;
5172          change=TRUE;
5173        }
5174      }
5175    }
5176 // lp(i),lp(j) -> lp(i+j)
5177    if(((*iv)[1]==ringorder_lp)
5178    && (h->next!=NULL))
5179    {
5180      intvec *iv2 = (intvec *)(h->next->data);
5181      if ((*iv2)[1]==ringorder_lp)
5182      {
5183        leftv hh=h->next;
5184        h->next=hh->next;
5185        hh->next=NULL;
5186        if ((*iv2)[0]==1)
5187          (*iv)[2] += 1; // last block unspecified, at least 1
5188        else
5189          (*iv)[2] += (*iv2)[2];
5190        hh->CleanUp();
5191        omFree(hh);
5192        change=TRUE;
5193      }
5194    }
5195   // -------------------
5196    if (!change) h=h->next;
5197 }
5198 return ord;
5199}
5200
5201
5202BOOLEAN rSleftvOrdering2Ordering(sleftv *ord, ring R)
5203{
5204  int last = 0, o=0, n = 1, i=0, typ = 1, j;
5205  ord=rOptimizeOrdAsSleftv(ord);
5206  sleftv *sl = ord;
5207
5208  // determine nBlocks
5209  while (sl!=NULL)
5210  {
5211    intvec *iv = (intvec *)(sl->data);
5212    if (((*iv)[1]==ringorder_c)||((*iv)[1]==ringorder_C))
5213      i++;
5214    else if ((*iv)[1]==ringorder_L)
5215    {
5216      R->bitmask=(*iv)[2];
5217      n--;
5218    }
5219    else if (((*iv)[1]!=ringorder_a)
5220    && ((*iv)[1]!=ringorder_a64)
5221    && ((*iv)[1]!=ringorder_am))
5222      o++;
5223    n++;
5224    sl=sl->next;
5225  }
5226  // check whether at least one real ordering
5227  if (o==0)
5228  {
5229    WerrorS("invalid combination of orderings");
5230    return TRUE;
5231  }
5232  // if no c/C ordering is given, increment n
5233  if (i==0) n++;
5234  else if (i != 1)
5235  {
5236    // throw error if more than one is given
5237    WerrorS("more than one ordering c/C specified");
5238    return TRUE;
5239  }
5240
5241  // initialize fields of R
5242  R->order=(int *)omAlloc0(n*sizeof(int));
5243  R->block0=(int *)omAlloc0(n*sizeof(int));
5244  R->block1=(int *)omAlloc0(n*sizeof(int));
5245  R->wvhdl=(int**)omAlloc0(n*sizeof(int_ptr));
5246
5247  int *weights=(int*)omAlloc0((R->N+1)*sizeof(int));
5248
5249  // init order, so that rBlocks works correctly
5250  for (j=0; j < n-1; j++)
5251    R->order[j] = (int) ringorder_unspec;
5252  // set last _C order, if no c/C order was given
5253  if (i == 0) R->order[n-2] = ringorder_C;
5254
5255  /* init orders */
5256  sl=ord;
5257  n=-1;
5258  while (sl!=NULL)
5259  {
5260    intvec *iv;
5261    iv = (intvec *)(sl->data);
5262    if ((*iv)[1]!=ringorder_L)
5263    {
5264      n++;
5265
5266      /* the format of an ordering:
5267       *  iv[0]: factor
5268       *  iv[1]: ordering
5269       *  iv[2..end]: weights
5270       */
5271      R->order[n] = (*iv)[1];
5272      typ=1;
5273      switch ((*iv)[1])
5274      {
5275          case ringorder_ws:
5276          case ringorder_Ws:
5277            typ=-1;
5278          case ringorder_wp:
5279          case ringorder_Wp:
5280            R->wvhdl[n]=(int*)omAlloc((iv->length()-1)*sizeof(int));
5281            R->block0[n] = last+1;
5282            for (i=2; i<iv->length(); i++)
5283            {
5284              R->wvhdl[n][i-2] = (*iv)[i];
5285              last++;
5286              if (weights[last]==0) weights[last]=(*iv)[i]*typ;
5287            }
5288            R->block1[n] = si_min(last,R->N);
5289            break;
5290          case ringorder_ls:
5291          case ringorder_ds:
5292          case ringorder_Ds:
5293          case ringorder_rs:
5294            typ=-1;
5295          case ringorder_lp:
5296          case ringorder_dp:
5297          case ringorder_Dp:
5298          case ringorder_rp:
5299            R->block0[n] = last+1;
5300            if (iv->length() == 3) last+=(*iv)[2];
5301            else last += (*iv)[0];
5302            R->block1[n] = si_min(last,R->N);
5303            if (rCheckIV(iv)) return TRUE;
5304            for(i=si_min(rVar(R),R->block1[n]);i>=R->block0[n];i--)
5305            {
5306              if (weights[i]==0) weights[i]=typ;
5307            }
5308            break;
5309
5310          case ringorder_s: // no 'rank' params!
5311          {
5312
5313            if(iv->length() > 3)
5314              return TRUE;
5315
5316            if(iv->length() == 3)
5317            {
5318              const int s = (*iv)[2];
5319              R->block0[n] = s;
5320              R->block1[n] = s;
5321            }
5322            break;
5323          }
5324          case ringorder_IS:
5325          {
5326            if(iv->length() != 3) return TRUE;
5327
5328            const int s = (*iv)[2];
5329
5330            if( 1 < s || s < -1 ) return TRUE;
5331
5332            R->block0[n] = s;
5333            R->block1[n] = s;
5334            break;
5335          }
5336          case ringorder_S:
5337          case ringorder_c:
5338          case ringorder_C:
5339          {
5340            if (rCheckIV(iv)) return TRUE;
5341            break;
5342          }
5343          case ringorder_aa:
5344          case ringorder_a:
5345          {
5346            R->block0[n] = last+1;
5347            R->block1[n] = si_min(last+iv->length()-2 , R->N);
5348            R->wvhdl[n] = (int*)omAlloc((iv->length()-1)*sizeof(int));
5349            for (i=2; i<iv->length(); i++)
5350            {
5351              R->wvhdl[n][i-2]=(*iv)[i];
5352              last++;
5353              if (weights[last]==0) weights[last]=(*iv)[i]*typ;
5354            }
5355            last=R->block0[n]-1;
5356            break;
5357          }
5358          case ringorder_am:
5359          {
5360            R->block0[n] = last+1;
5361            R->block1[n] = si_min(last+iv->length()-2 , R->N);
5362            R->wvhdl[n] = (int*)omAlloc(iv->length()*sizeof(int));
5363            if (R->block1[n]- R->block0[n]+2>=iv->length())
5364               WarnS("missing module weights");
5365            for (i=2; i<=(R->block1[n]-R->block0[n]+2); i++)
5366            {
5367              R->wvhdl[n][i-2]=(*iv)[i];
5368              last++;
5369              if (weights[last]==0) weights[last]=(*iv)[i]*typ;
5370            }
5371            R->wvhdl[n][i-2]=iv->length() -3 -(R->block1[n]- R->block0[n]);
5372            for (; i<iv->length(); i++)
5373            {
5374              R->wvhdl[n][i-1]=(*iv)[i];
5375            }
5376            last=R->block0[n]-1;
5377            break;
5378          }
5379          case ringorder_a64:
5380          {
5381            R->block0[n] = last+1;
5382            R->block1[n] = si_min(last+iv->length()-2 , R->N);
5383            R->wvhdl[n] = (int*)omAlloc((iv->length()-1)*sizeof(int64));
5384            int64 *w=(int64 *)R->wvhdl[n];
5385            for (i=2; i<iv->length(); i++)
5386            {
5387              w[i-2]=(*iv)[i];
5388              last++;
5389              if (weights[last]==0) weights[last]=(*iv)[i]*typ;
5390            }
5391            last=R->block0[n]-1;
5392            break;
5393          }
5394          case ringorder_M:
5395          {
5396            int Mtyp=rTypeOfMatrixOrder(iv);
5397            if (Mtyp==0) return TRUE;
5398            if (Mtyp==-1) typ = -1;
5399
5400            R->wvhdl[n] =( int *)omAlloc((iv->length()-1)*sizeof(int));
5401            for (i=2; i<iv->length();i++)
5402              R->wvhdl[n][i-2]=(*iv)[i];
5403
5404            R->block0[n] = last+1;
5405            last += (int)sqrt((double)(iv->length()-2));
5406            R->block1[n] = si_min(last,R->N);
5407            for(i=R->block1[n];i>=R->block0[n];i--)
5408            {
5409              if (weights[i]==0) weights[i]=typ;
5410            }
5411            break;
5412          }
5413
5414          case ringorder_no:
5415            R->order[n] = ringorder_unspec;
5416            return TRUE;
5417
5418          default:
5419            Werror("Internal Error: Unknown ordering %d", (*iv)[1]);
5420            R->order[n] = ringorder_unspec;
5421            return TRUE;
5422      }
5423    }
5424    if (last>R->N)
5425    {
5426      Werror("mismatch of number of vars (%d) and ordering (>=%d vars)",
5427             R->N,last);
5428      return TRUE;
5429    }
5430    sl=sl->next;
5431  }
5432  // find OrdSgn:
5433  R->OrdSgn = 1;
5434  for(i=1;i<=R->N;i++)
5435  { if (weights[i]<0) { R->OrdSgn=-1;break; }}
5436  omFree(weights);
5437
5438  // check for complete coverage
5439  while ( n >= 0 && (
5440          (R->order[n]==ringorder_c)
5441      ||  (R->order[n]==ringorder_C)
5442      ||  (R->order[n]==ringorder_s)
5443      ||  (R->order[n]==ringorder_S)
5444      ||  (R->order[n]==ringorder_IS)
5445                    )) n--;
5446
5447  assume( n >= 0 );
5448
5449  if (R->block1[n] != R->N)
5450  {
5451    if (((R->order[n]==ringorder_dp) ||
5452         (R->order[n]==ringorder_ds) ||
5453         (R->order[n]==ringorder_Dp) ||
5454         (R->order[n]==ringorder_Ds) ||
5455         (R->order[n]==ringorder_rp) ||
5456         (R->order[n]==ringorder_rs) ||
5457         (R->order[n]==ringorder_lp) ||
5458         (R->order[n]==ringorder_ls))
5459        &&
5460        R->block0[n] <= R->N)
5461    {
5462      R->block1[n] = R->N;
5463    }
5464    else
5465    {
5466      Werror("mismatch of number of vars (%d) and ordering (%d vars)",
5467             R->N,R->block1[n]);
5468      return TRUE;
5469    }
5470  }
5471  return FALSE;
5472}
5473
5474static BOOLEAN rSleftvList2StringArray(leftv sl, char** p)
5475{
5476
5477  while(sl!=NULL)
5478  {
5479    if ((sl->rtyp == IDHDL)||(sl->rtyp==ALIAS_CMD))
5480    {
5481      *p = omStrDup(sl->Name());
5482    }
5483    else if (sl->name!=NULL)
5484    {
5485      *p = (char*)sl->name;
5486      sl->name=NULL;
5487    }
5488    else if (sl->rtyp==POLY_CMD)
5489    {
5490      sleftv s_sl;
5491      iiConvert(POLY_CMD,ANY_TYPE,-1,sl,&s_sl);
5492      if (s_sl.name != NULL)
5493      {
5494        *p = (char*)s_sl.name; s_sl.name=NULL;
5495      }
5496      else
5497        *p = NULL;
5498      sl->next = s_sl.next;
5499      s_sl.next = NULL;
5500      s_sl.CleanUp();
5501      if (*p == NULL) return TRUE;
5502    }
5503    else return TRUE;
5504    p++;
5505    sl=sl->next;
5506  }
5507  return FALSE;
5508}
5509
5510const short MAX_SHORT = 32767; // (1 << (sizeof(short)*8)) - 1;
5511
5512////////////////////
5513//
5514// rInit itself:
5515//
5516// INPUT:  pn: ch & parameter (names), rv: variable (names)
5517//         ord: ordering (all !=NULL)
5518// RETURN: currRingHdl on success
5519//         NULL        on error
5520// NOTE:   * makes new ring to current ring, on success
5521//         * considers input sleftv's as read-only
5522ring rInit(leftv pn, leftv rv, leftv ord)
5523{
5524#ifdef HAVE_RINGS
5525  //unsigned int ringtype = 0;
5526  mpz_ptr modBase = NULL;
5527  unsigned int modExponent = 1;
5528#endif
5529  int float_len=0;
5530  int float_len2=0;
5531  ring R = NULL;
5532  //BOOLEAN ffChar=FALSE;
5533
5534  /* ch -------------------------------------------------------*/
5535  // get ch of ground field
5536
5537  // allocated ring
5538  R = (ring) omAlloc0Bin(sip_sring_bin);
5539
5540  coeffs cf = NULL;
5541
5542  assume( pn != NULL );
5543  const int P = pn->listLength();
5544
5545  #ifdef SINGULAR_4_1
5546  if (pn->Typ()==CRING_CMD)
5547  {
5548    cf=(coeffs)pn->CopyD();
5549    leftv pnn=pn;
5550    if(P>1) /*parameter*/
5551    {
5552      pnn = pnn->next;
5553      const int pars = pnn->listLength();
5554      assume( pars > 0 );
5555      char ** names = (char**)omAlloc0(pars * sizeof(char_ptr));
5556
5557      if (rSleftvList2StringArray(pnn, names))
5558      {
5559        WerrorS("parameter expected");
5560        goto rInitError;
5561      }
5562
5563      TransExtInfo extParam;
5564
5565      extParam.r = rDefault( cf, pars, names); // Q/Zp [ p_1, ... p_pars ]
5566      for(int i=pars-1; i>=0;i--)
5567      {
5568        omFree(names[i]);
5569      }
5570      omFree(names);
5571
5572      cf = nInitChar(n_transExt, &extParam);
5573    }
5574    assume( cf != NULL );
5575  }
5576  else
5577  #endif
5578  if (pn->Typ()==INT_CMD)
5579  {
5580    int ch = (int)(long)pn->Data();
5581    leftv pnn=pn;
5582
5583    /* parameter? -------------------------------------------------------*/
5584    pnn = pnn->next;
5585
5586    if (pnn == NULL) // no params!?
5587    {
5588      if (ch!=0)
5589      {
5590        int ch2=IsPrime(ch);
5591        if ((ch<2)||(ch!=ch2))
5592        {
5593          Warn("%d is invalid as characteristic of the ground field. 32003 is used.", ch);
5594          ch=32003;
5595        }
5596        cf = nInitChar(n_Zp, (void*)(long)ch);
5597      }
5598      else
5599        cf = nInitChar(n_Q, (void*)(long)ch);
5600    }
5601    else
5602    {
5603      const int pars = pnn->listLength();
5604
5605      assume( pars > 0 );
5606
5607      // predefined finite field: (p^k, a)
5608      if ((ch!=0) && (ch!=IsPrime(ch)) && (pars == 1))
5609      {
5610        GFInfo param;
5611
5612        param.GFChar = ch;
5613        param.GFDegree = 1;
5614        param.GFPar_name = pnn->name;
5615
5616        cf = nInitChar(n_GF, &param);
5617      }
5618      else // (0/p, a, b, ..., z)
5619      {
5620        if ((ch!=0) && (ch!=IsPrime(ch)))
5621        {
5622          WerrorS("too many parameters");
5623          goto rInitError;
5624        }
5625
5626        char ** names = (char**)omAlloc0(pars * sizeof(char_ptr));
5627
5628        if (rSleftvList2StringArray(pnn, names))
5629        {
5630          WerrorS("parameter expected");
5631          goto rInitError;
5632        }
5633
5634        TransExtInfo extParam;
5635
5636        extParam.r = rDefault( ch, pars, names); // Q/Zp [ p_1, ... p_pars ]
5637        for(int i=pars-1; i>=0;i--)
5638        {
5639          omFree(names[i]);
5640        }
5641        omFree(names);
5642
5643        cf = nInitChar(n_transExt, &extParam);
5644      }
5645    }
5646
5647//    if (cf==NULL) goto rInitError;
5648    assume( cf != NULL );
5649  }
5650  else if ((pn->name != NULL)
5651  && ((strcmp(pn->name,"real")==0) || (strcmp(pn->name,"complex")==0)))
5652  {
5653    leftv pnn=pn->next;
5654    BOOLEAN complex_flag=(strcmp(pn->name,"complex")==0);
5655    if ((pnn!=NULL) && (pnn->Typ()==INT_CMD))
5656    {
5657      float_len=(int)(long)pnn->Data();
5658      float_len2=float_len;
5659      pnn=pnn->next;
5660      if ((pnn!=NULL) && (pnn->Typ()==INT_CMD))
5661      {
5662        float_len2=(int)(long)pnn->Data();
5663        pnn=pnn->next;
5664      }
5665    }
5666
5667    if (!complex_flag)
5668      complex_flag= (pnn!=NULL) && (pnn->name!=NULL);
5669    if( !complex_flag && (float_len2 <= (short)SHORT_REAL_LENGTH))
5670       cf=nInitChar(n_R, NULL);
5671    else // longR or longC?
5672    {
5673       LongComplexInfo param;
5674
5675       param.float_len = si_min (float_len, 32767);
5676       param.float_len2 = si_min (float_len2, 32767);
5677
5678       // set the parameter name
5679       if (complex_flag)
5680       {
5681         if (param.float_len < SHORT_REAL_LENGTH)
5682         {
5683           param.float_len= SHORT_REAL_LENGTH;
5684           param.float_len2= SHORT_REAL_LENGTH;
5685         }
5686         if ((pnn == NULL) || (pnn->name == NULL))
5687           param.par_name=(const char*)"i"; //default to i
5688         else
5689           param.par_name = (const char*)pnn->name;
5690       }
5691
5692       cf = nInitChar(complex_flag ? n_long_C: n_long_R, (void*)&param);
5693    }
5694    assume( cf != NULL );
5695  }
5696#ifdef HAVE_RINGS
5697  else if ((pn->name != NULL) && (strcmp(pn->name, "integer") == 0))
5698  {
5699    // TODO: change to use coeffs_BIGINT!?
5700    modBase = (mpz_ptr) omAlloc(sizeof(mpz_t));
5701    mpz_init_set_si(modBase, 0);
5702    if (pn->next!=NULL)
5703    {
5704      leftv pnn=pn;
5705      if (pnn->next->Typ()==INT_CMD)
5706      {
5707        pnn=pnn->next;
5708        mpz_set_ui(modBase, (int)(long) pnn->Data());
5709        if ((pnn->next!=NULL) && (pnn->next->Typ()==INT_CMD))
5710        {
5711          pnn=pnn->next;
5712          modExponent = (long) pnn->Data();
5713        }
5714        while ((pnn->next!=NULL) && (pnn->next->Typ()==INT_CMD))
5715        {
5716          pnn=pnn->next;
5717          mpz_mul_ui(modBase, modBase, (int)(long) pnn->Data());
5718        }
5719      }
5720      else if (pnn->next->Typ()==BIGINT_CMD)
5721      {
5722        number p=(number)pnn->next->CopyD();
5723        nlGMP(p,(number)modBase,coeffs_BIGINT); // TODO? // extern void   nlGMP(number &i, number n, const coeffs r); // FIXME: n_MPZ( modBase, p, coeffs_BIGINT); ?
5724        n_Delete(&p,coeffs_BIGINT);
5725      }
5726    }
5727    else
5728      cf=nInitChar(n_Z,NULL);
5729
5730    if ((mpz_cmp_ui(modBase, 1) == 0) && (mpz_cmp_ui(modBase, 0) < 0))
5731    {
5732      Werror("Wrong ground ring specification (module is 1)");
5733      goto rInitError;
5734    }
5735    if (modExponent < 1)
5736    {
5737      Werror("Wrong ground ring specification (exponent smaller than 1");
5738      goto rInitError;
5739    }
5740    // module is 0 ---> integers ringtype = 4;
5741    // we have an exponent
5742    if (modExponent > 1 && cf == NULL)
5743    {
5744      if ((mpz_cmp_ui(modBase, 2) == 0) && (modExponent <= 8*sizeof(unsigned long)))
5745      {
5746        /* this branch should be active for modExponent = 2..32 resp. 2..64,
5747           depending on the size of a long on the respective platform */
5748        //ringtype = 1;       // Use Z/2^ch
5749        cf=nInitChar(n_Z2m,(void*)(long)modExponent);
5750        mpz_clear(modBase);
5751        omFreeSize (modBase, sizeof (mpz_t));
5752      }
5753      else
5754      {
5755        if (mpz_cmp_ui(modBase,0)==0)
5756        {
5757          WerrorS("modulus must not be 0 or parameter not allowed");
5758          goto rInitError;
5759        }
5760        //ringtype = 3;
5761        ZnmInfo info;
5762        info.base= modBase;
5763        info.exp= modExponent;
5764        cf=nInitChar(n_Znm,(void*) &info); //exponent is missing
5765      }
5766    }
5767    // just a module m > 1
5768    else if (cf == NULL)
5769    {
5770      if (mpz_cmp_ui(modBase,0)==0)
5771      {
5772        WerrorS("modulus must not be 0 or parameter not allowed");
5773        goto rInitError;
5774      }
5775      //ringtype = 2;
5776      ZnmInfo info;
5777      info.base= modBase;
5778      info.exp= modExponent;
5779      cf=nInitChar(n_Zn,(void*) &info);
5780    }
5781    assume( cf != NULL );
5782  }
5783#endif
5784  // ring NEW = OLD, (), (); where OLD is a polynomial ring...
5785  else if ((pn->Typ()==RING_CMD) && (P == 1))
5786  {
5787    TransExtInfo extParam;
5788    extParam.r = (ring)pn->Data();
5789    cf = nInitChar(n_transExt, &extParam);
5790  }
5791  else if ((pn->Typ()==QRING_CMD) && (P == 1)) // same for qrings - which should be fields!?
5792  {
5793    AlgExtInfo extParam;
5794    extParam.r = (ring)pn->Data();
5795
5796    cf = nInitChar(n_algExt, &extParam);   // Q[a]/<minideal>
5797  }
5798  else
5799  {
5800    Werror("Wrong or unknown ground field specification");
5801#ifndef SING_NDEBUG
5802    sleftv* p = pn;
5803    while (p != NULL)
5804    {
5805      Print( "pn[%p]: type: %d [%s]: %p, name: %s", (void*)p, p->Typ(), Tok2Cmdname(p->Typ()), p->Data(), (p->name == NULL? "NULL" : p->name) );
5806      PrintLn();
5807      p = p->next;
5808    }
5809#endif
5810    goto rInitError;
5811  }
5812
5813  /*every entry in the new ring is initialized to 0*/
5814
5815  /* characteristic -----------------------------------------------*/
5816  /* input: 0 ch=0 : Q     parameter=NULL    ffChar=FALSE   float_len
5817   *         0    1 : Q(a,...)        *names         FALSE
5818   *         0   -1 : R               NULL           FALSE  0
5819   *         0   -1 : R               NULL           FALSE  prec. >6
5820   *         0   -1 : C               *names         FALSE  prec. 0..?
5821   *         p    p : Fp              NULL           FALSE
5822   *         p   -p : Fp(a)           *names         FALSE
5823   *         q    q : GF(q=p^n)       *names         TRUE
5824  */
5825  if (cf==NULL)
5826  {
5827    Werror("Invalid ground field specification");
5828    goto rInitError;
5829//    const int ch=32003;
5830//    cf=nInitChar(n_Zp, (void*)(long)ch);
5831  }
5832
5833  assume( R != NULL );
5834
5835  R->cf = cf;
5836
5837  /* names and number of variables-------------------------------------*/
5838  {
5839    int l=rv->listLength();
5840
5841    if (l>MAX_SHORT)
5842    {
5843      Werror("too many ring variables(%d), max is %d",l,MAX_SHORT);
5844       goto rInitError;
5845    }
5846    R->N = l; /*rv->listLength();*/
5847  }
5848  R->names   = (char **)omAlloc0(R->N * sizeof(char_ptr));
5849  if (rSleftvList2StringArray(rv, R->names))
5850  {
5851    WerrorS("name of ring variable expected");
5852    goto rInitError;
5853  }
5854
5855  /* check names and parameters for conflicts ------------------------- */
5856  rRenameVars(R); // conflicting variables will be renamed
5857  /* ordering -------------------------------------------------------------*/
5858  if (rSleftvOrdering2Ordering(ord, R))
5859    goto rInitError;
5860
5861  // Complete the initialization
5862  if (rComplete(R,1))
5863    goto rInitError;
5864
5865/*#ifdef HAVE_RINGS
5866// currently, coefficients which are ring elements require a global ordering:
5867  if (rField_is_Ring(R) && (R->OrdSgn==-1))
5868  {
5869    WerrorS("global ordering required for these coefficients");
5870    goto rInitError;
5871  }
5872#endif*/
5873
5874  rTest(R);
5875
5876  // try to enter the ring into the name list
5877  // need to clean up sleftv here, before this ring can be set to
5878  // new currRing or currRing can be killed beacuse new ring has
5879  // same name
5880  pn->CleanUp();
5881  rv->CleanUp();
5882  ord->CleanUp();
5883  //if ((tmp = enterid(s, myynest, RING_CMD, &IDROOT))==NULL)
5884  //  goto rInitError;
5885
5886  //memcpy(IDRING(tmp),R,sizeof(*R));
5887  // set current ring
5888  //omFreeBin(R,  ip_sring_bin);
5889  //return tmp;
5890  return R;
5891
5892  // error case:
5893  rInitError:
5894  if  ((R != NULL)&&(R->cf!=NULL)) rDelete(R);
5895  pn->CleanUp();
5896  rv->CleanUp();
5897  ord->CleanUp();
5898  return NULL;
5899}
5900
5901ring rSubring(ring org_ring, sleftv* rv)
5902{
5903  ring R = rCopy0(org_ring);
5904  int *perm=(int *)omAlloc0((org_ring->N+1)*sizeof(int));
5905  int n = rBlocks(org_ring), i=0, j;
5906
5907  /* names and number of variables-------------------------------------*/
5908  {
5909    int l=rv->listLength();
5910    if (l>MAX_SHORT)
5911    {
5912      Werror("too many ring variables(%d), max is %d",l,MAX_SHORT);
5913       goto rInitError;
5914    }
5915    R->N = l; /*rv->listLength();*/
5916  }
5917  omFree(R->names);
5918  R->names   = (char **)omAlloc0(R->N * sizeof(char_ptr));
5919  if (rSleftvList2StringArray(rv, R->names))
5920  {
5921    WerrorS("name of ring variable expected");
5922    goto rInitError;
5923  }
5924
5925  /* check names for subring in org_ring ------------------------- */
5926  {
5927    i=0;
5928
5929    for(j=0;j<R->N;j++)
5930    {
5931      for(;i<org_ring->N;i++)
5932      {
5933        if (strcmp(org_ring->names[i],R->names[j])==0)
5934        {
5935          perm[i+1]=j+1;
5936          break;
5937        }
5938      }
5939      if (i>org_ring->N)
5940      {
5941        Werror("variable %d (%s) not in basering",j+1,R->names[j]);
5942        break;
5943      }
5944    }
5945  }
5946  //Print("perm=");
5947  //for(i=1;i<org_ring->N;i++) Print("v%d -> v%d\n",i,perm[i]);
5948  /* ordering -------------------------------------------------------------*/
5949
5950  for(i=0;i<n;i++)
5951  {
5952    int min_var=-1;
5953    int max_var=-1;
5954    for(j=R->block0[i];j<=R->block1[i];j++)
5955    {
5956      if (perm[j]>0)
5957      {
5958        if (min_var==-1) min_var=perm[j];
5959        max_var=perm[j];
5960      }
5961    }
5962    if (min_var!=-1)
5963    {
5964      //Print("block %d: old %d..%d, now:%d..%d\n",
5965      //      i,R->block0[i],R->block1[i],min_var,max_var);
5966      R->block0[i]=min_var;
5967      R->block1[i]=max_var;
5968      if (R->wvhdl[i]!=NULL)
5969      {
5970        omFree(R->wvhdl[i]);
5971        R->wvhdl[i]=(int*)omAlloc0((max_var-min_var+1)*sizeof(int));
5972        for(j=org_ring->block0[i];j<=org_ring->block1[i];j++)
5973        {
5974          if (perm[j]>0)
5975          {
5976            R->wvhdl[i][perm[j]-R->block0[i]]=
5977                org_ring->wvhdl[i][j-org_ring->block0[i]];
5978            //Print("w%d=%d (orig_w%d)\n",perm[j],R->wvhdl[i][perm[j]-R->block0[i]],j);
5979          }
5980        }
5981      }
5982    }
5983    else
5984    {
5985      if(R->block0[i]>0)
5986      {
5987        //Print("skip block %d\n",i);
5988        R->order[i]=ringorder_unspec;
5989        if (R->wvhdl[i] !=NULL) omFree(R->wvhdl[i]);
5990        R->wvhdl[i]=NULL;
5991      }
5992      //else Print("keep block %d\n",i);
5993    }
5994  }
5995  i=n-1;
5996  while(i>0)
5997  {
5998    // removed unneded blocks
5999    if(R->order[i-1]==ringorder_unspec)
6000    {
6001      for(j=i;j<=n;j++)
6002      {
6003        R->order[j-1]=R->order[j];
6004        R->block0[j-1]=R->block0[j];
6005        R->block1[j-1]=R->block1[j];
6006        if (R->wvhdl[j-1] !=NULL) omFree(R->wvhdl[j-1]);
6007        R->wvhdl[j-1]=R->wvhdl[j];
6008      }
6009      R->order[n]=ringorder_unspec;
6010      n--;
6011    }
6012    i--;
6013  }
6014  n=rBlocks(org_ring)-1;
6015  while (R->order[n]==0)  n--;
6016  while (R->order[n]==ringorder_unspec)  n--;
6017  if ((R->order[n]==ringorder_c) ||  (R->order[n]==ringorder_C)) n--;
6018  if (R->block1[n] != R->N)
6019  {
6020    if (((R->order[n]==ringorder_dp) ||
6021         (R->order[n]==ringorder_ds) ||
6022         (R->order[n]==ringorder_Dp) ||
6023         (R->order[n]==ringorder_Ds) ||
6024         (R->order[n]==ringorder_rp) ||
6025         (R->order[n]==ringorder_rs) ||
6026         (R->order[n]==ringorder_lp) ||
6027         (R->order[n]==ringorder_ls))
6028        &&
6029        R->block0[n] <= R->N)
6030    {
6031      R->block1[n] = R->N;
6032    }
6033    else
6034    {
6035      Werror("mismatch of number of vars (%d) and ordering (%d vars) in block %d",
6036             R->N,R->block1[n],n);
6037      return NULL;
6038    }
6039  }
6040  omFree(perm);
6041  // find OrdSgn:
6042  R->OrdSgn = org_ring->OrdSgn; // IMPROVE!
6043  //for(i=1;i<=R->N;i++)
6044  //{ if (weights[i]<0) { R->OrdSgn=-1;break; }}
6045  //omFree(weights);
6046  // Complete the initialization
6047  if (rComplete(R,1))
6048    goto rInitError;
6049
6050  rTest(R);
6051
6052  if (rv != NULL) rv->CleanUp();
6053
6054  return R;
6055
6056  // error case:
6057  rInitError:
6058  if  (R != NULL) rDelete(R);
6059  if (rv != NULL) rv->CleanUp();
6060  return NULL;
6061}
6062
6063void rKill(ring r)
6064{
6065  if ((r->ref<=0)&&(r->order!=NULL))
6066  {
6067#ifdef RDEBUG
6068    if (traceit &TRACE_SHOW_RINGS) Print("kill ring %lx\n",(long)r);
6069#endif
6070    if (r->qideal!=NULL)
6071    {
6072      id_Delete(&r->qideal, r);
6073      r->qideal = NULL;
6074    }
6075    int j;
6076#ifdef USE_IILOCALRING
6077    for (j=0;j<myynest;j++)
6078    {
6079      if (iiLocalRing[j]==r)
6080      {
6081        if (j+1==myynest) Warn("killing the basering for level %d",j);
6082        iiLocalRing[j]=NULL;
6083      }
6084    }
6085#else /* USE_IILOCALRING */
6086//#endif /* USE_IILOCALRING */
6087    {
6088      proclevel * nshdl = procstack;
6089      int lev=myynest-1;
6090
6091      for(; nshdl != NULL; nshdl = nshdl->next)
6092      {
6093        if (nshdl->cRing==r)
6094        {
6095          Warn("killing the basering for level %d",lev);
6096          nshdl->cRing=NULL;
6097          nshdl->cRingHdl=NULL;
6098        }
6099      }
6100    }
6101#endif /* USE_IILOCALRING */
6102// any variables depending on r ?
6103    while (r->idroot!=NULL)
6104    {
6105      r->idroot->lev=myynest; // avoid warning about kill global objects
6106      killhdl2(r->idroot,&(r->idroot),r);
6107    }
6108    if (r==currRing)
6109    {
6110      // all dependend stuff is done, clean global vars:
6111      if ((currRing->ppNoether)!=NULL) pDelete(&(currRing->ppNoether));
6112      if (sLastPrinted.RingDependend())
6113      {
6114        sLastPrinted.CleanUp();
6115      }
6116      //if ((myynest>0) && (iiRETURNEXPR.RingDependend()))
6117      //{
6118      //  WerrorS("return value depends on local ring variable (export missing ?)");
6119      //  iiRETURNEXPR.CleanUp();
6120      //}
6121      currRing=NULL;
6122      currRingHdl=NULL;
6123    }
6124
6125    /* nKillChar(r); will be called from inside of rDelete */
6126    rDelete(r);
6127    return;
6128  }
6129  r->ref--;
6130}
6131
6132void rKill(idhdl h)
6133{
6134  ring r = IDRING(h);
6135  int ref=0;
6136  if (r!=NULL)
6137  {
6138    ref=r->ref;
6139    rKill(r);
6140  }
6141  if (h==currRingHdl)
6142  {
6143    if (ref<=0) { currRing=NULL; currRingHdl=NULL;}
6144    else
6145    {
6146      currRingHdl=rFindHdl(r,currRingHdl);
6147    }
6148  }
6149}
6150
6151idhdl rSimpleFindHdl(ring r, idhdl root, idhdl n)
6152{
6153  //idhdl next_best=NULL;
6154  idhdl h=root;
6155  while (h!=NULL)
6156  {
6157    if (((IDTYP(h)==RING_CMD)||(IDTYP(h)==QRING_CMD))
6158    && (h!=n)
6159    && (IDRING(h)==r)
6160    )
6161    {
6162   //   if (IDLEV(h)==myynest)
6163   //     return h;
6164   //   if ((IDLEV(h)==0) || (next_best==NULL))
6165   //     next_best=h;
6166   //   else if (IDLEV(next_best)<IDLEV(h))
6167   //     next_best=h;
6168      return h;
6169    }
6170    h=IDNEXT(h);
6171  }
6172  //return next_best;
6173  return NULL;
6174}
6175
6176extern BOOLEAN jjPROC(leftv res, leftv u, leftv v);
6177ideal kGroebner(ideal F, ideal Q)
6178{
6179  //test|=Sy_bit(OPT_PROT);
6180  idhdl save_ringhdl=currRingHdl;
6181  ideal resid;
6182  idhdl new_ring=NULL;
6183  if ((currRingHdl==NULL) || (IDRING(currRingHdl)!=currRing))
6184  {
6185    currRingHdl=enterid(omStrDup(" GROEBNERring"),0,RING_CMD,&IDROOT,FALSE);
6186    new_ring=currRingHdl;
6187    IDRING(currRingHdl)=currRing;
6188  }
6189  sleftv v; memset(&v,0,sizeof(v)); v.rtyp=IDEAL_CMD; v.data=(char *) F;
6190  idhdl h=ggetid("groebner");
6191  sleftv u; memset(&u,0,sizeof(u)); u.rtyp=IDHDL; u.data=(char *) h;
6192            u.name=IDID(h);
6193
6194  sleftv res; memset(&res,0,sizeof(res));
6195  if(jjPROC(&res,&u,&v))
6196  {
6197    resid=kStd(F,Q,testHomog,NULL);
6198  }
6199  else
6200  {
6201    //printf("typ:%d\n",res.rtyp);
6202    resid=(ideal)(res.data);
6203  }
6204  // cleanup GROEBNERring, save_ringhdl, u,v,(res )
6205  if (new_ring!=NULL)
6206  {
6207    idhdl h=IDROOT;
6208    if (h==new_ring) IDROOT=h->next;
6209    else
6210    {
6211      while ((h!=NULL) &&(h->next!=new_ring)) h=h->next;
6212      if (h!=NULL) h->next=h->next->next;
6213    }
6214    if (h!=NULL) omFreeSize(h,sizeof(*h));
6215  }
6216  currRingHdl=save_ringhdl;
6217  u.CleanUp();
6218  v.CleanUp();
6219  return resid;
6220}
6221
6222static void jjINT_S_TO_ID(int n,int *e, leftv res)
6223{
6224  if (n==0) n=1;
6225  ideal l=idInit(n,1);
6226  int i;
6227  poly p;
6228  for(i=rVar(currRing);i>0;i--)
6229  {
6230    if (e[i]>0)
6231    {
6232      n--;
6233      p=pOne();
6234      pSetExp(p,i,1);
6235      pSetm(p);
6236      l->m[n]=p;
6237      if (n==0) break;
6238    }
6239  }
6240  res->data=(char*)l;
6241  setFlag(res,FLAG_STD);
6242  omFreeSize((ADDRESS)e,(rVar(currRing)+1)*sizeof(int));
6243}
6244BOOLEAN jjVARIABLES_P(leftv res, leftv u)
6245{
6246  int *e=(int *)omAlloc0((rVar(currRing)+1)*sizeof(int));
6247  int n=pGetVariables((poly)u->Data(),e);
6248  jjINT_S_TO_ID(n,e,res);
6249  return FALSE;
6250}
6251
6252BOOLEAN jjVARIABLES_ID(leftv res, leftv u)
6253{
6254  int *e=(int *)omAlloc0((rVar(currRing)+1)*sizeof(int));
6255  ideal I=(ideal)u->Data();
6256  int i;
6257  int n=0;
6258  for(i=I->nrows*I->ncols-1;i>=0;i--)
6259  {
6260    int n0=pGetVariables(I->m[i],e);
6261    if (n0>n) n=n0;
6262  }
6263  jjINT_S_TO_ID(n,e,res);
6264  return FALSE;
6265}
6266
6267void paPrint(const char *n,package p)
6268{
6269  Print(" %s (",n);
6270  switch (p->language)
6271  {
6272    case LANG_SINGULAR: PrintS("S"); break;
6273    case LANG_C:        PrintS("C"); break;
6274    case LANG_TOP:      PrintS("T"); break;
6275    case LANG_NONE:     PrintS("N"); break;
6276    default:            PrintS("U");
6277  }
6278  if(p->libname!=NULL)
6279  Print(",%s", p->libname);
6280  PrintS(")");
6281}
6282
6283BOOLEAN iiApplyINTVEC(leftv res, leftv a, int op, leftv proc)
6284{
6285  intvec *aa=(intvec*)a->Data();
6286  sleftv tmp_out;
6287  sleftv tmp_in;
6288  leftv curr=res;
6289  BOOLEAN bo=FALSE;
6290  for(int i=0;i<aa->length(); i++)
6291  {
6292    memset(&tmp_in,0,sizeof(tmp_in));
6293    tmp_in.rtyp=INT_CMD;
6294    tmp_in.data=(void*)(long)(*aa)[i];
6295    if (proc==NULL)
6296      bo=iiExprArith1(&tmp_out,&tmp_in,op);
6297    else
6298      bo=jjPROC(&tmp_out,proc,&tmp_in);
6299    if (bo)
6300    {
6301      res->CleanUp(currRing);
6302      Werror("apply fails at index %d",i+1);
6303      return TRUE;
6304    }
6305    if (i==0) { memcpy(res,&tmp_out,sizeof(tmp_out)); }
6306    else
6307    {
6308      curr->next=(leftv)omAllocBin(sleftv_bin);
6309      curr=curr->next;
6310      memcpy(curr,&tmp_out,sizeof(tmp_out));
6311    }
6312  }
6313  return FALSE;
6314}
6315BOOLEAN iiApplyBIGINTMAT(leftv, leftv, int, leftv)
6316{
6317  WerrorS("not implemented");
6318  return TRUE;
6319}
6320BOOLEAN iiApplyIDEAL(leftv, leftv, int, leftv)
6321{
6322  WerrorS("not implemented");
6323  return TRUE;
6324}
6325BOOLEAN iiApplyLIST(leftv res, leftv a, int op, leftv proc)
6326{
6327  lists aa=(lists)a->Data();
6328  sleftv tmp_out;
6329  sleftv tmp_in;
6330  leftv curr=res;
6331  BOOLEAN bo=FALSE;
6332  for(int i=0;i<=aa->nr; i++)
6333  {
6334    memset(&tmp_in,0,sizeof(tmp_in));
6335    tmp_in.Copy(&(aa->m[i]));
6336    if (proc==NULL)
6337      bo=iiExprArith1(&tmp_out,&tmp_in,op);
6338    else
6339      bo=jjPROC(&tmp_out,proc,&tmp_in);
6340    tmp_in.CleanUp();
6341    if (bo)
6342    {
6343      res->CleanUp(currRing);
6344      Werror("apply fails at index %d",i+1);
6345      return TRUE;
6346    }
6347    if (i==0) { memcpy(res,&tmp_out,sizeof(tmp_out)); }
6348    else
6349    {
6350      curr->next=(leftv)omAllocBin(sleftv_bin);
6351      curr=curr->next;
6352      memcpy(curr,&tmp_out,sizeof(tmp_out));
6353    }
6354  }
6355  return FALSE;
6356}
6357BOOLEAN iiApply(leftv res, leftv a, int op, leftv proc)
6358{
6359  memset(res,0,sizeof(sleftv));
6360  res->rtyp=a->Typ();
6361  switch (res->rtyp /*a->Typ()*/)
6362  {
6363    case INTVEC_CMD:
6364    case INTMAT_CMD:
6365        return iiApplyINTVEC(res,a,op,proc);
6366    case BIGINTMAT_CMD:
6367        return iiApplyBIGINTMAT(res,a,op,proc);
6368    case IDEAL_CMD:
6369    case MODUL_CMD:
6370    case MATRIX_CMD:
6371        return iiApplyIDEAL(res,a,op,proc);
6372    case LIST_CMD:
6373        return iiApplyLIST(res,a,op,proc);
6374  }
6375  WerrorS("first argument to `apply` must allow an index");
6376  return TRUE;
6377}
6378
6379BOOLEAN iiTestAssume(leftv a, leftv b)
6380{
6381  // assume a: level
6382  if ((a->Typ()==INT_CMD)&&((long)a->Data()>=0))
6383  {
6384    if ((TEST_V_ALLWARN) && (myynest==0)) WarnS("ASSUME at top level is of no use: see documentation");
6385    char       assume_yylinebuf[80];
6386    strncpy(assume_yylinebuf,my_yylinebuf,79);
6387    int lev=(long)a->Data();
6388    int startlev=0;
6389    idhdl h=ggetid("assumeLevel");
6390    if ((h!=NULL)&&(IDTYP(h)==INT_CMD)) startlev=(long)IDINT(h);
6391    if(lev <=startlev)
6392    {
6393      BOOLEAN bo=b->Eval();
6394      if (bo) { WerrorS("syntax error in ASSUME");return TRUE;}
6395      if (b->Typ()!=INT_CMD) { WerrorS("ASUMME(<level>,<int expr>)");return TRUE; }
6396      if (b->Data()==NULL) { Werror("ASSUME failed:%s",assume_yylinebuf);return TRUE;}
6397    }
6398  }
6399  b->CleanUp();
6400  a->CleanUp();
6401  return FALSE;
6402}
6403
6404#include "libparse.h"
6405
6406BOOLEAN iiARROW(leftv r, char* a, char *s)
6407{
6408  char *ss=(char*)omAlloc(strlen(a)+strlen(s)+30); /* max. 27 currently */
6409  // find end of s:
6410  int end_s=strlen(s);
6411  while ((end_s>0) && ((s[end_s]<=' ')||(s[end_s]==';'))) end_s--;
6412  s[end_s+1]='\0';
6413  char *name=(char *)omAlloc(strlen(a)+strlen(s)+30);
6414  sprintf(name,"%s->%s",a,s);
6415  // find start of last expression
6416  int start_s=end_s-1;
6417  while ((start_s>=0) && (s[start_s]!=';')) start_s--;
6418  if (start_s<0) // ';' not found
6419  {
6420    sprintf(ss,"parameter def %s;return(%s);\n",a,s);
6421  }
6422  else // s[start_s] is ';'
6423  {
6424    s[start_s]='\0';
6425    sprintf(ss,"parameter def %s;%s;return(%s);\n",a,s,s+start_s+1);
6426  }
6427  memset(r,0,sizeof(*r));
6428  // now produce procinfo for PROC_CMD:
6429  r->data = (void *)omAlloc0Bin(procinfo_bin);
6430  ((procinfo *)(r->data))->language=LANG_NONE;
6431  iiInitSingularProcinfo((procinfo *)r->data,"",name,0,0);
6432  ((procinfo *)r->data)->data.s.body=ss;
6433  omFree(name);
6434  r->rtyp=PROC_CMD;
6435  //r->rtyp=STRING_CMD;
6436  //r->data=ss;
6437  return FALSE;
6438}
6439
6440BOOLEAN iiAssignCR(leftv r, leftv arg)
6441{
6442  int t=arg->Typ();
6443  char* ring_name=omStrDup((char*)r->Name());
6444  if ((t==RING_CMD) ||(t==QRING_CMD))
6445  {
6446    sleftv tmp;
6447    memset(&tmp,0,sizeof(tmp));
6448    tmp.rtyp=IDHDL;
6449    tmp.data=(char*)rDefault(ring_name);
6450    if (tmp.data!=NULL)
6451    {
6452      BOOLEAN b=iiAssign(&tmp,arg);
6453      if (b) return TRUE;
6454      rSetHdl(ggetid(ring_name));
6455      omFree(ring_name);
6456      return FALSE;
6457    }
6458    else
6459      return TRUE;
6460  }
6461  #ifdef SINGULAR_4_1
6462  else if (t==CRING_CMD)
6463  {
6464    sleftv tmp;
6465    sleftv n;
6466    memset(&n,0,sizeof(n));
6467    n.name=ring_name;
6468    if (iiDeclCommand(&tmp,&n,myynest,CRING_CMD,&IDROOT)) return TRUE;
6469    if (iiAssign(&tmp,arg)) return TRUE;
6470    //Print("create %s\n",r->Name());
6471    //Print("from %s(%d)\n",Tok2Cmdname(arg->Typ()),arg->Typ());
6472    return FALSE;
6473  }
6474  #endif
6475  //Print("create %s\n",r->Name());
6476  //Print("from %s(%d)\n",Tok2Cmdname(arg->Typ()),arg->Typ());
6477  return TRUE;// not handled -> error for now
6478}
6479
6480static void iiReportTypes(int nr,int t,const short *T)
6481{
6482  char *buf=(char*)omAlloc(250);
6483  buf[0]='\0';
6484  if (nr==0)
6485    sprintf(buf,"wrong length of parameters(%d), expected ",t);
6486  else
6487    sprintf(buf,"par. %d is of type `%s`, expected ",nr,Tok2Cmdname(t));
6488  for(int i=1;i<=T[0];i++)
6489  {
6490    strcat(buf,"`");
6491    strcat(buf,Tok2Cmdname(T[i]));
6492    strcat(buf,"`");
6493    if (i<T[0]) strcat(buf,",");
6494  }
6495  WerrorS(buf);
6496}
6497
6498BOOLEAN iiCheckTypes(leftv args, const short *type_list, int report)
6499{
6500  if (args==NULL)
6501  {
6502    if (type_list[0]==0) return TRUE;
6503    else
6504    {
6505      if (report) WerrorS("no arguments expected");
6506      return FALSE;
6507    }
6508  }
6509  int l=args->listLength();
6510  if (l!=(int)type_list[0])
6511  {
6512    if (report) iiReportTypes(0,l,type_list);
6513    return FALSE;
6514  }
6515  for(int i=1;i<=l;i++,args=args->next)
6516  {
6517    short t=type_list[i];
6518    if (t!=ANY_TYPE)
6519    {
6520      if (((t==IDHDL)&&(args->rtyp!=IDHDL))
6521      || (t!=args->Typ()))
6522      {
6523        if (report) iiReportTypes(i,args->Typ(),type_list);
6524        return FALSE;
6525      }
6526    }
6527  }
6528  return TRUE;
6529}
Note: See TracBrowser for help on using the repository browser.