source: git/Singular/ipshell.cc @ 8d1432e

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