source: git/Singular/ipshell.cc @ 23f2f57

spielwiese
Last change on this file since 23f2f57 was 23f2f57, checked in by Hans Schönemann <hannes@…>, 19 years ago
*hannes; local rings again git-svn-id: file:///usr/local/Singular/svn/trunk@7909 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 105.6 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/* $Id: ipshell.cc,v 1.104 2005-04-26 18:21:26 Singular Exp $ */
5/*
6* ABSTRACT:
7*/
8
9//#include <stdlib.h>
10#include <stdio.h>
11#include <string.h>
12#include <ctype.h>
13#include <math.h>
14
15#include "mod2.h"
16#include "tok.h"
17#include "ipid.h"
18#include "intvec.h"
19#include "omalloc.h"
20#include "febase.h"
21#include "polys.h"
22#include "ideals.h"
23#include "matpol.h"
24#include "kstd1.h"
25#include "ring.h"
26#include "subexpr.h"
27#include "maps.h"
28#include "syz.h"
29#include "numbers.h"
30#include "lists.h"
31#include "attrib.h"
32#include "ipconv.h"
33#include "silink.h"
34#include "stairc.h"
35#include "weight.h"
36#include "semic.h"
37#include "splist.h"
38#include "spectrum.h"
39#include "gnumpfl.h"
40#include "mpr_base.h"
41#include "ffields.h"
42#include "clapsing.h"
43#include "hutil.h"
44#include "ipshell.h"
45#ifdef HAVE_FACTORY
46#define SI_DONT_HAVE_GLOBAL_VARS
47#include <factory.h>
48#endif
49
50// define this if you want to use the fast_map routine for mapping ideals
51#define FAST_MAP
52
53#ifdef FAST_MAP
54#include "fast_maps.h"
55#endif
56
57leftv iiCurrArgs=NULL;
58int  traceit = 0;
59char *lastreserved=NULL;
60
61int  myynest = -1;
62
63static BOOLEAN iiNoKeepRing=TRUE;
64
65/*0 implementation*/
66
67char * Tok2Cmdname(int tok)
68{
69  int i = 0;
70  if (tok < 0)
71  {
72    return cmds[0].name;
73  }
74  if (tok==ANY_TYPE) return "any_type";
75  if (tok==NONE) return "nothing";
76  //if (tok==IFBREAK) return "if_break";
77  //if (tok==VECTOR_FROM_POLYS) return "vector_from_polys";
78  //if (tok==ORDER_VECTOR) return "ordering";
79  //if (tok==REF_VAR) return "ref";
80  //if (tok==OBJECT) return "object";
81  //if (tok==PRINT_EXPR) return "print_expr";
82  if (tok==IDHDL) return "identifier";
83  while (cmds[i].tokval!=0)
84  {
85    if ((cmds[i].tokval == tok)&&(cmds[i].alias==0))
86    {
87      return cmds[i].name;
88    }
89    i++;
90  }
91  return cmds[0].name;
92}
93
94char * iiTwoOps(int t)
95{
96  if (t<127)
97  {
98    static char ch[2];
99    switch (t)
100    {
101      case '&':
102        return "and";
103      case '|':
104        return "or";
105      default:
106        ch[0]=t;
107        ch[1]='\0';
108        return ch;
109    }
110  }
111  switch (t)
112  {
113    case COLONCOLON:  return "::";
114    case DOTDOT:      return "..";
115    //case PLUSEQUAL:   return "+=";
116    //case MINUSEQUAL:  return "-=";
117    case MINUSMINUS:  return "--";
118    case PLUSPLUS:    return "++";
119    case EQUAL_EQUAL: return "==";
120    case LE:          return "<=";
121    case GE:          return ">=";
122    case NOTEQUAL:    return "<>";
123    default:          return Tok2Cmdname(t);
124  }
125}
126
127static void list1(char* s, idhdl h,BOOLEAN c, BOOLEAN fullname)
128{
129  char buffer[22];
130  int l;
131  char buf2[128];
132
133  if(fullname) sprintf(buf2, "%s::%s", "", IDID(h));
134  else sprintf(buf2, "%s", IDID(h));
135
136  Print("%s%-20.20s [%d]  ",s,buf2,IDLEV(h));
137  if (h == currRingHdl) PrintS("*");
138  PrintS(Tok2Cmdname((int)IDTYP(h)));
139
140  ipListFlag(h);
141  switch(IDTYP(h))
142  {
143    case INT_CMD:   Print(" %d",IDINT(h)); break;
144    case INTVEC_CMD:Print(" (%d)",IDINTVEC(h)->length()); break;
145    case INTMAT_CMD:Print(" %d x %d",IDINTVEC(h)->rows(),IDINTVEC(h)->cols());
146                    break;
147    case POLY_CMD:
148    case VECTOR_CMD:if (c)
149                    {
150                      PrintS(" ");wrp(IDPOLY(h));
151                      if(IDPOLY(h) != NULL)
152                      {
153                        Print(", %d monomial(s)",pLength(IDPOLY(h)));
154                      }
155                    }
156                    break;
157    case MODUL_CMD: Print(", rk %d", IDIDEAL(h)->rank);
158    case IDEAL_CMD: Print(", %u generator(s)",
159                    IDELEMS(IDIDEAL(h)),IDIDEAL(h)->rank); break;
160    case MAP_CMD:
161                    Print(" from %s",IDMAP(h)->preimage); break;
162    case MATRIX_CMD:Print(" %u x %u"
163                      ,MATROWS(IDMATRIX(h))
164                      ,MATCOLS(IDMATRIX(h))
165                    );
166                    break;
167    case PACKAGE_CMD:
168                    PrintS(" (");
169                    switch (IDPACKAGE(h)->language)
170                    {
171                        case LANG_SINGULAR: PrintS("S"); break;
172                        case LANG_C:        PrintS("C"); break;
173                        case LANG_TOP:      PrintS("T"); break;
174                        case LANG_NONE:     PrintS("N"); break;
175                        default:            PrintS("U");
176                    }
177                    if(IDPACKAGE(h)->libname!=NULL)
178                      Print(",%s", IDPACKAGE(h)->libname);
179                    PrintS(")");
180                    break;
181    case PROC_CMD: if(strlen(IDPROC(h)->libname)>0)
182                     Print(" from %s",IDPROC(h)->libname);
183                   if(IDPROC(h)->is_static)
184                     PrintS(" (static)");
185                   break;
186    case STRING_CMD:
187                   {
188                     char *s;
189                     l=strlen(IDSTRING(h));
190                     memset(buffer,0,22);
191                     strncpy(buffer,IDSTRING(h),min(l,20));
192                     if ((s=strchr(buffer,'\n'))!=NULL)
193                     {
194                       *s='\0';
195                     }
196                     PrintS(" ");
197                     PrintS(buffer);
198                     if((s!=NULL) ||(l>20))
199                     {
200                       Print("..., %d char(s)",l);
201                     }
202                     break;
203                   }
204    case LIST_CMD: Print(", size: %d",IDLIST(h)->nr+1);
205                   break;
206    case QRING_CMD:
207    case RING_CMD:
208                   if ((IDRING(h)==currRing) && (currRingHdl!=h))
209                     PrintS("(*)"); /* this is an alias to currRing */
210#ifdef RDEBUG
211                   if (traceit &TRACE_SHOW_RINGS)
212                     Print(" <%x>",IDRING(h));
213#endif
214                   break;
215    /*default:     break;*/
216  }
217  PrintLn();
218}
219
220void type_cmd(idhdl h)
221{
222  BOOLEAN oldShortOut = FALSE;
223
224  if (currRing != NULL)
225  {
226    oldShortOut = currRing->ShortOut;
227    currRing->ShortOut = 1;
228  }
229  list1("// ",h,FALSE,FALSE);
230  if (IDTYP(h)!=INT_CMD)
231  {
232    sleftv expr;
233    memset(&expr,0,sizeof(expr));
234    expr.rtyp=IDHDL;
235    expr.name=IDID(h);
236    expr.data=(void *)h;
237    expr.Print();
238  }
239  if (currRing != NULL)
240    currRing->ShortOut = oldShortOut;
241}
242
243static void killlocals0(int v, idhdl * localhdl, const ring r)
244{
245  idhdl h = *localhdl;
246  while (h!=NULL)
247  {
248    int vv;
249    //Print("consider %s, lev: %d:",IDID(h),IDLEV(h));
250    if ((vv=IDLEV(h))>0)
251    {
252      if (vv < v)
253      {
254        if (iiNoKeepRing)
255        {
256          //PrintS(" break\n");
257          return;
258        }
259        h = IDNEXT(h);
260        //PrintLn();
261      }
262      else if (vv >= v)
263      {
264        idhdl nexth = IDNEXT(h);
265        killhdl2(h,localhdl,r);
266        h = nexth;
267        //PrintS("kill\n");
268      }
269    }
270    else
271    {
272      h = IDNEXT(h);
273      //PrintLn();
274    }
275  }
276}
277#ifndef HAVE_NS
278void killlocals(int v)
279{
280  killlocals0(v,&IDROOT,currRing);
281
282  if ((iiRETURNEXPR_len > myynest)
283  && ((iiRETURNEXPR[myynest].Typ()==RING_CMD)
284    || (iiRETURNEXPR[myynest].Typ()==QRING_CMD)))
285  {
286    leftv h=&iiRETURNEXPR[myynest];
287    killlocals0(v,&(((ring)h->data)->idroot),(ring)h->data);
288  }
289
290  idhdl sh=currRingHdl;
291  ring sr=currRing;
292  BOOLEAN changed=FALSE;
293  idhdl h = IDROOT;
294
295//  Print("killlocals in %s\n",IDID(currPackHdl));
296  while (h!=NULL)
297  {
298    if (((IDTYP(h)==QRING_CMD) || (IDTYP(h) == RING_CMD))
299    && (IDRING(h)->idroot!=NULL))
300    {
301      if (IDRING(h)!=currRing) {changed=TRUE;rSetHdl(h);}
302      killlocals0(v,&(IDRING(h)->idroot),IDRING(h));
303    }
304    else if (IDTYP(h) == PACKAGE_CMD)
305    {
306      killlocals0(v,&(IDPACKAGE(h)->idroot),IDRING(h));
307    }
308    h = IDNEXT(h);
309  }
310  if (changed)
311  {
312    currRing=NULL;
313    currRingHdl=NULL;
314    if (sh!=NULL) rSetHdl(sh);
315    else if (sr!=NULL)
316    {
317      sh=rFindHdl(sr,NULL,NULL);
318      rSetHdl(sh);
319    }
320  }
321
322  if (myynest<=1) iiNoKeepRing=TRUE;
323  //Print("end killlocals  >= %d\n",v);
324  //listall();
325}
326#endif
327#ifdef HAVE_NS
328void killlocals_rec(idhdl *root,int v, ring r)
329{
330  idhdl h=*root;
331  while (h!=NULL)
332  {
333    if (IDLEV(h)>=v)
334    {
335//      Print("kill %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
336      idhdl n=IDNEXT(h);
337      killhdl2(h,root,r);
338      h=n;
339    }
340    else if (IDTYP(h)==PACKAGE_CMD)
341    {
342 //     Print("into pack %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
343      if (IDPACKAGE(h)!=basePack)
344        killlocals_rec(&(IDRING(h)->idroot),v,r);
345      h=IDNEXT(h);
346    }
347    else if ((IDTYP(h)==RING_CMD)
348    ||(IDTYP(h)==QRING_CMD))
349    {
350      if (IDRING(h)->idroot!=NULL)
351      {
352  //    Print("into ring %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
353        killlocals_rec(&(IDRING(h)->idroot),v,IDRING(h));
354      }
355      h=IDNEXT(h);
356    }
357    else
358    {
359//      Print("skip %s lev %d for lev %d\n",IDID(h),IDLEV(h),v);
360      h=IDNEXT(h);
361    }
362  }
363}
364BOOLEAN killlocals_list(int v, lists L)
365{
366  if (L==NULL) return FALSE;
367  BOOLEAN changed=FALSE;
368  int n=L->nr;
369  for(;n>=0;n--)
370  {
371    leftv h=&(L->m[n]);
372    void *d=h->data;
373    if (((h->rtyp==RING_CMD) || (h->rtyp==QRING_CMD))
374    && (((ring)d)->idroot!=NULL))
375    {
376      if (d!=currRing) {changed=TRUE;rChangeCurrRing((ring)d);}
377      killlocals0(v,&(((ring)h->data)->idroot),(ring)h->data);
378    }
379    else if (h->rtyp==LIST_CMD)
380      changed|=killlocals_list(v,(lists)d);
381  }
382  return changed;
383}
384void killlocals(int v)
385{
386  BOOLEAN changed=FALSE;
387  idhdl sh=currRingHdl;
388  ring cr=currRing;
389  if (sh!=NULL) changed=((IDLEV(sh)<v) || (IDRING(sh)->ref>0));
390  //if (changed) Print("currRing=%s(%x), lev=%d,ref=%d\n",IDID(sh),IDRING(sh),IDLEV(sh),IDRING(sh)->ref);
391
392  killlocals_rec(&(basePack->idroot),v,currRing);
393
394  if (iiRETURNEXPR_len > myynest)
395  {
396    int t=iiRETURNEXPR[myynest].Typ();
397    if ((/*iiRETURNEXPR[myynest].Typ()*/ t==RING_CMD)
398    || (/*iiRETURNEXPR[myynest].Typ()*/ t==QRING_CMD))
399    {
400      leftv h=&iiRETURNEXPR[myynest];
401      if (((ring)h->data)->idroot!=NULL)
402        killlocals0(v,&(((ring)h->data)->idroot),(ring)h->data);
403    }
404    else if (/*iiRETURNEXPR[myynest].Typ()*/ t==LIST_CMD)
405    {
406      leftv h=&iiRETURNEXPR[myynest];
407      changed |=killlocals_list(v,(lists)h->data);
408    }
409  }
410  if (changed)
411  {
412    currRingHdl=rFindHdl(cr,NULL,NULL);
413    if (currRingHdl==NULL) 
414      currRing=NULL;
415    else
416      rChangeCurrRing(cr);
417  }
418
419  if (myynest<=1) iiNoKeepRing=TRUE;
420  //Print("end killlocals  >= %d\n",v);
421  //listall();
422}
423#endif
424
425void list_cmd(int typ, const char* what, char *prefix,BOOLEAN iterate, BOOLEAN fullname)
426{
427  idhdl h,start;
428  BOOLEAN all = typ<0;
429  BOOLEAN really_all=FALSE;
430  BOOLEAN do_packages=FALSE;
431
432  if ( typ == -1 ) do_packages=TRUE;
433  if ( typ==0 )
434  {
435    if (strcmp(what,"all")==0)
436    {
437      really_all=TRUE;
438#ifdef HAVE_NS
439      h=basePack->idroot;
440#else
441      h=IDROOT;
442#endif
443    }
444    else
445    {
446      h = ggetid(what);
447      if (h!=NULL)
448      {
449        if (iterate) list1(prefix,h,TRUE,fullname);
450        if ((IDTYP(h)==RING_CMD)
451            || (IDTYP(h)==QRING_CMD)
452#ifdef HAVE_NS
453            //|| (IDTYP(h)==PACKE_CMD)
454#endif
455        )
456        {
457          h=IDRING(h)->idroot;
458        }
459        else if((IDTYP(h)==PACKAGE_CMD) || (IDTYP(h)==POINTER_CMD))
460        {
461          Print("list_cmd:package or pointer\n");
462          all=TRUE;typ=PROC_CMD;fullname=TRUE;really_all=TRUE;
463          h=IDPACKAGE(h)->idroot;
464        }
465        else
466          return;
467      }
468      else
469      {
470        Werror("%s is undefined",what);
471        return;
472      }
473    }
474    all=TRUE;
475  }
476  else if (RingDependend(typ))
477  {
478    h = currRing->idroot;
479  }
480  else
481    h = IDROOT;
482  start=h;
483  while (h!=NULL)
484  {
485    if ((all && (IDTYP(h)!=PROC_CMD) &&(IDTYP(h)!=PACKAGE_CMD))
486    || (typ == IDTYP(h))
487    || ((IDTYP(h)==QRING_CMD) && (typ==RING_CMD)))
488    {
489      list1(prefix,h,start==currRingHdl, fullname);
490      if (((IDTYP(h)==RING_CMD)||(IDTYP(h)==QRING_CMD))
491        && (really_all || (all && (h==currRingHdl)))
492        && ((IDLEV(h)==0)||(IDLEV(h)==myynest)))
493      {
494        list_cmd(0,IDID(h),"//      ",FALSE);
495      }
496#ifdef HAVE_NS
497      if (IDTYP(h)==PACKAGE_CMD && really_all)
498      {
499        list_cmd(0,IDID(h),"//      ",FALSE);
500      }
501#endif /* HAVE_NS */
502    }
503    h = IDNEXT(h);
504  }
505}
506
507void test_cmd(int i)
508{
509  int ii=(char)i;
510
511  if (i == (-32))
512  {
513    test = 0;
514  }
515  else
516  {
517    if (i<0)
518    {
519      ii= -i;
520      if (Sy_bit(ii) & kOptions)
521      {
522        Warn("Gerhard, use the option command");
523        test &= ~Sy_bit(ii);
524      }
525      else if (Sy_bit(ii) & validOpts)
526        test &= ~Sy_bit(ii);
527    }
528    else if (i<32)
529    {
530      if (Sy_bit(ii) & kOptions)
531      {
532        Warn("Gerhard, use the option command");
533        test |= Sy_bit(ii);
534      }
535      else if (Sy_bit(ii) & validOpts)
536        test |= Sy_bit(ii);
537    }
538  }
539}
540
541int exprlist_length(leftv v)
542{
543  int rc = 0;
544  while (v!=NULL)
545  {
546    switch (v->Typ())
547    {
548      case INT_CMD:
549      case POLY_CMD:
550      case VECTOR_CMD:
551      case NUMBER_CMD:
552        rc++;
553        break;
554      case INTVEC_CMD:
555      case INTMAT_CMD:
556        rc += ((intvec *)(v->Data()))->length();
557        break;
558      case MATRIX_CMD:
559      case IDEAL_CMD:
560      case MODUL_CMD:
561        {
562          matrix mm = (matrix)(v->Data());
563          rc += mm->rows() * mm->cols();
564        }
565        break;
566      case LIST_CMD:
567        rc+=((lists)v->Data())->nr+1;
568        break;
569      default:
570        rc++;
571    }
572    v = v->next;
573  }
574  return rc;
575}
576
577int IsPrime(int p)  /* brute force !!!! */
578{
579  int i,j;
580  if      (p == 0)    return 0;
581  else if (p == 1)    return 1/*1*/;
582  else if (p == 2)    return p;
583  else if (p < 0)     return (-IsPrime(-p));
584  else if (!(p & 1)) return IsPrime(p-1);
585#ifdef HAVE_FACTORY
586  else if (p<=32749) // max. small prime in factory
587  {
588    int a=0;
589    int e=cf_getNumSmallPrimes()-1;
590    i=e/2;
591    do
592    {
593      if (p==(j=cf_getSmallPrime(i))) return p;
594      if (p<j) e=i-1;
595      else     a=i+1;
596      i=a+(e-a)/2;
597    } while ( a<= e);
598    if (p>j) return j;
599    else     return cf_getSmallPrime(i-1);
600  }
601#endif
602#ifdef HAVE_FACTORY
603  int end_i=cf_getNumSmallPrimes()-1;
604#else
605  int end_i=p/2;
606#endif 
607  int end_p=(int)sqrt((double)p);
608restart:
609  for (i=0; i<end_i; i++)
610  {
611#ifdef HAVE_FACTORY
612    j=cf_getSmallPrime(i);
613#else
614    if (i==0) j=2;
615    else j=2*i-1;
616#endif   
617    if ((p%j) == 0)
618    {
619    #ifdef HAVE_FACTORY
620      if (p<=32751) return IsPrime(p-2);
621    #endif 
622      p-=2;
623      goto restart;
624    }
625    if (j > end_p) return p;
626  }
627  return p;
628}
629
630BOOLEAN iiWRITE(leftv res,leftv v)
631{
632  sleftv vf;
633  if (iiConvert(v->Typ(),LINK_CMD,iiTestConvert(v->Typ(),LINK_CMD),v,&vf))
634  {
635    WerrorS("link expected");
636    return TRUE;
637  }
638  si_link l=(si_link)vf.Data();
639  if (vf.next == NULL)
640  {
641    WerrorS("write: need at least two arguments");
642    return TRUE;
643  }
644
645  BOOLEAN b=slWrite(l,vf.next); /* iiConvert preserves next */
646  if (b)
647  {
648    const char *s;
649    if ((l!=NULL)&&(l->name!=NULL)) s=l->name;
650    else                            s=sNoName;
651    Werror("cannot write to %s",s);
652  }
653  vf.CleanUp();
654  return b;
655}
656
657leftv iiMap(map theMap, char * what)
658{
659  idhdl w,r;
660  leftv v;
661  int i;
662  nMapFunc nMap;
663
664  r=IDROOT->get(theMap->preimage,myynest);
665#ifdef HAVE_NS
666  if ((currPack!=basePack)
667  &&((r==NULL) || ((r->typ != RING_CMD) && (r->typ != QRING_CMD))))
668    r=basePack->idroot->get(theMap->preimage,myynest);
669  if ((r==NULL) && (currRingHdl!=NULL)
670  && (strcmp(theMap->preimage,IDID(currRingHdl))==0))
671  {
672    r=currRingHdl;
673  }
674#endif /* HAVE_NS */
675  if ((r!=NULL) && ((r->typ == RING_CMD) || (r->typ== QRING_CMD)))
676  {
677    //if ((nMap=nSetMap(rInternalChar(IDRING(r)),
678    //             IDRING(r)->parameter,
679    //             rPar(IDRING(r)),
680    //             IDRING(r)->minpoly)))
681    if ((nMap=nSetMap(IDRING(r)))==NULL)
682    {
683      if (rEqual(IDRING(r),currRing))
684      {
685        nMap=nCopy;
686      }
687      else
688      {
689        Werror("can not map from ground field of %s to current ground field",
690          theMap->preimage);
691        return NULL;
692      }
693    }
694    if (IDELEMS(theMap)<IDRING(r)->N)
695    {
696      theMap->m=(polyset)omReallocSize((ADDRESS)theMap->m,
697                                 IDELEMS(theMap)*sizeof(poly),
698                                 (IDRING(r)->N)*sizeof(poly));
699      for(i=IDELEMS(theMap);i<IDRING(r)->N;i++)
700        theMap->m[i]=NULL;
701      IDELEMS(theMap)=IDRING(r)->N;
702    }
703    if (what==NULL)
704    {
705      WerrorS("argument of a map must have a name");
706    }
707    else if ((w=IDRING(r)->idroot->get(what,myynest))!=NULL)
708    {
709      v=(leftv)omAlloc0Bin(sleftv_bin);
710      sleftv tmpW;
711      memset(&tmpW,0,sizeof(sleftv));
712      tmpW.rtyp=IDTYP(w);
713      tmpW.data=IDDATA(w);
714      #ifdef FAST_MAP
715      if ((tmpW.rtyp==IDEAL_CMD) && (nMap==nCopy)
716      #ifdef HAVE_PLURAL
717      && (!rIsPluralRing(currRing))
718      #endif
719      )
720      {
721        v->rtyp=IDEAL_CMD;
722        v->data=fast_map(IDIDEAL(w), IDRING(r), (ideal)theMap, currRing);
723      }
724      else
725      #endif
726      if (maApplyFetch(MAP_CMD,theMap,v,&tmpW,IDRING(r),NULL,NULL,0,nMap))
727      {
728        Werror("cannot map %s(%d)",Tok2Cmdname(w->typ),w->typ);
729        omFreeBin((ADDRESS)v, sleftv_bin);
730        return NULL;
731      }
732      return v;
733    }
734    else
735    {
736      Werror("%s undefined in %s",what,theMap->preimage);
737    }
738  }
739  else
740  {
741    Werror("cannot find preimage %s",theMap->preimage);
742  }
743  return NULL;
744}
745
746#ifdef OLD_RES
747void  iiMakeResolv(resolvente r, int length, int rlen, char * name, int typ0,
748                   intvec ** weights)
749{
750  lists L=liMakeResolv(r,length,rlen,typ0,weights);
751  int i=0;
752  idhdl h;
753  char * s=(char *)omAlloc(strlen(name)+5);
754
755  while (i<=L->nr)
756  {
757    sprintf(s,"%s(%d)",name,i+1);
758    if (i==0)
759      h=enterid(s,myynest,typ0,&(currRing->idroot), FALSE);
760    else
761      h=enterid(s,myynest,MODUL_CMD,&(currRing->idroot), FALSE);
762    if (h!=NULL)
763    {
764      h->data.uideal=(ideal)L->m[i].data;
765      h->attribute=L->m[i].attribute;
766      if (BVERBOSE(V_DEF_RES))
767        Print("//defining: %s as %d-th syzygy module\n",s,i+1);
768    }
769    else
770    {
771      idDelete((ideal *)&(L->m[i].data));
772      Warn("cannot define %s",s);
773    }
774    //L->m[i].data=NULL;
775    //L->m[i].rtyp=0;
776    //L->m[i].attribute=NULL;
777    i++;
778  }
779  omFreeSize((ADDRESS)L->m,(L->nr+1)*sizeof(sleftv));
780  omFreeBin((ADDRESS)L, slists_bin);
781  omFreeSize((ADDRESS)s,strlen(name)+5);
782}
783#endif
784
785//resolvente iiFindRes(char * name, int * len, int *typ0)
786//{
787//  char *s=(char *)omAlloc(strlen(name)+5);
788//  int i=-1;
789//  resolvente r;
790//  idhdl h;
791//
792//  do
793//  {
794//    i++;
795//    sprintf(s,"%s(%d)",name,i+1);
796//    h=currRing->idroot->get(s,myynest);
797//  } while (h!=NULL);
798//  *len=i-1;
799//  if (*len<=0)
800//  {
801//    Werror("no objects %s(1),.. found",name);
802//    omFreeSize((ADDRESS)s,strlen(name)+5);
803//    return NULL;
804//  }
805//  r=(ideal *)omAlloc(/*(len+1)*/ i*sizeof(ideal));
806//  memset(r,0,(*len)*sizeof(ideal));
807//  i=-1;
808//  *typ0=MODUL_CMD;
809//  while (i<(*len))
810//  {
811//    i++;
812//    sprintf(s,"%s(%d)",name,i+1);
813//    h=currRing->idroot->get(s,myynest);
814//    if (h->typ != MODUL_CMD)
815//    {
816//      if ((i!=0) || (h->typ!=IDEAL_CMD))
817//      {
818//        Werror("%s is not of type module",s);
819//        omFreeSize((ADDRESS)r,(*len)*sizeof(ideal));
820//        omFreeSize((ADDRESS)s,strlen(name)+5);
821//        return NULL;
822//      }
823//      *typ0=IDEAL_CMD;
824//    }
825//    if ((i>0) && (idIs0(r[i-1])))
826//    {
827//      *len=i-1;
828//      break;
829//    }
830//    r[i]=IDIDEAL(h);
831//  }
832//  omFreeSize((ADDRESS)s,strlen(name)+5);
833//  return r;
834//}
835
836static resolvente iiCopyRes(resolvente r, int l)
837{
838  int i;
839  resolvente res=(ideal *)omAlloc0((l+1)*sizeof(ideal));
840
841  for (i=0; i<l; i++)
842    res[i]=idCopy(r[i]);
843  return res;
844}
845
846BOOLEAN jjMINRES(leftv res, leftv v)
847{
848  int len=0;
849  int typ0;
850  resolvente rr=liFindRes((lists)v->Data(),&len,&typ0);
851  if (rr==NULL) return TRUE;
852  resolvente r=iiCopyRes(rr,len);
853
854  syMinimizeResolvente(r,len,0);
855  omFreeSize((ADDRESS)rr,len*sizeof(ideal));
856  len++;
857  res->data=(char *)liMakeResolv(r,len,-1,typ0,NULL);
858  return FALSE;
859}
860
861BOOLEAN jjBETTI(leftv res, leftv u)
862{
863  sleftv tmp;
864  memset(&tmp,0,sizeof(tmp));
865  tmp.rtyp=INT_CMD;
866  tmp.data=(void *)1; 
867  if ((u->Typ()==IDEAL_CMD)
868  || (u->Typ()==MODUL_CMD))
869    return jjBETTI2_ID(res,u,&tmp);
870  else
871    return jjBETTI2(res,u,&tmp);
872}
873
874BOOLEAN jjBETTI2_ID(leftv res, leftv u, leftv v)
875{
876  lists l=(lists) omAllocBin(slists_bin);
877  l->Init(1);
878  l->m[0].rtyp=u->Typ();
879  l->m[0].data=u->Data();
880  l->m[0].attribute=u->attribute;
881  sleftv tmp2;
882  memset(&tmp2,0,sizeof(tmp2));
883  tmp2.rtyp=LIST_CMD;
884  tmp2.data=(void *)l;
885  BOOLEAN r=jjBETTI2(res,&tmp2,v);
886  l->m[0].data=NULL;
887  l->m[0].attribute=NULL;
888  l->m[0].rtyp=DEF_CMD;
889  l->Clean();
890  return r;
891}
892
893BOOLEAN jjBETTI2(leftv res, leftv u, leftv v)
894{
895  resolvente r;
896  int len;
897  int reg,typ0;
898  lists l=(lists)u->Data();
899
900  intvec *weights=NULL;
901  int add_row_shift=0;
902  intvec *ww=(intvec *)atGet(&(l->m[0]),"isHomog",INTVEC_CMD); 
903  if (ww!=NULL)
904  {
905     weights=ivCopy(ww);
906     add_row_shift = ww->min_in();
907     (*weights) -= add_row_shift;
908  }
909  //Print("attr:%x\n",weights);
910
911  r=liFindRes(l,&len,&typ0);
912  if (r==NULL) return TRUE;
913  res->data=(char *)syBetti(r,len,&reg,weights,(int)v->Data());
914  omFreeSize((ADDRESS)r,(len)*sizeof(ideal));
915  atSet(res,omStrDup("rowShift"),(void*)add_row_shift,INT_CMD);
916  if (weights!=NULL) delete weights;
917  return FALSE;
918}
919
920int iiRegularity(lists L)
921{
922  int len,reg,typ0;
923
924  resolvente r=liFindRes(L,&len,&typ0);
925
926  if (r==NULL)
927    return -2;
928  intvec * dummy=syBetti(r,len,&reg);
929  omFreeSize((ADDRESS)r,len*sizeof(ideal));
930  delete dummy;
931  return reg+1;
932}
933
934BOOLEAN iiDebugMarker=TRUE;
935#define BREAK_LINE_LENGTH 80
936void iiDebug()
937{
938  Print("\n-- break point in %s --\n",VoiceName());
939  if (iiDebugMarker) VoiceBackTrack();
940  char * s;
941  iiDebugMarker=FALSE;
942  s = (char *)omAlloc(BREAK_LINE_LENGTH+4);
943  loop
944  {
945    memset(s,0,80);
946    fe_fgets_stdin("",s,BREAK_LINE_LENGTH);
947    if (s[BREAK_LINE_LENGTH-1]!='\0')
948    {
949      Print("line too long, max is %d chars\n",BREAK_LINE_LENGTH);
950    }
951    else
952      break;
953  }
954  if (*s=='\n')
955  {
956    iiDebugMarker=TRUE;
957  }
958#if MDEBUG
959  else if(strncmp(s,"cont;",5)==0)
960  {
961    iiDebugMarker=TRUE;
962  }
963#endif /* MDEBUG */
964  else
965  {
966    strcat( s, "\n;~\n");
967    newBuffer(s,BT_execute);
968  }
969}
970
971lists scIndIndset(ideal S, BOOLEAN all, ideal Q)
972{
973  int i;
974  indset save;
975  lists res=(lists)omAlloc0Bin(slists_bin);
976
977  hexist = hInit(S, Q, &hNexist);
978  if ((hNexist == 0) || (hisModule!=0))
979  {
980    res->Init(0);
981    return res;
982  }
983  save = ISet = (indset)omAlloc0Bin(indlist_bin);
984  hMu = 0;
985  hwork = (scfmon)omAlloc(hNexist * sizeof(scmon));
986  hvar = (varset)omAlloc((pVariables + 1) * sizeof(int));
987  hpure = (scmon)omAlloc((1 + (pVariables * pVariables)) * sizeof(Exponent_t));
988  hrad = hexist;
989  hNrad = hNexist;
990  radmem = hCreate(pVariables - 1);
991  hCo = pVariables + 1;
992  hNvar = pVariables;
993  hRadical(hrad, &hNrad, hNvar);
994  hSupp(hrad, hNrad, hvar, &hNvar);
995  if (hNvar)
996  {
997    hCo = hNvar;
998    memset(hpure, 0, (pVariables + 1) * sizeof(Exponent_t));
999    hPure(hrad, 0, &hNrad, hvar, hNvar, hpure, &hNpure);
1000    hLexR(hrad, hNrad, hvar, hNvar);
1001    hDimSolve(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1002  }
1003  if (hCo && (hCo < pVariables))
1004  {
1005    hIndMult(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1006  }
1007  if (hMu!=0)
1008  {
1009    ISet = save;
1010    hMu2 = 0;
1011    if (all && (hCo+1 < pVariables))
1012    {
1013      JSet = (indset)omAlloc0Bin(indlist_bin);
1014      hIndAllMult(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1015      i=hMu+hMu2;
1016      res->Init(i);
1017      if (hMu2 == 0)
1018      {
1019        omFreeBin((ADDRESS)JSet, indlist_bin);
1020      }
1021    }
1022    else
1023    {
1024      res->Init(hMu);
1025    }
1026    for (i=0;i<hMu;i++)
1027    {
1028      res->m[i].data = (void *)save->set;
1029      res->m[i].rtyp = INTVEC_CMD;
1030      ISet = save;
1031      save = save->nx;
1032      omFreeBin((ADDRESS)ISet, indlist_bin);
1033    }
1034    omFreeBin((ADDRESS)save, indlist_bin);
1035    if (hMu2 != 0)
1036    {
1037      save = JSet;
1038      for (i=hMu;i<hMu+hMu2;i++)
1039      {
1040        res->m[i].data = (void *)save->set;
1041        res->m[i].rtyp = INTVEC_CMD;
1042        JSet = save;
1043        save = save->nx;
1044        omFreeBin((ADDRESS)JSet, indlist_bin);
1045      }
1046      omFreeBin((ADDRESS)save, indlist_bin);
1047    }
1048  }
1049  else
1050  {
1051    res->Init(0);
1052    omFreeBin((ADDRESS)ISet,  indlist_bin);
1053  }
1054  hKill(radmem, pVariables - 1);
1055  omFreeSize((ADDRESS)hpure, (1 + (pVariables * pVariables)) * sizeof(Exponent_t));
1056  omFreeSize((ADDRESS)hvar, (pVariables + 1) * sizeof(int));
1057  omFreeSize((ADDRESS)hwork, hNexist * sizeof(scmon));
1058  hDelete(hexist, hNexist);
1059  return res;
1060}
1061
1062int iiDeclCommand(leftv sy, leftv name, int lev,int t, idhdl* root,BOOLEAN isring, BOOLEAN init_b)
1063{
1064  BOOLEAN res=FALSE;
1065  char *id = name->name;
1066
1067  memset(sy,0,sizeof(sleftv));
1068  if ((name->name==NULL)||(isdigit(name->name[0])))
1069  {
1070    WerrorS("object to declare is not a name");
1071    res=TRUE;
1072  }
1073  else
1074  {
1075    //if (name->rtyp!=0)
1076    //{
1077    //  Warn("`%s` is already in use",name->name);
1078    //}
1079    {
1080      sy->data = (char *)enterid(id,lev,t,root,init_b);
1081    }
1082    if (sy->data!=NULL)
1083    {
1084      sy->rtyp=IDHDL;
1085      currid=sy->name=IDID((idhdl)sy->data);
1086      // name->name=NULL; /* used in enterid */
1087      //sy->e = NULL;
1088      if (name->next!=NULL)
1089      {
1090        sy->next=(leftv)omAllocBin(sleftv_bin);
1091        res=iiDeclCommand(sy->next,name->next,lev,t,root, isring);
1092      }
1093    }
1094    else res=TRUE;
1095  }
1096  name->CleanUp();
1097  return res;
1098}
1099
1100BOOLEAN iiParameter(leftv p)
1101{
1102  if (iiCurrArgs==NULL)
1103  {
1104    if (strcmp(p->name,"#")==0) return FALSE;
1105    Werror("not enough arguments for proc %s",VoiceName());
1106    p->CleanUp();
1107    return TRUE;
1108  }
1109  leftv h=iiCurrArgs;
1110  if (strcmp(p->name,"#")==0)
1111  {
1112    iiCurrArgs=NULL;
1113  }
1114  else
1115  {
1116    iiCurrArgs=h->next;
1117    h->next=NULL;
1118  }
1119  BOOLEAN res=iiAssign(p,h);
1120  omFreeBin((ADDRESS)h, sleftv_bin);
1121  return res;
1122}
1123
1124static BOOLEAN iiInternalExport (leftv v, int toLev)
1125{
1126  idhdl h=(idhdl)v->data;
1127  //Print("iiInternalExport('%s',%d)%s\n", v->name, toLev,"");
1128  if (IDLEV(h)==0) Warn("`%s` is already global",IDID(h));
1129  else
1130  {
1131    h=IDROOT->get(v->name,toLev);
1132    idhdl *root=&IDROOT;
1133    if ((h==NULL)&&(currRing!=NULL))
1134    {
1135      h=currRing->idroot->get(v->name,toLev);
1136      root=&currRing->idroot;
1137    }
1138    BOOLEAN keepring=FALSE;
1139    if ((h!=NULL)&&(IDLEV(h)==toLev))
1140    {
1141      if (IDTYP(h)==v->Typ())
1142      {
1143        if (((IDTYP(h)==RING_CMD)||(IDTYP(h)==QRING_CMD))
1144        && (v->Data()==IDDATA(h)))
1145        {
1146          IDRING(h)->ref++;
1147          keepring=TRUE;
1148          IDLEV(h)=toLev;
1149          //WarnS("keepring");
1150          return FALSE;
1151        }
1152        if (BVERBOSE(V_REDEFINE))
1153        {
1154          Warn("redefining %s",IDID(h));
1155        }
1156#ifdef USE_IILOCALRING
1157        if (iiLocalRing[0]==IDRING(h) && (!keepring)) iiLocalRing[0]=NULL;
1158#else
1159        proclevel *p=procstack;
1160        while (p->next!=NULL) p=p->next;
1161        if ((p->cRing==IDRING(h)) && (!keepring))
1162        {
1163          p->cRing=NULL;
1164          p->cRingHdl=NULL;
1165        }
1166#endif
1167        killhdl2(h,root,currRing);
1168      }
1169      else
1170      {
1171        return TRUE;
1172      }
1173    }
1174    h=(idhdl)v->data;
1175    IDLEV(h)=toLev;
1176    if (keepring) IDRING(h)->ref--;
1177    iiNoKeepRing=FALSE;
1178    //Print("export %s\n",IDID(h));
1179  }
1180  return FALSE;
1181}
1182
1183#ifdef HAVE_NS
1184BOOLEAN iiInternalExport (leftv v, int toLev, idhdl roothdl)
1185{
1186  idhdl h=(idhdl)v->data;
1187  if(h==NULL)
1188  {
1189    Warn("'%s': no such identifier\n", v->name);
1190    return FALSE;
1191  }
1192  package frompack=v->req_packhdl; 
1193  if (frompack==NULL) frompack=currPack;
1194  package rootpack = IDPACKAGE(roothdl);
1195//  Print("iiInternalExport('%s',%d,%s->%s) typ:%d\n", v->name, toLev, IDID(currPackHdl),IDID(roothdl),v->Typ());
1196  if ((RingDependend(IDTYP(h)))
1197  || ((IDTYP(h)==LIST_CMD)
1198     && (lRingDependend(IDLIST(h)))
1199     )
1200  )
1201  {
1202    //Print("// ==> Ringdependent set nesting to 0\n");
1203    return (iiInternalExport(v, toLev));
1204  }
1205  else
1206  {
1207    IDLEV(h)=toLev;
1208    v->req_packhdl=rootpack;
1209    if (h==frompack->idroot)
1210    {
1211      frompack->idroot=h->next;
1212    }
1213    else
1214    {
1215      idhdl hh=frompack->idroot;
1216      while ((hh!=NULL) && (hh->next!=h))
1217        hh=hh->next;
1218      if ((hh!=NULL) && (hh->next==h))
1219        hh->next=h->next;
1220      else
1221      {
1222        Werror("`%s` not found",v->Name());
1223        return TRUE;
1224      }
1225    }
1226    h->next=rootpack->idroot;
1227    rootpack->idroot=h;
1228  }
1229  return FALSE;
1230}
1231#endif /* HAVE_NS */
1232
1233BOOLEAN iiExport (leftv v, int toLev)
1234{
1235#ifdef HAVE_NS
1236  checkall();
1237#endif
1238  BOOLEAN nok=FALSE;
1239  leftv r=v;
1240  while (v!=NULL)
1241  {
1242    if ((v->name==NULL)||(v->rtyp==0)||(v->e!=NULL))
1243    {
1244      WerrorS("cannot export");
1245      nok=TRUE;
1246    }
1247    else
1248    {
1249      if(iiInternalExport(v, toLev))
1250      {
1251        r->CleanUp();
1252        return TRUE;
1253      }
1254    }
1255    v=v->next;
1256  }
1257  r->CleanUp();
1258#ifdef HAVE_NS
1259  checkall();
1260#endif
1261  return nok;
1262}
1263
1264/*assume root!=idroot*/
1265#ifdef HAVE_NS
1266BOOLEAN iiExport (leftv v, int toLev, idhdl root)
1267{
1268  checkall();
1269  //  Print("iiExport1: pack=%s\n",IDID(root));
1270  BOOLEAN nok=FALSE;
1271  leftv rv=v;
1272  while (v!=NULL)
1273  {
1274    if ((v->name==NULL)||(v->rtyp==0)||(v->e!=NULL)
1275    )
1276    {
1277      WerrorS("cannot export");
1278      nok=TRUE;
1279    }
1280    else
1281    {
1282      idhdl old=root->get(v->name,toLev);
1283      if (old!=NULL)
1284      {
1285        if (IDTYP(old)==v->Typ())
1286        {
1287          if (BVERBOSE(V_REDEFINE))
1288          {
1289            Warn("redefining %s",IDID(old));
1290          }
1291          killhdl2(old,&root,currRing);
1292        }
1293        else
1294        {
1295          rv->CleanUp();
1296          return TRUE;
1297        }
1298      }
1299      //Print("iiExport: pack=%s\n",IDID(root));
1300      if(iiInternalExport(v, toLev, root))
1301      {
1302        rv->CleanUp();
1303        return TRUE;
1304      }
1305    }
1306    v=v->next;
1307  }
1308  rv->CleanUp();
1309  checkall();
1310  return nok;
1311}
1312#endif
1313
1314BOOLEAN iiCheckRing(int i)
1315{
1316  if (currRingHdl==NULL)
1317  {
1318    #ifdef SIQ
1319    if (siq<=0)
1320    {
1321    #endif
1322      if (RingDependend(i))
1323      {
1324        WerrorS("no ring active");
1325        return TRUE;
1326      }
1327    #ifdef SIQ
1328    }
1329    #endif
1330  }
1331  return FALSE;
1332}
1333
1334poly    iiHighCorner(ideal I, int ak)
1335{
1336  int i;
1337  if(!idIsZeroDim(I)) return NULL; // not zero-dim.
1338  poly po=NULL;
1339  if (currRing->OrdSgn== -1)
1340  {
1341    scComputeHC(I,currQuotient,ak,po);
1342    if (po!=NULL)
1343    {
1344      pGetCoeff(po)=nInit(1);
1345      for (i=pVariables; i>0; i--)
1346      {
1347        if (pGetExp(po, i) > 0) pDecrExp(po,i);
1348      }
1349      pSetComp(po,ak);
1350      pSetm(po);
1351    }
1352  }
1353  else
1354    po=pOne();
1355  return po;
1356}
1357
1358#ifdef HAVE_NS
1359void iiCheckPack(package &p)
1360{
1361  if (p==basePack) return;
1362
1363  idhdl t=basePack->idroot;
1364
1365  while ((t!=NULL) && (IDTYP(t)!=PACKAGE_CMD) && (IDPACKAGE(t)!=p)) t=t->next;
1366
1367  if (t==NULL)
1368  {
1369    WarnS("package not found\n");
1370    p=basePack;
1371  }
1372  return;
1373}
1374#endif
1375
1376idhdl rDefault(char *s)
1377{
1378  idhdl tmp=NULL;
1379
1380  if (s!=NULL) tmp = enterid(s, myynest, RING_CMD, &IDROOT);
1381  if (tmp==NULL) return NULL;
1382
1383  if (ppNoether!=NULL) pDelete(&ppNoether);
1384  if (sLastPrinted.RingDependend())
1385  {
1386    sLastPrinted.CleanUp();
1387    memset(&sLastPrinted,0,sizeof(sleftv));
1388  }
1389
1390  ring r = IDRING(tmp);
1391
1392  r->ch    = 32003;
1393  r->N     = 3;
1394  /*r->P     = 0; Alloc0 in idhdl::set, ipid.cc*/
1395  /*names*/
1396  r->names = (char **) omAlloc0(3 * sizeof(char_ptr));
1397  r->names[0]  = omStrDup("x");
1398  r->names[1]  = omStrDup("y");
1399  r->names[2]  = omStrDup("z");
1400  /*weights: entries for 3 blocks: NULL*/
1401  r->wvhdl = (int **)omAlloc0(3 * sizeof(int_ptr));
1402  /*order: dp,C,0*/
1403  r->order = (int *) omAlloc(3 * sizeof(int *));
1404  r->block0 = (int *)omAlloc0(3 * sizeof(int *));
1405  r->block1 = (int *)omAlloc0(3 * sizeof(int *));
1406  /* ringorder dp for the first block: var 1..3 */
1407  r->order[0]  = ringorder_dp;
1408  r->block0[0] = 1;
1409  r->block1[0] = 3;
1410  /* ringorder C for the second block: no vars */
1411  r->order[1]  = ringorder_C;
1412  /* the last block: everything is 0 */
1413  r->order[2]  = 0;
1414  /*polynomial ring*/
1415  r->OrdSgn    = 1;
1416
1417  /* complete ring intializations */
1418  rComplete(r);
1419  rSetHdl(tmp);
1420  return currRingHdl;
1421}
1422
1423idhdl rFindHdl(ring r, idhdl n, idhdl w)
1424{
1425  idhdl h=rSimpleFindHdl(r,IDROOT,n);
1426  if (h!=NULL)  return h;
1427#ifdef HAVE_NS
1428  if (IDROOT!=basePack->idroot) h=rSimpleFindHdl(r,basePack->idroot,n);
1429  if (h!=NULL)  return h;
1430  proclevel *p=procstack;
1431  while(p!=NULL)
1432  {
1433    if ((p->cPack!=basePack)
1434    && (p->cPack!=currPack))
1435      h=rSimpleFindHdl(r,p->cPack->idroot,n);
1436    if (h!=NULL)  return h;
1437    p=p->next;
1438  }
1439  idhdl tmp=basePack->idroot;
1440  while (tmp!=NULL)
1441  {
1442    if (IDTYP(tmp)==PACKAGE_CMD)
1443      h=rSimpleFindHdl(r,IDPACKAGE(tmp)->idroot,n);
1444    if (h!=NULL)  return h;
1445    tmp=IDNEXT(tmp);
1446  }
1447#endif
1448  return NULL;
1449}
1450
1451void rDecomposeCF(leftv h,const ring r,const ring R)
1452{
1453  lists L=(lists)omAlloc0Bin(slists_bin);
1454  L->Init(4);
1455  h->rtyp=LIST_CMD;
1456  h->data=(void *)L;
1457  // 0: char/ cf - ring
1458  // 1: list (var)
1459  // 2: list (ord)
1460  // 3: qideal
1461  // ----------------------------------------
1462  // 0: char/ cf - ring
1463  L->m[0].rtyp=INT_CMD;
1464  L->m[0].data=(void *)r->ch;
1465  // ----------------------------------------
1466  // 1: list (var)
1467  lists LL=(lists)omAlloc0Bin(slists_bin);
1468  LL->Init(r->N);
1469  int i;
1470  for(i=0; i<r->N; i++)
1471  {
1472    LL->m[i].rtyp=STRING_CMD;
1473    LL->m[i].data=(void *)omStrDup(r->names[i]);
1474  }
1475  L->m[1].rtyp=LIST_CMD;
1476  L->m[1].data=(void *)LL;
1477  // ----------------------------------------
1478  // 2: list (ord)
1479  LL=(lists)omAlloc0Bin(slists_bin);
1480  i=rBlocks(r)-1;
1481  LL->Init(i);
1482  i--;
1483  lists LLL;
1484  for(; i>=0; i--)
1485  {
1486    intvec *iv;
1487    int j;
1488    LL->m[i].rtyp=LIST_CMD;
1489    LLL=(lists)omAlloc0Bin(slists_bin);
1490    LLL->Init(2);
1491    LLL->m[0].rtyp=STRING_CMD;
1492    LLL->m[0].data=(void *)omStrDup(rSimpleOrdStr(r->order[i]));
1493    if (r->block1[i]-r->block0[i] >=0 )
1494    {
1495      j=r->block1[i]-r->block0[i];
1496      iv=new intvec(j+1);
1497      if ((r->wvhdl!=NULL) && (r->wvhdl[i]!=NULL))
1498      {
1499        for(;j>=0; j--) (*iv)[j]=r->wvhdl[i][j];
1500      }
1501      else switch (r->order[i])
1502      {
1503        case ringorder_dp:
1504        case ringorder_Dp:
1505        case ringorder_ds:
1506        case ringorder_Ds:
1507        case ringorder_lp:
1508          for(;j>=0; j--) (*iv)[j]=1;
1509          break;
1510        default: /* do nothing */;
1511      }
1512    }
1513    else
1514    {
1515      iv=new intvec(1);
1516    }
1517    LLL->m[1].rtyp=INTVEC_CMD;
1518    LLL->m[1].data=(void *)iv;
1519    LL->m[i].data=(void *)LLL;
1520  }
1521  L->m[2].rtyp=LIST_CMD;
1522  L->m[2].data=(void *)LL;
1523  // ----------------------------------------
1524  // 3: qideal
1525  L->m[3].rtyp=IDEAL_CMD;
1526  if (R->minpoly==NULL)
1527    L->m[3].data=(void *)idInit(1,1);
1528  else
1529  {
1530    ideal I=idInit(1,1);
1531    L->m[3].data=(void *)I;
1532    I->m[0]=pOne();
1533    pSetCoeff(I->m[0],R->minpoly);
1534  }
1535  // ----------------------------------------
1536}
1537void rDecomposeC(leftv h,const ring R)
1538{
1539  lists L=(lists)omAlloc0Bin(slists_bin);
1540  L->Init(4);
1541  h->rtyp=LIST_CMD;
1542  h->data=(void *)L;
1543  // 0: char/ cf - ring
1544  // 1: list (var)
1545  // 2: list (ord)
1546  // 3: qideal
1547  // ----------------------------------------
1548  // 0: char/ cf - ring
1549  L->m[0].rtyp=INT_CMD;
1550  L->m[0].data=(void *)R->ch;
1551  // ----------------------------------------
1552  // 1: list (var)
1553  lists LL=(lists)omAlloc0Bin(slists_bin);
1554  LL->Init(1);
1555    LL->m[0].rtyp=STRING_CMD;
1556    LL->m[0].data=(void *)omStrDup(R->parameter[0]);
1557  L->m[1].rtyp=LIST_CMD;
1558  L->m[1].data=(void *)LL;
1559  // ----------------------------------------
1560  // 2: list (ord)
1561  LL=(lists)omAlloc0Bin(slists_bin);
1562  LL->Init(1);
1563  lists LLL;
1564  {
1565    intvec *iv;
1566    int j;
1567    LL->m[0].rtyp=LIST_CMD;
1568    LLL=(lists)omAlloc0Bin(slists_bin);
1569    LLL->Init(2);
1570    LLL->m[0].rtyp=STRING_CMD;
1571    LLL->m[0].data=(void *)omStrDup("lp");
1572    {
1573      iv=new intvec(1);
1574      for(;j>=0; j--) (*iv)[0]=1;
1575    }
1576    LLL->m[1].rtyp=INTVEC_CMD;
1577    LLL->m[1].data=(void *)iv;
1578    LL->m[0].data=(void *)LLL;
1579  }
1580  L->m[2].rtyp=LIST_CMD;
1581  L->m[2].data=(void *)LL;
1582  // ----------------------------------------
1583  // 3: qideal
1584  L->m[3].rtyp=IDEAL_CMD;
1585  {
1586    ideal I=idInit(1,1);
1587    L->m[3].data=(void *)I;
1588    //I->m[0]=pOne();
1589    //pSetCoeff(I->m[0],R->minpoly);
1590  }
1591  // ----------------------------------------
1592}
1593
1594lists rDecompose(ring r)
1595{
1596  // 0: char/ cf - ring
1597  // 1: list (var)
1598  // 2: list (ord)
1599  // 3: qideal
1600  // possibly:
1601  // 4: C
1602  // 5: D
1603  lists L=(lists)omAlloc0Bin(slists_bin);
1604  if (rIsPluralRing(r))
1605    L->Init(6);
1606  else
1607    L->Init(4);
1608  // ----------------------------------------
1609  // 0: char/ cf - ring
1610  #if 1 /* TODO */
1611  if (rIsExtension(r))
1612  {
1613    if (rField_is_long_C(r))
1614      rDecomposeC(&(L->m[0]),r);
1615    else
1616      rDecomposeCF(&(L->m[0]),r->algring,r);
1617  }
1618  else
1619  #endif
1620  {
1621    L->m[0].rtyp=INT_CMD;
1622    L->m[0].data=(void *)r->ch;
1623  }
1624  // ----------------------------------------
1625  // 1: list (var)
1626  lists LL=(lists)omAlloc0Bin(slists_bin);
1627  LL->Init(r->N);
1628  int i;
1629  for(i=0; i<r->N; i++)
1630  {
1631    LL->m[i].rtyp=STRING_CMD;
1632    LL->m[i].data=(void *)omStrDup(r->names[i]);
1633  }
1634  L->m[1].rtyp=LIST_CMD;
1635  L->m[1].data=(void *)LL;
1636  // ----------------------------------------
1637  // 2: list (ord)
1638  LL=(lists)omAlloc0Bin(slists_bin);
1639  i=rBlocks(r)-1;
1640  LL->Init(i);
1641  i--;
1642  lists LLL;
1643  for(; i>=0; i--)
1644  {
1645    intvec *iv;
1646    int j;
1647    LL->m[i].rtyp=LIST_CMD;
1648    LLL=(lists)omAlloc0Bin(slists_bin);
1649    LLL->Init(2);
1650    LLL->m[0].rtyp=STRING_CMD;
1651    LLL->m[0].data=(void *)omStrDup(rSimpleOrdStr(r->order[i]));
1652    if (r->block1[i]-r->block0[i] >=0 )
1653    {
1654      j=r->block1[i]-r->block0[i];
1655      iv=new intvec(j+1);
1656      if ((r->wvhdl!=NULL) && (r->wvhdl[i]!=NULL))
1657      {
1658        for(;j>=0; j--) (*iv)[j]=r->wvhdl[i][j];
1659      }
1660      else switch (r->order[i])
1661      {
1662        case ringorder_dp:
1663        case ringorder_Dp:
1664        case ringorder_ds:
1665        case ringorder_Ds:
1666        case ringorder_lp:
1667          for(;j>=0; j--) (*iv)[j]=1;
1668          break;
1669        default: /* do nothing */;
1670      }
1671    }
1672    else
1673    {
1674      iv=new intvec(1);
1675    }
1676    LLL->m[1].rtyp=INTVEC_CMD;
1677    LLL->m[1].data=(void *)iv;
1678    LL->m[i].data=(void *)LLL;
1679  }
1680  L->m[2].rtyp=LIST_CMD;
1681  L->m[2].data=(void *)LL;
1682  // ----------------------------------------
1683  // 3: qideal
1684  L->m[3].rtyp=IDEAL_CMD;
1685  if (r->qideal==NULL)
1686    L->m[3].data=(void *)idInit(1,1);
1687  else
1688    L->m[3].data=(void *)idCopy(r->qideal);
1689  // ----------------------------------------
1690  #ifdef HAVE_PLURAL
1691  if (rIsPluralRing(r))
1692  {
1693    L->m[4].rtyp=MATRIX_CMD;
1694    L->m[4].data=(void *)mpCopy(r->nc->C);
1695    L->m[5].rtyp=MATRIX_CMD;
1696    L->m[5].data=(void *)mpCopy(r->nc->D);
1697  }
1698  #endif
1699  return L;
1700}
1701
1702ring rCompose(lists  L)
1703{
1704  if ((L->nr!=3)
1705#ifdef HAVE_PLURAL
1706  &&(L->nr!=5)
1707#endif
1708  )
1709    return NULL;
1710  // 0: char/ cf - ring
1711  // 1: list (var)
1712  // 2: list (ord)
1713  // 3: qideal
1714  // possibly:
1715  // 4: C
1716  // 5: D
1717  ring R=(ring) omAlloc0Bin(sip_sring_bin);
1718  if (L->m[0].Typ()==INT_CMD)
1719  {
1720    R->ch=(int)L->m[0].Data();
1721  }
1722  else if (L->m[0].Typ()==LIST_CMD)
1723  {
1724    R->algring=rCompose((lists)L->m[0].Data());
1725    if (R->algring==NULL)
1726    {
1727      WerrorS("could not create rational function coefficient field");
1728      goto rCompose_err;
1729    }
1730    if (R->algring->ch>0)
1731       R->ch= -R->algring->ch;
1732    else
1733       R->ch=1;
1734    R->P=R->algring->N;
1735    R->parameter=(char**)omAlloc0(R->P*sizeof(char_ptr));
1736    int i;
1737    for(i=R->P-1;i>=0;i--)
1738      R->parameter[i]=omStrDup(R->algring->names[i]);
1739    if (R->algring->qideal!=NULL)
1740    {
1741      R->minpoly=pGetCoeff(R->algring->qideal->m[0]);
1742    }
1743  }
1744  else
1745  {
1746    WerrorS("coefficient field must be described by `int` or `list`");
1747    goto rCompose_err;
1748  }
1749  // ------------------------- VARS ---------------------------
1750  if (L->m[1].Typ()==LIST_CMD)
1751  {
1752    lists v=(lists)L->m[1].Data();
1753    R->N = v->nr+1;
1754    R->names   = (char **)omAlloc0(R->N * sizeof(char_ptr));
1755    int i;
1756    for(i=0;i<R->N;i++)
1757    {
1758      if (v->m[i].Typ()==STRING_CMD)
1759        R->names[i]=omStrDup((char *)v->m[i].Data());
1760      else if (v->m[i].Typ()==POLY_CMD)
1761      {
1762        poly p=(poly)v->m[i].Data();
1763        int nr=pIsPurePower(p);
1764        if (nr>0)
1765          R->names[i]=omStrDup(currRing->names[nr-1]);
1766        else
1767        {
1768          Werror("var name %d must be a string or a ring variable",i+1);
1769          goto rCompose_err;
1770        }
1771      }
1772      else
1773      {
1774        Werror("var name %d must be `string`",i+1);
1775        goto rCompose_err;
1776      }
1777    }
1778  }
1779  else
1780  {
1781    WerrorS("variable must be given as `list`");
1782    goto rCompose_err;
1783  }
1784  // ------------------------ ORDER ------------------------------
1785  if (L->m[2].Typ()==LIST_CMD)
1786  {
1787    lists v=(lists)L->m[2].Data();
1788    int n= v->nr+2;
1789    int j;
1790    // initialize fields of R
1791    R->order=(int *)omAlloc0(n*sizeof(int));
1792    R->block0=(int *)omAlloc0(n*sizeof(int));
1793    R->block1=(int *)omAlloc0(n*sizeof(int));
1794    R->wvhdl=(int**)omAlloc0(n*sizeof(int_ptr));
1795    // init order, so that rBlocks works correctly
1796    for (j=0; j < n-2; j++)
1797      R->order[j] = (int) ringorder_unspec;
1798    // orderings
1799    R->OrdSgn=1;
1800    for(j=0;j<n-1;j++)
1801    {
1802    // todo: a(..), M
1803      if (v->m[j].Typ()!=LIST_CMD)
1804      {
1805        WerrorS("ordering must be list of lists");
1806        goto rCompose_err;
1807      }
1808      lists vv=(lists)v->m[j].Data();
1809      if ((vv->nr!=1)
1810      || (vv->m[0].Typ()!=STRING_CMD)
1811      || ((vv->m[1].Typ()!=INTVEC_CMD) && (vv->m[1].Typ()!=INT_CMD)))
1812      {
1813        WerrorS("ordering name must be a (string,intvec)");
1814        goto rCompose_err;
1815      }
1816      R->order[j]=rOrderName(omStrDup((char*)vv->m[0].Data())); // assume STRING
1817      if (j==0) R->block0[0]=1;
1818      else      R->block0[j]=R->block1[j-1]+1;
1819      intvec *iv;
1820      if (vv->m[1].Typ()==INT_CMD)
1821        iv=new intvec((int)vv->m[1].Data(),(int)vv->m[1].Data());
1822      else
1823        iv=ivCopy((intvec*)vv->m[1].Data()); //assume INTVEC
1824      R->block1[j]=max(R->block0[j],R->block0[j]+iv->length()-1);
1825      int i;
1826      switch (R->order[j])
1827      {
1828         case ringorder_ws:
1829         case ringorder_Ws:
1830            R->OrdSgn=-1;
1831         case ringorder_wp:
1832         case ringorder_Wp:
1833           R->wvhdl[j] =( int *)omAlloc((iv->length())*sizeof(int));
1834           for (i=0; i<iv->length();i++) R->wvhdl[j][i]=(*iv)[i];
1835           break;
1836         case ringorder_ls:
1837         case ringorder_ds:
1838         case ringorder_Ds:
1839           R->OrdSgn=-1;
1840         case ringorder_lp:
1841         case ringorder_dp:
1842         case ringorder_Dp:
1843         case ringorder_rp:
1844           break;
1845         case ringorder_S:
1846           break;
1847         case ringorder_c:
1848         case ringorder_C:
1849           R->block1[j]=R->block0[j]-1;
1850           break;
1851         case ringorder_aa:
1852         case ringorder_a:
1853           R->wvhdl[j] =( int *)omAlloc((iv->length())*sizeof(int));
1854           for (i=1; i<iv->length();i++) R->wvhdl[n][i-1]=(*iv)[i];
1855         // todo
1856           break;
1857         case ringorder_M:
1858         // todo
1859           break;
1860      }
1861    }
1862    // sanity check
1863    j=n-2;
1864    if ((R->order[j]==ringorder_c)
1865    || (R->order[j]==ringorder_C)) j--;
1866    if (R->block1[j] != R->N)
1867    {
1868      if (((R->order[j]==ringorder_dp) ||
1869           (R->order[j]==ringorder_ds) ||
1870           (R->order[j]==ringorder_Dp) ||
1871           (R->order[j]==ringorder_Ds) ||
1872           (R->order[j]==ringorder_rp) ||
1873           (R->order[j]==ringorder_lp) ||
1874           (R->order[j]==ringorder_ls))
1875          &&
1876            R->block0[j] <= R->N)
1877      {
1878        R->block1[j] = R->N;
1879      }
1880      else
1881      {
1882        Werror("ordering incomplete: size (%d) should be %d",R->block1[j],R->N);
1883        goto rCompose_err;
1884      }
1885    }
1886  }
1887  else
1888  {
1889    WerrorS("ordering must be given as `list`");
1890    goto rCompose_err;
1891  }
1892  // ------------------------ Q-IDEAL ------------------------
1893  if (L->m[3].Typ()==IDEAL_CMD)
1894  {
1895    ideal q=(ideal)L->m[3].Data();
1896    if (q->m[0]!=NULL)
1897      R->qideal=idCopy(q);
1898  }
1899  else
1900  {
1901    WerrorS("q-ideal must be given as `ideal`");
1902    goto rCompose_err;
1903  }
1904
1905  // todo
1906  rComplete(R);
1907
1908  // ---------------------------------------------------------------
1909  #ifdef HAVE_PLURAL
1910  if (L->nr==5)
1911  {
1912    if (nc_CallPlural((matrix)L->m[4].Data(),(matrix)L->m[5].Data(),
1913        NULL,NULL,R))
1914      goto rCompose_err;
1915  }
1916  #endif
1917  return R;
1918
1919rCompose_err:
1920  if (R->N>0)
1921  {
1922    int i;
1923    if (R->names!=NULL)
1924    {
1925      i=R->N-1;
1926      while (i>=0) { if (R->names[i]!=NULL) omFree(R->names[i]); i--; }
1927      omFree(R->names);
1928    }
1929  }
1930  if (R->order!=NULL) omFree(R->order);
1931  if (R->block0!=NULL) omFree(R->block0);
1932  if (R->block1!=NULL) omFree(R->block1);
1933  if (R->wvhdl!=NULL) omFree(R->wvhdl);
1934  omFree(R);
1935  return NULL;
1936}
1937
1938// from matpol.cc
1939
1940/*2
1941* compute the jacobi matrix of an ideal
1942*/
1943BOOLEAN mpJacobi(leftv res,leftv a)
1944{
1945  int     i,j;
1946  matrix result;
1947  ideal id=(ideal)a->Data();
1948
1949  result =mpNew(IDELEMS(id),pVariables);
1950  for (i=1; i<=IDELEMS(id); i++)
1951  {
1952    for (j=1; j<=pVariables; j++)
1953    {
1954      MATELEM(result,i,j) = pDiff(id->m[i-1],j);
1955    }
1956  }
1957  res->data=(char *)result;
1958  return FALSE;
1959}
1960
1961/*2
1962* returns the Koszul-matrix of degree d of a vectorspace with dimension n
1963* uses the first n entrees of id, if id <> NULL
1964*/
1965BOOLEAN mpKoszul(leftv res,leftv c/*ip*/, leftv b/*in*/, leftv id)
1966{
1967  int n=(int)b->Data();
1968  int d=(int)c->Data();
1969  int     k,l,sign,row,col;
1970  matrix  result;
1971  ideal temp;
1972  BOOLEAN bo;
1973  poly    p;
1974
1975  if ((d>n) || (d<1) || (n<1))
1976  {
1977    res->data=(char *)mpNew(1,1);
1978    return FALSE;
1979  }
1980  int *choise = (int*)omAlloc(d*sizeof(int));
1981  if (id==NULL)
1982    temp=idMaxIdeal(1);
1983  else
1984    temp=(ideal)id->Data();
1985
1986  k = binom(n,d);
1987  l = k*d;
1988  l /= n-d+1;
1989  result =mpNew(l,k);
1990  col = 1;
1991  idInitChoise(d,1,n,&bo,choise);
1992  while (!bo)
1993  {
1994    sign = 1;
1995    for (l=1;l<=d;l++)
1996    {
1997      if (choise[l-1]<=IDELEMS(temp))
1998      {
1999        p = pCopy(temp->m[choise[l-1]-1]);
2000        if (sign == -1) p = pNeg(p);
2001        sign *= -1;
2002        row = idGetNumberOfChoise(l-1,d,1,n,choise);
2003        MATELEM(result,row,col) = p;
2004      }
2005    }
2006    col++;
2007    idGetNextChoise(d,n,&bo,choise);
2008  }
2009  if (id==NULL) idDelete(&temp);
2010
2011  res->data=(char *)result;
2012  return FALSE;
2013}
2014
2015// from syz1.cc
2016/*2
2017* read out the Betti numbers from resolution
2018* (interpreter interface)
2019*/
2020BOOLEAN syBetti2(leftv res, leftv u, leftv w)
2021{
2022  syStrategy syzstr=(syStrategy)u->Data();
2023  BOOLEAN minim=(int)w->Data();
2024  int row_shift=0;
2025  int add_row_shift=0;
2026  intvec *weights=NULL;
2027  intvec *ww=(intvec *)atGet(u,"isHomog",INTVEC_CMD);
2028  if (ww!=NULL)
2029  {
2030     weights=ivCopy(ww);
2031     add_row_shift = ww->min_in();
2032     (*weights) -= add_row_shift;
2033  }
2034  res->data=(void *)syBettiOfComputation(syzstr,minim,&row_shift,weights);
2035  //row_shift += add_row_shift;
2036  //Print("row_shift=%d, add_row_shift=%d\n",row_shift,add_row_shift);
2037  atSet(res,omStrDup("rowShift"),(void*)add_row_shift,INT_CMD);
2038  return FALSE;
2039}
2040BOOLEAN syBetti1(leftv res, leftv u)
2041{
2042  sleftv tmp;
2043  memset(&tmp,0,sizeof(tmp));
2044  tmp.rtyp=INT_CMD;
2045  tmp.data=(void *)1; 
2046  return syBetti2(res,u,&tmp);
2047}
2048
2049/*3
2050* converts a resolution into a list of modules
2051*/
2052lists syConvRes(syStrategy syzstr,BOOLEAN toDel)
2053{
2054  if ((syzstr->fullres==NULL) && (syzstr->minres==NULL))
2055  {
2056    if (syzstr->hilb_coeffs==NULL)
2057    {
2058      syzstr->fullres = syReorder(syzstr->res,syzstr->length,syzstr);
2059    }
2060    else
2061    {
2062      syzstr->minres = syReorder(syzstr->orderedRes,syzstr->length,syzstr);
2063      syKillEmptyEntres(syzstr->minres,syzstr->length);
2064    }
2065  }
2066  resolvente tr;
2067  int typ0=IDEAL_CMD;
2068  if (syzstr->minres!=NULL)
2069    tr = syzstr->minres;
2070  else
2071    tr = syzstr->fullres;
2072  resolvente trueres=NULL;
2073  intvec ** w=NULL;
2074  if (syzstr->length>0)
2075  {
2076    trueres=(resolvente)omAlloc0((syzstr->length)*sizeof(ideal));
2077    for (int i=(syzstr->length)-1;i>=0;i--)
2078    {
2079      if (tr[i]!=NULL)
2080      {
2081        trueres[i] = idCopy(tr[i]);
2082      }
2083    }
2084    if (idRankFreeModule(trueres[0]) > 0)
2085      typ0 = MODUL_CMD;
2086    if (syzstr->weights!=NULL)
2087    {
2088      w = (intvec**)omAlloc0((syzstr->length)*sizeof(intvec*));
2089      for (int i=(syzstr->length)-1;i>=0;i--)
2090      {
2091        if (syzstr->weights[i]!=NULL) w[i] = ivCopy(syzstr->weights[i]);
2092      }
2093    }
2094  }
2095  lists li = liMakeResolv(trueres,syzstr->length,syzstr->list_length,typ0,w);
2096  if (w != NULL) omFreeSize(w, (syzstr->length)*sizeof(intvec*));
2097  if (toDel) syKillComputation(syzstr);
2098  return li;
2099}
2100
2101/*3
2102* converts a list of modules into a resolution
2103*/
2104syStrategy syConvList(lists li,BOOLEAN toDel)
2105{
2106  int typ0;
2107  syStrategy result=(syStrategy)omAlloc0(sizeof(ssyStrategy));
2108
2109  resolvente fr = liFindRes(li,&(result->length),&typ0,&(result->weights));
2110  if (fr != NULL)
2111  {
2112
2113    result->fullres = (resolvente)omAlloc0((result->length+1)*sizeof(ideal));
2114    for (int i=result->length-1;i>=0;i--)
2115    {
2116      if (fr[i]!=NULL)
2117        result->fullres[i] = idCopy(fr[i]);
2118    }
2119    result->list_length=result->length;
2120    omFreeSize((ADDRESS)fr,(result->length)*sizeof(ideal));
2121  }
2122  else
2123  {
2124    omFreeSize(result, sizeof(ssyStrategy));
2125    result = NULL;
2126  }
2127  if (toDel) li->Clean();
2128  return result;
2129}
2130
2131/*3
2132* converts a list of modules into a minimal resolution
2133*/
2134syStrategy syForceMin(lists li)
2135{
2136  int typ0;
2137  syStrategy result=(syStrategy)omAlloc0(sizeof(ssyStrategy));
2138
2139  resolvente fr = liFindRes(li,&(result->length),&typ0);
2140  result->minres = (resolvente)omAlloc0((result->length+1)*sizeof(ideal));
2141  for (int i=result->length-1;i>=0;i--)
2142  {
2143    if (fr[i]!=NULL)
2144      result->minres[i] = idCopy(fr[i]);
2145  }
2146  omFreeSize((ADDRESS)fr,(result->length)*sizeof(ideal));
2147  return result;
2148}
2149// from weight.cc
2150BOOLEAN kWeight(leftv res,leftv id)
2151{
2152  ideal F=(ideal)id->Data();
2153  intvec * iv = new intvec(pVariables);
2154  polyset s;
2155  int  sl, n, i;
2156  int  *x;
2157
2158  res->data=(char *)iv;
2159  s = F->m;
2160  sl = IDELEMS(F) - 1;
2161  n = pVariables;
2162  wNsqr = (double)2.0 / (double)n;
2163  wFunctional = wFunctionalBuch;
2164  x = (int * )omAlloc(2 * (n + 1) * sizeof(int));
2165  wCall(s, sl, x);
2166  for (i = n; i!=0; i--)
2167    (*iv)[i-1] = x[i + n + 1];
2168  omFreeSize((ADDRESS)x, 2 * (n + 1) * sizeof(int));
2169  return FALSE;
2170}
2171
2172BOOLEAN kQHWeight(leftv res,leftv v)
2173{
2174  res->data=(char *)idQHomWeight((ideal)v->Data());
2175  if (res->data==NULL)
2176    res->data=(char *)new intvec(pVariables);
2177  return FALSE;
2178}
2179/*==============================================================*/
2180// from clapsing.cc
2181#if 0
2182BOOLEAN jjIS_SQR_FREE(leftv res, leftv u)
2183{
2184  BOOLEAN b=singclap_factorize((poly)(u->Data()), &v, 0);
2185  res->data=(void *)b;
2186}
2187#endif
2188
2189BOOLEAN jjRESULTANT(leftv res, leftv u, leftv v, leftv w)
2190{
2191  res->data=singclap_resultant((poly)u->Data(),(poly)v->Data(), (poly)w->Data());
2192  return errorreported;
2193}
2194BOOLEAN jjCHARSERIES(leftv res, leftv u)
2195{
2196  res->data=singclap_irrCharSeries((ideal)u->Data());
2197  return (res->data==NULL);
2198}
2199
2200// from semic.cc
2201#ifdef HAVE_SPECTRUM
2202
2203// ----------------------------------------------------------------------------
2204//  Initialize a  spectrum  deep from another  spectrum
2205// ----------------------------------------------------------------------------
2206
2207void spectrum::copy_deep( const spectrum &spec )
2208{
2209    mu = spec.mu;
2210    pg = spec.pg;
2211    n  = spec.n;
2212
2213    copy_new( n );
2214
2215    for( int i=0; i<n; i++ )
2216    {
2217        s[i] = spec.s[i];
2218        w[i] = spec.w[i];
2219    }
2220}
2221
2222// ----------------------------------------------------------------------------
2223//  Initialize a  spectrum  deep from a  singular  lists
2224// ----------------------------------------------------------------------------
2225
2226void spectrum::copy_deep( lists l )
2227{
2228    mu = (int)(l->m[0].Data( ));
2229    pg = (int)(l->m[1].Data( ));
2230    n  = (int)(l->m[2].Data( ));
2231
2232    copy_new( n );
2233
2234    intvec  *num = (intvec*)l->m[3].Data( );
2235    intvec  *den = (intvec*)l->m[4].Data( );
2236    intvec  *mul = (intvec*)l->m[5].Data( );
2237
2238    for( int i=0; i<n; i++ )
2239    {
2240        s[i] = (Rational)((*num)[i])/(Rational)((*den)[i]);
2241        w[i] = (*mul)[i];
2242    }
2243}
2244
2245// ----------------------------------------------------------------------------
2246//  singular lists  constructor for  spectrum
2247// ----------------------------------------------------------------------------
2248
2249spectrum::spectrum( lists l )
2250{
2251    copy_deep( l );
2252}
2253
2254// ----------------------------------------------------------------------------
2255//  generate a Singular  lists  from a spectrum
2256// ----------------------------------------------------------------------------
2257
2258lists   spectrum::thelist( void )
2259{
2260    lists   L  = (lists)omAllocBin( slists_bin);
2261
2262    L->Init( 6 );
2263
2264    intvec            *num  = new intvec( n );
2265    intvec            *den  = new intvec( n );
2266    intvec            *mult = new intvec( n );
2267
2268    for( int i=0; i<n; i++ )
2269    {
2270        (*num) [i] = s[i].get_num_si( );
2271        (*den) [i] = s[i].get_den_si( );
2272        (*mult)[i] = w[i];
2273    }
2274
2275    L->m[0].rtyp = INT_CMD;    //  milnor number
2276    L->m[1].rtyp = INT_CMD;    //  geometrical genus
2277    L->m[2].rtyp = INT_CMD;    //  # of spectrum numbers
2278    L->m[3].rtyp = INTVEC_CMD; //  numerators
2279    L->m[4].rtyp = INTVEC_CMD; //  denomiantors
2280    L->m[5].rtyp = INTVEC_CMD; //  multiplicities
2281
2282    L->m[0].data = (void*)mu;
2283    L->m[1].data = (void*)pg;
2284    L->m[2].data = (void*)n;
2285    L->m[3].data = (void*)num;
2286    L->m[4].data = (void*)den;
2287    L->m[5].data = (void*)mult;
2288
2289    return  L;
2290}
2291// from spectrum.cc
2292// ----------------------------------------------------------------------------
2293//  print out an error message for a spectrum list
2294// ----------------------------------------------------------------------------
2295
2296void    list_error( semicState state )
2297{
2298    switch( state )
2299    {
2300        case semicListTooShort:
2301            WerrorS( "the list is too short" );
2302            break;
2303        case semicListTooLong:
2304            WerrorS( "the list is too long" );
2305            break;
2306
2307        case semicListFirstElementWrongType:
2308            WerrorS( "first element of the list should be int" );
2309            break;
2310        case semicListSecondElementWrongType:
2311            WerrorS( "second element of the list should be int" );
2312            break;
2313        case semicListThirdElementWrongType:
2314            WerrorS( "third element of the list should be int" );
2315            break;
2316        case semicListFourthElementWrongType:
2317            WerrorS( "fourth element of the list should be intvec" );
2318            break;
2319        case semicListFifthElementWrongType:
2320            WerrorS( "fifth element of the list should be intvec" );
2321            break;
2322        case semicListSixthElementWrongType:
2323            WerrorS( "sixth element of the list should be intvec" );
2324            break;
2325
2326        case semicListNNegative:
2327            WerrorS( "first element of the list should be positive" );
2328            break;
2329        case semicListWrongNumberOfNumerators:
2330            WerrorS( "wrong number of numerators" );
2331            break;
2332        case semicListWrongNumberOfDenominators:
2333            WerrorS( "wrong number of denominators" );
2334            break;
2335        case semicListWrongNumberOfMultiplicities:
2336            WerrorS( "wrong number of multiplicities" );
2337            break;
2338
2339        case semicListMuNegative:
2340            WerrorS( "the Milnor number should be positive" );
2341            break;
2342        case semicListPgNegative:
2343            WerrorS( "the geometrical genus should be nonnegative" );
2344            break;
2345        case semicListNumNegative:
2346            WerrorS( "all numerators should be positive" );
2347            break;
2348        case semicListDenNegative:
2349            WerrorS( "all denominators should be positive" );
2350            break;
2351        case semicListMulNegative:
2352            WerrorS( "all multiplicities should be positive" );
2353            break;
2354
2355        case semicListNotSymmetric:
2356            WerrorS( "it is not symmetric" );
2357            break;
2358        case semicListNotMonotonous:
2359            WerrorS( "it is not monotonous" );
2360            break;
2361
2362        case semicListMilnorWrong:
2363            WerrorS( "the Milnor number is wrong" );
2364            break;
2365        case semicListPGWrong:
2366            WerrorS( "the geometrical genus is wrong" );
2367            break;
2368
2369        default:
2370            WerrorS( "unspecific error" );
2371            break;
2372    }
2373}
2374// ----------------------------------------------------------------------------
2375//  this is the main spectrum computation function
2376// ----------------------------------------------------------------------------
2377
2378spectrumState   spectrumCompute( poly h,lists *L,int fast )
2379{
2380  int i,j;
2381
2382  #ifdef SPECTRUM_DEBUG
2383  #ifdef SPECTRUM_PRINT
2384  #ifdef SPECTRUM_IOSTREAM
2385    cout << "spectrumCompute\n";
2386    if( fast==0 ) cout << "    no optimization" << endl;
2387    if( fast==1 ) cout << "    weight optimization" << endl;
2388    if( fast==2 ) cout << "    symmetry optimization" << endl;
2389  #else
2390    fprintf( stdout,"spectrumCompute\n" );
2391    if( fast==0 ) fprintf( stdout,"    no optimization\n" );
2392    if( fast==1 ) fprintf( stdout,"    weight optimization\n" );
2393    if( fast==2 ) fprintf( stdout,"    symmetry optimization\n" );
2394  #endif
2395  #endif
2396  #endif
2397
2398  // ----------------------
2399  //  check if  h  is zero
2400  // ----------------------
2401
2402  if( h==(poly)NULL )
2403  {
2404    return  spectrumZero;
2405  }
2406
2407  // ----------------------------------
2408  //  check if  h  has a constant term
2409  // ----------------------------------
2410
2411  if( hasConstTerm( h ) )
2412  {
2413    return  spectrumBadPoly;
2414  }
2415
2416  // --------------------------------
2417  //  check if  h  has a linear term
2418  // --------------------------------
2419
2420  if( hasLinearTerm( h ) )
2421  {
2422    *L = (lists)omAllocBin( slists_bin);
2423    (*L)->Init( 1 );
2424    (*L)->m[0].rtyp = INT_CMD;    //  milnor number
2425    /* (*L)->m[0].data = (void*)0;a  -- done by Init */
2426
2427    return  spectrumNoSingularity;
2428  }
2429
2430  // ----------------------------------
2431  //  compute the jacobi ideal of  (h)
2432  // ----------------------------------
2433
2434  ideal J = NULL;
2435  J = idInit( pVariables,1 );
2436
2437  #ifdef SPECTRUM_DEBUG
2438  #ifdef SPECTRUM_PRINT
2439  #ifdef SPECTRUM_IOSTREAM
2440    cout << "\n   computing the Jacobi ideal...\n";
2441  #else
2442    fprintf( stdout,"\n   computing the Jacobi ideal...\n" );
2443  #endif
2444  #endif
2445  #endif
2446
2447  for( i=0; i<pVariables; i++ )
2448  {
2449    J->m[i] = pDiff( h,i+1); //j );
2450
2451    #ifdef SPECTRUM_DEBUG
2452    #ifdef SPECTRUM_PRINT
2453    #ifdef SPECTRUM_IOSTREAM
2454      cout << "        ";
2455    #else
2456      fprintf( stdout,"        " );
2457    #endif
2458      pWrite( J->m[i] );
2459    #endif
2460    #endif
2461  }
2462
2463  // --------------------------------------------
2464  //  compute a standard basis  stdJ  of  jac(h)
2465  // --------------------------------------------
2466
2467  #ifdef SPECTRUM_DEBUG
2468  #ifdef SPECTRUM_PRINT
2469  #ifdef SPECTRUM_IOSTREAM
2470    cout << endl;
2471    cout << "    computing a standard basis..." << endl;
2472  #else
2473    fprintf( stdout,"\n" );
2474    fprintf( stdout,"    computing a standard basis...\n" );
2475  #endif
2476  #endif
2477  #endif
2478
2479  ideal stdJ = kStd(J,currQuotient,isNotHomog,NULL);
2480  idSkipZeroes( stdJ );
2481
2482  #ifdef SPECTRUM_DEBUG
2483  #ifdef SPECTRUM_PRINT
2484    for( i=0; i<IDELEMS(stdJ); i++ )
2485    {
2486      #ifdef SPECTRUM_IOSTREAM
2487        cout << "        ";
2488      #else
2489        fprintf( stdout,"        " );
2490      #endif
2491
2492      pWrite( stdJ->m[i] );
2493    }
2494  #endif
2495  #endif
2496
2497  idDelete( &J );
2498
2499  // ------------------------------------------
2500  //  check if the  h  has a singularity
2501  // ------------------------------------------
2502
2503  if( hasOne( stdJ ) )
2504  {
2505    // -------------------------------
2506    //  h is smooth in the origin
2507    //  return only the Milnor number
2508    // -------------------------------
2509
2510    *L = (lists)omAllocBin( slists_bin);
2511    (*L)->Init( 1 );
2512    (*L)->m[0].rtyp = INT_CMD;    //  milnor number
2513    /* (*L)->m[0].data = (void*)0;a  -- done by Init */
2514
2515    return  spectrumNoSingularity;
2516  }
2517
2518  // ------------------------------------------
2519  //  check if the singularity  h  is isolated
2520  // ------------------------------------------
2521
2522  for( i=pVariables; i>0; i-- )
2523  {
2524    if( hasAxis( stdJ,i )==FALSE )
2525    {
2526      return  spectrumNotIsolated;
2527    }
2528  }
2529
2530  // ------------------------------------------
2531  //  compute the highest corner  hc  of  stdJ
2532  // ------------------------------------------
2533
2534  #ifdef SPECTRUM_DEBUG
2535  #ifdef SPECTRUM_PRINT
2536  #ifdef SPECTRUM_IOSTREAM
2537    cout << "\n    computing the highest corner...\n";
2538  #else
2539    fprintf( stdout,"\n    computing the highest corner...\n" );
2540  #endif
2541  #endif
2542  #endif
2543
2544  poly hc = (poly)NULL;
2545
2546  scComputeHC( stdJ,currQuotient, 0,hc );
2547
2548  if( hc!=(poly)NULL )
2549  {
2550    pGetCoeff(hc) = nInit(1);
2551
2552    for( i=pVariables; i>0; i-- )
2553    {
2554      if( pGetExp( hc,i )>0 ) pDecrExp( hc,i );
2555    }
2556    pSetm( hc );
2557  }
2558  else
2559  {
2560    return  spectrumNoHC;
2561  }
2562
2563  #ifdef SPECTRUM_DEBUG
2564  #ifdef SPECTRUM_PRINT
2565  #ifdef SPECTRUM_IOSTREAM
2566    cout << "       ";
2567  #else
2568    fprintf( stdout,"       " );
2569  #endif
2570    pWrite( hc );
2571  #endif
2572  #endif
2573
2574  // ----------------------------------------
2575  //  compute the Newton polygon  nph  of  h
2576  // ----------------------------------------
2577
2578  #ifdef SPECTRUM_DEBUG
2579  #ifdef SPECTRUM_PRINT
2580  #ifdef SPECTRUM_IOSTREAM
2581    cout << "\n    computing the newton polygon...\n";
2582  #else
2583    fprintf( stdout,"\n    computing the newton polygon...\n" );
2584  #endif
2585  #endif
2586  #endif
2587
2588  newtonPolygon nph( h );
2589
2590  #ifdef SPECTRUM_DEBUG
2591  #ifdef SPECTRUM_PRINT
2592    cout << nph;
2593  #endif
2594  #endif
2595
2596  // -----------------------------------------------
2597  //  compute the weight corner  wc  of  (stdj,nph)
2598  // -----------------------------------------------
2599
2600  #ifdef SPECTRUM_DEBUG
2601  #ifdef SPECTRUM_PRINT
2602  #ifdef SPECTRUM_IOSTREAM
2603    cout << "\n    computing the weight corner...\n";
2604  #else
2605    fprintf( stdout,"\n    computing the weight corner...\n" );
2606  #endif
2607  #endif
2608  #endif
2609
2610  poly    wc = ( fast==0 ? pCopy( hc ) :
2611               ( fast==1 ? computeWC( nph,(Rational)pVariables ) :
2612              /* fast==2 */computeWC( nph,((Rational)pVariables)/(Rational)2 ) ) );
2613
2614  #ifdef SPECTRUM_DEBUG
2615  #ifdef SPECTRUM_PRINT
2616  #ifdef SPECTRUM_IOSTREAM
2617    cout << "        ";
2618  #else
2619    fprintf( stdout,"        " );
2620  #endif
2621    pWrite( wc );
2622  #endif
2623  #endif
2624
2625  // -------------
2626  //  compute  NF
2627  // -------------
2628
2629  #ifdef SPECTRUM_DEBUG
2630  #ifdef SPECTRUM_PRINT
2631  #ifdef SPECTRUM_IOSTREAM
2632    cout << "\n    computing NF...\n" << endl;
2633  #else
2634    fprintf( stdout,"\n    computing NF...\n" );
2635  #endif
2636  #endif
2637  #endif
2638
2639  spectrumPolyList NF( &nph );
2640
2641  computeNF( stdJ,hc,wc,&NF );
2642
2643  #ifdef SPECTRUM_DEBUG
2644  #ifdef SPECTRUM_PRINT
2645    cout << NF;
2646  #ifdef SPECTRUM_IOSTREAM
2647    cout << endl;
2648  #else
2649    fprintf( stdout,"\n" );
2650  #endif
2651  #endif
2652  #endif
2653
2654  // ----------------------------
2655  //  compute the spectrum of  h
2656  // ----------------------------
2657
2658  return  NF.spectrum( L,fast );
2659}
2660
2661// ----------------------------------------------------------------------------
2662//  this procedure is called from the interpreter
2663// ----------------------------------------------------------------------------
2664//  first  = polynomial
2665//  result = list of spectrum numbers
2666// ----------------------------------------------------------------------------
2667
2668BOOLEAN spectrumProc( leftv result,leftv first )
2669{
2670  spectrumState state = spectrumOK;
2671
2672  // -------------------
2673  //  check consistency
2674  // -------------------
2675
2676  //  check for a local ring
2677
2678  if( !ringIsLocal( ) )
2679  {
2680    WerrorS( "only works for local orderings" );
2681    state = spectrumWrongRing;
2682  }
2683
2684  //  no quotient rings are allowed
2685
2686  else if( currRing->qideal != NULL )
2687  {
2688    WerrorS( "does not work in quotient rings" );
2689    state = spectrumWrongRing;
2690  }
2691  else
2692  {
2693    lists   L    = (lists)NULL;
2694    int     flag = 1; // weight corner optimization is safe
2695
2696    state = spectrumCompute( (poly)first->Data( ),&L,flag );
2697
2698    if( state==spectrumOK )
2699    {
2700      result->rtyp = LIST_CMD;
2701      result->data = (char*)L;
2702    }
2703    else
2704    {
2705      spectrumPrintError(state);
2706    }
2707  }
2708
2709  return  (state!=spectrumOK);
2710}
2711
2712// ----------------------------------------------------------------------------
2713//  this procedure is called from the interpreter
2714// ----------------------------------------------------------------------------
2715//  first  = polynomial
2716//  result = list of spectrum numbers
2717// ----------------------------------------------------------------------------
2718
2719BOOLEAN spectrumfProc( leftv result,leftv first )
2720{
2721  spectrumState state = spectrumOK;
2722
2723  // -------------------
2724  //  check consistency
2725  // -------------------
2726
2727  //  check for a local polynomial ring
2728
2729  if( currRing->OrdSgn != -1 )
2730  // ?? HS: the test above is also true for k[x][[y]], k[[x]][y]
2731  // or should we use:
2732  //if( !ringIsLocal( ) )
2733  {
2734    WerrorS( "only works for local orderings" );
2735    state = spectrumWrongRing;
2736  }
2737  else if( currRing->qideal != NULL )
2738  {
2739    WerrorS( "does not work in quotient rings" );
2740    state = spectrumWrongRing;
2741  }
2742  else
2743  {
2744    lists   L    = (lists)NULL;
2745    int     flag = 2; // symmetric optimization
2746
2747    state = spectrumCompute( (poly)first->Data( ),&L,flag );
2748
2749    if( state==spectrumOK )
2750    {
2751      result->rtyp = LIST_CMD;
2752      result->data = (char*)L;
2753    }
2754    else
2755    {
2756      spectrumPrintError(state);
2757    }
2758  }
2759
2760  return  (state!=spectrumOK);
2761}
2762
2763// ----------------------------------------------------------------------------
2764//  check if a list is a spectrum
2765//  check for:
2766//      list has 6 elements
2767//      1st element is int (mu=Milnor number)
2768//      2nd element is int (pg=geometrical genus)
2769//      3rd element is int (n =number of different spectrum numbers)
2770//      4th element is intvec (num=numerators)
2771//      5th element is intvec (den=denomiantors)
2772//      6th element is intvec (mul=multiplicities)
2773//      exactly n numerators
2774//      exactly n denominators
2775//      exactly n multiplicities
2776//      mu>0
2777//      pg>=0
2778//      n>0
2779//      num>0
2780//      den>0
2781//      mul>0
2782//      symmetriy with respect to numberofvariables/2
2783//      monotony
2784//      mu = sum of all multiplicities
2785//      pg = sum of all multiplicities where num/den<=1
2786// ----------------------------------------------------------------------------
2787
2788semicState  list_is_spectrum( lists l )
2789{
2790    // -------------------
2791    //  check list length
2792    // -------------------
2793
2794    if( l->nr < 5 )
2795    {
2796        return  semicListTooShort;
2797    }
2798    else if( l->nr > 5 )
2799    {
2800        return  semicListTooLong;
2801    }
2802
2803    // -------------
2804    //  check types
2805    // -------------
2806
2807    if( l->m[0].rtyp != INT_CMD )
2808    {
2809        return  semicListFirstElementWrongType;
2810    }
2811    else if( l->m[1].rtyp != INT_CMD )
2812    {
2813        return  semicListSecondElementWrongType;
2814    }
2815    else if( l->m[2].rtyp != INT_CMD )
2816    {
2817        return  semicListThirdElementWrongType;
2818    }
2819    else if( l->m[3].rtyp != INTVEC_CMD )
2820    {
2821        return  semicListFourthElementWrongType;
2822    }
2823    else if( l->m[4].rtyp != INTVEC_CMD )
2824    {
2825        return  semicListFifthElementWrongType;
2826    }
2827    else if( l->m[5].rtyp != INTVEC_CMD )
2828    {
2829        return  semicListSixthElementWrongType;
2830    }
2831
2832    // -------------------------
2833    //  check number of entries
2834    // -------------------------
2835
2836    int     mu = (int)(l->m[0].Data( ));
2837    int     pg = (int)(l->m[1].Data( ));
2838    int     n  = (int)(l->m[2].Data( ));
2839
2840    if( n <= 0 )
2841    {
2842        return  semicListNNegative;
2843    }
2844
2845    intvec  *num = (intvec*)l->m[3].Data( );
2846    intvec  *den = (intvec*)l->m[4].Data( );
2847    intvec  *mul = (intvec*)l->m[5].Data( );
2848
2849    if( n != num->length( ) )
2850    {
2851        return  semicListWrongNumberOfNumerators;
2852    }
2853    else if( n != den->length( ) )
2854    {
2855        return  semicListWrongNumberOfDenominators;
2856    }
2857    else if( n != mul->length( ) )
2858    {
2859        return  semicListWrongNumberOfMultiplicities;
2860    }
2861
2862    // --------
2863    //  values
2864    // --------
2865
2866    if( mu <= 0 )
2867    {
2868        return  semicListMuNegative;
2869    }
2870    if( pg < 0 )
2871    {
2872        return  semicListPgNegative;
2873    }
2874
2875    int i;
2876
2877    for( i=0; i<n; i++ )
2878    {
2879        if( (*num)[i] <= 0 )
2880        {
2881            return  semicListNumNegative;
2882        }
2883        if( (*den)[i] <= 0 )
2884        {
2885            return  semicListDenNegative;
2886        }
2887        if( (*mul)[i] <= 0 )
2888        {
2889            return  semicListMulNegative;
2890        }
2891    }
2892
2893    // ----------------
2894    //  check symmetry
2895    // ----------------
2896
2897    int     j;
2898
2899    for( i=0, j=n-1; i<=j; i++,j-- )
2900    {
2901        if( (*num)[i] != pVariables*((*den)[i]) - (*num)[j] ||
2902            (*den)[i] != (*den)[j] ||
2903            (*mul)[i] != (*mul)[j] )
2904        {
2905            return  semicListNotSymmetric;
2906        }
2907    }
2908
2909    // ----------------
2910    //  check monotony
2911    // ----------------
2912
2913    for( i=0, j=1; i<n/2; i++,j++ )
2914    {
2915        if( (*num)[i]*(*den)[j] >= (*num)[j]*(*den)[i] )
2916        {
2917            return  semicListNotMonotonous;
2918        }
2919    }
2920
2921    // ---------------------
2922    //  check Milnor number
2923    // ---------------------
2924
2925    for( mu=0, i=0; i<n; i++ )
2926    {
2927        mu += (*mul)[i];
2928    }
2929
2930    if( mu != (int)(l->m[0].Data( )) )
2931    {
2932        return  semicListMilnorWrong;
2933    }
2934
2935    // -------------------------
2936    //  check geometrical genus
2937    // -------------------------
2938
2939    for( pg=0, i=0; i<n; i++ )
2940    {
2941        if( (*num)[i]<=(*den)[i] )
2942        {
2943            pg += (*mul)[i];
2944        }
2945    }
2946
2947    if( pg != (int)(l->m[1].Data( )) )
2948    {
2949        return  semicListPGWrong;
2950    }
2951
2952    return  semicOK;
2953}
2954
2955// ----------------------------------------------------------------------------
2956//  this procedure is called from the interpreter
2957// ----------------------------------------------------------------------------
2958//  first  = list of spectrum numbers
2959//  second = list of spectrum numbers
2960//  result = sum of the two lists
2961// ----------------------------------------------------------------------------
2962
2963BOOLEAN spaddProc( leftv result,leftv first,leftv second )
2964{
2965    semicState  state;
2966
2967    // -----------------
2968    //  check arguments
2969    // -----------------
2970
2971    lists l1 = (lists)first->Data( );
2972    lists l2 = (lists)second->Data( );
2973
2974    if( (state=list_is_spectrum( l1 )) != semicOK )
2975    {
2976        WerrorS( "first argument is not a spectrum:" );
2977        list_error( state );
2978    }
2979    else if( (state=list_is_spectrum( l2 )) != semicOK )
2980    {
2981        WerrorS( "second argument is not a spectrum:" );
2982        list_error( state );
2983    }
2984    else
2985    {
2986        spectrum s1( l1 );
2987        spectrum s2( l2 );
2988        spectrum sum( s1+s2 );
2989
2990        result->rtyp = LIST_CMD;
2991        result->data = (char*)(sum.thelist( ));
2992    }
2993
2994    return  (state!=semicOK);
2995}
2996
2997// ----------------------------------------------------------------------------
2998//  this procedure is called from the interpreter
2999// ----------------------------------------------------------------------------
3000//  first  = list of spectrum numbers
3001//  second = integer
3002//  result = the multiple of the first list by the second factor
3003// ----------------------------------------------------------------------------
3004
3005BOOLEAN spmulProc( leftv result,leftv first,leftv second )
3006{
3007    semicState  state;
3008
3009    // -----------------
3010    //  check arguments
3011    // -----------------
3012
3013    lists   l = (lists)first->Data( );
3014    int     k = (int)second->Data( );
3015
3016    if( (state=list_is_spectrum( l ))!=semicOK )
3017    {
3018        WerrorS( "first argument is not a spectrum" );
3019        list_error( state );
3020    }
3021    else if( k < 0 )
3022    {
3023        WerrorS( "second argument should be positive" );
3024        state = semicMulNegative;
3025    }
3026    else
3027    {
3028        spectrum s( l );
3029        spectrum product( k*s );
3030
3031        result->rtyp = LIST_CMD;
3032        result->data = (char*)product.thelist( );
3033    }
3034
3035    return  (state!=semicOK);
3036}
3037
3038// ----------------------------------------------------------------------------
3039//  this procedure is called from the interpreter
3040// ----------------------------------------------------------------------------
3041//  first  = list of spectrum numbers
3042//  second = list of spectrum numbers
3043//  result = semicontinuity index
3044// ----------------------------------------------------------------------------
3045
3046BOOLEAN    semicProc3   ( leftv res,leftv u,leftv v,leftv w )
3047{
3048  semicState  state;
3049  BOOLEAN qh=(((int)w->Data())==1);
3050
3051  // -----------------
3052  //  check arguments
3053  // -----------------
3054
3055  lists l1 = (lists)u->Data( );
3056  lists l2 = (lists)v->Data( );
3057
3058  if( (state=list_is_spectrum( l1 ))!=semicOK )
3059  {
3060    WerrorS( "first argument is not a spectrum" );
3061    list_error( state );
3062  }
3063  else if( (state=list_is_spectrum( l2 ))!=semicOK )
3064  {
3065    WerrorS( "second argument is not a spectrum" );
3066    list_error( state );
3067  }
3068  else
3069  {
3070    spectrum s1( l1 );
3071    spectrum s2( l2 );
3072
3073    res->rtyp = INT_CMD;
3074    if (qh)
3075      res->data = (void*)(s1.mult_spectrumh( s2 ));
3076    else
3077      res->data = (void*)(s1.mult_spectrum( s2 ));
3078  }
3079
3080  // -----------------
3081  //  check status
3082  // -----------------
3083
3084  return  (state!=semicOK);
3085}
3086BOOLEAN    semicProc   ( leftv res,leftv u,leftv v )
3087{
3088  sleftv tmp;
3089  memset(&tmp,0,sizeof(tmp));
3090  tmp.rtyp=INT_CMD;
3091  /* tmp.data = (void *)0;  -- done by memset */
3092
3093  return  semicProc3(res,u,v,&tmp);
3094}
3095// from splist.cc
3096// ----------------------------------------------------------------------------
3097//  Compute the spectrum of a  spectrumPolyList
3098// ----------------------------------------------------------------------------
3099
3100spectrumState   spectrumPolyList::spectrum( lists *L,int fast )
3101{
3102    spectrumPolyNode  **node = &root;
3103    spectrumPolyNode  *search;
3104
3105    poly              f,tmp;
3106    int               found,cmp;
3107
3108    Rational smax( ( fast==0 ? 0 : pVariables ),
3109                   ( fast==2 ? 2 : 1 ) );
3110
3111    Rational weight_prev( 0,1 );
3112
3113    int     mu = 0;          // the milnor number
3114    int     pg = 0;          // the geometrical genus
3115    int     n  = 0;          // number of different spectral numbers
3116    int     z  = 0;          // number of spectral number equal to smax
3117
3118    int     k = 0;
3119
3120    while( (*node)!=(spectrumPolyNode*)NULL &&
3121           ( fast==0 || (*node)->weight<=smax ) )
3122    {
3123        // ---------------------------------------
3124        //  determine the first normal form which
3125        //  contains the monomial  node->mon
3126        // ---------------------------------------
3127
3128        found  = FALSE;
3129        search = *node;
3130
3131        while( search!=(spectrumPolyNode*)NULL && found==FALSE )
3132        {
3133            if( search->nf!=(poly)NULL )
3134            {
3135                f = search->nf;
3136
3137                do
3138                {
3139                    // --------------------------------
3140                    //  look for  (*node)->mon  in   f
3141                    // --------------------------------
3142
3143                    cmp = pCmp( (*node)->mon,f );
3144
3145                    if( cmp<0 )
3146                    {
3147                        f = pNext( f );
3148                    }
3149                    else if( cmp==0 )
3150                    {
3151                        // -----------------------------
3152                        //  we have found a normal form
3153                        // -----------------------------
3154
3155                        found = TRUE;
3156
3157                        //  normalize coefficient
3158
3159                        number inv = nInvers( pGetCoeff( f ) );
3160                        pMult_nn( search->nf,inv );
3161                        nDelete( &inv );
3162
3163                        //  exchange  normal forms
3164
3165                        tmp         = (*node)->nf;
3166                        (*node)->nf = search->nf;
3167                        search->nf  = tmp;
3168                    }
3169                }
3170                while( cmp<0 && f!=(poly)NULL );
3171            }
3172            search = search->next;
3173        }
3174
3175        if( found==FALSE )
3176        {
3177            // ------------------------------------------------
3178            //  the weight of  node->mon  is a spectrum number
3179            // ------------------------------------------------
3180
3181            mu++;
3182
3183            if( (*node)->weight<=(Rational)1 )              pg++;
3184            if( (*node)->weight==smax )           z++;
3185            if( (*node)->weight>weight_prev )     n++;
3186
3187            weight_prev = (*node)->weight;
3188            node = &((*node)->next);
3189        }
3190        else
3191        {
3192            // -----------------------------------------------
3193            //  determine all other normal form which contain
3194            //  the monomial  node->mon
3195            //  replace for  node->mon  its normal form
3196            // -----------------------------------------------
3197
3198            while( search!=(spectrumPolyNode*)NULL )
3199            {
3200                    if( search->nf!=(poly)NULL )
3201                {
3202                    f = search->nf;
3203
3204                    do
3205                    {
3206                        // --------------------------------
3207                        //  look for  (*node)->mon  in   f
3208                        // --------------------------------
3209
3210                        cmp = pCmp( (*node)->mon,f );
3211
3212                        if( cmp<0 )
3213                        {
3214                            f = pNext( f );
3215                        }
3216                        else if( cmp==0 )
3217                        {
3218                            search->nf = pSub( search->nf,
3219                                ppMult_nn( (*node)->nf,pGetCoeff( f ) ) );
3220                            pNorm( search->nf );
3221                        }
3222                    }
3223                    while( cmp<0 && f!=(poly)NULL );
3224                }
3225                search = search->next;
3226            }
3227            delete_node( node );
3228        }
3229
3230    }
3231
3232    // --------------------------------------------------------
3233    //  fast computation exploits the symmetry of the spectrum
3234    // --------------------------------------------------------
3235
3236    if( fast==2 )
3237    {
3238        mu = 2*mu - z;
3239        n  = ( z > 0 ? 2*n - 1 : 2*n );
3240    }
3241
3242    // --------------------------------------------------------
3243    //  compute the spectrum numbers with their multiplicities
3244    // --------------------------------------------------------
3245
3246    intvec            *nom  = new intvec( n );
3247    intvec            *den  = new intvec( n );
3248    intvec            *mult = new intvec( n );
3249
3250    int count         = 0;
3251    int multiplicity  = 1;
3252
3253    for( search=root; search!=(spectrumPolyNode*)NULL &&
3254                     ( fast==0 || search->weight<=smax );
3255                     search=search->next )
3256    {
3257        if( search->next==(spectrumPolyNode*)NULL ||
3258            search->weight<search->next->weight )
3259        {
3260            (*nom) [count] = search->weight.get_num_si( );
3261            (*den) [count] = search->weight.get_den_si( );
3262            (*mult)[count] = multiplicity;
3263
3264            multiplicity=1;
3265            count++;
3266        }
3267        else
3268        {
3269            multiplicity++;
3270        }
3271    }
3272
3273    // --------------------------------------------------------
3274    //  fast computation exploits the symmetry of the spectrum
3275    // --------------------------------------------------------
3276
3277    if( fast==2 )
3278    {
3279        int n1,n2;
3280        for( n1=0, n2=n-1; n1<n2; n1++, n2-- )
3281        {
3282            (*nom) [n2] = pVariables*(*den)[n1]-(*nom)[n1];
3283            (*den) [n2] = (*den)[n1];
3284            (*mult)[n2] = (*mult)[n1];
3285        }
3286    }
3287
3288    // -----------------------------------
3289    //  test if the spectrum is symmetric
3290    // -----------------------------------
3291
3292    if( fast==0 || fast==1 )
3293    {
3294        int symmetric=TRUE;
3295
3296        for( int n1=0, n2=n-1 ; n1<n2 && symmetric==TRUE; n1++, n2-- )
3297        {
3298            if( (*mult)[n1]!=(*mult)[n2] ||
3299                (*den) [n1]!= (*den)[n2] ||
3300                (*nom)[n1]+(*nom)[n2]!=pVariables*(*den) [n1] )
3301            {
3302                symmetric = FALSE;
3303            }
3304        }
3305
3306        if( symmetric==FALSE )
3307        {
3308            // ---------------------------------------------
3309            //  the spectrum is not symmetric => degenerate
3310            //  principal part
3311            // ---------------------------------------------
3312
3313            *L = (lists)omAllocBin( slists_bin);
3314            (*L)->Init( 1 );
3315            (*L)->m[0].rtyp = INT_CMD;    //  milnor number
3316            (*L)->m[0].data = (void*)mu;
3317
3318            return spectrumDegenerate;
3319        }
3320    }
3321
3322    *L = (lists)omAllocBin( slists_bin);
3323
3324    (*L)->Init( 6 );
3325
3326    (*L)->m[0].rtyp = INT_CMD;    //  milnor number
3327    (*L)->m[1].rtyp = INT_CMD;    //  geometrical genus
3328    (*L)->m[2].rtyp = INT_CMD;    //  number of spectrum values
3329    (*L)->m[3].rtyp = INTVEC_CMD; //  nominators
3330    (*L)->m[4].rtyp = INTVEC_CMD; //  denomiantors
3331    (*L)->m[5].rtyp = INTVEC_CMD; //  multiplicities
3332
3333    (*L)->m[0].data = (void*)mu;
3334    (*L)->m[1].data = (void*)pg;
3335    (*L)->m[2].data = (void*)n;
3336    (*L)->m[3].data = (void*)nom;
3337    (*L)->m[4].data = (void*)den;
3338    (*L)->m[5].data = (void*)mult;
3339
3340    return  spectrumOK;
3341}
3342
3343#endif
3344
3345//from mpr_inout.cc
3346extern void nPrint(number n);
3347
3348BOOLEAN loNewtonP( leftv res, leftv arg1 )
3349{
3350  res->data= (void*)loNewtonPolytope( (ideal)arg1->Data() );
3351  return FALSE;
3352}
3353
3354BOOLEAN loSimplex( leftv res, leftv args )
3355{
3356  if ( !(rField_is_long_R()) )
3357  {
3358    WerrorS("Ground field not implemented!");
3359    return TRUE;
3360  }
3361
3362  simplex * LP;
3363  matrix m;
3364
3365  leftv v= args;
3366  if ( v->Typ() != MATRIX_CMD ) // 1: matrix
3367    return TRUE;
3368  else
3369    m= (matrix)(v->CopyD());
3370
3371  LP = new simplex(MATROWS(m),MATCOLS(m));
3372  LP->mapFromMatrix(m);
3373
3374  v= v->next;
3375  if ( v->Typ() != INT_CMD )    // 2: m = number of constraints
3376    return TRUE;
3377  else
3378    LP->m= (int)(v->Data());
3379
3380  v= v->next;
3381  if ( v->Typ() != INT_CMD )    // 3: n = number of variables
3382    return TRUE;
3383  else
3384    LP->n= (int)(v->Data());
3385
3386  v= v->next;
3387  if ( v->Typ() != INT_CMD )    // 4: m1 = number of <= constraints
3388    return TRUE;
3389  else
3390    LP->m1= (int)(v->Data());
3391
3392  v= v->next;
3393  if ( v->Typ() != INT_CMD )    // 5: m2 = number of >= constraints
3394    return TRUE;
3395  else
3396    LP->m2= (int)(v->Data());
3397
3398  v= v->next;
3399  if ( v->Typ() != INT_CMD )    // 6: m3 = number of == constraints
3400    return TRUE;
3401  else
3402    LP->m3= (int)(v->Data());
3403
3404#ifdef mprDEBUG_PROT
3405  Print("m (constraints) %d\n",LP->m);
3406  Print("n (columns) %d\n",LP->n);
3407  Print("m1 (<=) %d\n",LP->m1);
3408  Print("m2 (>=) %d\n",LP->m2);
3409  Print("m3 (==) %d\n",LP->m3);
3410#endif
3411
3412  LP->compute();
3413
3414  lists lres= (lists)omAlloc( sizeof(slists) );
3415  lres->Init( 6 );
3416
3417  lres->m[0].rtyp= MATRIX_CMD; // output matrix
3418  lres->m[0].data=(void*)LP->mapToMatrix(m);
3419
3420  lres->m[1].rtyp= INT_CMD;   // found a solution?
3421  lres->m[1].data=(void*)LP->icase;
3422
3423  lres->m[2].rtyp= INTVEC_CMD;
3424  lres->m[2].data=(void*)LP->posvToIV();
3425
3426  lres->m[3].rtyp= INTVEC_CMD;
3427  lres->m[3].data=(void*)LP->zrovToIV();
3428
3429  lres->m[4].rtyp= INT_CMD;
3430  lres->m[4].data=(void*)LP->m;
3431
3432  lres->m[5].rtyp= INT_CMD;
3433  lres->m[5].data=(void*)LP->n;
3434
3435  res->data= (void*)lres;
3436
3437  return FALSE;
3438}
3439
3440BOOLEAN nuMPResMat( leftv res, leftv arg1, leftv arg2 )
3441{
3442  ideal gls = (ideal)(arg1->Data());
3443  int imtype= (int)arg2->Data();
3444
3445  uResultant::resMatType mtype= determineMType( imtype );
3446
3447  // check input ideal ( = polynomial system )
3448  if ( mprIdealCheck( gls, arg1->Name(), mtype, true ) != mprOk )
3449  {
3450    return TRUE;
3451  }
3452
3453  uResultant *resMat= new uResultant( gls, mtype, false );
3454
3455  res->rtyp = MODUL_CMD;
3456  res->data= (void*)resMat->accessResMat()->getMatrix();
3457
3458  delete resMat;
3459
3460  return FALSE;
3461}
3462
3463BOOLEAN nuLagSolve( leftv res, leftv arg1, leftv arg2, leftv arg3 )
3464{
3465
3466  poly gls;
3467  gls= (poly)(arg1->Data());
3468  int howclean= (int)arg3->Data();
3469
3470  if ( !(rField_is_R() ||
3471         rField_is_Q() ||
3472         rField_is_long_R() ||
3473         rField_is_long_C()) )
3474  {
3475    WerrorS("Ground field not implemented!");
3476    return TRUE;
3477  }
3478
3479  if ( !(rField_is_R()||rField_is_long_R()||rField_is_long_C()) )
3480  {
3481    unsigned long int ii = (unsigned long int)arg2->Data();
3482    setGMPFloatDigits( ii, ii );
3483  }
3484
3485  if ( gls == NULL || pIsConstant( gls ) )
3486  {
3487    WerrorS("Input polynomial is constant!");
3488    return TRUE;
3489  }
3490
3491  int ldummy;
3492  int deg= pLDeg( gls, &ldummy, currRing );
3493  //  int deg= pDeg( gls );
3494  int len= pLength( gls );
3495  int i,vpos;
3496  poly piter;
3497  lists elist;
3498  lists rlist;
3499
3500  elist= (lists)omAlloc( sizeof(slists) );
3501  elist->Init( 0 );
3502
3503  if ( pVariables > 1 )
3504  {
3505    piter= gls;
3506    for ( i= 1; i <= pVariables; i++ )
3507      if ( pGetExp( piter, i ) )
3508      {
3509        vpos= i;
3510        break;
3511      }
3512    while ( piter )
3513    {
3514      for ( i= 1; i <= pVariables; i++ )
3515        if ( (vpos != i) && (pGetExp( piter, i ) != 0) )
3516        {
3517          WerrorS("The input polynomial must be univariate!");
3518          return TRUE;
3519        }
3520      pIter( piter );
3521    }
3522  }
3523
3524  rootContainer * roots= new rootContainer();
3525  number * pcoeffs= (number *)omAlloc( (deg+1) * sizeof( number ) );
3526  piter= gls;
3527  for ( i= deg; i >= 0; i-- )
3528  {
3529    //if ( piter ) Print("deg %d, pDeg(piter) %d\n",i,pTotaldegree(piter));
3530    if ( piter && pTotaldegree(piter) == i )
3531    {
3532      pcoeffs[i]= nCopy( pGetCoeff( piter ) );
3533      //nPrint( pcoeffs[i] );PrintS("  ");
3534      pIter( piter );
3535    }
3536    else
3537    {
3538      pcoeffs[i]= nInit(0);
3539    }
3540  }
3541
3542#ifdef mprDEBUG_PROT
3543  for (i=deg; i >= 0; i--)
3544  {
3545    nPrint( pcoeffs[i] );PrintS("  ");
3546  }
3547  PrintLn();
3548#endif
3549
3550  roots->fillContainer( pcoeffs, NULL, 1, deg, rootContainer::onepoly, 1 );
3551  roots->solver( howclean );
3552
3553  int elem= roots->getAnzRoots();
3554  char *out;
3555  char *dummy;
3556  int j;
3557
3558  rlist= (lists)omAlloc( sizeof(slists) );
3559  rlist->Init( elem );
3560
3561  if (rField_is_long_C())
3562  {
3563    for ( j= 0; j < elem; j++ )
3564    {
3565      rlist->m[j].rtyp=NUMBER_CMD;
3566      rlist->m[j].data=(void *)nCopy((number)(roots->getRoot(j)));
3567      //rlist->m[j].data=(void *)(number)(roots->getRoot(j));
3568    }
3569  }
3570  else
3571  {
3572    for ( j= 0; j < elem; j++ )
3573    {
3574      dummy = complexToStr( (*roots)[j], gmp_output_digits );
3575      rlist->m[j].rtyp=STRING_CMD;
3576      rlist->m[j].data=(void *)dummy;
3577    }
3578  }
3579
3580  elist->Clean();
3581  //omFreeSize( (ADDRESS) elist, sizeof(slists) );
3582
3583  // this is (via fillContainer) the same data as in root
3584  //for ( i= deg; i >= 0; i-- ) nDelete( &pcoeffs[i] );
3585  //omFreeSize( (ADDRESS) pcoeffs, (deg+1) * sizeof( number ) );
3586
3587  delete roots;
3588
3589  res->rtyp= LIST_CMD;
3590  res->data= (void*)rlist;
3591
3592  return FALSE;
3593}
3594
3595BOOLEAN nuVanderSys( leftv res, leftv arg1, leftv arg2, leftv arg3)
3596{
3597  int i;
3598  ideal p,w;
3599  p= (ideal)arg1->Data();
3600  w= (ideal)arg2->Data();
3601
3602  // w[0] = f(p^0)
3603  // w[1] = f(p^1)
3604  // ...
3605  // p can be a vector of numbers (multivariate polynom)
3606  //   or one number (univariate polynom)
3607  // tdg = deg(f)
3608
3609  int n= IDELEMS( p );
3610  int m= IDELEMS( w );
3611  int tdg= (int)arg3->Data();
3612
3613  res->data= (void*)NULL;
3614
3615  // check the input
3616  if ( tdg < 1 )
3617  {
3618    WerrorS("Last input parameter must be > 0!");
3619    return TRUE;
3620  }
3621  if ( n != pVariables )
3622  {
3623    Werror("Size of first input ideal must be equal to %d!",pVariables);
3624    return TRUE;
3625  }
3626  if ( m != (int)pow((double)tdg+1,(double)n) )
3627  {
3628    Werror("Size of second input ideal must be equal to %d!",
3629      (int)pow((double)tdg+1,(double)n));
3630    return TRUE;
3631  }
3632  if ( !(rField_is_Q() /* ||
3633         rField_is_R() || rField_is_long_R() ||
3634         rField_is_long_C()*/ ) )
3635         {
3636    WerrorS("Ground field not implemented!");
3637    return TRUE;
3638  }
3639
3640  number tmp;
3641  number *pevpoint= (number *)omAlloc( n * sizeof( number ) );
3642  for ( i= 0; i < n; i++ )
3643  {
3644    pevpoint[i]=nInit(0);
3645    if (  (p->m)[i] )
3646    {
3647      tmp = pGetCoeff( (p->m)[i] );
3648      if ( nIsZero(tmp) || nIsOne(tmp) || nIsMOne(tmp) )
3649      {
3650        omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
3651        WerrorS("Elements of first input ideal must not be equal to -1, 0, 1!");
3652        return TRUE;
3653      }
3654    } else tmp= NULL;
3655    if ( !nIsZero(tmp) )
3656    {
3657      if ( !pIsConstant((p->m)[i]))
3658      {
3659        omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
3660        WerrorS("Elements of first input ideal must be numbers!");
3661        return TRUE;
3662      }
3663      pevpoint[i]= nCopy( tmp );
3664    }
3665  }
3666
3667  number *wresults= (number *)omAlloc( m * sizeof( number ) );
3668  for ( i= 0; i < m; i++ )
3669  {
3670    wresults[i]= nInit(0);
3671    if ( (w->m)[i] && !nIsZero(pGetCoeff((w->m)[i])) )
3672    {
3673      if ( !pIsConstant((w->m)[i]))
3674      {
3675        omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
3676        omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
3677        WerrorS("Elements of second input ideal must be numbers!");
3678        return TRUE;
3679      }
3680      wresults[i]= nCopy(pGetCoeff((w->m)[i]));
3681    }
3682  }
3683
3684  vandermonde vm( m, n, tdg, pevpoint, FALSE );
3685  number *ncpoly= vm.interpolateDense( wresults );
3686  // do not free ncpoly[]!!
3687  poly rpoly= vm.numvec2poly( ncpoly );
3688
3689  omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
3690  omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
3691
3692  res->data= (void*)rpoly;
3693  return FALSE;
3694}
3695
3696BOOLEAN nuUResSolve( leftv res, leftv args )
3697{
3698  leftv v= args;
3699
3700  ideal gls;
3701  int imtype;
3702  int howclean;
3703
3704  // get ideal
3705  if ( v->Typ() != IDEAL_CMD )
3706    return TRUE;
3707  else gls= (ideal)(v->Data());
3708  v= v->next;
3709
3710  // get resultant matrix type to use (0,1)
3711  if ( v->Typ() != INT_CMD )
3712    return TRUE;
3713  else imtype= (int)v->Data();
3714  v= v->next;
3715
3716  // get and set precision in digits ( > 0 )
3717  if ( v->Typ() != INT_CMD )
3718    return TRUE;
3719  else if ( !(rField_is_R()||rField_is_long_R()||rField_is_long_C()) )
3720  {
3721    unsigned long int ii=(unsigned long int)v->Data();
3722    setGMPFloatDigits( ii, ii );
3723  }
3724  v= v->next;
3725
3726  // get interpolation steps (0,1,2)
3727  if ( v->Typ() != INT_CMD )
3728    return TRUE;
3729  else howclean= (int)v->Data();
3730
3731  uResultant::resMatType mtype= determineMType( imtype );
3732  int i,c,count;
3733  lists listofroots= NULL;
3734  lists emptylist;
3735  number smv= NULL;
3736  BOOLEAN interpolate_det= (mtype==uResultant::denseResMat)?TRUE:FALSE;
3737
3738  //emptylist= (lists)omAlloc( sizeof(slists) );
3739  //emptylist->Init( 0 );
3740
3741  //res->rtyp = LIST_CMD;
3742  //res->data= (void *)emptylist;
3743
3744  // check input ideal ( = polynomial system )
3745  if ( mprIdealCheck( gls, args->Name(), mtype ) != mprOk )
3746  {
3747    return TRUE;
3748  }
3749
3750  uResultant * ures;
3751  rootContainer ** iproots;
3752  rootContainer ** muiproots;
3753  rootArranger * arranger;
3754
3755  // main task 1: setup of resultant matrix
3756  ures= new uResultant( gls, mtype );
3757  if ( ures->accessResMat()->initState() != resMatrixBase::ready )
3758  {
3759    WerrorS("Error occurred during matrix setup!");
3760    return TRUE;
3761  }
3762
3763  // if dense resultant, check if minor nonsingular
3764  if ( mtype == uResultant::denseResMat )
3765  {
3766    smv= ures->accessResMat()->getSubDet();
3767#ifdef mprDEBUG_PROT
3768    PrintS("// Determinant of submatrix: ");nPrint(smv);PrintLn();
3769#endif
3770    if ( nIsZero(smv) )
3771    {
3772      WerrorS("Unsuitable input ideal: Minor of resultant matrix is singular!");
3773      return TRUE;
3774    }
3775  }
3776
3777  // main task 2: Interpolate specialized resultant polynomials
3778  if ( interpolate_det )
3779    iproots= ures->interpolateDenseSP( false, smv );
3780  else
3781    iproots= ures->specializeInU( false, smv );
3782
3783  // main task 3: Interpolate specialized resultant polynomials
3784  if ( interpolate_det )
3785    muiproots= ures->interpolateDenseSP( true, smv );
3786  else
3787    muiproots= ures->specializeInU( true, smv );
3788
3789#ifdef mprDEBUG_PROT
3790  c= iproots[0]->getAnzElems();
3791  for (i=0; i < c; i++) pWrite(iproots[i]->getPoly());
3792  c= muiproots[0]->getAnzElems();
3793  for (i=0; i < c; i++) pWrite(muiproots[i]->getPoly());
3794#endif
3795
3796  // main task 4: Compute roots of specialized polys and match them up
3797  arranger= new rootArranger( iproots, muiproots, howclean );
3798  arranger->solve_all();
3799
3800  // get list of roots
3801  if ( arranger->success() )
3802  {
3803    arranger->arrange();
3804    listofroots= arranger->listOfRoots( gmp_output_digits );
3805  }
3806  else
3807  {
3808    WerrorS("Solver was unable to find any roots!");
3809    return TRUE;
3810  }
3811
3812  // free everything
3813  count= iproots[0]->getAnzElems();
3814  for (i=0; i < count; i++) delete iproots[i];
3815  omFreeSize( (ADDRESS) iproots, count * sizeof(rootContainer*) );
3816  count= muiproots[0]->getAnzElems();
3817  for (i=0; i < count; i++) delete muiproots[i];
3818  omFreeSize( (ADDRESS) muiproots, count * sizeof(rootContainer*) );
3819
3820  delete ures;
3821  delete arranger;
3822  nDelete( &smv );
3823
3824  res->data= (void *)listofroots;
3825
3826  //emptylist->Clean();
3827  //  omFreeSize( (ADDRESS) emptylist, sizeof(slists) );
3828
3829  return FALSE;
3830}
3831
3832// from mpr_numeric.cc
3833lists rootArranger::listOfRoots( const unsigned int oprec )
3834{
3835  int i,j,tr;
3836  int count= roots[0]->getAnzRoots(); // number of roots
3837  int elem= roots[0]->getAnzElems();  // number of koordinates per root
3838
3839  lists listofroots= (lists)omAlloc( sizeof(slists) ); // must be done this way!
3840
3841  if ( found_roots )
3842  {
3843    listofroots->Init( count );
3844
3845    for (i=0; i < count; i++)
3846    {
3847      lists onepoint= (lists)omAlloc(sizeof(slists)); // must be done this way!
3848      onepoint->Init(elem);
3849      for ( j= 0; j < elem; j++ )
3850      {
3851        if ( !rField_is_long_C() )
3852        {
3853          onepoint->m[j].rtyp=STRING_CMD;
3854          onepoint->m[j].data=(void *)complexToStr((*roots[j])[i],oprec);
3855        }
3856        else
3857        {
3858          onepoint->m[j].rtyp=NUMBER_CMD;
3859          onepoint->m[j].data=(void *)nCopy((number)(roots[j]->getRoot(i)));
3860        }
3861        onepoint->m[j].next= NULL;
3862        onepoint->m[j].name= NULL;
3863      }
3864      listofroots->m[i].rtyp=LIST_CMD;
3865      listofroots->m[i].data=(void *)onepoint;
3866      listofroots->m[j].next= NULL;
3867      listofroots->m[j].name= NULL;
3868    }
3869
3870  }
3871  else
3872  {
3873    listofroots->Init( 0 );
3874  }
3875
3876  return listofroots;
3877}
3878
3879#ifdef PDEBUG
3880
3881#if (OM_TRACK > 2) && defined(OM_TRACK_CUSTOM)
3882
3883void p_SetRingOfPoly(poly p, ring r)
3884{
3885  while (p != NULL)
3886  {
3887    p_SetRingOfLm(p, r);
3888    pIter(p);
3889  }
3890}
3891
3892void p_SetRingOfIdeal(ideal id, ring r)
3893{
3894  if (id == NULL) return;
3895
3896  int i, n = id->ncols*id->nrows;
3897
3898  for (i=0; i<n; i++)
3899  {
3900    p_SetRingOfPoly(id->m[i], r);
3901  }
3902}
3903
3904void p_SetRingOfList(lists L, ring r)
3905{
3906  int i;
3907  for (i=0; i<L->nr; i++)
3908  {
3909    p_SetRingOfLeftv(&(L->m[i]), r);
3910  }
3911}
3912
3913void p_SetRingOfCommand(command cmd, ring r)
3914{
3915  if (cmd->op == PROC_CMD && cmd->argc == 2)
3916    p_SetRingOfLeftv(&(cmd->arg2), r);
3917  else if (cmd->argc > 0)
3918  {
3919    p_SetRingOfLeftv(&(cmd->arg1), r);
3920    if (cmd->argc > 1)
3921    {
3922      p_SetRingOfLeftv(&(cmd->arg2), r);
3923      if (cmd->argc > 2)
3924        p_SetRingOfLeftv(&(cmd->arg3), r);
3925    }
3926  }
3927}
3928
3929void p_SetRingOfLeftv(leftv l, ring r)
3930{
3931  while (l != NULL)
3932  {
3933    switch(l->rtyp)
3934    {
3935        case POLY_CMD:
3936        case VECTOR_CMD:
3937          p_SetRingOfPoly((poly) l->data, r);
3938      break;
3939
3940      case IDEAL_CMD:
3941      case MODUL_CMD:
3942      case MATRIX_CMD:
3943      case MAP_CMD:
3944        p_SetRingOfIdeal((ideal) l->data, r);
3945        break;
3946
3947        case LIST_CMD:
3948          p_SetRingOfList((lists) l->data, r);
3949          break;
3950
3951        case COMMAND:
3952          p_SetRingOfCommand((command)l->data, r);
3953        default:
3954          break;
3955    }
3956    l = l->next;
3957  }
3958}
3959#endif // (OM_TRACK > 2) && defined(OM_TRACK_CUSTOM)
3960
3961#endif // PDEBUG
3962
3963// from ring.cc
3964void rSetHdl(idhdl h)
3965{
3966  int i;
3967  ring rg = NULL;
3968  if (h!=NULL)
3969  {
3970//   Print(" new ring:%s (l:%d)\n",IDID(h),IDLEV(h));
3971    rg = IDRING(h);
3972    omCheckAddrSize((ADDRESS)h,sizeof(idrec));
3973    if (IDID(h))  // OB: ????
3974      omCheckAddr((ADDRESS)IDID(h));
3975    rTest(rg);
3976  }
3977
3978  // clean up history
3979  if (sLastPrinted.RingDependend())
3980  {
3981    sLastPrinted.CleanUp();
3982    memset(&sLastPrinted,0,sizeof(sleftv));
3983  }
3984
3985   /*------------ change the global ring -----------------------*/
3986  rChangeCurrRing(rg);
3987  currRingHdl = h;
3988}
3989
3990BOOLEAN rSleftvOrdering2Ordering(sleftv *ord, ring R)
3991{
3992  int last = 0, o=0, n = 1, i=0, typ = 1, j;
3993  sleftv *sl = ord;
3994
3995  // determine nBlocks
3996  while (sl!=NULL)
3997  {
3998    intvec *iv = (intvec *)(sl->data);
3999    if (((*iv)[1]==ringorder_c)||((*iv)[1]==ringorder_C)) i++;
4000    else if ((*iv)[1]==ringorder_L)
4001    {
4002      R->bitmask=(*iv)[2];
4003      n--;
4004    }
4005    else if ((*iv)[1]!=ringorder_a) o++;
4006    n++;
4007    sl=sl->next;
4008  }
4009  // check whether at least one real ordering
4010  if (o==0)
4011  {
4012    WerrorS("invalid combination of orderings");
4013    return TRUE;
4014  }
4015  // if no c/C ordering is given, increment n
4016  if (i==0) n++;
4017  else if (i != 1)
4018  {
4019    // throw error if more than one is given
4020    WerrorS("more than one ordering c/C specified");
4021    return TRUE;
4022  }
4023
4024  // initialize fields of R
4025  R->order=(int *)omAlloc0(n*sizeof(int));
4026  R->block0=(int *)omAlloc0(n*sizeof(int));
4027  R->block1=(int *)omAlloc0(n*sizeof(int));
4028  R->wvhdl=(int**)omAlloc0(n*sizeof(int_ptr));
4029
4030  // init order, so that rBlocks works correctly
4031  for (j=0; j < n-1; j++)
4032    R->order[j] = (int) ringorder_unspec;
4033  // set last _C order, if no c/C order was given
4034  if (i == 0) R->order[n-2] = ringorder_C;
4035
4036  /* init orders */
4037  sl=ord;
4038  n=-1;
4039  while (sl!=NULL)
4040  {
4041    intvec *iv;
4042    iv = (intvec *)(sl->data);
4043    if ((*iv)[1]!=ringorder_L)
4044    {
4045      n++;
4046
4047      /* the format of an ordering:
4048       *  iv[0]: factor
4049       *  iv[1]: ordering
4050       *  iv[2..end]: weights
4051       */
4052      R->order[n] = (*iv)[1];
4053      switch ((*iv)[1])
4054      {
4055          case ringorder_ws:
4056          case ringorder_Ws:
4057            typ=-1;
4058          case ringorder_wp:
4059          case ringorder_Wp:
4060            R->wvhdl[n]=(int*)omAlloc((iv->length()-1)*sizeof(int));
4061            for (i=2; i<iv->length(); i++)
4062              R->wvhdl[n][i-2] = (*iv)[i];
4063            R->block0[n] = last+1;
4064            last += iv->length()-2;
4065            R->block1[n] = last;
4066            break;
4067          case ringorder_ls:
4068          case ringorder_ds:
4069          case ringorder_Ds:
4070            typ=-1;
4071          case ringorder_lp:
4072          case ringorder_dp:
4073          case ringorder_Dp:
4074          case ringorder_rp:
4075            R->block0[n] = last+1;
4076            if (iv->length() == 3) last+=(*iv)[2];
4077            else last += (*iv)[0];
4078            R->block1[n] = last;
4079            if (rCheckIV(iv)) return TRUE;
4080            break;
4081          case ringorder_S:
4082          case ringorder_c:
4083          case ringorder_C:
4084            if (rCheckIV(iv)) return TRUE;
4085            break;
4086          case ringorder_aa:
4087          case ringorder_a:
4088            R->block0[n] = last+1;
4089            R->block1[n] = si_min(last+iv->length()-2 , R->N);
4090            R->wvhdl[n] = (int*)omAlloc((iv->length()-1)*sizeof(int));
4091            for (i=2; i<iv->length(); i++)
4092            {
4093              R->wvhdl[n][i-2]=(*iv)[i];
4094              if ((*iv)[i]<0) typ=-1;
4095            }
4096            break;
4097          case ringorder_M:
4098          {
4099            int Mtyp=rTypeOfMatrixOrder(iv);
4100            if (Mtyp==0) return TRUE;
4101            if (Mtyp==-1) typ = -1;
4102
4103            R->wvhdl[n] =( int *)omAlloc((iv->length()-1)*sizeof(int));
4104            for (i=2; i<iv->length();i++)
4105              R->wvhdl[n][i-2]=(*iv)[i];
4106
4107            R->block0[n] = last+1;
4108            last += (int)sqrt((double)(iv->length()-2));
4109            R->block1[n] = last;
4110            break;
4111          }
4112
4113          case ringorder_no:
4114            R->order[n] = ringorder_unspec;
4115            return TRUE;
4116
4117          default:
4118            Werror("Internal Error: Unknown ordering %d", (*iv)[1]);
4119            R->order[n] = ringorder_unspec;
4120            return TRUE;
4121      }
4122    }
4123    sl=sl->next;
4124  }
4125
4126  // check for complete coverage
4127  if ((R->order[n]==ringorder_c) ||  (R->order[n]==ringorder_C)) n--;
4128  if (R->block1[n] != R->N)
4129  {
4130    if (((R->order[n]==ringorder_dp) ||
4131         (R->order[n]==ringorder_ds) ||
4132         (R->order[n]==ringorder_Dp) ||
4133         (R->order[n]==ringorder_Ds) ||
4134         (R->order[n]==ringorder_rp) ||
4135         (R->order[n]==ringorder_lp) ||
4136         (R->order[n]==ringorder_ls))
4137        &&
4138        R->block0[n] <= R->N)
4139    {
4140      R->block1[n] = R->N;
4141    }
4142    else
4143    {
4144      Werror("mismatch of number of vars (%d) and ordering (%d vars)",
4145             R->N,R->block1[n]);
4146      return TRUE;
4147    }
4148  }
4149  R->OrdSgn = typ;
4150  return FALSE;
4151}
4152
4153BOOLEAN rSleftvList2StringArray(sleftv* sl, char** p)
4154{
4155
4156  while(sl!=NULL)
4157  {
4158    if (sl->Name() == sNoName)
4159    {
4160      if (sl->Typ()==POLY_CMD)
4161      {
4162        sleftv s_sl;
4163        iiConvert(POLY_CMD,ANY_TYPE,-1,sl,&s_sl);
4164        if (s_sl.Name() != sNoName)
4165          *p = omStrDup(s_sl.Name());
4166        else
4167          *p = NULL;
4168        sl->next = s_sl.next;
4169        s_sl.next = NULL;
4170        s_sl.CleanUp();
4171        if (*p == NULL) return TRUE;
4172      }
4173      else
4174        return TRUE;
4175    }
4176    else
4177      *p = omStrDup(sl->Name());
4178    p++;
4179    sl=sl->next;
4180  }
4181  return FALSE;
4182}
4183
4184////////////////////
4185//
4186// rInit itself:
4187//
4188// INPUT:  s: name, pn: ch & parameter (names), rv: variable (names)
4189//         ord: ordering
4190// RETURN: currRingHdl on success
4191//         NULL        on error
4192// NOTE:   * makes new ring to current ring, on success
4193//         * considers input sleftv's as read-only
4194//idhdl rInit(char *s, sleftv* pn, sleftv* rv, sleftv* ord)
4195ring rInit(sleftv* pn, sleftv* rv, sleftv* ord)
4196{
4197  int ch;
4198  int float_len=0;
4199  int float_len2=0;
4200  ring R = NULL;
4201  idhdl tmp = NULL;
4202  BOOLEAN ffChar=FALSE;
4203  int typ = 1;
4204
4205  /* ch -------------------------------------------------------*/
4206  // get ch of ground field
4207  int numberOfAllocatedBlocks;
4208
4209  if (pn->Typ()==INT_CMD)
4210  {
4211    ch=(int)pn->Data();
4212  }
4213  else if ((pn->name != NULL)
4214  && ((strcmp(pn->name,"real")==0) || (strcmp(pn->name,"complex")==0)))
4215  {
4216    BOOLEAN complex_flag=(strcmp(pn->name,"complex")==0);
4217    ch=-1;
4218    if ((pn->next!=NULL) && (pn->next->Typ()==INT_CMD))
4219    {
4220      float_len=(int)pn->next->Data();
4221      float_len2=float_len;
4222      pn=pn->next;
4223      if ((pn->next!=NULL) && (pn->next->Typ()==INT_CMD))
4224      {
4225        float_len2=(int)pn->next->Data();
4226        pn=pn->next;
4227      }
4228    }
4229    if ((pn->next==NULL) && complex_flag)
4230    {
4231      pn->next=(leftv)omAlloc0Bin(sleftv_bin);
4232      pn->next->name=omStrDup("i");
4233    }
4234  }
4235  else
4236  {
4237    Werror("Wrong ground field specification");
4238    goto rInitError;
4239  }
4240  pn=pn->next;
4241
4242  int l, last;
4243  sleftv * sl;
4244  /*every entry in the new ring is initialized to 0*/
4245
4246  /* characteristic -----------------------------------------------*/
4247  /* input: 0 ch=0 : Q     parameter=NULL    ffChar=FALSE   float_len
4248   *         0    1 : Q(a,...)        *names         FALSE
4249   *         0   -1 : R               NULL           FALSE  0
4250   *         0   -1 : R               NULL           FALSE  prec. >6
4251   *         0   -1 : C               *names         FALSE  prec. 0..?
4252   *         p    p : Fp              NULL           FALSE
4253   *         p   -p : Fp(a)           *names         FALSE
4254   *         q    q : GF(q=p^n)       *names         TRUE
4255  */
4256  if (ch!=-1)
4257  {
4258    int l = 0;
4259
4260    if (ch!=0 && (ch<2)
4261    #ifndef NV_OPS
4262    || (ch > 32003)
4263    #endif
4264    )
4265    {
4266      Warn("%d is invalid characteristic of ground field. 32003 is used.", ch);
4267      ch=32003;
4268    }
4269    // load fftable, if necessary
4270    if (pn!=NULL)
4271    {
4272      while ((ch!=fftable[l]) && (fftable[l])) l++;
4273      if (fftable[l]==0) ch = IsPrime(ch);
4274      else
4275      {
4276        char *m[1]={(char *)sNoName};
4277        nfSetChar(ch,m);
4278        if (errorreported) goto rInitError;
4279        else ffChar=TRUE;
4280      }
4281    }
4282    else
4283      ch = IsPrime(ch);
4284  }
4285  // allocated ring and set ch
4286  R = (ring) omAlloc0Bin(sip_sring_bin);
4287  R->ch = ch;
4288  if (ch == -1)
4289  {
4290    R->float_len= si_min(float_len,32767);
4291    R->float_len2= si_min(float_len2,32767);
4292  }
4293
4294  /* parameter -------------------------------------------------------*/
4295  if (pn!=NULL)
4296  {
4297    R->P=pn->listLength();
4298    //if ((ffChar|| (ch == 1)) && (R->P > 1))
4299    if ((R->P > 1) && (ffChar || (ch == -1)))
4300    {
4301      WerrorS("too many parameters");
4302      goto rInitError;
4303    }
4304    R->parameter=(char**)omAlloc0(R->P*sizeof(char_ptr));
4305    if (rSleftvList2StringArray(pn, R->parameter))
4306    {
4307      WerrorS("parameter expected");
4308      goto rInitError;
4309    }
4310    if (ch>1 && !ffChar) R->ch=-ch;
4311    else if (ch==0) R->ch=1;
4312  }
4313  else if (ffChar)
4314  {
4315    WerrorS("need one parameter");
4316    goto rInitError;
4317  }
4318  /* post-processing of field description */
4319  // we have short reals, but no short complex
4320  if ((R->ch == - 1)
4321  && (R->parameter !=NULL)
4322  && (R->float_len < SHORT_REAL_LENGTH))
4323  {
4324    R->float_len = SHORT_REAL_LENGTH;
4325    R->float_len2 = SHORT_REAL_LENGTH;
4326  }
4327
4328  /* names and number of variables-------------------------------------*/
4329  {
4330    int l=rv->listLength();
4331#if SIZEOF_SHORT == 2
4332#define MAX_SHORT 0x7fff
4333#endif
4334    if (l>MAX_SHORT)
4335    {
4336      Werror("too many ring variables(%d), max is %d",l,MAX_SHORT);
4337       goto rInitError;
4338    }
4339    R->N = l; /*rv->listLength();*/
4340  }
4341  R->names   = (char **)omAlloc0(R->N * sizeof(char_ptr));
4342  if (rSleftvList2StringArray(rv, R->names))
4343  {
4344    WerrorS("name of ring variable expected");
4345    goto rInitError;
4346  }
4347
4348  /* check names and parameters for conflicts ------------------------- */
4349  {
4350    int i,j;
4351    for(i=0;i<R->P; i++)
4352    {
4353      for(j=0;j<R->N;j++)
4354      {
4355        if (strcmp(R->parameter[i],R->names[j])==0)
4356        {
4357          Werror("parameter %d conflicts with variable %d",i+1,j+1);
4358          goto rInitError;
4359        }
4360      }
4361    }
4362  }
4363  /* ordering -------------------------------------------------------------*/
4364  if (rSleftvOrdering2Ordering(ord, R))
4365    goto rInitError;
4366
4367  // Complete the initialization
4368  if (rComplete(R,1))
4369    goto rInitError;
4370
4371  rTest(R);
4372
4373  // try to enter the ring into the name list
4374  // need to clean up sleftv here, before this ring can be set to
4375  // new currRing or currRing can be killed beacuse new ring has
4376  // same name
4377  if (pn != NULL) pn->CleanUp();
4378  if (rv != NULL) rv->CleanUp();
4379  if (ord != NULL) ord->CleanUp();
4380  //if ((tmp = enterid(s, myynest, RING_CMD, &IDROOT))==NULL)
4381  //  goto rInitError;
4382
4383  //memcpy(IDRING(tmp),R,sizeof(*R));
4384  // set current ring
4385  //omFreeBin(R,  ip_sring_bin);
4386  //return tmp;
4387  return R;
4388
4389  // error case:
4390  rInitError:
4391  if  (R != NULL) rDelete(R);
4392  if (pn != NULL) pn->CleanUp();
4393  if (rv != NULL) rv->CleanUp();
4394  if (ord != NULL) ord->CleanUp();
4395  return NULL;
4396}
4397
4398void rKill(ring r)
4399{
4400  if ((r->ref<=0)&&(r->order!=NULL))
4401  {
4402#ifdef RDEBUG
4403    if (traceit &TRACE_SHOW_RINGS) Print("kill ring %x\n",r);
4404#endif
4405    if (r==currRing)
4406    {
4407      if (r->qideal!=NULL)
4408      {
4409        idDelete(&r->qideal);
4410        r->qideal=NULL;
4411        currQuotient=NULL;
4412      }
4413      if (ppNoether!=NULL) pDelete(&ppNoether);
4414      if (sLastPrinted.RingDependend())
4415      {
4416        sLastPrinted.CleanUp();
4417      }
4418      if ((myynest>0) && (iiRETURNEXPR[myynest].RingDependend()))
4419      {
4420        WerrorS("return value depends on local ring variable (export missing ?)");
4421        iiRETURNEXPR[myynest].CleanUp();
4422      }
4423      currRing=NULL;
4424      currRingHdl=NULL;
4425    }
4426    else if (r->qideal!=NULL)
4427    {
4428      id_Delete(&r->qideal, r);
4429      r->qideal = NULL;
4430    }
4431    int i=1;
4432    int j;
4433    int *pi=r->order;
4434#ifdef USE_IILOCALRING
4435    for (j=0;j<iiRETURNEXPR_len;j++)
4436    {
4437      if (iiLocalRing[j]==r)
4438      {
4439        if (j<myynest) Warn("killing the basering for level %d",j);
4440        iiLocalRing[j]=NULL;
4441      }
4442    }
4443#else /* USE_IILOCALRING */
4444//#endif /* USE_IILOCALRING */
4445    {
4446      proclevel * nshdl = procstack;
4447      int lev=myynest-1;
4448
4449      for(; nshdl != NULL; nshdl = nshdl->next)
4450      {
4451        if (nshdl->cRing==r)
4452        {
4453          Warn("killing the basering for level %d",lev);
4454          nshdl->cRing=NULL;
4455          nshdl->cRingHdl=NULL;
4456        }
4457      }
4458    }
4459#endif /* USE_IILOCALRING */
4460
4461    /* nKillChar(r); will be called from inside of rDelete */
4462    rDelete(r);
4463    return;
4464  }
4465  r->ref--;
4466}
4467
4468void rKill(idhdl h)
4469{
4470  ring r = IDRING(h);
4471  int ref=0;
4472  if (r!=NULL)
4473  {
4474    ref=r->ref;
4475    rKill(r);
4476  }
4477  if (h==currRingHdl)
4478  {
4479    if (ref<=0) { currRing=NULL; currRingHdl=NULL;}
4480    else
4481    {
4482      currRingHdl=rFindHdl(r,currRingHdl,NULL);
4483    }
4484  }
4485}
4486
4487idhdl rSimpleFindHdl(ring r, idhdl root, idhdl n=NULL)
4488{
4489  //idhdl next_best=NULL;
4490  idhdl h=root;
4491  while (h!=NULL)
4492  {
4493    if (((IDTYP(h)==RING_CMD)||(IDTYP(h)==QRING_CMD))
4494    && (h!=n)
4495    && (IDRING(h)==r)
4496    )
4497    {
4498   //   if (IDLEV(h)==myynest)
4499   //     return h;
4500   //   if ((IDLEV(h)==0) || (next_best==NULL))
4501   //     next_best=h;
4502   //   else if (IDLEV(next_best)<IDLEV(h))
4503   //     next_best=h;
4504      return h;
4505    }
4506    h=IDNEXT(h);
4507  }
4508  //return next_best;
4509  return NULL;
4510}
4511
Note: See TracBrowser for help on using the repository browser.