source: git/Singular/ipshell.cc @ 417a91a

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