source: git/libpolys/polys/monomials/ring.cc @ 1ddd04

spielwiese
Last change on this file since 1ddd04 was 1ddd04, checked in by Hans Schoenemann <hannes@…>, 13 years ago
rDefault via coeffs
  • Property mode set to 100644
File size: 141.6 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/* $Id$ */
5
6/*
7* ABSTRACT - the interpreter related ring operations
8*/
9
10/* includes */
11#include <math.h>
12
13// #include <polys/options.h>
14#include <omalloc/omalloc.h>
15#include <misc/options.h>
16#include <polys/monomials/ring.h>
17#include <polys/monomials/p_polys.h>
18#include <polys/simpleideals.h>
19#include <coeffs/numbers.h>
20// #include <???/febase.h>
21// #include <???/intvec.h>
22#include <polys/ext_fields/longalg.h>
23#include <polys/ext_fields/longtrans.h>
24#include <coeffs/ffields.h>
25#include <polys/monomials/ring.h>
26#include <polys/monomials/maps.h>
27#include <polys/prCopy.h>
28// #include "../Singular/ipshell.h"
29#include <polys/templates/p_Procs.h>
30
31#include <polys/matpol.h>
32
33#ifdef HAVE_PLURAL
34#include <polys/nc/nc.h>
35#include <polys/nc/sca.h>
36#endif
37// #include <???/maps.h>
38// #include <???/matpol.h>
39#ifdef HAVE_FACTORY
40#define SI_DONT_HAVE_GLOBAL_VARS
41#  include <factory/factory.h>
42#endif
43
44#define BITS_PER_LONG 8*SIZEOF_LONG
45
46omBin sip_sring_bin = omGetSpecBin(sizeof(ip_sring));
47
48static const char * const ringorder_name[] =
49{
50  " ?", ///< ringorder_no = 0,
51  "a", ///< ringorder_a,
52  "A", ///< ringorder_a64,
53  "c", ///< ringorder_c,
54  "C", ///< ringorder_C,
55  "M", ///< ringorder_M,
56  "S", ///< ringorder_S,
57  "s", ///< ringorder_s,
58  "lp", ///< ringorder_lp,
59  "dp", ///< ringorder_dp,
60  "rp", ///< ringorder_rp,
61  "Dp", ///< ringorder_Dp,
62  "wp", ///< ringorder_wp,
63  "Wp", ///< ringorder_Wp,
64  "ls", ///< ringorder_ls,
65  "ds", ///< ringorder_ds,
66  "Ds", ///< ringorder_Ds,
67  "ws", ///< ringorder_ws,
68  "Ws", ///< ringorder_Ws,
69  "L", ///< ringorder_L,
70  "aa", ///< ringorder_aa
71  "rs", ///< ringorder_rs,
72  "IS", ///<  ringorder_IS
73  " _" ///< ringorder_unspec
74};
75
76#define  rInternalChar(r) ((r)->ch)
77
78const char * rSimpleOrdStr(int ord)
79{
80  return ringorder_name[ord];
81}
82
83/// unconditionally deletes fields in r
84void rDelete(ring r);
85/// set r->VarL_Size, r->VarL_Offset, r->VarL_LowIndex
86static void rSetVarL(ring r);
87/// get r->divmask depending on bits per exponent
88static unsigned long rGetDivMask(int bits);
89/// right-adjust r->VarOffset
90static void rRightAdjustVarOffset(ring r);
91static void rOptimizeLDeg(ring r);
92
93/*0 implementation*/
94//BOOLEAN rField_is_R(ring r)
95//{
96//  if (r->ch== -1)
97//  {
98//    if (r->float_len==(short)0) return TRUE;
99//  }
100//  return FALSE;
101//}
102
103/*
104
105/// internally changes the gloabl ring and resets the relevant
106/// global variables:
107/// SHOULD BE DEPRECATED NOW...?
108void rChangeCurrRing(ring r)
109{
110 // if ((currRing!=NULL) && (currRing->minpoly!=NULL))
111 // {
112 //   omCheckAddr(currRing->minpoly);
113 // }
114  //------------ set global ring vars --------------------------------
115  //currRing = r;
116  //currQuotient=NULL;
117  if (r != NULL)
118  {
119    rTest(r);
120    //------------ set global ring vars --------------------------------
121    //currQuotient=r->qideal;
122
123    //------------ global variables related to coefficients ------------
124    nSetChar(r->cf);
125
126    //------------ global variables related to polys -------------------
127    p_SetGlobals(r);
128    //------------ global variables related to factory -----------------
129#ifdef HAVE_FACTORY
130    //int c=ABS(nGetChar());
131    //if (c==1) c=0;
132    //setCharacteristic( c );
133#endif
134  }
135}
136*/
137
138ring rDefault(coeffs cf, int N, char **n,int ord_size, int *ord, int *block0, int *block1)
139{
140  ring r=(ring) omAlloc0Bin(sip_sring_bin);
141  r->ch    = n_GetChar(cf);
142  r->N     = N;
143  r->cf = cf;
144  /*r->P     = 0; Alloc0 */
145  /*names*/
146  r->names = (char **) omAlloc0(N * sizeof(char *));
147  int i;
148  for(i=0;i<N;i++)
149  {
150    r->names[i]  = omStrDup(n[i]);
151  }
152  /*weights: entries for 2 blocks: NULL*/
153  r->wvhdl = (int **)omAlloc0((ord_size+1) * sizeof(int *));
154  r->order = ord;
155  r->block0 = block0;
156  r->block1 = block1;
157  /*polynomial ring*/
158  r->OrdSgn    = 1;
159
160  /* complete ring intializations */
161  rComplete(r);
162  return r;
163}
164ring rDefault(int ch, int N, char **n,int ord_size, int *ord, int *block0, int *block1)
165{
166  coeffs cf;
167  if (ch==0) cf=nInitChar(n_Q,NULL);
168  else       cf=nInitChar(n_Zp,(void*)(long)ch);
169  return rDefault(cf,N,n,ord_size,ord,block0,block1);
170}
171ring   rDefault(coeffs cf, int N, char **n)
172{ 
173  /*order: lp,0*/
174  int *order = (int *) omAlloc(2* sizeof(int));
175  int *block0 = (int *)omAlloc0(2 * sizeof(int));
176  int *block1 = (int *)omAlloc0(2 * sizeof(int));
177  /* ringorder dp for the first block: var 1..N */
178  order[0]  = ringorder_lp;
179  block0[0] = 1;
180  block1[0] = N;
181  /* the last block: everything is 0 */
182  order[1]  = 0;
183
184  return rDefault(cf,N,n,2,order,block0,block1);
185}
186
187ring rDefault(int ch, int N, char **n)
188{
189  coeffs cf;
190  if (ch==0) cf=nInitChar(n_Q,NULL);
191  else       cf=nInitChar(n_Zp,(void*)(long)ch);
192  return rDefault(cf,N,n);
193}
194
195///////////////////////////////////////////////////////////////////////////
196//
197// rInit: define a new ring from sleftv's
198//
199//-> ipshell.cc
200
201/////////////////////////////
202// Auxillary functions
203//
204
205// check intvec, describing the ordering
206BOOLEAN rCheckIV(intvec *iv)
207{
208  if ((iv->length()!=2)&&(iv->length()!=3))
209  {
210    WerrorS("weights only for orderings wp,ws,Wp,Ws,a,M");
211    return TRUE;
212  }
213  return FALSE;
214}
215
216int rTypeOfMatrixOrder(intvec * order)
217{
218  int i=0,j,typ=1;
219  int sz = (int)sqrt((double)(order->length()-2));
220  if ((sz*sz)!=(order->length()-2))
221  {
222    WerrorS("Matrix order is not a square matrix");
223    typ=0;
224  }
225  while ((i<sz) && (typ==1))
226  {
227    j=0;
228    while ((j<sz) && ((*order)[j*sz+i+2]==0)) j++;
229    if (j>=sz)
230    {
231      typ = 0;
232      WerrorS("Matrix order not complete");
233    }
234    else if ((*order)[j*sz+i+2]<0)
235      typ = -1;
236    else
237      i++;
238  }
239  return typ;
240}
241
242/*2
243 * set a new ring from the data:
244 s: name, chr: ch, varnames: rv, ordering: ord, typ: typ
245 */
246
247int r_IsRingVar(const char *n, ring r)
248{
249  if ((r!=NULL) && (r->names!=NULL))
250  {
251    for (int i=0; i<r->N; i++)
252    {
253      if (r->names[i]==NULL) return -1;
254      if (strcmp(n,r->names[i]) == 0) return (int)i;
255    }
256  }
257  return -1;
258}
259
260
261void rWrite(ring r)
262{
263  if ((r==NULL)||(r->order==NULL))
264    return; /*to avoid printing after errors....*/
265
266  int nblocks=rBlocks(r);
267
268  // omCheckAddrSize(r,sizeof(ip_sring));
269  omCheckAddrSize(r->order,nblocks*sizeof(int));
270  omCheckAddrSize(r->block0,nblocks*sizeof(int));
271  omCheckAddrSize(r->block1,nblocks*sizeof(int));
272  omCheckAddrSize(r->wvhdl,nblocks*sizeof(int *));
273  omCheckAddrSize(r->names,r->N*sizeof(char *));
274
275  nblocks--;
276
277
278  if (rField_is_GF(r))
279  {
280    Print("//   # ground field : %d\n",r->ch);
281    Print("//   primitive element : %s\n", r->parameter[0]);
282    StringSetS("//   minpoly        : ");
283    nfShowMipo(r->cf);PrintS(StringAppendS("\n"));
284  }
285#ifdef HAVE_RINGS
286  else if (rField_is_Ring(r))
287  {
288    PrintS("//   coeff. ring is : ");
289    if (rField_is_Ring_Z(r)) PrintS("Integers\n");
290    long l = (long)mpz_sizeinbase(r->ringflaga, 10) + 2;
291    char* s = (char*) omAlloc(l);
292    mpz_get_str(s,10,r->ringflaga);
293    if (rField_is_Ring_ModN(r)) Print("Z/%s\n", s);
294    if (rField_is_Ring_2toM(r)) Print("Z/2^%lu\n", r->ringflagb);
295    if (rField_is_Ring_PtoM(r)) Print("Z/%s^%lu\n", s, r->ringflagb);
296    omFreeSize((ADDRESS)s, l);
297  }
298#endif
299  else
300  {
301    PrintS("//   characteristic : ");
302    if ( rField_is_R(r) )             PrintS("0 (real)\n");  /* R */
303    else if ( rField_is_long_R(r) )
304      Print("0 (real:%d digits, additional %d digits)\n",
305             r->float_len,r->float_len2);  /* long R */
306    else if ( rField_is_long_C(r) )
307      Print("0 (complex:%d digits, additional %d digits)\n",
308             r->float_len, r->float_len2);  /* long C */
309    else
310      Print ("%d\n",rChar(r)); /* Fp(a) */
311    if (r->parameter!=NULL)
312    {
313      Print ("//   %d parameter    : ",rPar(r));
314      char **sp=r->parameter;
315      int nop=0;
316      while (nop<rPar(r))
317      {
318        PrintS(*sp);
319        PrintS(" ");
320        sp++; nop++;
321      }
322      PrintS("\n//   minpoly        : ");
323      if ( rField_is_long_C(r) )
324      {
325        // i^2+1:
326        Print("(%s^2+1)\n",r->parameter[0]);
327      }
328      else if (r->minpoly==NULL)
329      {
330        PrintS("0\n");
331      }
332      else
333      {
334        StringSetS(""); n_Write(r->minpoly,r->cf); PrintS(StringAppendS("\n"));
335      }
336      //if (r->minideal!=NULL)
337      //{
338      //  iiWriteMatrix((matrix)r->minideal,"//   minpolys",1,0);
339      //  PrintLn();
340      //}
341    }
342  }
343  Print("//   number of vars : %d",r->N);
344
345  //for (nblocks=0; r->order[nblocks]; nblocks++);
346  nblocks=rBlocks(r)-1;
347
348  for (int l=0, nlen=0 ; l<nblocks; l++)
349  {
350    int i;
351    Print("\n//        block %3d : ",l+1);
352
353    Print("ordering %s", rSimpleOrdStr(r->order[l]));
354
355
356    if (r->order[l] == ringorder_s)
357    {
358      assume( l == 0 );
359#ifndef NDEBUG
360      Print("  syzcomp at %d",r->typ[l].data.syz.limit);
361#endif
362      continue;
363    }
364    else if (r->order[l] == ringorder_IS)
365    {
366      assume( r->block0[l] == r->block1[l] );
367      const int s = r->block0[l];
368      assume( (-2 < s) && (s < 2) );
369      Print("(%d)", s); // 0 => prefix! +/-1 => suffix!
370      continue;
371    }
372    else if (
373    (  (r->order[l] >= ringorder_lp)
374    ||(r->order[l] == ringorder_M)
375    ||(r->order[l] == ringorder_a)
376    ||(r->order[l] == ringorder_a64)
377    ||(r->order[l] == ringorder_aa) ) && (r->order[l] < ringorder_IS) )
378    {
379      PrintS("\n//                  : names   ");
380      for (i = r->block0[l]-1; i<r->block1[l]; i++)
381      {
382        nlen = strlen(r->names[i]);
383        Print(" %s",r->names[i]);
384      }
385    }
386
387    if (r->wvhdl[l]!=NULL)
388    {
389      for (int j= 0;
390           j<(r->block1[l]-r->block0[l]+1)*(r->block1[l]-r->block0[l]+1);
391           j+=i)
392      {
393        PrintS("\n//                  : weights ");
394        for (i = 0; i<=r->block1[l]-r->block0[l]; i++)
395        {
396          if (r->order[l] == ringorder_a64)
397          {
398            int64 *w=(int64 *)r->wvhdl[l];
399            #if SIZEOF_LONG == 4
400            Print("%*lld " ,nlen,w[i+j]);
401            #else
402            Print(" %*ld"  ,nlen,w[i+j]);
403            #endif
404          }
405          else
406            Print(" %*d" ,nlen,r->wvhdl[l][i+j]);
407        }
408        if (r->order[l]!=ringorder_M) break;
409      }
410    }
411  }
412#ifdef HAVE_PLURAL
413  if(rIsPluralRing(r))
414  {
415    PrintS("\n//   noncommutative relations:");
416    poly pl=NULL;
417    int nl;
418    int i,j;
419    for (i = 1; i<r->N; i++)
420    {
421      for (j = i+1; j<=r->N; j++)
422      {
423        nl = n_IsOne(p_GetCoeff(MATELEM(r->GetNC()->C,i,j),r->GetNC()->basering), r->cf);
424        if ( (MATELEM(r->GetNC()->D,i,j)!=NULL) || (!nl) )
425        {
426          Print("\n//    %s%s=",r->names[j-1],r->names[i-1]);
427          pl = MATELEM(r->GetNC()->MT[UPMATELEM(i,j,r->N)],1,1);
428          p_Write0(pl, r, r);
429        }
430      }
431    }
432#if 0  /*Singularg should not differ from Singular except in error case*/
433    Print("\n//   noncommutative type:%d", (int)ncRingType(r));
434    Print("\n//      is skew constant:%d",r->GetNC()->IsSkewConstant);
435    if( rIsSCA(r) )
436    {
437      Print("\n//   alternating variables: [%d, %d]", scaFirstAltVar(r), scaLastAltVar(r));
438      const ideal Q = SCAQuotient(r); // resides within r!
439      PrintS("\n//   quotient of sca by ideal");
440
441      if (Q!=NULL)
442      {
443//        if (r==currRing)
444//        {
445//          PrintLn();
446//          iiWriteMatrix((matrix)Q,"scaQ",1);
447//        }
448//        else
449            PrintS(" ...");
450      }
451      else
452        PrintS(" (NULL)");
453    }
454#endif
455  }
456#endif
457  //if (r->qideal!=NULL)
458  //{
459  //  PrintS("\n// quotient ring from ideal\n");
460  //  iiWriteMatrix((matrix)r->qideal,"_",1);
461  //}
462}
463
464void rDelete(ring r)
465{
466  int i, j;
467
468  if (r == NULL) return;
469
470#ifdef HAVE_PLURAL
471  if (rIsPluralRing(r))
472    nc_rKill(r);
473#endif
474
475  nKillChar(r->cf);
476  rUnComplete(r);
477  // delete order stuff
478  if (r->order != NULL)
479  {
480    i=rBlocks(r);
481    assume(r->block0 != NULL && r->block1 != NULL && r->wvhdl != NULL);
482    // delete order
483    omFreeSize((ADDRESS)r->order,i*sizeof(int));
484    omFreeSize((ADDRESS)r->block0,i*sizeof(int));
485    omFreeSize((ADDRESS)r->block1,i*sizeof(int));
486    // delete weights
487    for (j=0; j<i; j++)
488    {
489      if (r->wvhdl[j]!=NULL)
490        omFree(r->wvhdl[j]);
491    }
492    omFreeSize((ADDRESS)r->wvhdl,i*sizeof(int *));
493  }
494  else
495  {
496    assume(r->block0 == NULL && r->block1 == NULL && r->wvhdl == NULL);
497  }
498
499  // delete varnames
500  if(r->names!=NULL)
501  {
502    for (i=0; i<r->N; i++)
503    {
504      if (r->names[i] != NULL) omFree((ADDRESS)r->names[i]);
505    }
506    omFreeSize((ADDRESS)r->names,r->N*sizeof(char *));
507  }
508
509  // delete parameter
510  if (r->parameter!=NULL)
511  {
512    char **s=r->parameter;
513    j = 0;
514    while (j < rPar(r))
515    {
516      if (*s != NULL) omFree((ADDRESS)*s);
517      s++;
518      j++;
519    }
520    omFreeSize((ADDRESS)r->parameter,rPar(r)*sizeof(char *));
521  }
522#ifdef HAVE_RINGS
523  if (r->ringflaga != NULL)
524  {
525    mpz_clear(r->ringflaga);
526    omFree((ADDRESS)r->ringflaga);
527  }
528  if (r->nrnModul != NULL)
529  {
530    mpz_clear(r->nrnModul);
531    omFree((ADDRESS)r->nrnModul);
532  }
533#endif
534  omFreeBin(r, sip_sring_bin);
535}
536
537int rOrderName(char * ordername)
538{
539  int order=ringorder_unspec;
540  while (order!= 0)
541  {
542    if (strcmp(ordername,rSimpleOrdStr(order))==0)
543      break;
544    order--;
545  }
546  if (order==0) Werror("wrong ring order `%s`",ordername);
547  omFree((ADDRESS)ordername);
548  return order;
549}
550
551char * rOrdStr(ring r)
552{
553  if ((r==NULL)||(r->order==NULL)) return omStrDup("");
554  int nblocks,l,i;
555
556  for (nblocks=0; r->order[nblocks]; nblocks++);
557  nblocks--;
558
559  StringSetS("");
560  for (l=0; ; l++)
561  {
562    StringAppendS((char *)rSimpleOrdStr(r->order[l]));
563    if (
564           (r->order[l] != ringorder_c)
565        && (r->order[l] != ringorder_C)
566        && (r->order[l] != ringorder_s)
567        && (r->order[l] != ringorder_S)
568        && (r->order[l] != ringorder_IS)
569       )
570    {
571      if (r->wvhdl[l]!=NULL)
572      {
573        StringAppendS("(");
574        for (int j= 0;
575             j<(r->block1[l]-r->block0[l]+1)*(r->block1[l]-r->block0[l]+1);
576             j+=i+1)
577        {
578          char c=',';
579          if(r->order[l]==ringorder_a64)
580          {
581            int64 * w=(int64 *)r->wvhdl[l];
582            for (i = 0; i<r->block1[l]-r->block0[l]; i++)
583            {
584              StringAppend("%lld," ,w[i]);
585            }
586            StringAppend("%lld)" ,w[i]);
587            break;
588          }
589          else
590          {
591            for (i = 0; i<r->block1[l]-r->block0[l]; i++)
592            {
593              StringAppend("%d," ,r->wvhdl[l][i+j]);
594            }
595          }
596          if (r->order[l]!=ringorder_M)
597          {
598            StringAppend("%d)" ,r->wvhdl[l][i+j]);
599            break;
600          }
601          if (j+i+1==(r->block1[l]-r->block0[l]+1)*(r->block1[l]-r->block0[l]+1))
602            c=')';
603          StringAppend("%d%c" ,r->wvhdl[l][i+j],c);
604        }
605      }
606      else
607        StringAppend("(%d)",r->block1[l]-r->block0[l]+1);
608    }
609    else if (r->order[l] == ringorder_IS)
610    {
611      assume( r->block0[l] == r->block1[l] );
612      const int s = r->block0[l];
613      assume( (-2 < s) && (s < 2) );
614
615      StringAppend("(%d)", s);
616    }
617
618    if (l==nblocks) return omStrDup(StringAppendS(""));
619    StringAppendS(",");
620  }
621}
622
623char * rVarStr(ring r)
624{
625  if ((r==NULL)||(r->names==NULL)) return omStrDup("");
626  int i;
627  int l=2;
628  char *s;
629
630  for (i=0; i<r->N; i++)
631  {
632    l+=strlen(r->names[i])+1;
633  }
634  s=(char *)omAlloc((long)l);
635  s[0]='\0';
636  for (i=0; i<r->N-1; i++)
637  {
638    strcat(s,r->names[i]);
639    strcat(s,",");
640  }
641  strcat(s,r->names[i]);
642  return s;
643}
644
645char * rCharStr(ring r)
646{
647  char *s;
648  int i;
649
650#ifdef HAVE_RINGS
651  if (rField_is_Ring_Z(r))
652  {
653    s=omStrDup("integer");                   /* Z */
654    return s;
655  }
656  if(rField_is_Ring_2toM(r))
657  {
658    char* s = (char*) omAlloc(7+10+2);
659    sprintf(s,"integer,%lu",r->ringflagb);
660    return s;
661  }
662  if(rField_is_Ring_ModN(r))
663  {
664    long l = (long)mpz_sizeinbase(r->ringflaga, 10) + 2+7;
665    char* s = (char*) omAlloc(l);
666    gmp_sprintf(s,"integer,%Zd",r->ringflaga);
667    return s;
668  }
669  if(rField_is_Ring_PtoM(r))
670  {
671    long l = (long)mpz_sizeinbase(r->ringflaga, 10) + 2+7+10;
672    char* s = (char*) omAlloc(l);
673    gmp_sprintf(s,"integer,%Zd^%lu",r->ringflaga,r->ringflagb);
674    return s;
675  }
676#endif
677  if (r->parameter==NULL)
678  {
679    i=r->ch;
680    if(i==-1)
681      s=omStrDup("real");                    /* R */
682    else
683    {
684      s=(char *)omAlloc(MAX_INT_LEN+1);
685      sprintf(s,"%d",i);                   /* Q, Z/p */
686    }
687    return s;
688  }
689  if (rField_is_long_C(r))
690  {
691    s=(char *)omAlloc(21+strlen(r->parameter[0]));
692    sprintf(s,"complex,%d,%s",r->float_len,r->parameter[0]);   /* C */
693    return s;
694  }
695  int l=0;
696  for(i=0; i<rPar(r);i++)
697  {
698    l+=(strlen(r->parameter[i])+1);
699  }
700  s=(char *)omAlloc((long)(l+MAX_INT_LEN+1));
701  s[0]='\0';
702  if (r->ch<0)       sprintf(s,"%d",-r->ch); /* Fp(a) */
703  else if (r->ch==1) sprintf(s,"0");         /* Q(a)  */
704  else
705  {
706    sprintf(s,"%d,%s",r->ch,r->parameter[0]); /* GF(q)  */
707    return s;
708  }
709  char tt[2];
710  tt[0]=',';
711  tt[1]='\0';
712  for(i=0; i<rPar(r);i++)
713  {
714    strcat(s,tt);
715    strcat(s,r->parameter[i]);
716  }
717  return s;
718}
719
720char * rParStr(ring r)
721{
722  if ((r==NULL)||(r->parameter==NULL)) return omStrDup("");
723
724  int i;
725  int l=2;
726
727  for (i=0; i<rPar(r); i++)
728  {
729    l+=strlen(r->parameter[i])+1;
730  }
731  char *s=(char *)omAlloc((long)l);
732  s[0]='\0';
733  for (i=0; i<rPar(r)-1; i++)
734  {
735    strcat(s,r->parameter[i]);
736    strcat(s,",");
737  }
738  strcat(s,r->parameter[i]);
739  return s;
740}
741
742char * rString(ring r)
743{
744  char *ch=rCharStr(r);
745  char *var=rVarStr(r);
746  char *ord=rOrdStr(r);
747  char *res=(char *)omAlloc(strlen(ch)+strlen(var)+strlen(ord)+9);
748  sprintf(res,"(%s),(%s),(%s)",ch,var,ord);
749  omFree((ADDRESS)ch);
750  omFree((ADDRESS)var);
751  omFree((ADDRESS)ord);
752  return res;
753}
754
755int  rIsExtension(const ring r)
756{
757  return (r->parameter!=NULL); /* R, Q, Fp: FALSE */
758}
759
760static int binaryPower (const int a, const int b)
761{
762  /* computes a^b according to the binary representation of b,
763     i.e., a^7 = a^4 * a^2 * a^1. This saves some multiplications. */
764  int result = 1;
765  int factor = a;
766  int bb = b;
767  while (bb != 0)
768  {
769    if (bb % 2 != 0) result = result * factor;
770    bb = bb / 2;
771    factor = factor * factor;
772  }
773  return result;
774}
775
776int rChar(ring r)
777{
778#ifdef HAVE_RINGS
779  if (rField_is_Ring_2toM(r))
780    return binaryPower(2, (int)(unsigned long)r->ringflagb);
781  if (rField_is_Ring_ModN(r))
782    return (int)mpz_get_ui(r->ringflaga);
783  if (rField_is_Ring_PtoM(r))
784    return binaryPower((int)mpz_get_ui(r->ringflaga),
785                       (int)(unsigned long)r->ringflagb);
786#endif
787  if (rField_is_numeric(r))
788    return 0;
789  if (!rIsExtension(r)) /* Q, Fp */
790    return r->ch;
791  if (rField_is_Zp_a(r))  /* Fp(a)  */
792    return -r->ch;
793  if (rField_is_Q_a(r))   /* Q(a)  */
794    return 0;
795  /*else*/               /* GF(p,n) */
796  {
797    if ((r->ch & 1)==0) return 2;
798    int i=3;
799    while ((r->ch % i)!=0) i+=2;
800    return i;
801  }
802}
803
804typedef char *             char_ptr;
805omBin char_ptr_bin = omGetSpecBin(sizeof(char_ptr)); // deallocation?
806
807
808// creates a commutative nc extension; "converts" comm.ring to a Plural ring
809#ifdef HAVE_PLURAL
810ring nc_rCreateNCcomm_rCopy(ring r)
811{
812  r = rCopy(r);
813  if (rIsPluralRing(r)) 
814    return r;
815
816  matrix C = mpNew(r->N,r->N); // ring-independent!?!
817  matrix D = mpNew(r->N,r->N);
818
819  for(int i=1; i<r->N; i++)
820    for(int j=i+1; j<=r->N; j++)
821      MATELEM(C,i,j) = p_One( r);
822
823  if (nc_CallPlural(C, D, NULL, NULL, r, false, true, false, r/*??currRing??*/, TRUE)) // TODO: what about quotient ideal?
824    WarnS("Error initializing multiplication!"); // No reaction!???
825
826  return r;
827}
828#endif
829
830
831/*2
832 *returns -1 for not compatible, (sum is undefined)
833 *         1 for compatible (and sum)
834 */
835/* vartest: test for variable/paramter names
836* dp_dp: for comm. rings: use block order dp + dp/ds/wp
837*/
838int rSumInternal(ring r1, ring r2, ring &sum, BOOLEAN vartest, BOOLEAN dp_dp)
839{
840  ip_sring tmpR;
841  memset(&tmpR,0,sizeof(tmpR));
842  /* check coeff. field =====================================================*/
843  if ((rFieldType(r1)==rFieldType(r2))
844  && (rInternalChar(r1)==rInternalChar(r2)))
845  {
846    tmpR.ch=r1->ch;
847    if (rField_is_Q(r1)||rField_is_Zp(r1)||rField_is_GF(r1)) /*Q, Z/p, GF(p,n)*/
848    {
849      if (r1->parameter!=NULL)
850      {
851        if (!vartest || (strcmp(r1->parameter[0],r2->parameter[0])==0)) /* 1 par */
852        {
853          //tmpR.parameter=(char **)omAllocBin(char *);
854          //tmpR.parameter[0]=omStrDup(r1->parameter[0]);
855          //tmpR.P=1;
856        }
857        else
858        {
859          WerrorS("GF(p,n)+GF(p,n)");
860          return -1;
861        }
862      }
863    }
864    else if (rField_is_Extension(r1)) /* Q(a),Z/p(a) */
865    {
866      if (r1->minpoly!=NULL)
867      {
868        if (r2->minpoly!=NULL)
869        {
870          // HANNES: TODO: delete nSetChar
871         
872          // TODO: eliminate rChangeCurrRing
873//        rChangeCurrRing(r1);
874           
875          if ((!vartest || (strcmp(r1->parameter[0],r2->parameter[0])==0)) /* 1 par */
876              && n_Equal(r1->minpoly,r2->minpoly, r1->cf))
877          {
878            //tmpR.parameter=(char **)omAllocBin(char *);
879            //tmpR.parameter[0]=omStrDup(r1->parameter[0]);
880            //tmpR.minpoly=n_Copy(r1->minpoly, r1->cf);
881            //tmpR.P=1;
882            // HANNES: TODO: delete nSetChar
883          }
884          else
885          {
886            // HANNES: TODO: delete nSetChar
887            WerrorS("different minpolys");
888            return -1;
889          }
890        }
891        else
892        {
893          if ((!vartest || (strcmp(r1->parameter[0],r2->parameter[0])==0)) /* 1 par */
894              && (rPar(r2)==1))
895          {
896            //tmpR.parameter=(char **)omAllocBin(char *);
897            //tmpR.parameter[0]=omStrDup(r1->parameter[0]);
898            //tmpR.P=1;
899            //tmpR.minpoly=n_Copy(r1->minpoly, r1->cf);
900          }
901          else
902          {
903            WerrorS("different parameters and minpoly!=0");
904            return -1;
905          }
906        }
907      }
908      else /* r1->minpoly==NULL */
909      {
910        if (r2->minpoly!=NULL)
911        {
912          if ((!vartest || (strcmp(r1->parameter[0],r2->parameter[0])==0)) /* 1 par */
913              && (rPar(r1)==1))
914          {
915            //tmpR.parameter=(char **)omAllocBin(char *);
916            //tmpR.parameter[0]=omStrDup(r1->parameter[0]);
917            //tmpR.P=1;
918            //tmpR.minpoly=n_Copy(r2->minpoly, r2->cf);
919          }
920          else
921          {
922            WerrorS("different parameters and minpoly!=0");
923            return -1;
924          }
925        }
926        else
927        {
928          int len=rPar(r1)+rPar(r2);
929          tmpR.parameter=(char **)omAlloc0(len*sizeof(char *));
930          int i;
931          for (i=0;i<rPar(r1);i++)
932          {
933            tmpR.parameter[i]=omStrDup(r1->parameter[i]);
934          }
935          int j,l;
936          for(j=0;j<rPar(r2);j++)
937          {
938            if (vartest)
939            {
940              for(l=0;l<i;l++)
941              {
942                if(strcmp(tmpR.parameter[l],r2->parameter[j])==0)
943                  break;
944              }
945            }
946            else
947              l=i;
948            if (l==i)
949            {
950              tmpR.parameter[i]=omStrDup(r2->parameter[j]);
951              i++;
952            }
953          }
954          if (i!=len)
955          {
956            tmpR.parameter=(char**)omReallocSize(tmpR.parameter,
957                                                 len*sizeof(char *),
958                                                 i*sizeof(char *));
959          }
960          tmpR.P=i;
961        }
962      }
963    }
964    #ifdef HAVE_RINGS
965    else if (rField_is_Ring(r1)||rField_is_Ring(r2))
966    {
967      if (r1->ringtype != r2->ringtype)
968      {
969        Werror("rSumInternal not yet implemented for %s",
970               "different coefficient rings");
971        return -1;
972      }
973      else
974      {
975        tmpR.ch        = rInternalChar(r1);
976        tmpR.ringtype  = r1->ringtype;
977        if (r1->ringflaga != NULL)
978        {
979          omBin tmpBin = omGetSpecBin(sizeof(mpz_t));
980          tmpR.ringflaga = (int_number)omAllocBin(tmpBin);
981          mpz_init_set(tmpR.ringflaga, (int_number)r1->ringflaga);
982        }
983        tmpR.ringflagb = r1->ringflagb;
984        tmpR.nr2mModul = r1->nr2mModul;
985        if (r1->nrnModul != NULL)
986        {
987          omBin tmpBin = omGetSpecBin(sizeof(mpz_t));
988          tmpR.nrnModul = (int_number)omAllocBin(tmpBin);
989          mpz_init_set(tmpR.nrnModul, (int_number)r1->nrnModul);
990        }
991      }
992    }
993    #endif
994  }
995  else /* r1->ch!=r2->ch */
996  {
997    if (r1->ch<-1) /* Z/p(a) */
998    {
999      if ((r2->ch==0) /* Q */
1000          || (r2->ch==-r1->ch)) /* Z/p */
1001      {
1002        tmpR.ch=r1->ch;
1003        tmpR.P=rPar(r1);
1004        tmpR.parameter=(char **)omAlloc(rPar(r1)*sizeof(char *));
1005        int i;
1006        for (i=0;i<rPar(r1);i++)
1007        {
1008          tmpR.parameter[i]=omStrDup(r1->parameter[i]);
1009        }
1010        if (r1->minpoly!=NULL)
1011        {
1012          tmpR.minpoly=n_Copy(r1->minpoly, r1->cf);
1013        }
1014      }
1015      else  /* R, Q(a),Z/q,Z/p(a),GF(p,n) */
1016      {
1017        WerrorS("Z/p(a)+(R,Q(a),Z/q(a),GF(q,n))");
1018        return -1;
1019      }
1020    }
1021    else if (r1->ch==-1) /* R */
1022    {
1023      WerrorS("R+..");
1024      return -1;
1025    }
1026    else if (r1->ch==0) /* Q */
1027    {
1028      if ((r2->ch<-1)||(r2->ch==1)) /* Z/p(a),Q(a) */
1029      {
1030        tmpR.ch=r2->ch;
1031        tmpR.P=rPar(r2);
1032        tmpR.parameter=(char **)omAlloc(rPar(r2)*sizeof(char *));
1033        int i;
1034        for (i=0;i<rPar(r2);i++)
1035        {
1036          tmpR.parameter[i]=omStrDup(r2->parameter[i]);
1037        }
1038        if (r2->minpoly!=NULL)
1039        {
1040          tmpR.minpoly=n_Copy(r2->minpoly, r2->cf);
1041        }
1042      }
1043      else if (r2->ch>1) /* Z/p,GF(p,n) */
1044      {
1045        tmpR.ch=r2->ch;
1046        if (r2->parameter!=NULL)
1047        {
1048          tmpR.parameter=(char_ptr *)omAllocBin(char_ptr_bin);
1049          tmpR.P=1;
1050          tmpR.parameter[0]=omStrDup(r2->parameter[0]);
1051        }
1052      }
1053      else
1054      {
1055        WerrorS("Q+R");
1056        return -1; /* R */
1057      }
1058    }
1059    else if (r1->ch==1) /* Q(a) */
1060    {
1061      if (r2->ch==0) /* Q */
1062      {
1063        tmpR.ch=r1->ch;
1064        tmpR.P=rPar(r1);
1065        tmpR.parameter=(char **)omAlloc(rPar(r1)*sizeof(char *));
1066        int i;
1067        for(i=0;i<rPar(r1);i++)
1068        {
1069          tmpR.parameter[i]=omStrDup(r1->parameter[i]);
1070        }
1071        if (r1->minpoly!=NULL)
1072        {
1073          tmpR.minpoly=n_Copy(r1->minpoly, r1->cf);
1074        }
1075      }
1076      else  /* R, Z/p,GF(p,n) */
1077      {
1078        WerrorS("Q(a)+(R,Z/p,GF(p,n))");
1079        return -1;
1080      }
1081    }
1082    else /* r1->ch >=2 , Z/p */
1083    {
1084      if (r2->ch==0) /* Q */
1085      {
1086        tmpR.ch=r1->ch;
1087      }
1088      else if (r2->ch==-r1->ch) /* Z/p(a) */
1089      {
1090        tmpR.ch=r2->ch;
1091        tmpR.P=rPar(r2);
1092        tmpR.parameter=(char **)omAlloc(rPar(r2)*sizeof(char *));
1093        int i;
1094        for(i=0;i<rPar(r2);i++)
1095        {
1096          tmpR.parameter[i]=omStrDup(r2->parameter[i]);
1097        }
1098        if (r2->minpoly!=NULL)
1099        {
1100          tmpR.minpoly=n_Copy(r2->minpoly, r2->cf);
1101        }
1102      }
1103      else
1104      {
1105        WerrorS("Z/p+(GF(q,n),Z/q(a),R,Q(a))");
1106        return -1; /* GF(p,n),Z/q(a),R,Q(a) */
1107      }
1108    }
1109  }
1110  /* variable names ========================================================*/
1111  int i,j,k;
1112  int l=r1->N+r2->N;
1113  char **names=(char **)omAlloc0(l*sizeof(char *));
1114  k=0;
1115
1116  // collect all varnames from r1, except those which are parameters
1117  // of r2, or those which are the empty string
1118  for (i=0;i<r1->N;i++)
1119  {
1120    BOOLEAN b=TRUE;
1121
1122    if (*(r1->names[i]) == '\0')
1123      b = FALSE;
1124    else if ((r2->parameter!=NULL) && (strlen(r1->names[i])==1))
1125    {
1126      if (vartest)
1127      {
1128        for(j=0;j<rPar(r2);j++)
1129        {
1130          if (strcmp(r1->names[i],r2->parameter[j])==0)
1131          {
1132            b=FALSE;
1133            break;
1134          }
1135        }
1136      }
1137    }
1138
1139    if (b)
1140    {
1141      //Print("name : %d: %s\n",k,r1->names[i]);
1142      names[k]=omStrDup(r1->names[i]);
1143      k++;
1144    }
1145    //else
1146    //  Print("no name (par1) %s\n",r1->names[i]);
1147  }
1148  // Add variables from r2, except those which are parameters of r1
1149  // those which are empty strings, and those which equal a var of r1
1150  for(i=0;i<r2->N;i++)
1151  {
1152    BOOLEAN b=TRUE;
1153
1154    if (*(r2->names[i]) == '\0')
1155      b = FALSE;
1156    else if ((r1->parameter!=NULL) && (strlen(r2->names[i])==1))
1157    {
1158      if (vartest)
1159      {
1160        for(j=0;j<rPar(r1);j++)
1161        {
1162          if (strcmp(r2->names[i],r1->parameter[j])==0)
1163          {
1164            b=FALSE;
1165            break;
1166          }
1167        }
1168      }
1169    }
1170
1171    if (b)
1172    {
1173      if (vartest)
1174      {
1175        for(j=0;j<r1->N;j++)
1176        {
1177          if (strcmp(r1->names[j],r2->names[i])==0)
1178          {
1179            b=FALSE;
1180            break;
1181          }
1182        }
1183      }
1184      if (b)
1185      {
1186        //Print("name : %d : %s\n",k,r2->names[i]);
1187        names[k]=omStrDup(r2->names[i]);
1188        k++;
1189      }
1190      //else
1191      //  Print("no name (var): %s\n",r2->names[i]);
1192    }
1193    //else
1194    //  Print("no name (par): %s\n",r2->names[i]);
1195  }
1196  // check whether we found any vars at all
1197  if (k == 0)
1198  {
1199    names[k]=omStrDup("");
1200    k=1;
1201  }
1202  tmpR.N=k;
1203  tmpR.names=names;
1204  /* ordering *======================================================== */
1205  tmpR.OrdSgn=1;
1206  if (dp_dp
1207#ifdef HAVE_PLURAL
1208      && !rIsPluralRing(r1) && !rIsPluralRing(r2)
1209#endif
1210     )
1211  {
1212    tmpR.order=(int*)omAlloc(4*sizeof(int));
1213    tmpR.block0=(int*)omAlloc0(4*sizeof(int));
1214    tmpR.block1=(int*)omAlloc0(4*sizeof(int));
1215    tmpR.wvhdl=(int**)omAlloc0(4*sizeof(int *));
1216    tmpR.order[0]=ringorder_dp;
1217    tmpR.block0[0]=1;
1218    tmpR.block1[0]=rVar(r1);
1219    if (r2->OrdSgn==1)
1220    {
1221      if ((r2->block0[0]==1)
1222      && (r2->block1[0]==rVar(r2))
1223      && ((r2->order[0]==ringorder_wp)
1224        || (r2->order[0]==ringorder_Wp)
1225        || (r2->order[0]==ringorder_Dp))
1226     )
1227     {
1228       tmpR.order[1]=r2->order[0];
1229       if (r2->wvhdl[0]!=NULL)
1230         tmpR.wvhdl[1]=(int *)omMemDup(r2->wvhdl[0]);
1231     }
1232     else
1233        tmpR.order[1]=ringorder_dp;
1234    }
1235    else
1236    {
1237      tmpR.order[1]=ringorder_ds;
1238      tmpR.OrdSgn=-1;
1239    }
1240    tmpR.block0[1]=rVar(r1)+1;
1241    tmpR.block1[1]=rVar(r1)+rVar(r2);
1242    tmpR.order[2]=ringorder_C;
1243    tmpR.order[3]=0;
1244  }
1245  else
1246  {
1247    if ((r1->order[0]==ringorder_unspec)
1248        && (r2->order[0]==ringorder_unspec))
1249    {
1250      tmpR.order=(int*)omAlloc(3*sizeof(int));
1251      tmpR.block0=(int*)omAlloc(3*sizeof(int));
1252      tmpR.block1=(int*)omAlloc(3*sizeof(int));
1253      tmpR.wvhdl=(int**)omAlloc0(3*sizeof(int *));
1254      tmpR.order[0]=ringorder_unspec;
1255      tmpR.order[1]=ringorder_C;
1256      tmpR.order[2]=0;
1257      tmpR.block0[0]=1;
1258      tmpR.block1[0]=tmpR.N;
1259    }
1260    else if (l==k) /* r3=r1+r2 */
1261    {
1262      int b;
1263      ring rb;
1264      if (r1->order[0]==ringorder_unspec)
1265      {
1266        /* extend order of r2 to r3 */
1267        b=rBlocks(r2);
1268        rb=r2;
1269        tmpR.OrdSgn=r2->OrdSgn;
1270      }
1271      else if (r2->order[0]==ringorder_unspec)
1272      {
1273        /* extend order of r1 to r3 */
1274        b=rBlocks(r1);
1275        rb=r1;
1276        tmpR.OrdSgn=r1->OrdSgn;
1277      }
1278      else
1279      {
1280        b=rBlocks(r1)+rBlocks(r2)-2; /* for only one order C, only one 0 */
1281        rb=NULL;
1282      }
1283      tmpR.order=(int*)omAlloc0(b*sizeof(int));
1284      tmpR.block0=(int*)omAlloc0(b*sizeof(int));
1285      tmpR.block1=(int*)omAlloc0(b*sizeof(int));
1286      tmpR.wvhdl=(int**)omAlloc0(b*sizeof(int *));
1287      /* weights not implemented yet ...*/
1288      if (rb!=NULL)
1289      {
1290        for (i=0;i<b;i++)
1291        {
1292          tmpR.order[i]=rb->order[i];
1293          tmpR.block0[i]=rb->block0[i];
1294          tmpR.block1[i]=rb->block1[i];
1295          if (rb->wvhdl[i]!=NULL)
1296            WarnS("rSum: weights not implemented");
1297        }
1298        tmpR.block0[0]=1;
1299      }
1300      else /* ring sum for complete rings */
1301      {
1302        for (i=0;r1->order[i]!=0;i++)
1303        {
1304          tmpR.order[i]=r1->order[i];
1305          tmpR.block0[i]=r1->block0[i];
1306          tmpR.block1[i]=r1->block1[i];
1307          if (r1->wvhdl[i]!=NULL)
1308            tmpR.wvhdl[i] = (int*) omMemDup(r1->wvhdl[i]);
1309        }
1310        j=i;
1311        i--;
1312        if ((r1->order[i]==ringorder_c)
1313            ||(r1->order[i]==ringorder_C))
1314        {
1315          j--;
1316          tmpR.order[b-2]=r1->order[i];
1317        }
1318        for (i=0;r2->order[i]!=0;i++)
1319        {
1320          if ((r2->order[i]!=ringorder_c)
1321              &&(r2->order[i]!=ringorder_C))
1322          {
1323            tmpR.order[j]=r2->order[i];
1324            tmpR.block0[j]=r2->block0[i]+rVar(r1);
1325            tmpR.block1[j]=r2->block1[i]+rVar(r1);
1326            if (r2->wvhdl[i]!=NULL)
1327            {
1328              tmpR.wvhdl[j] = (int*) omMemDup(r2->wvhdl[i]);
1329            }
1330            j++;
1331          }
1332        }
1333        if((r1->OrdSgn==-1)||(r2->OrdSgn==-1))
1334          tmpR.OrdSgn=-1;
1335      }
1336    }
1337    else if ((k==rVar(r1)) && (k==rVar(r2))) /* r1 and r2 are "quite"
1338                                                the same ring */
1339      /* copy r1, because we have the variables from r1 */
1340    {
1341      int b=rBlocks(r1);
1342
1343      tmpR.order=(int*)omAlloc0(b*sizeof(int));
1344      tmpR.block0=(int*)omAlloc0(b*sizeof(int));
1345      tmpR.block1=(int*)omAlloc0(b*sizeof(int));
1346      tmpR.wvhdl=(int**)omAlloc0(b*sizeof(int *));
1347      /* weights not implemented yet ...*/
1348      for (i=0;i<b;i++)
1349      {
1350        tmpR.order[i]=r1->order[i];
1351        tmpR.block0[i]=r1->block0[i];
1352        tmpR.block1[i]=r1->block1[i];
1353        if (r1->wvhdl[i]!=NULL)
1354        {
1355          tmpR.wvhdl[i] = (int*) omMemDup(r1->wvhdl[i]);
1356        }
1357      }
1358      tmpR.OrdSgn=r1->OrdSgn;
1359    }
1360    else
1361    {
1362      for(i=0;i<k;i++) omFree((ADDRESS)tmpR.names[i]);
1363      omFreeSize((ADDRESS)names,tmpR.N*sizeof(char *));
1364      Werror("difficulties with variables: %d,%d -> %d",rVar(r1),rVar(r2),k);
1365      return -1;
1366    }
1367  }
1368  sum=(ring)omAllocBin(sip_sring_bin);
1369  memcpy(sum,&tmpR,sizeof(ip_sring));
1370  rComplete(sum);
1371
1372//#ifdef RDEBUG
1373//  rDebugPrint(sum);
1374//#endif
1375
1376
1377
1378#ifdef HAVE_PLURAL
1379  if(1)
1380  {
1381//    ring old_ring = currRing;
1382
1383    BOOLEAN R1_is_nc = rIsPluralRing(r1);
1384    BOOLEAN R2_is_nc = rIsPluralRing(r2);
1385
1386    if ( (R1_is_nc) || (R2_is_nc))
1387    {
1388      ring R1 = nc_rCreateNCcomm_rCopy(r1);
1389      assume( rIsPluralRing(R1) );
1390
1391#if 0
1392#ifdef RDEBUG
1393      rWrite(R1);
1394      rDebugPrint(R1);
1395#endif
1396#endif
1397      ring R2 = nc_rCreateNCcomm_rCopy(r2);
1398#if 0
1399#ifdef RDEBUG
1400      rWrite(R2);
1401      rDebugPrint(R2);
1402#endif
1403#endif
1404
1405//      rChangeCurrRing(sum); // ?
1406
1407      // Projections from R_i into Sum:
1408      /* multiplication matrices business: */
1409      /* find permutations of vars and pars */
1410      int *perm1 = (int *)omAlloc0((rVar(R1)+1)*sizeof(int));
1411      int *par_perm1 = NULL;
1412      if (rPar(R1)!=0) par_perm1=(int *)omAlloc0((rPar(R1)+1)*sizeof(int));
1413
1414      int *perm2 = (int *)omAlloc0((rVar(R2)+1)*sizeof(int));
1415      int *par_perm2 = NULL;
1416      if (rPar(R2)!=0) par_perm2=(int *)omAlloc0((rPar(R2)+1)*sizeof(int));
1417
1418      maFindPerm(R1->names,  rVar(R1),  R1->parameter,  rPar(R1),
1419                 sum->names, rVar(sum), sum->parameter, rPar(sum),
1420                 perm1, par_perm1, sum->ch);
1421
1422      maFindPerm(R2->names,  rVar(R2),  R2->parameter,  rPar(R2),
1423                 sum->names, rVar(sum), sum->parameter, rPar(sum),
1424                 perm2, par_perm2, sum->ch);
1425
1426
1427      matrix C1 = R1->GetNC()->C, C2 = R2->GetNC()->C;
1428      matrix D1 = R1->GetNC()->D, D2 = R2->GetNC()->D;
1429
1430      // !!!! BUG? C1 and C2 might live in different baserings!!!
1431
1432      int l = rVar(R1) + rVar(R2);
1433
1434      matrix C  = mpNew(l,l);
1435      matrix D  = mpNew(l,l);
1436
1437      for (i = 1; i <= rVar(R1); i++)
1438        for (j= rVar(R1)+1; j <= l; j++)
1439          MATELEM(C,i,j) = p_One(sum); // in 'sum'
1440
1441      id_Test((ideal)C, sum);
1442
1443      nMapFunc nMap1 = n_SetMap(R1->cf,sum->cf); /* can change something global: not usable
1444                                                    after the next nSetMap call :( */
1445      // Create blocked C and D matrices:
1446      for (i=1; i<= rVar(R1); i++)
1447        for (j=i+1; j<=rVar(R1); j++)
1448        {
1449          assume(MATELEM(C1,i,j) != NULL);
1450          MATELEM(C,i,j) = p_PermPoly(MATELEM(C1,i,j), perm1, R1, sum, nMap1, par_perm1, rPar(R1)); // need ADD + CMP ops.
1451
1452          if (MATELEM(D1,i,j) != NULL)
1453            MATELEM(D,i,j) = p_PermPoly(MATELEM(D1,i,j), perm1, R1, sum, nMap1, par_perm1, rPar(R1));
1454        }
1455
1456      id_Test((ideal)C, sum);
1457      id_Test((ideal)D, sum);
1458
1459
1460      nMapFunc nMap2 = n_SetMap(R2->cf,sum->cf); /* can change something global: not usable
1461                                                    after the next nSetMap call :( */
1462      for (i=1; i<= rVar(R2); i++)
1463        for (j=i+1; j<=rVar(R2); j++)
1464        {
1465          assume(MATELEM(C2,i,j) != NULL);
1466          MATELEM(C,rVar(R1)+i,rVar(R1)+j) = p_PermPoly(MATELEM(C2,i,j),perm2,R2,sum, nMap2,par_perm2,rPar(R2));
1467
1468          if (MATELEM(D2,i,j) != NULL)
1469            MATELEM(D,rVar(R1)+i,rVar(R1)+j) = p_PermPoly(MATELEM(D2,i,j),perm2,R2,sum, nMap2,par_perm2,rPar(R2));
1470        }
1471
1472      id_Test((ideal)C, sum);
1473      id_Test((ideal)D, sum);
1474
1475      // Now sum is non-commutative with blocked structure constants!
1476      if (nc_CallPlural(C, D, NULL, NULL, sum, false, false, true, sum))
1477        WarnS("Error initializing non-commutative multiplication!");
1478
1479      /* delete R1, R2*/
1480
1481#if 0
1482#ifdef RDEBUG
1483      rWrite(sum);
1484      rDebugPrint(sum);
1485
1486      Print("\nRefs: R1: %d, R2: %d\n", R1->GetNC()->ref, R2->GetNC()->ref);
1487
1488#endif
1489#endif
1490
1491
1492      rDelete(R1);
1493      rDelete(R2);
1494
1495      /* delete perm arrays */
1496      if (perm1!=NULL) omFree((ADDRESS)perm1);
1497      if (perm2!=NULL) omFree((ADDRESS)perm2);
1498      if (par_perm1!=NULL) omFree((ADDRESS)par_perm1);
1499      if (par_perm2!=NULL) omFree((ADDRESS)par_perm2);
1500
1501//      rChangeCurrRing(old_ring);
1502    }
1503
1504  }
1505#endif
1506
1507  ideal Q=NULL;
1508  ideal Q1=NULL, Q2=NULL;
1509  if (r1->qideal!=NULL)
1510  {
1511//    rChangeCurrRing(sum);
1512//     if (r2->qideal!=NULL)
1513//     {
1514//       WerrorS("todo: qring+qring");
1515//       return -1;
1516//     }
1517//     else
1518//     {}
1519    /* these were defined in the Plural Part above... */
1520    int *perm1 = (int *)omAlloc0((rVar(r1)+1)*sizeof(int));
1521    int *par_perm1 = NULL;
1522    if (rPar(r1)!=0) par_perm1=(int *)omAlloc0((rPar(r1)+1)*sizeof(int));
1523    maFindPerm(r1->names,  rVar(r1),  r1->parameter,  rPar(r1),
1524               sum->names, rVar(sum), sum->parameter, rPar(sum),
1525               perm1, par_perm1, sum->ch);
1526    nMapFunc nMap1 = n_SetMap(r1->cf,sum->cf);
1527    Q1 = idInit(IDELEMS(r1->qideal),1);
1528
1529    for (int for_i=0;for_i<IDELEMS(r1->qideal);for_i++)
1530      Q1->m[for_i] = p_PermPoly(
1531                                r1->qideal->m[for_i], perm1,
1532                                r1, sum,
1533                                nMap1, 
1534                                par_perm1, rPar(r1));
1535   
1536    omFree((ADDRESS)perm1);
1537  }
1538
1539  if (r2->qideal!=NULL)
1540  {
1541    //if (currRing!=sum)
1542    //  rChangeCurrRing(sum);
1543    int *perm2 = (int *)omAlloc0((rVar(r2)+1)*sizeof(int));
1544    int *par_perm2 = NULL;
1545    if (rPar(r2)!=0) par_perm2=(int *)omAlloc0((rPar(r2)+1)*sizeof(int));
1546    maFindPerm(r2->names,  rVar(r2),  r2->parameter,  rPar(r2),
1547               sum->names, rVar(sum), sum->parameter, rPar(sum),
1548               perm2, par_perm2, sum->ch);
1549    nMapFunc nMap2 = n_SetMap(r2->cf,sum->cf);
1550    Q2 = idInit(IDELEMS(r2->qideal),1);
1551     
1552    for (int for_i=0;for_i<IDELEMS(r2->qideal);for_i++)
1553      Q2->m[for_i] = p_PermPoly( 
1554                  r2->qideal->m[for_i], perm2, 
1555                  r2, sum, 
1556                  nMap2, 
1557                  par_perm2, rPar(r2));
1558   
1559    omFree((ADDRESS)perm2);
1560  }
1561  if ( (Q1!=NULL) || ( Q2!=NULL))
1562  {
1563    Q = id_SimpleAdd(Q1,Q2,sum);
1564  }
1565  sum->qideal = Q;
1566
1567#ifdef HAVE_PLURAL
1568  if( rIsPluralRing(sum) )
1569    nc_SetupQuotient( sum );
1570#endif
1571  return 1;
1572}
1573
1574/*2
1575 *returns -1 for not compatible, (sum is undefined)
1576 *         0 for equal, (and sum)
1577 *         1 for compatible (and sum)
1578 */
1579int rSum(ring r1, ring r2, ring &sum)
1580{
1581  if (r1==r2)
1582  {
1583    sum=r1;
1584    r1->ref++;
1585    return 0;
1586  }
1587  return rSumInternal(r1,r2,sum,TRUE,FALSE);
1588}
1589
1590/*2
1591 * create a copy of the ring r, which must be equivalent to currRing
1592 * used for qring definition,..
1593 * (i.e.: normal rings: same nCopy as currRing;
1594 *        qring:        same nCopy, same idCopy as currRing)
1595 * DOES NOT CALL rComplete
1596 */
1597ring rCopy0(const ring r, BOOLEAN copy_qideal, BOOLEAN copy_ordering)
1598{
1599  if (r == NULL) return NULL;
1600  int i,j;
1601  ring res=(ring)omAllocBin(sip_sring_bin);
1602  memset(res,0,sizeof(ip_sring));
1603  //memcpy(res,r,sizeof(ip_sring));
1604  //memset: res->idroot=NULL; /* local objects */
1605  //ideal      minideal;
1606  res->options=r->options; /* ring dependent options */
1607
1608  //memset: res->ordsgn=NULL;
1609  //memset: res->typ=NULL;
1610  //memset: res->VarOffset=NULL;
1611  //memset: res->firstwv=NULL;
1612
1613  //struct omBin   PolyBin; /* Bin from where monoms are allocated */
1614  //memset: res->PolyBin=NULL; // rComplete
1615  res->ch=r->ch;     /* characteristic */
1616#ifdef HAVE_RINGS
1617  res->ringtype=r->ringtype;  /* cring = 0 => coefficient field, cring = 1 => coeffs from Z/2^m */
1618  if (r->ringflaga!=NULL)
1619  {
1620    res->ringflaga = (int_number) omAlloc(sizeof(mpz_t));
1621    mpz_init_set(res->ringflaga,r->ringflaga);
1622  }
1623  res->ringflagb=r->ringflagb;
1624  if (r->nrnModul!=NULL)
1625  {
1626    res->nrnModul = (int_number) omAlloc(sizeof(mpz_t));
1627    mpz_init_set(res->nrnModul,r->nrnModul);
1628  }
1629#endif
1630  //memset: res->ref=0; /* reference counter to the ring */
1631
1632  res->float_len=r->float_len; /* additional char-flags */
1633  res->float_len2=r->float_len2; /* additional char-flags */
1634
1635  res->N=r->N;      /* number of vars */
1636  res->P=r->P;      /* number of pars */
1637  res->OrdSgn=r->OrdSgn; /* 1 for polynomial rings, -1 otherwise */
1638
1639  res->firstBlockEnds=r->firstBlockEnds;
1640#ifdef HAVE_PLURAL
1641  res->real_var_start=r->real_var_start;
1642  res->real_var_end=r->real_var_end;
1643#endif
1644
1645#ifdef HAVE_SHIFTBBA
1646  res->isLPring=r->isLPring; /* 0 for non-letterplace rings, otherwise the number of LP blocks, at least 1, known also as lV */
1647#endif
1648
1649  res->VectorOut=r->VectorOut;
1650  res->ShortOut=r->ShortOut;
1651  res->CanShortOut=r->CanShortOut;
1652  res->LexOrder=r->LexOrder; // TRUE if the monomial ordering has polynomial and power series blocks
1653  res->MixedOrder=r->MixedOrder; // ?? 1 for lex ordering (except ls), -1 otherwise
1654  res->ComponentOrder=r->ComponentOrder;
1655
1656  //memset: res->ExpL_Size=0;
1657  //memset: res->CmpL_Size=0;
1658  //memset: res->VarL_Size=0;
1659  //memset: res->pCompIndex=0;
1660  //memset: res->pOrdIndex=0;
1661  //memset: res->OrdSize=0;
1662  //memset: res->VarL_LowIndex=0;
1663  //memset: res->MinExpPerLong=0;
1664  //memset: res->NegWeightL_Size=0;
1665  //memset: res->NegWeightL_Offset=NULL;
1666  //memset: res->VarL_Offset=NULL;
1667
1668  // the following are set by rComplete unless predefined
1669  // therefore, we copy these values: maybe they are non-standard
1670  /* mask for getting single exponents */
1671  res->bitmask=r->bitmask;
1672  res->divmask=r->divmask;
1673  res->BitsPerExp = r->BitsPerExp;
1674  res->ExpPerLong =  r->ExpPerLong;
1675
1676  //memset: res->p_Procs=NULL;
1677  //memset: res->pFDeg=NULL;
1678  //memset: res->pLDeg=NULL;
1679  //memset: res->pFDegOrig=NULL;
1680  //memset: res->pLDegOrig=NULL;
1681  //memset: res->p_Setm=NULL;
1682  //memset: res->cf=NULL;
1683  res->options=r->options;
1684  #ifdef HAVE_RINGS
1685  res->ringtype=r->ringtype;
1686  #endif
1687  //
1688  if (r->algring!=NULL)
1689    r->algring->ref++;
1690  res->algring=r->algring;
1691  //memset: res->minideal=NULL;
1692  if (r->parameter!=NULL)
1693  {
1694    if (r->minpoly!=NULL) res->minpoly=n_Copy(r->minpoly,r->cf);
1695    int l=rPar(r);
1696    res->parameter=(char **)omAlloc(l*sizeof(char *));
1697    int i;
1698    for(i=0;i<rPar(r);i++)
1699    {
1700      res->parameter[i]=omStrDup(r->parameter[i]);
1701    }
1702    if (r->minideal!=NULL)
1703    {
1704      res->minideal=id_Copy(r->minideal,r->algring);
1705    }
1706  }
1707  if (copy_ordering == TRUE)
1708  {
1709    i=rBlocks(r);
1710    res->wvhdl   = (int **)omAlloc(i * sizeof(int *));
1711    res->order   = (int *) omAlloc(i * sizeof(int));
1712    res->block0  = (int *) omAlloc(i * sizeof(int));
1713    res->block1  = (int *) omAlloc(i * sizeof(int));
1714    for (j=0; j<i; j++)
1715    {
1716      if (r->wvhdl[j]!=NULL)
1717      {
1718        res->wvhdl[j] = (int*) omMemDup(r->wvhdl[j]);
1719      }
1720      else
1721        res->wvhdl[j]=NULL;
1722    }
1723    memcpy(res->order,r->order,i * sizeof(int));
1724    memcpy(res->block0,r->block0,i * sizeof(int));
1725    memcpy(res->block1,r->block1,i * sizeof(int));
1726  }
1727  //memset: else
1728  //memset: {
1729  //memset:   res->wvhdl = NULL;
1730  //memset:   res->order = NULL;
1731  //memset:   res->block0 = NULL;
1732  //memset:   res->block1 = NULL;
1733  //memset: }
1734
1735  res->names   = (char **)omAlloc0(rVar(r) * sizeof(char *));
1736  for (i=0; i<rVar(res); i++)
1737  {
1738    res->names[i] = omStrDup(r->names[i]);
1739  }
1740  if (r->qideal!=NULL)
1741  {
1742    if (copy_qideal)
1743    {
1744      #ifndef NDEBUG
1745      if (!copy_ordering)
1746        WerrorS("internal error: rCopy0(Q,TRUE,FALSE)");
1747      else
1748      #endif
1749      {
1750      #ifndef NDEBUG
1751        WarnS("internal bad stuff: rCopy0(Q,TRUE,TRUE)");
1752      #endif
1753        rComplete(res);
1754        res->qideal= idrCopyR_NoSort(r->qideal, r, res);
1755        rUnComplete(res);
1756      }
1757    }
1758    //memset: else res->qideal = NULL;
1759  }
1760  //memset: else res->qideal = NULL;
1761  //memset: res->GetNC() = NULL; // copy is purely commutative!!!
1762  return res;
1763}
1764
1765/*2
1766 * create a copy of the ring r, which must be equivalent to currRing
1767 * used for qring definition,..
1768 * (i.e.: normal rings: same nCopy as currRing;
1769 *        qring:        same nCopy, same idCopy as currRing)
1770 */
1771ring rCopy(ring r)
1772{
1773  if (r == NULL) return NULL;
1774  ring res=rCopy0(r,FALSE,TRUE);
1775  rComplete(res, 1); // res is purely commutative so far
1776  if (r->qideal!=NULL) res->qideal=idrCopyR_NoSort(r->qideal, r, res);
1777
1778#ifdef HAVE_PLURAL
1779  if (rIsPluralRing(r))
1780    if( nc_rCopy(res, r, true) );
1781#endif
1782
1783  return res;
1784}
1785
1786// returns TRUE, if r1 equals r2 FALSE, otherwise Equality is
1787// determined componentwise, if qr == 1, then qrideal equality is
1788// tested, as well
1789BOOLEAN rEqual(ring r1, ring r2, BOOLEAN qr)
1790{
1791  int i, j;
1792
1793  if (r1 == r2) return TRUE;
1794
1795  if (r1 == NULL || r2 == NULL) return FALSE;
1796
1797  if ((r1->cf->type != r2->cf->type)
1798  || (r1->float_len != r2->float_len)
1799  || (r1->float_len2 != r2->float_len2)
1800  || (rVar(r1) != rVar(r2))
1801  || (r1->OrdSgn != r2->OrdSgn)
1802  || (rPar(r1) != rPar(r2)))
1803    return FALSE;
1804
1805  for (i=0; i<rVar(r1); i++)
1806  {
1807    if (r1->names[i] != NULL && r2->names[i] != NULL)
1808    {
1809      if (strcmp(r1->names[i], r2->names[i])) return FALSE;
1810    }
1811    else if ((r1->names[i] != NULL) ^ (r2->names[i] != NULL))
1812    {
1813      return FALSE;
1814    }
1815  }
1816
1817  i=0;
1818  while (r1->order[i] != 0)
1819  {
1820    if (r2->order[i] == 0) return FALSE;
1821    if ((r1->order[i] != r2->order[i])
1822    || (r1->block0[i] != r2->block0[i])
1823    || (r1->block1[i] != r2->block1[i]))
1824      return FALSE;
1825    if (r1->wvhdl[i] != NULL)
1826    {
1827      if (r2->wvhdl[i] == NULL)
1828        return FALSE;
1829      for (j=0; j<r1->block1[i]-r1->block0[i]+1; j++)
1830        if (r2->wvhdl[i][j] != r1->wvhdl[i][j])
1831          return FALSE;
1832    }
1833    else if (r2->wvhdl[i] != NULL) return FALSE;
1834    i++;
1835  }
1836  if (r2->order[i] != 0) return FALSE;
1837
1838  for (i=0; i<rPar(r1);i++)
1839  {
1840      if (strcmp(r1->parameter[i], r2->parameter[i])!=0)
1841        return FALSE;
1842  }
1843
1844  if (r1->minpoly != NULL)
1845  {
1846    if (r2->minpoly == NULL) return FALSE;
1847    if (! n_Equal(r1->minpoly, r2->minpoly, r1->cf)) return FALSE;
1848  }
1849  else if (r2->minpoly != NULL) return FALSE;
1850
1851  if (qr)
1852  {
1853    if (r1->qideal != NULL)
1854    {
1855      ideal id1 = r1->qideal, id2 = r2->qideal;
1856      int i, n;
1857      poly *m1, *m2;
1858
1859      if (id2 == NULL) return FALSE;
1860      if ((n = IDELEMS(id1)) != IDELEMS(id2)) return FALSE;
1861
1862      {
1863        m1 = id1->m;
1864        m2 = id2->m;
1865        for (i=0; i<n; i++)
1866          if (! p_EqualPolys(m1[i],m2[i],r1)) return FALSE;
1867      }
1868    }
1869    else if (r2->qideal != NULL) return FALSE;
1870  }
1871
1872  return TRUE;
1873}
1874
1875// returns TRUE, if r1 and r2 represents the monomials in the same way
1876// FALSE, otherwise
1877// this is an analogue to rEqual but not so strict
1878BOOLEAN rSamePolyRep(ring r1, ring r2)
1879{
1880  int i, j;
1881
1882  if (r1 == r2) return TRUE;
1883
1884  if (r1 == NULL || r2 == NULL) return FALSE;
1885
1886  if ((r1->cf->type != r2->cf->type)
1887  || (r1->float_len != r2->float_len)
1888  || (r1->float_len2 != r2->float_len2)
1889  || (rVar(r1) != rVar(r2))
1890  || (r1->OrdSgn != r2->OrdSgn)
1891  || (rPar(r1) != rPar(r2)))
1892    return FALSE;
1893
1894  if (rVar(r1)!=rVar(r2)) return FALSE;
1895  if (rPar(r1)!=rPar(r2)) return FALSE;
1896
1897  i=0;
1898  while (r1->order[i] != 0)
1899  {
1900    if (r2->order[i] == 0) return FALSE;
1901    if ((r1->order[i] != r2->order[i])
1902    || (r1->block0[i] != r2->block0[i])
1903    || (r1->block1[i] != r2->block1[i]))
1904      return FALSE;
1905    if (r1->wvhdl[i] != NULL)
1906    {
1907      if (r2->wvhdl[i] == NULL)
1908        return FALSE;
1909      for (j=0; j<r1->block1[i]-r1->block0[i]+1; j++)
1910        if (r2->wvhdl[i][j] != r1->wvhdl[i][j])
1911          return FALSE;
1912    }
1913    else if (r2->wvhdl[i] != NULL) return FALSE;
1914    i++;
1915  }
1916  if (r2->order[i] != 0) return FALSE;
1917
1918  // we do not check minpoly
1919  // we do not check qideal
1920
1921  return TRUE;
1922}
1923
1924rOrderType_t rGetOrderType(ring r)
1925{
1926  // check for simple ordering
1927  if (rHasSimpleOrder(r))
1928  {
1929    if ((r->order[1] == ringorder_c)
1930    || (r->order[1] == ringorder_C))
1931    {
1932      switch(r->order[0])
1933      {
1934          case ringorder_dp:
1935          case ringorder_wp:
1936          case ringorder_ds:
1937          case ringorder_ws:
1938          case ringorder_ls:
1939          case ringorder_unspec:
1940            if (r->order[1] == ringorder_C
1941            ||  r->order[0] == ringorder_unspec)
1942              return rOrderType_ExpComp;
1943            return rOrderType_Exp;
1944
1945          default:
1946            assume(r->order[0] == ringorder_lp ||
1947                   r->order[0] == ringorder_rs ||
1948                   r->order[0] == ringorder_Dp ||
1949                   r->order[0] == ringorder_Wp ||
1950                   r->order[0] == ringorder_Ds ||
1951                   r->order[0] == ringorder_Ws);
1952
1953            if (r->order[1] == ringorder_c) return rOrderType_ExpComp;
1954            return rOrderType_Exp;
1955      }
1956    }
1957    else
1958    {
1959      assume((r->order[0]==ringorder_c)||(r->order[0]==ringorder_C));
1960      return rOrderType_CompExp;
1961    }
1962  }
1963  else
1964    return rOrderType_General;
1965}
1966
1967BOOLEAN rHasSimpleOrder(const ring r)
1968{
1969  if (r->order[0] == ringorder_unspec) return TRUE;
1970  int blocks = rBlocks(r) - 1;
1971  assume(blocks >= 1);
1972  if (blocks == 1) return TRUE;
1973
1974  int s = 0; 
1975  while( (s < blocks) && (r->order[s] == ringorder_IS) && (r->order[blocks-1] == ringorder_IS) )
1976  {
1977    s++;
1978    blocks--;
1979  }
1980
1981  if ((blocks - s) > 2)  return FALSE;
1982
1983  assume( blocks == s + 2 );
1984 
1985  if (
1986     (r->order[s] != ringorder_c)
1987  && (r->order[s] != ringorder_C)
1988  && (r->order[s+1] != ringorder_c)
1989  && (r->order[s+1] != ringorder_C)
1990     )
1991    return FALSE;
1992  if ((r->order[s+1] == ringorder_M)
1993  || (r->order[s] == ringorder_M))
1994    return FALSE;
1995  return TRUE;
1996}
1997
1998// returns TRUE, if simple lp or ls ordering
1999BOOLEAN rHasSimpleLexOrder(const ring r)
2000{
2001  return rHasSimpleOrder(r) &&
2002    (r->order[0] == ringorder_ls ||
2003     r->order[0] == ringorder_lp ||
2004     r->order[1] == ringorder_ls ||
2005     r->order[1] == ringorder_lp);
2006}
2007
2008BOOLEAN rOrder_is_DegOrdering(const rRingOrder_t order)
2009{
2010  switch(order)
2011  {
2012      case ringorder_dp:
2013      case ringorder_Dp:
2014      case ringorder_ds:
2015      case ringorder_Ds:
2016      case ringorder_Ws:
2017      case ringorder_Wp:
2018      case ringorder_ws:
2019      case ringorder_wp:
2020        return TRUE;
2021
2022      default:
2023        return FALSE;
2024  }
2025}
2026
2027BOOLEAN rOrder_is_WeightedOrdering(rRingOrder_t order)
2028{
2029  switch(order)
2030  {
2031      case ringorder_Ws:
2032      case ringorder_Wp:
2033      case ringorder_ws:
2034      case ringorder_wp:
2035        return TRUE;
2036
2037      default:
2038        return FALSE;
2039  }
2040}
2041
2042BOOLEAN rHasSimpleOrderAA(ring r)
2043{
2044  if (r->order[0] == ringorder_unspec) return TRUE;
2045  int blocks = rBlocks(r) - 1;
2046  assume(blocks >= 1);
2047  if (blocks == 1) return TRUE;
2048
2049  int s = 0; 
2050  while( (s < blocks) && (r->order[s] == ringorder_IS) && (r->order[blocks-1] == ringorder_IS) )
2051  {
2052    s++;
2053    blocks--;
2054  }
2055
2056  if ((blocks - s) > 3)  return FALSE;
2057 
2058//  if ((blocks > 3) || (blocks < 2)) return FALSE;
2059  if ((blocks - s) == 3)
2060  {
2061    return (((r->order[s] == ringorder_aa) && (r->order[s+1] != ringorder_M) &&
2062             ((r->order[s+2] == ringorder_c) || (r->order[s+2] == ringorder_C))) ||
2063            (((r->order[s] == ringorder_c) || (r->order[s] == ringorder_C)) &&
2064             (r->order[s+1] == ringorder_aa) && (r->order[s+2] != ringorder_M)));
2065  }
2066  else
2067  {
2068    return ((r->order[s] == ringorder_aa) && (r->order[s+1] != ringorder_M));
2069  }
2070}
2071
2072// return TRUE if p_SetComp requires p_Setm
2073BOOLEAN rOrd_SetCompRequiresSetm(ring r)
2074{
2075  if (r->typ != NULL)
2076  {
2077    int pos;
2078    for (pos=0;pos<r->OrdSize;pos++)
2079    {
2080      sro_ord* o=&(r->typ[pos]);
2081      if ((o->ord_typ == ro_syzcomp) || (o->ord_typ == ro_syz) || (o->ord_typ == ro_is) || (o->ord_typ == ro_isTemp)) return TRUE;
2082    }
2083  }
2084  return FALSE;
2085}
2086
2087// return TRUE if p->exp[r->pOrdIndex] holds total degree of p */
2088BOOLEAN rOrd_is_Totaldegree_Ordering(ring r)
2089{
2090  // Hmm.... what about Syz orderings?
2091  return (rVar(r) > 1 &&
2092          ((rHasSimpleOrder(r) &&
2093           (rOrder_is_DegOrdering((rRingOrder_t)r->order[0]) ||
2094            rOrder_is_DegOrdering(( rRingOrder_t)r->order[1]))) ||
2095           (rHasSimpleOrderAA(r) &&
2096            (rOrder_is_DegOrdering((rRingOrder_t)r->order[1]) ||
2097             rOrder_is_DegOrdering((rRingOrder_t)r->order[2])))));
2098}
2099
2100// return TRUE if p->exp[r->pOrdIndex] holds a weighted degree of p */
2101BOOLEAN rOrd_is_WeightedDegree_Ordering(ring r )
2102{
2103  // Hmm.... what about Syz orderings?
2104  return ((rVar(r) > 1) &&
2105          rHasSimpleOrder(r) &&
2106          (rOrder_is_WeightedOrdering((rRingOrder_t)r->order[0]) ||
2107           rOrder_is_WeightedOrdering(( rRingOrder_t)r->order[1])));
2108}
2109
2110BOOLEAN rIsPolyVar(int v, ring r)
2111{
2112  int  i=0;
2113  while(r->order[i]!=0)
2114  {
2115    if((r->block0[i]<=v)
2116    && (r->block1[i]>=v))
2117    {
2118      switch(r->order[i])
2119      {
2120        case ringorder_a:
2121          return (r->wvhdl[i][v-r->block0[i]]>0);
2122        case ringorder_M:
2123          return 2; /*don't know*/
2124        case ringorder_a64: /* assume: all weight are non-negative!*/
2125        case ringorder_lp:
2126        case ringorder_rs:
2127        case ringorder_dp:
2128        case ringorder_Dp:
2129        case ringorder_wp:
2130        case ringorder_Wp:
2131          return TRUE;
2132        case ringorder_ls:
2133        case ringorder_ds:
2134        case ringorder_Ds:
2135        case ringorder_ws:
2136        case ringorder_Ws:
2137          return FALSE;
2138        default:
2139          break;
2140      }
2141    }
2142    i++;
2143  }
2144  return 3; /* could not find var v*/
2145}
2146
2147#ifdef RDEBUG
2148// This should eventually become a full-fledge ring check, like pTest
2149BOOLEAN rDBTest(ring r, const char* fn, const int l)
2150{
2151  int i,j;
2152
2153  if (r == NULL)
2154  {
2155    dReportError("Null ring in %s:%d", fn, l);
2156    return FALSE;
2157  }
2158
2159
2160  if (r->N == 0) return TRUE;
2161
2162//  omCheckAddrSize(r,sizeof(ip_sring));
2163#if OM_CHECK > 0
2164  i=rBlocks(r);
2165  omCheckAddrSize(r->order,i*sizeof(int));
2166  omCheckAddrSize(r->block0,i*sizeof(int));
2167  omCheckAddrSize(r->block1,i*sizeof(int));
2168  if (r->wvhdl!=NULL)
2169  {
2170    omCheckAddrSize(r->wvhdl,i*sizeof(int *));
2171    for (j=0;j<i; j++)
2172    {
2173      if (r->wvhdl[j] != NULL) omCheckAddr(r->wvhdl[j]);
2174    }
2175  }
2176#endif
2177  if (r->VarOffset == NULL)
2178  {
2179    dReportError("Null ring VarOffset -- no rComplete (?) in n %s:%d", fn, l);
2180    return FALSE;
2181  }
2182  omCheckAddrSize(r->VarOffset,(r->N+1)*sizeof(int));
2183
2184  if ((r->OrdSize==0)!=(r->typ==NULL))
2185  {
2186    dReportError("mismatch OrdSize and typ-pointer in %s:%d");
2187    return FALSE;
2188  }
2189  omcheckAddrSize(r->typ,r->OrdSize*sizeof(*(r->typ)));
2190  omCheckAddrSize(r->VarOffset,(r->N+1)*sizeof(*(r->VarOffset)));
2191  // test assumptions:
2192  for(i=0;i<=r->N;i++) // for all variables (i = 0..N)
2193  {
2194    if(r->typ!=NULL)
2195    {
2196      for(j=0;j<r->OrdSize;j++) // for all ordering blocks (j =0..OrdSize-1)
2197      {
2198        if(r->typ[j].ord_typ == ro_isTemp)
2199        {
2200          const int p = r->typ[j].data.isTemp.suffixpos;
2201
2202          if(p <= j)
2203            dReportError("ordrec prefix %d is unmatched",j);
2204
2205          assume( p < r->OrdSize );
2206
2207          if(r->typ[p].ord_typ != ro_is)
2208            dReportError("ordrec prefix %d is unmatched (suffix: %d is wrong!!!)",j, p);
2209
2210          // Skip all intermediate blocks for undone variables:
2211          if(r->typ[j].data.isTemp.pVarOffset[i] != -1) // Check i^th variable
2212          {
2213            j = p - 1; // SKIP ALL INTERNAL BLOCKS...???
2214            continue; // To make for check OrdSize bound...
2215          }
2216        }
2217        else if (r->typ[j].ord_typ == ro_is)
2218        {
2219          // Skip all intermediate blocks for undone variables:
2220          if(r->typ[j].data.is.pVarOffset[i] != -1)
2221          {
2222            // ???
2223          }
2224
2225        }
2226        else
2227        {
2228          if (r->typ[j].ord_typ==ro_cp)
2229          {
2230            if(((short)r->VarOffset[i]) == r->typ[j].data.cp.place)
2231              dReportError("ordrec %d conflicts with var %d",j,i);
2232          }
2233          else
2234            if ((r->typ[j].ord_typ!=ro_syzcomp)
2235            && (r->VarOffset[i] == r->typ[j].data.dp.place))
2236              dReportError("ordrec %d conflicts with var %d",j,i);
2237        }
2238      }
2239    }
2240    int tmp;
2241      tmp=r->VarOffset[i] & 0xffffff;
2242      #if SIZEOF_LONG == 8
2243        if ((r->VarOffset[i] >> 24) >63)
2244      #else
2245        if ((r->VarOffset[i] >> 24) >31)
2246      #endif
2247          dReportError("bit_start out of range:%d",r->VarOffset[i] >> 24);
2248      if (i > 0 && ((tmp<0) ||(tmp>r->ExpL_Size-1)))
2249      {
2250        dReportError("varoffset out of range for var %d: %d",i,tmp);
2251      }
2252  }
2253  if(r->typ!=NULL)
2254  {
2255    for(j=0;j<r->OrdSize;j++)
2256    {
2257      if ((r->typ[j].ord_typ==ro_dp)
2258      || (r->typ[j].ord_typ==ro_wp)
2259      || (r->typ[j].ord_typ==ro_wp_neg))
2260      {
2261        if (r->typ[j].data.dp.start > r->typ[j].data.dp.end)
2262          dReportError("in ordrec %d: start(%d) > end(%d)",j,
2263            r->typ[j].data.dp.start, r->typ[j].data.dp.end);
2264        if ((r->typ[j].data.dp.start < 1)
2265        || (r->typ[j].data.dp.end > r->N))
2266          dReportError("in ordrec %d: start(%d)<1 or end(%d)>vars(%d)",j,
2267            r->typ[j].data.dp.start, r->typ[j].data.dp.end,r->N);
2268      }
2269    }
2270  }
2271  if (r->minpoly!=NULL)
2272  {
2273    omCheckAddr(r->minpoly);
2274  }
2275  //assume(r->cf!=NULL);
2276
2277  return TRUE;
2278}
2279#endif
2280
2281static void rO_Align(int &place, int &bitplace)
2282{
2283  // increment place to the next aligned one
2284  // (count as Exponent_t,align as longs)
2285  if (bitplace!=BITS_PER_LONG)
2286  {
2287    place++;
2288    bitplace=BITS_PER_LONG;
2289  }
2290}
2291
2292static void rO_TDegree(int &place, int &bitplace, int start, int end,
2293    long *o, sro_ord &ord_struct)
2294{
2295  // degree (aligned) of variables v_start..v_end, ordsgn 1
2296  rO_Align(place,bitplace);
2297  ord_struct.ord_typ=ro_dp;
2298  ord_struct.data.dp.start=start;
2299  ord_struct.data.dp.end=end;
2300  ord_struct.data.dp.place=place;
2301  o[place]=1;
2302  place++;
2303  rO_Align(place,bitplace);
2304}
2305
2306static void rO_TDegree_neg(int &place, int &bitplace, int start, int end,
2307    long *o, sro_ord &ord_struct)
2308{
2309  // degree (aligned) of variables v_start..v_end, ordsgn -1
2310  rO_Align(place,bitplace);
2311  ord_struct.ord_typ=ro_dp;
2312  ord_struct.data.dp.start=start;
2313  ord_struct.data.dp.end=end;
2314  ord_struct.data.dp.place=place;
2315  o[place]=-1;
2316  place++;
2317  rO_Align(place,bitplace);
2318}
2319
2320static void rO_WDegree(int &place, int &bitplace, int start, int end,
2321    long *o, sro_ord &ord_struct, int *weights)
2322{
2323  // weighted degree (aligned) of variables v_start..v_end, ordsgn 1
2324  while((start<end) && (weights[0]==0)) { start++; weights++; }
2325  while((start<end) && (weights[end-start]==0)) { end--; }
2326  int i;
2327  int pure_tdeg=1;
2328  for(i=start;i<=end;i++)
2329  {
2330    if(weights[i-start]!=1)
2331    {
2332      pure_tdeg=0;
2333      break;
2334    }
2335  }
2336  if (pure_tdeg)
2337  {
2338    rO_TDegree(place,bitplace,start,end,o,ord_struct);
2339    return;
2340  }
2341  rO_Align(place,bitplace);
2342  ord_struct.ord_typ=ro_wp;
2343  ord_struct.data.wp.start=start;
2344  ord_struct.data.wp.end=end;
2345  ord_struct.data.wp.place=place;
2346  ord_struct.data.wp.weights=weights;
2347  o[place]=1;
2348  place++;
2349  rO_Align(place,bitplace);
2350  for(i=start;i<=end;i++)
2351  {
2352    if(weights[i-start]<0)
2353    {
2354      ord_struct.ord_typ=ro_wp_neg;
2355      break;
2356    }
2357  }
2358}
2359
2360static void rO_WDegree64(int &place, int &bitplace, int start, int end,
2361    long *o, sro_ord &ord_struct, int64 *weights)
2362{
2363  // weighted degree (aligned) of variables v_start..v_end, ordsgn 1,
2364  // reserved 2 places
2365  rO_Align(place,bitplace);
2366  ord_struct.ord_typ=ro_wp64;
2367  ord_struct.data.wp64.start=start;
2368  ord_struct.data.wp64.end=end;
2369  ord_struct.data.wp64.place=place;
2370  ord_struct.data.wp64.weights64=weights;
2371  o[place]=1;
2372  place++;
2373  o[place]=1;
2374  place++;
2375  rO_Align(place,bitplace);
2376}
2377
2378static void rO_WDegree_neg(int &place, int &bitplace, int start, int end,
2379    long *o, sro_ord &ord_struct, int *weights)
2380{
2381  // weighted degree (aligned) of variables v_start..v_end, ordsgn -1
2382  while((start<end) && (weights[0]==0)) { start++; weights++; }
2383  while((start<end) && (weights[end-start]==0)) { end--; }
2384  rO_Align(place,bitplace);
2385  ord_struct.ord_typ=ro_wp;
2386  ord_struct.data.wp.start=start;
2387  ord_struct.data.wp.end=end;
2388  ord_struct.data.wp.place=place;
2389  ord_struct.data.wp.weights=weights;
2390  o[place]=-1;
2391  place++;
2392  rO_Align(place,bitplace);
2393  int i;
2394  for(i=start;i<=end;i++)
2395  {
2396    if(weights[i-start]<0)
2397    {
2398      ord_struct.ord_typ=ro_wp_neg;
2399      break;
2400    }
2401  }
2402}
2403
2404static void rO_LexVars(int &place, int &bitplace, int start, int end,
2405  int &prev_ord, long *o,int *v, int bits, int opt_var)
2406{
2407  // a block of variables v_start..v_end with lex order, ordsgn 1
2408  int k;
2409  int incr=1;
2410  if(prev_ord==-1) rO_Align(place,bitplace);
2411
2412  if (start>end)
2413  {
2414    incr=-1;
2415  }
2416  for(k=start;;k+=incr)
2417  {
2418    bitplace-=bits;
2419    if (bitplace < 0) { bitplace=BITS_PER_LONG-bits; place++; }
2420    o[place]=1;
2421    v[k]= place | (bitplace << 24);
2422    if (k==end) break;
2423  }
2424  prev_ord=1;
2425  if (opt_var!= -1)
2426  {
2427    assume((opt_var == end+1) ||(opt_var == end-1));
2428    if((opt_var != end+1) &&(opt_var != end-1)) WarnS("hier-2");
2429    int save_bitplace=bitplace;
2430    bitplace-=bits;
2431    if (bitplace < 0)
2432    {
2433      bitplace=save_bitplace;
2434      return;
2435    }
2436    // there is enough space for the optional var
2437    v[opt_var]=place | (bitplace << 24);
2438  }
2439}
2440
2441static void rO_LexVars_neg(int &place, int &bitplace, int start, int end,
2442  int &prev_ord, long *o,int *v, int bits, int opt_var)
2443{
2444  // a block of variables v_start..v_end with lex order, ordsgn -1
2445  int k;
2446  int incr=1;
2447  if(prev_ord==1) rO_Align(place,bitplace);
2448
2449  if (start>end)
2450  {
2451    incr=-1;
2452  }
2453  for(k=start;;k+=incr)
2454  {
2455    bitplace-=bits;
2456    if (bitplace < 0) { bitplace=BITS_PER_LONG-bits; place++; }
2457    o[place]=-1;
2458    v[k]=place | (bitplace << 24);
2459    if (k==end) break;
2460  }
2461  prev_ord=-1;
2462//  #if 0
2463  if (opt_var!= -1)
2464  {
2465    assume((opt_var == end+1) ||(opt_var == end-1));
2466    if((opt_var != end+1) &&(opt_var != end-1)) WarnS("hier-1");
2467    int save_bitplace=bitplace;
2468    bitplace-=bits;
2469    if (bitplace < 0)
2470    {
2471      bitplace=save_bitplace;
2472      return;
2473    }
2474    // there is enough space for the optional var
2475    v[opt_var]=place | (bitplace << 24);
2476  }
2477//  #endif
2478}
2479
2480static void rO_Syzcomp(int &place, int &bitplace, int &prev_ord,
2481    long *o, sro_ord &ord_struct)
2482{
2483  // ordering is derived from component number
2484  rO_Align(place,bitplace);
2485  ord_struct.ord_typ=ro_syzcomp;
2486  ord_struct.data.syzcomp.place=place;
2487  ord_struct.data.syzcomp.Components=NULL;
2488  ord_struct.data.syzcomp.ShiftedComponents=NULL;
2489  o[place]=1;
2490  prev_ord=1;
2491  place++;
2492  rO_Align(place,bitplace);
2493}
2494
2495static void rO_Syz(int &place, int &bitplace, int &prev_ord,
2496    long *o, sro_ord &ord_struct)
2497{
2498  // ordering is derived from component number
2499  // let's reserve one Exponent_t for it
2500  if ((prev_ord== 1) || (bitplace!=BITS_PER_LONG))
2501    rO_Align(place,bitplace);
2502  ord_struct.ord_typ=ro_syz;
2503  ord_struct.data.syz.place=place;
2504  ord_struct.data.syz.limit=0;
2505  ord_struct.data.syz.syz_index = NULL;
2506  ord_struct.data.syz.curr_index = 1;
2507  o[place]= -1;
2508  prev_ord=-1;
2509  place++;
2510}
2511
2512#ifndef NDEBUG
2513# define MYTEST 0
2514#else /* ifndef NDEBUG */
2515# define MYTEST 0
2516#endif /* ifndef NDEBUG */
2517
2518static void rO_ISPrefix(int &place, int &bitplace, int &prev_ord,
2519    long *o, int N, int *v, sro_ord &ord_struct)
2520{
2521  if ((prev_ord== 1) || (bitplace!=BITS_PER_LONG))
2522    rO_Align(place,bitplace);
2523  // since we add something afterwards - it's better to start with anew!?
2524
2525  ord_struct.ord_typ = ro_isTemp;
2526  ord_struct.data.isTemp.start = place;
2527  ord_struct.data.isTemp.pVarOffset = (int *)omMemDup(v);
2528  ord_struct.data.isTemp.suffixpos = -1;
2529
2530  // We will act as rO_Syz on our own!!!
2531  // Here we allocate an exponent as a level placeholder
2532  o[place]= -1;
2533  prev_ord=-1;
2534  place++;
2535}
2536static void rO_ISSuffix(int &place, int &bitplace, int &prev_ord, long *o,
2537  int N, int *v, sro_ord *tmp_typ, int &typ_i, int sgn)
2538{
2539
2540  // Let's find previous prefix:
2541  int typ_j = typ_i - 1;
2542  while(typ_j >= 0)
2543  {
2544    if( tmp_typ[typ_j].ord_typ == ro_isTemp)
2545      break;
2546    typ_j --;
2547  }
2548
2549  assume( typ_j >= 0 );
2550
2551  if( typ_j < 0 ) // Found NO prefix!!! :(
2552    return;
2553
2554  assume( tmp_typ[typ_j].ord_typ == ro_isTemp );
2555
2556  // Get saved state:
2557  const int start = tmp_typ[typ_j].data.isTemp.start;
2558  int *pVarOffset = tmp_typ[typ_j].data.isTemp.pVarOffset;
2559
2560/*
2561  // shift up all blocks
2562  while(typ_j < (typ_i-1))
2563  {
2564    tmp_typ[typ_j] = tmp_typ[typ_j+1];
2565    typ_j++;
2566  }
2567  typ_j = typ_i - 1; // No increment for typ_i
2568*/
2569  tmp_typ[typ_j].data.isTemp.suffixpos = typ_i;
2570
2571  // Let's keep that dummy for now...
2572  typ_j = typ_i; // the typ to change!
2573  typ_i++; // Just for now...
2574
2575
2576  for( int i = 0; i <= N; i++ ) // Note [0] == component !!! No Skip?
2577  {
2578    // Was i-th variable allocated inbetween?
2579    if( v[i] != pVarOffset[i] )
2580    {
2581      pVarOffset[i] = v[i]; // Save for later...
2582      v[i] = -1; // Undo!
2583      assume( pVarOffset[i] != -1 );
2584    }
2585    else
2586      pVarOffset[i] = -1; // No change here...
2587  }
2588
2589  if( pVarOffset[0] != -1 )
2590    pVarOffset[0] &= 0x0fff;
2591
2592  sro_ord &ord_struct = tmp_typ[typ_j];
2593
2594
2595  ord_struct.ord_typ = ro_is;
2596  ord_struct.data.is.start = start;
2597  ord_struct.data.is.end   = place;
2598  ord_struct.data.is.pVarOffset = pVarOffset;
2599
2600
2601  // What about component???
2602//   if( v[0] != -1 ) // There is a component already...???
2603//     if( o[ v[0] & 0x0fff ] == sgn )
2604//     {
2605//       pVarOffset[0] = -1; // NEVER USED Afterwards...
2606//       return;
2607//     }
2608
2609
2610  // Moreover: we need to allocate the module component (v[0]) here!
2611  if( v[0] == -1) // It's possible that there was module component v0 at the begining (before prefix)!
2612  {
2613    // Start with a whole long exponent
2614    if( bitplace != BITS_PER_LONG )
2615      rO_Align(place, bitplace);
2616
2617    assume( bitplace == BITS_PER_LONG );
2618    bitplace -= BITS_PER_LONG;
2619    assume(bitplace == 0);
2620    v[0] = place | (bitplace << 24); // Never mind whether pVarOffset[0] > 0!!!
2621    o[place] = sgn; // Singnum for component ordering
2622    prev_ord = sgn;
2623  }
2624}
2625
2626
2627static unsigned long rGetExpSize(unsigned long bitmask, int & bits)
2628{
2629  if (bitmask == 0)
2630  {
2631    bits=16; bitmask=0xffff;
2632  }
2633  else if (bitmask <= 1L)
2634  {
2635    bits=1; bitmask = 1L;
2636  }
2637  else if (bitmask <= 3L)
2638  {
2639    bits=2; bitmask = 3L;
2640  }
2641  else if (bitmask <= 7L)
2642  {
2643    bits=3; bitmask=7L;
2644  }
2645  else if (bitmask <= 0xfL)
2646  {
2647    bits=4; bitmask=0xfL;
2648  }
2649  else if (bitmask <= 0x1fL)
2650  {
2651    bits=5; bitmask=0x1fL;
2652  }
2653  else if (bitmask <= 0x3fL)
2654  {
2655    bits=6; bitmask=0x3fL;
2656  }
2657#if SIZEOF_LONG == 8
2658  else if (bitmask <= 0x7fL)
2659  {
2660    bits=7; bitmask=0x7fL; /* 64 bit longs only */
2661  }
2662#endif
2663  else if (bitmask <= 0xffL)
2664  {
2665    bits=8; bitmask=0xffL;
2666  }
2667#if SIZEOF_LONG == 8
2668  else if (bitmask <= 0x1ffL)
2669  {
2670    bits=9; bitmask=0x1ffL; /* 64 bit longs only */
2671  }
2672#endif
2673  else if (bitmask <= 0x3ffL)
2674  {
2675    bits=10; bitmask=0x3ffL;
2676  }
2677#if SIZEOF_LONG == 8
2678  else if (bitmask <= 0xfffL)
2679  {
2680    bits=12; bitmask=0xfff; /* 64 bit longs only */
2681  }
2682#endif
2683  else if (bitmask <= 0xffffL)
2684  {
2685    bits=16; bitmask=0xffffL;
2686  }
2687#if SIZEOF_LONG == 8
2688  else if (bitmask <= 0xfffffL)
2689  {
2690    bits=20; bitmask=0xfffffL; /* 64 bit longs only */
2691  }
2692  else if (bitmask <= 0xffffffffL)
2693  {
2694    bits=32; bitmask=0xffffffffL;
2695  }
2696  else if (bitmask <= 0x7fffffffffffffffL)
2697  {
2698    bits=63; bitmask=0x7fffffffffffffffL; /* for overflow tests*/
2699  }
2700  else
2701  {
2702    bits=63; bitmask=0x7fffffffffffffffL; /* for overflow tests*/
2703  }
2704#else
2705  else if (bitmask <= 0x7fffffff)
2706  {
2707    bits=31; bitmask=0x7fffffff; /* for overflow tests*/
2708  }
2709  else
2710  {
2711    bits=31; bitmask=0x7fffffffL; /* for overflow tests*/
2712  }
2713#endif
2714  return bitmask;
2715}
2716
2717/*2
2718* optimize rGetExpSize for a block of N variables, exp <=bitmask
2719*/
2720static unsigned long rGetExpSize(unsigned long bitmask, int & bits, int N)
2721{
2722  bitmask =rGetExpSize(bitmask, bits);
2723  int vars_per_long=BIT_SIZEOF_LONG/bits;
2724  int bits1;
2725  loop
2726  {
2727    if (bits == BIT_SIZEOF_LONG-1)
2728    {
2729      bits =  BIT_SIZEOF_LONG - 1;
2730      return LONG_MAX;
2731    }
2732    unsigned long bitmask1 =rGetExpSize(bitmask+1, bits1);
2733    int vars_per_long1=BIT_SIZEOF_LONG/bits1;
2734    if ((((N+vars_per_long-1)/vars_per_long) ==
2735         ((N+vars_per_long1-1)/vars_per_long1)))
2736    {
2737      vars_per_long=vars_per_long1;
2738      bits=bits1;
2739      bitmask=bitmask1;
2740    }
2741    else
2742    {
2743      return bitmask; /* and bits */
2744    }
2745  }
2746}
2747
2748
2749/*2
2750 * create a copy of the ring r, which must be equivalent to currRing
2751 * used for std computations
2752 * may share data structures with currRing
2753 * DOES CALL rComplete
2754 */
2755ring rModifyRing(ring r, BOOLEAN omit_degree,
2756                         BOOLEAN omit_comp,
2757                         unsigned long exp_limit)
2758{
2759  assume (r != NULL );
2760  assume (exp_limit > 1);
2761  BOOLEAN need_other_ring;
2762  BOOLEAN omitted_degree = FALSE;
2763
2764  int iNeedInducedOrderingSetup = 0; ///< How many induced ordering block do we have?
2765  int bits;
2766
2767  exp_limit=rGetExpSize(exp_limit, bits, r->N);
2768  need_other_ring = (exp_limit != r->bitmask);
2769
2770  int nblocks=rBlocks(r);
2771  int *order=(int*)omAlloc0((nblocks+1)*sizeof(int));
2772  int *block0=(int*)omAlloc0((nblocks+1)*sizeof(int));
2773  int *block1=(int*)omAlloc0((nblocks+1)*sizeof(int));
2774  int **wvhdl=(int**)omAlloc0((nblocks+1)*sizeof(int *));
2775
2776  int i=0;
2777  int j=0; /*  i index in r, j index in res */
2778
2779  for( int r_ord=r->order[i]; (r_ord != 0) && (i < nblocks); j++, r_ord=r->order[++i])
2780  {
2781    BOOLEAN copy_block_index=TRUE;
2782
2783    if (r->block0[i]==r->block1[i])
2784    {
2785      switch(r_ord)
2786      {
2787        case ringorder_wp:
2788        case ringorder_dp:
2789        case ringorder_Wp:
2790        case ringorder_Dp:
2791          r_ord=ringorder_lp;
2792          break;
2793        case ringorder_Ws:
2794        case ringorder_Ds:
2795        case ringorder_ws:
2796        case ringorder_ds:
2797          r_ord=ringorder_ls;
2798          break;
2799        default:
2800          break;
2801      }
2802    }
2803    switch(r_ord)
2804    {
2805      case ringorder_S:
2806      {
2807#ifndef NDEBUG
2808        Warn("Error: unhandled ordering in rModifyRing: ringorder_S = [%d]", r_ord);
2809#endif
2810        order[j]=r_ord; /*r->order[i];*/
2811        break;
2812      }
2813      case ringorder_C:
2814      case ringorder_c:
2815        if (!omit_comp)
2816        {
2817          order[j]=r_ord; /*r->order[i]*/;
2818        }
2819        else
2820        {
2821          j--;
2822          need_other_ring=TRUE;
2823          omit_comp=FALSE;
2824          copy_block_index=FALSE;
2825        }
2826        break;
2827      case ringorder_wp:
2828      case ringorder_dp:
2829      case ringorder_ws:
2830      case ringorder_ds:
2831        if(!omit_degree)
2832        {
2833          order[j]=r_ord; /*r->order[i]*/;
2834        }
2835        else
2836        {
2837          order[j]=ringorder_rs;
2838          need_other_ring=TRUE;
2839          omit_degree=FALSE;
2840          omitted_degree = TRUE;
2841        }
2842        break;
2843      case ringorder_Wp:
2844      case ringorder_Dp:
2845      case ringorder_Ws:
2846      case ringorder_Ds:
2847        if(!omit_degree)
2848        {
2849          order[j]=r_ord; /*r->order[i];*/
2850        }
2851        else
2852        {
2853          order[j]=ringorder_lp;
2854          need_other_ring=TRUE;
2855          omit_degree=FALSE;
2856          omitted_degree = TRUE;
2857        }
2858        break;
2859      case ringorder_IS:
2860      {
2861        if (omit_comp)
2862        {
2863#ifndef NDEBUG
2864          Warn("Error: WRONG USAGE of rModifyRing: cannot omit component due to the ordering block [%d]: %d (ringorder_IS)", i, r_ord);
2865#endif
2866          omit_comp = FALSE;
2867        }
2868        order[j]=r_ord; /*r->order[i];*/
2869        iNeedInducedOrderingSetup++;
2870        break;
2871      }
2872      case ringorder_s:
2873      {
2874        assume((i == 0) && (j == 0));
2875        if (omit_comp)
2876        {
2877#ifndef NDEBUG
2878          Warn("WRONG USAGE? of rModifyRing: omitting component due to the ordering block [%d]: %d (ringorder_s)", i, r_ord);
2879#endif
2880          omit_comp = FALSE;
2881        }
2882        order[j]=r_ord; /*r->order[i];*/
2883        break;
2884      }
2885      default:
2886        order[j]=r_ord; /*r->order[i];*/
2887        break;
2888    }
2889    if (copy_block_index)
2890    {
2891      block0[j]=r->block0[i];
2892      block1[j]=r->block1[i];
2893      wvhdl[j]=r->wvhdl[i];
2894    }
2895
2896    // order[j]=ringorder_no; //  done by omAlloc0
2897  }
2898  if(!need_other_ring)
2899  {
2900    omFreeSize(order,(nblocks+1)*sizeof(int));
2901    omFreeSize(block0,(nblocks+1)*sizeof(int));
2902    omFreeSize(block1,(nblocks+1)*sizeof(int));
2903    omFreeSize(wvhdl,(nblocks+1)*sizeof(int *));
2904    return r;
2905  }
2906  ring res=(ring)omAlloc0Bin(sip_sring_bin);
2907  *res = *r;
2908
2909#ifdef HAVE_PLURAL
2910  res->GetNC() = NULL;
2911#endif
2912
2913  // res->qideal, res->idroot ???
2914  res->wvhdl=wvhdl;
2915  res->order=order;
2916  res->block0=block0;
2917  res->block1=block1;
2918  res->bitmask=exp_limit;
2919  int tmpref=r->cf->ref;
2920  rComplete(res, 1);
2921  r->cf->ref=tmpref;
2922
2923  // adjust res->pFDeg: if it was changed globally, then
2924  // it must also be changed for new ring
2925  if (r->pFDegOrig != res->pFDegOrig &&
2926           rOrd_is_WeightedDegree_Ordering(r))
2927  {
2928    // still might need adjustment for weighted orderings
2929    // and omit_degree
2930    res->firstwv = r->firstwv;
2931    res->firstBlockEnds = r->firstBlockEnds;
2932    res->pFDeg = res->pFDegOrig = p_WFirstTotalDegree;
2933  }
2934  if (omitted_degree)
2935    res->pLDeg = res->pLDegOrig = r->pLDegOrig;
2936
2937  rOptimizeLDeg(res);
2938
2939  // set syzcomp
2940  if (res->typ != NULL)
2941  {
2942    if( res->typ[0].ord_typ == ro_syz) // "s" Always on [0] place!
2943    {
2944      res->typ[0] = r->typ[0]; // Copy struct!? + setup the same limit!
2945
2946      if (r->typ[0].data.syz.limit > 0)
2947      {
2948        res->typ[0].data.syz.syz_index
2949          = (int*) omAlloc((r->typ[0].data.syz.limit +1)*sizeof(int));
2950        memcpy(res->typ[0].data.syz.syz_index, r->typ[0].data.syz.syz_index,
2951              (r->typ[0].data.syz.limit +1)*sizeof(int));
2952      }
2953    }
2954
2955    if( iNeedInducedOrderingSetup > 0 )
2956    {
2957      for(j = 0, i = 0; (i < nblocks) && (iNeedInducedOrderingSetup > 0); i++)
2958        if( res->typ[i].ord_typ == ro_is ) // Search for suffixes!
2959        {
2960          ideal F = idrHeadR(r->typ[i].data.is.F, r, res); // Copy F from r into res!
2961          assume(
2962            rSetISReference( res, 
2963              F,  // WILL BE COPIED!
2964              r->typ[i].data.is.limit,
2965              j++,
2966              r->typ[i].data.is.componentWeights // WILL BE COPIED
2967              )
2968            );
2969          id_Delete(&F, res);
2970          iNeedInducedOrderingSetup--;
2971        }
2972    } // Process all induced Ordering blocks! ...
2973  }
2974  // the special case: homog (omit_degree) and 1 block rs: that is global:
2975  // it comes from dp
2976  res->OrdSgn=r->OrdSgn;
2977
2978
2979#ifdef HAVE_PLURAL
2980  if (rIsPluralRing(r))
2981  {
2982    if ( nc_rComplete(r, res, false) ) // no qideal!
2983    {
2984#ifndef NDEBUG
2985      WarnS("error in nc_rComplete");
2986#endif
2987      // cleanup?
2988
2989//      rDelete(res);
2990//      return r;
2991
2992      // just go on..
2993    }
2994
2995    if( rIsSCA(r) )
2996    {
2997      if( !sca_Force(res, scaFirstAltVar(r), scaLastAltVar(r)) )
2998      WarnS("error in sca_Force!");
2999    }
3000  }
3001#endif
3002
3003  return res;
3004}
3005
3006// construct Wp,C ring
3007ring rModifyRing_Wp(ring r, int* weights)
3008{
3009  ring res=(ring)omAlloc0Bin(sip_sring_bin);
3010  *res = *r;
3011#ifdef HAVE_PLURAL
3012  res->GetNC() = NULL;
3013#endif
3014
3015  /*weights: entries for 3 blocks: NULL*/
3016  res->wvhdl = (int **)omAlloc0(3 * sizeof(int *));
3017  /*order: Wp,C,0*/
3018  res->order = (int *) omAlloc(3 * sizeof(int *));
3019  res->block0 = (int *)omAlloc0(3 * sizeof(int *));
3020  res->block1 = (int *)omAlloc0(3 * sizeof(int *));
3021  /* ringorder Wp for the first block: var 1..r->N */
3022  res->order[0]  = ringorder_Wp;
3023  res->block0[0] = 1;
3024  res->block1[0] = r->N;
3025  res->wvhdl[0] = weights;
3026  /* ringorder C for the second block: no vars */
3027  res->order[1]  = ringorder_C;
3028  /* the last block: everything is 0 */
3029  res->order[2]  = 0;
3030  /*polynomial ring*/
3031  res->OrdSgn    = 1;
3032
3033  int tmpref=r->cf->ref;
3034  rComplete(res, 1);
3035  r->cf->ref=tmpref;
3036#ifdef HAVE_PLURAL
3037  if (rIsPluralRing(r))
3038  {
3039    if ( nc_rComplete(r, res, false) ) // no qideal!
3040    {
3041#ifndef NDEBUG
3042      WarnS("error in nc_rComplete");
3043#endif
3044      // cleanup?
3045
3046//      rDelete(res);
3047//      return r;
3048
3049      // just go on..
3050    }
3051  }
3052#endif
3053  return res;
3054}
3055
3056// construct lp, C ring with r->N variables, r->names vars....
3057ring rModifyRing_Simple(ring r, BOOLEAN ommit_degree, BOOLEAN ommit_comp, unsigned long exp_limit, BOOLEAN &simple)
3058{
3059  simple=TRUE;
3060  if (!rHasSimpleOrder(r))
3061  {
3062    simple=FALSE; // sorting needed
3063    assume (r != NULL );
3064    assume (exp_limit > 1);
3065    int bits;
3066
3067    exp_limit=rGetExpSize(exp_limit, bits, r->N);
3068
3069    int nblocks=1+(ommit_comp!=0);
3070    int *order=(int*)omAlloc0((nblocks+1)*sizeof(int));
3071    int *block0=(int*)omAlloc0((nblocks+1)*sizeof(int));
3072    int *block1=(int*)omAlloc0((nblocks+1)*sizeof(int));
3073    int **wvhdl=(int**)omAlloc0((nblocks+1)*sizeof(int *));
3074
3075    order[0]=ringorder_lp;
3076    block0[0]=1;
3077    block1[0]=r->N;
3078    if (!ommit_comp)
3079    {
3080      order[1]=ringorder_C;
3081    }
3082    ring res=(ring)omAlloc0Bin(sip_sring_bin);
3083    *res = *r;
3084#ifdef HAVE_PLURAL
3085    res->GetNC() = NULL;
3086#endif
3087    // res->qideal, res->idroot ???
3088    res->wvhdl=wvhdl;
3089    res->order=order;
3090    res->block0=block0;
3091    res->block1=block1;
3092    res->bitmask=exp_limit;
3093    int tmpref=r->cf->ref;
3094    rComplete(res, 1);
3095    r->cf->ref=tmpref;
3096
3097#ifdef HAVE_PLURAL
3098    if (rIsPluralRing(r))
3099    {
3100      if ( nc_rComplete(r, res, false) ) // no qideal!
3101      {
3102#ifndef NDEBUG
3103        WarnS("error in nc_rComplete");
3104#endif
3105        // cleanup?
3106
3107//      rDelete(res);
3108//      return r;
3109
3110      // just go on..
3111      }
3112    }
3113#endif
3114
3115    rOptimizeLDeg(res);
3116
3117    return res;
3118  }
3119  return rModifyRing(r, ommit_degree, ommit_comp, exp_limit);
3120}
3121
3122void rKillModifiedRing_Simple(ring r)
3123{
3124  rKillModifiedRing(r);
3125}
3126
3127
3128void rKillModifiedRing(ring r)
3129{
3130  rUnComplete(r);
3131  omFree(r->order);
3132  omFree(r->block0);
3133  omFree(r->block1);
3134  omFree(r->wvhdl);
3135  omFreeBin(r,sip_sring_bin);
3136}
3137
3138void rKillModified_Wp_Ring(ring r)
3139{
3140  rUnComplete(r);
3141  omFree(r->order);
3142  omFree(r->block0);
3143  omFree(r->block1);
3144  omFree(r->wvhdl[0]);
3145  omFree(r->wvhdl);
3146  omFreeBin(r,sip_sring_bin);
3147}
3148
3149static void rSetOutParams(ring r)
3150{
3151  r->VectorOut = (r->order[0] == ringorder_c);
3152  r->ShortOut = TRUE;
3153  {
3154    int i;
3155    if (r->parameter!=NULL)
3156    {
3157      for (i=0;i<rPar(r);i++)
3158      {
3159        if(strlen(r->parameter[i])>1)
3160        {
3161          r->ShortOut=FALSE;
3162          break;
3163        }
3164      }
3165    }
3166    if (r->ShortOut)
3167    {
3168      // Hmm... sometimes (e.g., from maGetPreimage) new variables
3169      // are introduced, but their names are never set
3170      // hence, we do the following awkward trick
3171      int N = omSizeWOfAddr(r->names);
3172      if (r->N < N) N = r->N;
3173
3174      for (i=(N-1);i>=0;i--)
3175      {
3176        if(r->names[i] != NULL && strlen(r->names[i])>1)
3177        {
3178          r->ShortOut=FALSE;
3179          break;
3180        }
3181      }
3182    }
3183  }
3184  r->CanShortOut = r->ShortOut;
3185}
3186
3187/*2
3188* sets r->MixedOrder and r->ComponentOrder for orderings with more than one block
3189* block of variables (ip is the block number, o_r the number of the ordering)
3190* o is the position of the orderingering in r
3191*/
3192static void rHighSet(ring r, int o_r, int o)
3193{
3194  switch(o_r)
3195  {
3196    case ringorder_lp:
3197    case ringorder_dp:
3198    case ringorder_Dp:
3199    case ringorder_wp:
3200    case ringorder_Wp:
3201    case ringorder_rp:
3202    case ringorder_a:
3203    case ringorder_aa:
3204    case ringorder_a64:
3205      if (r->OrdSgn==-1) r->MixedOrder=TRUE;
3206      break;
3207    case ringorder_ls:
3208    case ringorder_rs:
3209    case ringorder_ds:
3210    case ringorder_Ds:
3211    case ringorder_s:
3212      break;
3213    case ringorder_ws:
3214    case ringorder_Ws:
3215      if (r->wvhdl[o]!=NULL)
3216      {
3217        int i;
3218        for(i=r->block1[o]-r->block0[o];i>=0;i--)
3219          if (r->wvhdl[o][i]<0) { r->MixedOrder=TRUE; break; }
3220      }
3221      break;
3222    case ringorder_c:
3223      r->ComponentOrder=1;
3224      break;
3225    case ringorder_C:
3226    case ringorder_S:
3227      r->ComponentOrder=-1;
3228      break;
3229    case ringorder_M:
3230      r->LexOrder=TRUE;
3231      break;
3232    case ringorder_IS:
3233    { // TODO: What is r->ComponentOrder???
3234      r->MixedOrder=TRUE;
3235/*
3236      if( r->block0[o] != 0 ) // Suffix has the comonent
3237        r->ComponentOrder = r->block0[o];
3238      else // Prefix has level...
3239        r->ComponentOrder=-1;
3240*/
3241      break;
3242    }
3243
3244    default:
3245      dReportError("wrong internal ordering:%d at %s, l:%d\n",o_r,__FILE__,__LINE__);
3246  }
3247}
3248
3249static void rSetFirstWv(ring r, int i, int* order, int* block1, int** wvhdl)
3250{
3251  // cheat for ringorder_aa
3252  if (order[i] == ringorder_aa)
3253    i++;
3254  if(block1[i]!=r->N) r->LexOrder=TRUE;
3255  r->firstBlockEnds=block1[i];
3256  r->firstwv = wvhdl[i];
3257  if ((order[i]== ringorder_ws)
3258  || (order[i]==ringorder_Ws)
3259  || (order[i]== ringorder_wp)
3260  || (order[i]==ringorder_Wp)
3261  || (order[i]== ringorder_a)
3262   /*|| (order[i]==ringorder_A)*/)
3263  {
3264    int j;
3265    for(j=block1[i]-r->block0[i];j>=0;j--)
3266    {
3267      if (r->firstwv[j]<0) r->MixedOrder=TRUE;
3268      if (r->firstwv[j]==0) r->LexOrder=TRUE;
3269    }
3270  }
3271  else if (order[i]==ringorder_a64)
3272  {
3273    int j;
3274    int64 *w=rGetWeightVec(r);
3275    for(j=block1[i]-r->block0[i];j>=0;j--)
3276    {
3277      if (w[j]==0) r->LexOrder=TRUE;
3278    }
3279  }
3280}
3281
3282static void rOptimizeLDeg(ring r)
3283{
3284  if (r->pFDeg == p_Deg)
3285  {
3286    if (r->pLDeg == pLDeg1)
3287      r->pLDeg = pLDeg1_Deg;
3288    if (r->pLDeg == pLDeg1c)
3289      r->pLDeg = pLDeg1c_Deg;
3290  }
3291  else if (r->pFDeg == p_Totaldegree)
3292  {
3293    if (r->pLDeg == pLDeg1)
3294      r->pLDeg = pLDeg1_Totaldegree;
3295    if (r->pLDeg == pLDeg1c)
3296      r->pLDeg = pLDeg1c_Totaldegree;
3297  }
3298  else if (r->pFDeg == p_WFirstTotalDegree)
3299  {
3300    if (r->pLDeg == pLDeg1)
3301      r->pLDeg = pLDeg1_WFirstTotalDegree;
3302    if (r->pLDeg == pLDeg1c)
3303      r->pLDeg = pLDeg1c_WFirstTotalDegree;
3304  }
3305}
3306
3307// set pFDeg, pLDeg, MixOrder, ComponentOrder, etc
3308static void rSetDegStuff(ring r)
3309{
3310  int* order = r->order;
3311  int* block0 = r->block0;
3312  int* block1 = r->block1;
3313  int** wvhdl = r->wvhdl;
3314
3315  if (order[0]==ringorder_S ||order[0]==ringorder_s || order[0]==ringorder_IS)
3316  {
3317    order++;
3318    block0++;
3319    block1++;
3320    wvhdl++;
3321  }
3322  r->LexOrder = FALSE;
3323  r->MixedOrder = FALSE;
3324  r->ComponentOrder = 1;
3325  r->pFDeg = p_Totaldegree;
3326  r->pLDeg = (r->OrdSgn == 1 ? pLDegb : pLDeg0);
3327
3328  /*======== ordering type is (_,c) =========================*/
3329  if ((order[0]==ringorder_unspec) || (order[1] == 0)
3330      ||(
3331    ((order[1]==ringorder_c)||(order[1]==ringorder_C)
3332     ||(order[1]==ringorder_S) 
3333     ||(order[1]==ringorder_s))
3334    && (order[0]!=ringorder_M)
3335    && (order[2]==0))
3336    )
3337  {
3338    if ((order[0]!=ringorder_unspec)
3339    && ((order[1]==ringorder_C)||(order[1]==ringorder_S)||
3340        (order[1]==ringorder_s)))
3341      r->ComponentOrder=-1;
3342    if (r->OrdSgn == -1) r->pLDeg = pLDeg0c;
3343    if ((order[0] == ringorder_lp)
3344    || (order[0] == ringorder_ls)
3345    || (order[0] == ringorder_rp)
3346    || (order[0] == ringorder_rs))
3347    {
3348      r->LexOrder=TRUE;
3349      r->pLDeg = pLDeg1c;
3350      r->pFDeg = p_Totaldegree;
3351    }
3352    if ((order[0] == ringorder_a)
3353    || (order[0] == ringorder_wp)
3354    || (order[0] == ringorder_Wp)
3355    || (order[0] == ringorder_ws)
3356    || (order[0] == ringorder_Ws))
3357      r->pFDeg = p_WFirstTotalDegree;
3358    r->firstBlockEnds=block1[0];
3359    r->firstwv = wvhdl[0];
3360  }
3361  /*======== ordering type is (c,_) =========================*/
3362  else if (((order[0]==ringorder_c)
3363            ||(order[0]==ringorder_C)
3364            ||(order[0]==ringorder_S)
3365            ||(order[0]==ringorder_s))
3366  && (order[1]!=ringorder_M)
3367  &&  (order[2]==0))
3368  {
3369    if ((order[0]==ringorder_C)||(order[0]==ringorder_S)||
3370        order[0]==ringorder_s)
3371      r->ComponentOrder=-1;
3372    if ((order[1] == ringorder_lp)
3373    || (order[1] == ringorder_ls)
3374    || (order[1] == ringorder_rp)
3375    || order[1] == ringorder_rs)
3376    {
3377      r->LexOrder=TRUE;
3378      r->pLDeg = pLDeg1c;
3379      r->pFDeg = p_Totaldegree;
3380    }
3381    r->firstBlockEnds=block1[1];
3382    if (wvhdl!=NULL) r->firstwv = wvhdl[1];
3383    if ((order[1] == ringorder_a)
3384    || (order[1] == ringorder_wp)
3385    || (order[1] == ringorder_Wp)
3386    || (order[1] == ringorder_ws)
3387    || (order[1] == ringorder_Ws))
3388      r->pFDeg = p_WFirstTotalDegree;
3389  }
3390  /*------- more than one block ----------------------*/
3391  else
3392  {
3393    if ((r->VectorOut)||(order[0]==ringorder_C)||(order[0]==ringorder_S)||(order[0]==ringorder_s))
3394    {
3395      rSetFirstWv(r, 1, order, block1, wvhdl);
3396    }
3397    else
3398      rSetFirstWv(r, 0, order, block1, wvhdl);
3399
3400    /*the number of orderings:*/
3401    int i = 0;
3402    while (order[++i] != 0);
3403    do
3404    {
3405      i--;
3406      rHighSet(r, order[i],i);
3407    }
3408    while (i != 0);
3409
3410    if ((order[0]!=ringorder_c)
3411        && (order[0]!=ringorder_C)
3412        && (order[0]!=ringorder_S)
3413        && (order[0]!=ringorder_s))
3414    {
3415      r->pLDeg = pLDeg1c;
3416    }
3417    else
3418    {
3419      r->pLDeg = pLDeg1;
3420    }
3421    r->pFDeg = p_WTotaldegree; // may be improved: p_Totaldegree for lp/dp/ls/.. blocks
3422  }
3423 
3424  if (rOrd_is_Totaldegree_Ordering(r) || rOrd_is_WeightedDegree_Ordering(r))
3425    r->pFDeg = p_Deg;
3426
3427  r->pFDegOrig = r->pFDeg;
3428  r->pLDegOrig = r->pLDeg;
3429  rOptimizeLDeg(r);
3430}
3431
3432/*2
3433* set NegWeightL_Size, NegWeightL_Offset
3434*/
3435static void rSetNegWeight(ring r)
3436{
3437  int i,l;
3438  if (r->typ!=NULL)
3439  {
3440    l=0;
3441    for(i=0;i<r->OrdSize;i++)
3442    {
3443      if(r->typ[i].ord_typ==ro_wp_neg) l++;
3444    }
3445    if (l>0)
3446    {
3447      r->NegWeightL_Size=l;
3448      r->NegWeightL_Offset=(int *) omAlloc(l*sizeof(int));
3449      l=0;
3450      for(i=0;i<r->OrdSize;i++)
3451      {
3452        if(r->typ[i].ord_typ==ro_wp_neg)
3453        {
3454          r->NegWeightL_Offset[l]=r->typ[i].data.wp.place;
3455          l++;
3456        }
3457      }
3458      return;
3459    }
3460  }
3461  r->NegWeightL_Size = 0;
3462  r->NegWeightL_Offset = NULL;
3463}
3464
3465static void rSetOption(ring r)
3466{
3467  // set redthrough
3468  if (!TEST_OPT_OLDSTD && r->OrdSgn == 1 && ! r->LexOrder)
3469    r->options |= Sy_bit(OPT_REDTHROUGH);
3470  else
3471    r->options &= ~Sy_bit(OPT_REDTHROUGH);
3472
3473  // set intStrategy
3474#ifdef HAVE_RINGS
3475  if (rField_is_Extension(r) || rField_is_Q(r) || rField_is_Ring(r))
3476#else
3477  if (rField_is_Extension(r) || rField_is_Q(r))
3478#endif
3479    r->options |= Sy_bit(OPT_INTSTRATEGY);
3480  else
3481    r->options &= ~Sy_bit(OPT_INTSTRATEGY);
3482
3483  // set redTail
3484  if (r->LexOrder || r->OrdSgn == -1 || rField_is_Extension(r))
3485    r->options &= ~Sy_bit(OPT_REDTAIL);
3486  else
3487    r->options |= Sy_bit(OPT_REDTAIL);
3488}
3489
3490static void rCheckOrdSgn(ring r,int i/*current block*/);
3491
3492BOOLEAN rComplete(ring r, int force)
3493{
3494  if (r->VarOffset!=NULL && force == 0) return FALSE;
3495  rSetOutParams(r);
3496  int n=rBlocks(r)-1;
3497  int i;
3498  int bits;
3499  r->bitmask=rGetExpSize(r->bitmask,bits,r->N);
3500  r->BitsPerExp = bits;
3501  r->ExpPerLong = BIT_SIZEOF_LONG / bits;
3502  r->divmask=rGetDivMask(bits);
3503
3504  // will be used for ordsgn:
3505  long *tmp_ordsgn=(long *)omAlloc0(3*(n+r->N)*sizeof(long));
3506  // will be used for VarOffset:
3507  int *v=(int *)omAlloc((r->N+1)*sizeof(int));
3508  for(i=r->N; i>=0 ; i--)
3509  {
3510    v[i]=-1;
3511  }
3512  sro_ord *tmp_typ=(sro_ord *)omAlloc0(3*(n+r->N)*sizeof(sro_ord));
3513  int typ_i=0;
3514  int prev_ordsgn=0;
3515
3516  // fill in v, tmp_typ, tmp_ordsgn, determine typ_i (== ordSize)
3517  int j=0;
3518  int j_bits=BITS_PER_LONG;
3519
3520  BOOLEAN need_to_add_comp=FALSE; // Only for ringorder_s and ringorder_S!
3521
3522  for(i=0;i<n;i++)
3523  {
3524    tmp_typ[typ_i].order_index=i;
3525    switch (r->order[i])
3526    {
3527      case ringorder_a:
3528      case ringorder_aa:
3529        rO_WDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,tmp_typ[typ_i],
3530                   r->wvhdl[i]);
3531        typ_i++;
3532        break;
3533
3534      case ringorder_a64:
3535        rO_WDegree64(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3536                     tmp_typ[typ_i], (int64 *)(r->wvhdl[i]));
3537        typ_i++;
3538        break;
3539
3540      case ringorder_c:
3541        rO_Align(j, j_bits);
3542        rO_LexVars_neg(j, j_bits, 0,0, prev_ordsgn,tmp_ordsgn,v,BITS_PER_LONG, -1);
3543        break;
3544
3545      case ringorder_C:
3546        rO_Align(j, j_bits);
3547        rO_LexVars(j, j_bits, 0,0, prev_ordsgn,tmp_ordsgn,v,BITS_PER_LONG, -1);
3548        break;
3549
3550      case ringorder_M:
3551        {
3552          int k,l;
3553          k=r->block1[i]-r->block0[i]+1; // number of vars
3554          for(l=0;l<k;l++)
3555          {
3556            rO_WDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3557                       tmp_typ[typ_i],
3558                       r->wvhdl[i]+(r->block1[i]-r->block0[i]+1)*l);
3559            typ_i++;
3560          }
3561          break;
3562        }
3563
3564      case ringorder_lp:
3565        rO_LexVars(j, j_bits, r->block0[i],r->block1[i], prev_ordsgn,
3566                   tmp_ordsgn,v,bits, -1);
3567        break;
3568
3569      case ringorder_ls:
3570        rO_LexVars_neg(j, j_bits, r->block0[i],r->block1[i], prev_ordsgn,
3571                       tmp_ordsgn,v, bits, -1);
3572        rCheckOrdSgn(r,i);
3573        break;
3574
3575      case ringorder_rs:
3576        rO_LexVars_neg(j, j_bits, r->block1[i],r->block0[i], prev_ordsgn,
3577                       tmp_ordsgn,v, bits, -1);
3578        rCheckOrdSgn(r,i);
3579        break;
3580
3581      case ringorder_rp:
3582        rO_LexVars(j, j_bits, r->block1[i],r->block0[i], prev_ordsgn,
3583                       tmp_ordsgn,v, bits, -1);
3584        break;
3585
3586      case ringorder_dp:
3587        if (r->block0[i]==r->block1[i])
3588        {
3589          rO_LexVars(j, j_bits, r->block0[i],r->block0[i], prev_ordsgn,
3590                     tmp_ordsgn,v, bits, -1);
3591        }
3592        else
3593        {
3594          rO_TDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3595                     tmp_typ[typ_i]);
3596          typ_i++;
3597          rO_LexVars_neg(j, j_bits, r->block1[i],r->block0[i]+1,
3598                         prev_ordsgn,tmp_ordsgn,v,bits, r->block0[i]);
3599        }
3600        break;
3601
3602      case ringorder_Dp:
3603        if (r->block0[i]==r->block1[i])
3604        {
3605          rO_LexVars(j, j_bits, r->block0[i],r->block0[i], prev_ordsgn,
3606                     tmp_ordsgn,v, bits, -1);
3607        }
3608        else
3609        {
3610          rO_TDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3611                     tmp_typ[typ_i]);
3612          typ_i++;
3613          rO_LexVars(j, j_bits, r->block0[i],r->block1[i]-1, prev_ordsgn,
3614                     tmp_ordsgn,v, bits, r->block1[i]);
3615        }
3616        break;
3617
3618      case ringorder_ds:
3619        if (r->block0[i]==r->block1[i])
3620        {
3621          rO_LexVars_neg(j, j_bits,r->block0[i],r->block1[i],prev_ordsgn,
3622                         tmp_ordsgn,v,bits, -1);
3623        }
3624        else
3625        {
3626          rO_TDegree_neg(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3627                         tmp_typ[typ_i]);
3628          typ_i++;
3629          rO_LexVars_neg(j, j_bits, r->block1[i],r->block0[i]+1,
3630                         prev_ordsgn,tmp_ordsgn,v,bits, r->block0[i]);
3631        }
3632        rCheckOrdSgn(r,i);
3633        break;
3634
3635      case ringorder_Ds:
3636        if (r->block0[i]==r->block1[i])
3637        {
3638          rO_LexVars_neg(j, j_bits, r->block0[i],r->block0[i],prev_ordsgn,
3639                         tmp_ordsgn,v, bits, -1);
3640        }
3641        else
3642        {
3643          rO_TDegree_neg(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3644                         tmp_typ[typ_i]);
3645          typ_i++;
3646          rO_LexVars(j, j_bits, r->block0[i],r->block1[i]-1, prev_ordsgn,
3647                     tmp_ordsgn,v, bits, r->block1[i]);
3648        }
3649        rCheckOrdSgn(r,i);
3650        break;
3651
3652      case ringorder_wp:
3653        rO_WDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3654                   tmp_typ[typ_i], r->wvhdl[i]);
3655        typ_i++;
3656        { // check for weights <=0
3657          int jj;
3658          BOOLEAN have_bad_weights=FALSE;
3659          for(jj=r->block1[i]-r->block0[i];jj>=0; jj--)
3660          {
3661            if (r->wvhdl[i][jj]<=0) have_bad_weights=TRUE;
3662          }
3663          if (have_bad_weights)
3664          {
3665             rO_TDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3666                                     tmp_typ[typ_i]);
3667             typ_i++;
3668             rCheckOrdSgn(r,i);
3669          }
3670        }
3671        if (r->block1[i]!=r->block0[i])
3672        {
3673          rO_LexVars_neg(j, j_bits,r->block1[i],r->block0[i]+1, prev_ordsgn,
3674                         tmp_ordsgn, v,bits, r->block0[i]);
3675        }
3676        break;
3677
3678      case ringorder_Wp:
3679        rO_WDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3680                   tmp_typ[typ_i], r->wvhdl[i]);
3681        typ_i++;
3682        { // check for weights <=0
3683          int jj;
3684          BOOLEAN have_bad_weights=FALSE;
3685          for(jj=r->block1[i]-r->block0[i];jj>=0; jj--)
3686          {
3687            if (r->wvhdl[i][jj]<=0) have_bad_weights=TRUE;
3688          }
3689          if (have_bad_weights)
3690          {
3691             rO_TDegree(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3692                                     tmp_typ[typ_i]);
3693             typ_i++;
3694             rCheckOrdSgn(r,i);
3695          }
3696        }
3697        if (r->block1[i]!=r->block0[i])
3698        {
3699          rO_LexVars(j, j_bits,r->block0[i],r->block1[i]-1, prev_ordsgn,
3700                     tmp_ordsgn,v, bits, r->block1[i]);
3701        }
3702        break;
3703
3704      case ringorder_ws:
3705        rO_WDegree_neg(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3706                       tmp_typ[typ_i], r->wvhdl[i]);
3707        typ_i++;
3708        if (r->block1[i]!=r->block0[i])
3709        {
3710          rO_LexVars_neg(j, j_bits,r->block1[i],r->block0[i]+1, prev_ordsgn,
3711                         tmp_ordsgn, v,bits, r->block0[i]);
3712        }
3713        rCheckOrdSgn(r,i);
3714        break;
3715
3716      case ringorder_Ws:
3717        rO_WDegree_neg(j,j_bits,r->block0[i],r->block1[i],tmp_ordsgn,
3718                       tmp_typ[typ_i], r->wvhdl[i]);
3719        typ_i++;
3720        if (r->block1[i]!=r->block0[i])
3721        {
3722          rO_LexVars(j, j_bits,r->block0[i],r->block1[i]-1, prev_ordsgn,
3723                     tmp_ordsgn,v, bits, r->block1[i]);
3724        }
3725        rCheckOrdSgn(r,i);
3726        break;
3727
3728      case ringorder_S:
3729        assume(typ_i == 1); // For LaScala3 only: on the 2nd place ([1])!
3730        // TODO: for K[x]: it is 0...?!
3731        rO_Syzcomp(j, j_bits,prev_ordsgn, tmp_ordsgn,tmp_typ[typ_i]);
3732        need_to_add_comp=TRUE;
3733        typ_i++;
3734        break;
3735
3736      case ringorder_s:
3737        assume(typ_i == 0 && j == 0);
3738        rO_Syz(j, j_bits, prev_ordsgn, tmp_ordsgn, tmp_typ[typ_i]); // set syz-limit?
3739        need_to_add_comp=TRUE;
3740        typ_i++;
3741        break;
3742
3743      case ringorder_IS:
3744      {
3745
3746        assume( r->block0[i] == r->block1[i] );
3747        const int s = r->block0[i];
3748        assume( -2 < s && s < 2);
3749
3750        if(s == 0) // Prefix IS
3751          rO_ISPrefix(j, j_bits, prev_ordsgn, tmp_ordsgn, r->N, v, tmp_typ[typ_i++]); // What about prev_ordsgn?
3752        else // s = +1 or -1 // Note: typ_i might be incrimented here inside!
3753        {
3754          rO_ISSuffix(j, j_bits, prev_ordsgn, tmp_ordsgn, r->N, v, tmp_typ, typ_i, s); // Suffix.
3755          need_to_add_comp=FALSE;
3756        }
3757
3758        break;
3759      }
3760      case ringorder_unspec:
3761      case ringorder_no:
3762      default:
3763        dReportError("undef. ringorder used\n");
3764        break;
3765    }
3766  }
3767
3768  int j0=j; // save j
3769  int j_bits0=j_bits; // save jbits
3770  rO_Align(j,j_bits);
3771  r->CmpL_Size = j;
3772
3773  j_bits=j_bits0; j=j0;
3774
3775  // fill in some empty slots with variables not already covered
3776  // v0 is special, is therefore normally already covered
3777  // now we do have rings without comp...
3778  if((need_to_add_comp) && (v[0]== -1))
3779  {
3780    if (prev_ordsgn==1)
3781    {
3782      rO_Align(j, j_bits);
3783      rO_LexVars(j, j_bits, 0,0, prev_ordsgn,tmp_ordsgn,v,BITS_PER_LONG, -1);
3784    }
3785    else
3786    {
3787      rO_Align(j, j_bits);
3788      rO_LexVars_neg(j, j_bits, 0,0, prev_ordsgn,tmp_ordsgn,v,BITS_PER_LONG, -1);
3789    }
3790  }
3791  // the variables
3792  for(i=1 ; i<=r->N ; i++)
3793  {
3794    if(v[i]==(-1))
3795    {
3796      if (prev_ordsgn==1)
3797      {
3798        rO_LexVars(j, j_bits, i,i, prev_ordsgn,tmp_ordsgn,v,bits, -1);
3799      }
3800      else
3801      {
3802        rO_LexVars_neg(j,j_bits,i,i, prev_ordsgn,tmp_ordsgn,v,bits, -1);
3803      }
3804    }
3805  }
3806
3807  rO_Align(j,j_bits);
3808  // ----------------------------
3809  // finished with constructing the monomial, computing sizes:
3810
3811  r->ExpL_Size=j;
3812  r->PolyBin = omGetSpecBin(POLYSIZE + (r->ExpL_Size)*sizeof(long));
3813  assume(r->PolyBin != NULL);
3814
3815  // ----------------------------
3816  // indices and ordsgn vector for comparison
3817  //
3818  // r->pCompHighIndex already set
3819  r->ordsgn=(long *)omAlloc0(r->ExpL_Size*sizeof(long));
3820
3821  for(j=0;j<r->CmpL_Size;j++)
3822  {
3823    r->ordsgn[j] = tmp_ordsgn[j];
3824  }
3825
3826  omFreeSize((ADDRESS)tmp_ordsgn,(3*(n+r->N)*sizeof(long)));
3827
3828  // ----------------------------
3829  // description of orderings for setm:
3830  //
3831  r->OrdSize=typ_i;
3832  if (typ_i==0) r->typ=NULL;
3833  else
3834  {
3835    r->typ=(sro_ord*)omAlloc(typ_i*sizeof(sro_ord));
3836    memcpy(r->typ,tmp_typ,typ_i*sizeof(sro_ord));
3837  }
3838  omFreeSize((ADDRESS)tmp_typ,(3*(n+r->N)*sizeof(sro_ord)));
3839
3840  // ----------------------------
3841  // indices for (first copy of ) variable entries in exp.e vector (VarOffset):
3842  r->VarOffset=v;
3843
3844  // ----------------------------
3845  // other indicies
3846  r->pCompIndex=(r->VarOffset[0] & 0xffff); //r->VarOffset[0];
3847  i=0; // position
3848  j=0; // index in r->typ
3849  if (i==r->pCompIndex) i++; // IS???
3850  while ((j < r->OrdSize)
3851         && ((r->typ[j].ord_typ==ro_syzcomp) ||
3852             (r->typ[j].ord_typ==ro_syz) || (r->typ[j].ord_typ==ro_isTemp) || (r->typ[j].ord_typ==ro_is) ||
3853             (r->order[r->typ[j].order_index] == ringorder_aa)))
3854  {
3855    i++; j++;
3856  }
3857  // No use of j anymore!!!????
3858
3859  if (i==r->pCompIndex) i++;
3860  r->pOrdIndex=i; // How came it is "i" here???!!!! exp[r->pOrdIndex] is order of a poly... This may be wrong!!! IS
3861
3862  // ----------------------------
3863  rSetDegStuff(r);
3864  rSetOption(r);
3865  // ----------------------------
3866  // r->p_Setm
3867  r->p_Setm = p_GetSetmProc(r);
3868
3869  // ----------------------------
3870  // set VarL_*
3871  rSetVarL(r);
3872
3873  //  ----------------------------
3874  // right-adjust VarOffset
3875  rRightAdjustVarOffset(r);
3876
3877  // ----------------------------
3878  // set NegWeightL*
3879  rSetNegWeight(r);
3880
3881  // ----------------------------
3882  // p_Procs: call AFTER NegWeightL
3883  r->p_Procs = (p_Procs_s*)omAlloc(sizeof(p_Procs_s));
3884  p_ProcsSet(r, r->p_Procs);
3885  return FALSE;
3886}
3887
3888static void rCheckOrdSgn(ring r,int i/*current block*/)
3889{ // set r->OrdSgn
3890  if ( r->OrdSgn==1)
3891  {
3892    int oo=-1;
3893    int jj;
3894    for(jj=i-1;jj>=0;jj--)
3895    {
3896      if(((r->order[jj]==ringorder_a)
3897        ||(r->order[jj]==ringorder_aa)
3898        ||(r->order[jj]==ringorder_a64))
3899      &&(r->block0[jj]<=r->block0[i])
3900      &&(r->block1[jj]>=r->block1[i]))
3901      { oo=1; break;}
3902    }
3903    r->OrdSgn=oo;
3904  }
3905}
3906
3907
3908void rUnComplete(ring r)
3909{
3910  if (r == NULL) return;
3911  if (r->VarOffset != NULL)
3912  {
3913    if (r->OrdSize!=0 && r->typ != NULL)
3914    {
3915      for(int i = 0; i < r->OrdSize; i++)
3916        if( r->typ[i].ord_typ == ro_is) // Search for suffixes! (prefix have the same VarOffset)
3917        {
3918          id_Delete(&r->typ[i].data.is.F, r);
3919          r->typ[i].data.is.F = NULL; // ?
3920
3921          if( r->typ[i].data.is.componentWeights != NULL )
3922          {
3923            delete r->typ[i].data.is.componentWeights;
3924            r->typ[i].data.is.componentWeights = NULL; // ?
3925          }
3926
3927          if( r->typ[i].data.is.pVarOffset != NULL )
3928          {
3929            omFreeSize((ADDRESS)r->typ[i].data.is.pVarOffset, (r->N +1)*sizeof(int));
3930            r->typ[i].data.is.pVarOffset = NULL; // ?
3931          }
3932        }
3933        else if (r->typ[i].ord_typ == ro_syz)
3934        {
3935          if(r->typ[i].data.syz.limit > 0)
3936            omFreeSize(r->typ[i].data.syz.syz_index, ((r->typ[i].data.syz.limit) +1)*sizeof(int));
3937          r->typ[i].data.syz.syz_index = NULL;
3938        }
3939        else if (r->typ[i].ord_typ == ro_syzcomp)
3940        {
3941          assume( r->typ[i].data.syzcomp.ShiftedComponents == NULL );
3942          assume( r->typ[i].data.syzcomp.Components        == NULL );
3943//          WarnS( "rUnComplete : ord_typ == ro_syzcomp was unhandled!!! Possibly memory leak!!!"  );
3944#ifndef NDEBUG
3945//          assume(0);
3946#endif
3947        }
3948
3949      omFreeSize((ADDRESS)r->typ,r->OrdSize*sizeof(sro_ord)); r->typ = NULL;
3950    }
3951
3952    if (r->order != NULL)
3953    {
3954      // delete r->order!!!???
3955    }
3956
3957    if (r->PolyBin != NULL)
3958      omUnGetSpecBin(&(r->PolyBin));
3959
3960    omFreeSize((ADDRESS)r->VarOffset, (r->N +1)*sizeof(int));
3961
3962    if (r->ordsgn != NULL && r->CmpL_Size != 0)
3963      omFreeSize((ADDRESS)r->ordsgn,r->ExpL_Size*sizeof(long));
3964    if (r->p_Procs != NULL)
3965      omFreeSize(r->p_Procs, sizeof(p_Procs_s));
3966    omfreeSize(r->VarL_Offset, r->VarL_Size*sizeof(int));
3967  }
3968  if (r->NegWeightL_Offset!=NULL)
3969  {
3970    omFreeSize(r->NegWeightL_Offset, r->NegWeightL_Size*sizeof(int));
3971    r->NegWeightL_Offset=NULL;
3972  }
3973}
3974
3975// set r->VarL_Size, r->VarL_Offset, r->VarL_LowIndex
3976static void rSetVarL(ring r)
3977{
3978  int  min = INT_MAX, min_j = -1;
3979  int* VarL_Number = (int*) omAlloc0(r->ExpL_Size*sizeof(int));
3980
3981  int i,j;
3982
3983  // count how often a var long is occupied by an exponent
3984  for (i=1; i<=r->N; i++)
3985  {
3986    VarL_Number[r->VarOffset[i] & 0xffffff]++;
3987  }
3988
3989  // determine how many and min
3990  for (i=0, j=0; i<r->ExpL_Size; i++)
3991  {
3992    if (VarL_Number[i] != 0)
3993    {
3994      if (min > VarL_Number[i])
3995      {
3996        min = VarL_Number[i];
3997        min_j = j;
3998      }
3999      j++;
4000    }
4001  }
4002
4003  r->VarL_Size = j; // number of long with exp. entries in
4004                    //  in p->exp
4005  r->VarL_Offset = (int*) omAlloc(r->VarL_Size*sizeof(int));
4006  r->VarL_LowIndex = 0;
4007
4008  // set VarL_Offset
4009  for (i=0, j=0; i<r->ExpL_Size; i++)
4010  {
4011    if (VarL_Number[i] != 0)
4012    {
4013      r->VarL_Offset[j] = i;
4014      if (j > 0 && r->VarL_Offset[j-1] != r->VarL_Offset[j] - 1)
4015        r->VarL_LowIndex = -1;
4016      j++;
4017    }
4018  }
4019  if (r->VarL_LowIndex >= 0)
4020    r->VarL_LowIndex = r->VarL_Offset[0];
4021
4022  r->MinExpPerLong = min;
4023  if (min_j != 0)
4024  {
4025    j = r->VarL_Offset[min_j];
4026    r->VarL_Offset[min_j] = r->VarL_Offset[0];
4027    r->VarL_Offset[0] = j;
4028  }
4029  omFree(VarL_Number);
4030}
4031
4032static void rRightAdjustVarOffset(ring r)
4033{
4034  int* shifts = (int*) omAlloc(r->ExpL_Size*sizeof(int));
4035  int i;
4036  // initialize shifts
4037  for (i=0;i<r->ExpL_Size;i++)
4038    shifts[i] = BIT_SIZEOF_LONG;
4039
4040  // find minimal bit shift in each long exp entry
4041  for (i=1;i<=r->N;i++)
4042  {
4043    if (shifts[r->VarOffset[i] & 0xffffff] > r->VarOffset[i] >> 24)
4044      shifts[r->VarOffset[i] & 0xffffff] = r->VarOffset[i] >> 24;
4045  }
4046  // reset r->VarOffset: set the minimal shift to 0
4047  for (i=1;i<=r->N;i++)
4048  {
4049    if (shifts[r->VarOffset[i] & 0xffffff] != 0)
4050      r->VarOffset[i]
4051        = (r->VarOffset[i] & 0xffffff) |
4052        (((r->VarOffset[i] >> 24) - shifts[r->VarOffset[i] & 0xffffff]) << 24);
4053  }
4054  omFree(shifts);
4055}
4056
4057// get r->divmask depending on bits per exponent
4058static unsigned long rGetDivMask(int bits)
4059{
4060  unsigned long divmask = 1;
4061  int i = bits;
4062
4063  while (i < BIT_SIZEOF_LONG)
4064  {
4065    divmask |= (((unsigned long) 1) << (unsigned long) i);
4066    i += bits;
4067  }
4068  return divmask;
4069}
4070
4071#ifdef RDEBUG
4072void rDebugPrint(ring r)
4073{
4074  if (r==NULL)
4075  {
4076    PrintS("NULL ?\n");
4077    return;
4078  }
4079  // corresponds to ro_typ from ring.h:
4080  const char *TYP[]={"ro_dp","ro_wp","ro_wp64","ro_wp_neg","ro_cp",
4081                     "ro_syzcomp", "ro_syz", "ro_isTemp", "ro_is", "ro_none"};
4082  int i,j;
4083
4084  Print("ExpL_Size:%d ",r->ExpL_Size);
4085  Print("CmpL_Size:%d ",r->CmpL_Size);
4086  Print("VarL_Size:%d\n",r->VarL_Size);
4087  Print("bitmask=0x%lx (expbound=%ld) \n",r->bitmask, r->bitmask);
4088  Print("BitsPerExp=%d ExpPerLong=%d MinExpPerLong=%d at L[%d]\n", r->BitsPerExp, r->ExpPerLong, r->MinExpPerLong, r->VarL_Offset[0]);
4089  PrintS("varoffset:\n");
4090  if (r->VarOffset==NULL) PrintS(" NULL\n");
4091  else
4092    for(j=0;j<=r->N;j++)
4093      Print("  v%d at e-pos %d, bit %d\n",
4094            j,r->VarOffset[j] & 0xffffff, r->VarOffset[j] >>24);
4095  Print("divmask=%lx\n", r->divmask);
4096  PrintS("ordsgn:\n");
4097  for(j=0;j<r->CmpL_Size;j++)
4098    Print("  ordsgn %ld at pos %d\n",r->ordsgn[j],j);
4099  Print("OrdSgn:%d\n",r->OrdSgn);
4100  PrintS("ordrec:\n");
4101  for(j=0;j<r->OrdSize;j++)
4102  {
4103    Print("  typ %s", TYP[r->typ[j].ord_typ]);
4104
4105
4106    if (r->typ[j].ord_typ==ro_syz)
4107    {
4108      const short place = r->typ[j].data.syz.place;
4109      const int limit = r->typ[j].data.syz.limit;
4110      const int curr_index = r->typ[j].data.syz.curr_index;
4111      const int* syz_index = r->typ[j].data.syz.syz_index;
4112
4113      Print("  limit %d (place: %d, curr_index: %d), syz_index: ", limit, place, curr_index);
4114
4115      if( syz_index == NULL )
4116        PrintS("(NULL)");
4117      else
4118      {
4119        Print("{");
4120        for( i=0; i <= limit; i++ )
4121          Print("%d ", syz_index[i]);
4122        Print("}");
4123      }
4124
4125    }
4126    else if (r->typ[j].ord_typ==ro_isTemp)
4127    {
4128      Print("  start (level) %d, suffixpos: %d, VO: ",r->typ[j].data.isTemp.start, r->typ[j].data.isTemp.suffixpos);
4129
4130    }
4131    else if (r->typ[j].ord_typ==ro_is)
4132    {
4133      Print("  start %d, end: %d: ",r->typ[j].data.is.start, r->typ[j].data.is.end);
4134
4135//      for( int k = 0; k <= r->N; k++) if (r->typ[j].data.is.pVarOffset[k] != -1) Print("[%2d]: %04x; ", k, r->typ[j].data.is.pVarOffset[k]);
4136
4137      Print("  limit %d\n",r->typ[j].data.is.limit);
4138      #ifndef NDEBUG
4139      //PrintS("  F: ");idShow(r->typ[j].data.is.F, r, r, 1);
4140      #endif
4141
4142      PrintS("weights: ");
4143
4144      if( r->typ[j].data.is.componentWeights == NULL )
4145        PrintS("NULL == [0,...0]\n");
4146      else
4147      {
4148        (r->typ[j].data.is.componentWeights)->show(); PrintLn();
4149      }
4150    }
4151    else
4152    {
4153      Print("  place %d",r->typ[j].data.dp.place);
4154
4155      if (r->typ[j].ord_typ!=ro_syzcomp  && r->typ[j].ord_typ!=ro_syz)
4156      {
4157        Print("  start %d",r->typ[j].data.dp.start);
4158        Print("  end %d",r->typ[j].data.dp.end);
4159        if ((r->typ[j].ord_typ==ro_wp)
4160        || (r->typ[j].ord_typ==ro_wp_neg))
4161        {
4162          PrintS(" w:");
4163          for(int l=r->typ[j].data.wp.start;l<=r->typ[j].data.wp.end;l++)
4164            Print(" %d",r->typ[j].data.wp.weights[l-r->typ[j].data.wp.start]);
4165        }
4166        else if (r->typ[j].ord_typ==ro_wp64)
4167        {
4168          PrintS(" w64:");
4169          int l;
4170          for(l=r->typ[j].data.wp64.start;l<=r->typ[j].data.wp64.end;l++)
4171            Print(" %ld",(long)(((int64*)r->typ[j].data.wp64.weights64)+l-r->typ[j].data.wp64.start));
4172          }
4173        }
4174    }
4175    PrintLn();
4176  }
4177  Print("pOrdIndex:%d pCompIndex:%d\n", r->pOrdIndex, r->pCompIndex);
4178  Print("OrdSize:%d\n",r->OrdSize);
4179  PrintS("--------------------\n");
4180  for(j=0;j<r->ExpL_Size;j++)
4181  {
4182    Print("L[%d]: ",j);
4183    if (j< r->CmpL_Size)
4184      Print("ordsgn %ld ", r->ordsgn[j]);
4185    else
4186      PrintS("no comp ");
4187    i=1;
4188    for(;i<=r->N;i++)
4189    {
4190      if( (r->VarOffset[i] & 0xffffff) == j )
4191      {  Print("v%d at e[%d], bit %d; ", i,r->VarOffset[i] & 0xffffff,
4192                                         r->VarOffset[i] >>24 ); }
4193    }
4194    if( r->pCompIndex==j ) PrintS("v0; ");
4195    for(i=0;i<r->OrdSize;i++)
4196    {
4197      if (r->typ[i].data.dp.place == j)
4198      {
4199        Print("ordrec:%s (start:%d, end:%d) ",TYP[r->typ[i].ord_typ],
4200          r->typ[i].data.dp.start, r->typ[i].data.dp.end);
4201      }
4202    }
4203
4204    if (j==r->pOrdIndex)
4205      PrintS("pOrdIndex\n");
4206    else
4207      PrintLn();
4208  }
4209  Print("LexOrder:%d, MixedOrder:%d\n",r->LexOrder, r->MixedOrder);
4210
4211  // p_Procs stuff
4212  p_Procs_s proc_names;
4213  const char* field;
4214  const char* length;
4215  const char* ord;
4216  p_Debug_GetProcNames(r, &proc_names); // changes p_Procs!!!
4217  p_Debug_GetSpecNames(r, field, length, ord);
4218
4219  Print("p_Spec  : %s, %s, %s\n", field, length, ord);
4220  PrintS("p_Procs :\n");
4221  for (i=0; i<(int) (sizeof(p_Procs_s)/sizeof(void*)); i++)
4222  {
4223    Print(" %s,\n", ((char**) &proc_names)[i]);
4224  }
4225
4226  {
4227#define pFDeg_CASE(A) if(r->pFDeg == A) PrintS( "" #A "" )
4228    Print("\npFDeg   : ");
4229   
4230    pFDeg_CASE(p_Totaldegree); else
4231      pFDeg_CASE(p_WFirstTotalDegree); else
4232      pFDeg_CASE(p_WTotaldegree); else
4233      pFDeg_CASE(p_Deg); else
4234      Print("(%p)", r->pFDeg); // default case
4235   
4236    PrintS("\n");
4237#undef pFDeg_CASE
4238  }
4239   
4240}
4241
4242void p_DebugPrint(poly p, const ring r)
4243{
4244  int i,j;
4245  p_Write(p,r);
4246  j=2;
4247  while(p!=NULL)
4248  {
4249    Print("\nexp[0..%d]\n",r->ExpL_Size-1);
4250    for(i=0;i<r->ExpL_Size;i++)
4251      Print("%ld ",p->exp[i]);
4252    PrintLn();
4253    Print("v0:%ld ",p_GetComp(p, r));
4254    for(i=1;i<=r->N;i++) Print(" v%d:%ld",i,p_GetExp(p,i, r));
4255    PrintLn();
4256    pIter(p);
4257    j--;
4258    if (j==0) { PrintS("...\n"); break; }
4259  }
4260}
4261
4262#endif // RDEBUG
4263
4264/// debug-print monomial poly/vector p, assuming that it lives in the ring R
4265static inline void m_DebugPrint(const poly p, const ring R)
4266{
4267  Print("\nexp[0..%d]\n", R->ExpL_Size - 1);
4268  for(int i = 0; i < R->ExpL_Size; i++)
4269    Print("%09lx ", p->exp[i]);
4270  PrintLn();
4271  Print("v0:%9ld ", p_GetComp(p, R));
4272  for(int i = 1; i <= R->N; i++) Print(" v%d:%5ld",i, p_GetExp(p, i, R));
4273  PrintLn();
4274}
4275
4276
4277#ifndef NDEBUG
4278/// debug-print at most nTerms (2 by default) terms from poly/vector p,
4279/// assuming that lt(p) lives in lmRing and tail(p) lives in tailRing.
4280void p_DebugPrint(const poly p, const ring lmRing, const ring tailRing, const int nTerms)
4281{
4282  assume( nTerms >= 0 );
4283  if( p != NULL )
4284  {
4285    assume( p != NULL );
4286
4287    p_Write(p, lmRing, tailRing);
4288
4289    if( (p != NULL) && (nTerms > 0) )
4290    {
4291      assume( p != NULL );
4292      assume( nTerms > 0 );
4293
4294      // debug pring leading term
4295      m_DebugPrint(p, lmRing);
4296
4297      poly q = pNext(p); // q = tail(p)
4298
4299      // debug pring tail (at most nTerms-1 terms from it)
4300      for(int j = nTerms - 1; (q !=NULL) && (j > 0); pIter(q), --j)
4301        m_DebugPrint(q, tailRing);
4302
4303      if (q != NULL)
4304        PrintS("...\n");
4305    }
4306  }
4307  else
4308    PrintS("0\n");
4309}
4310#endif
4311
4312
4313//    F = system("ISUpdateComponents", F, V, MIN );
4314//    // replace gen(i) -> gen(MIN + V[i-MIN]) for all i > MIN in all terms from F!
4315void pISUpdateComponents(ideal F, const intvec *const V, const int MIN, const ring r )
4316{
4317  assume( V != NULL );
4318  assume( MIN >= 0 );
4319
4320  if( F == NULL )
4321    return;
4322
4323  for( int j = (F->ncols*F->nrows) - 1; j >= 0; j-- )
4324  {
4325#ifdef PDEBUG
4326    Print("F[%d]:", j);
4327    p_DebugPrint(F->m[j], r, r, 0);
4328#endif
4329
4330    for( poly p = F->m[j]; p != NULL; pIter(p) )
4331    {
4332      int c = p_GetComp(p, r);
4333
4334      if( c > MIN )
4335      {
4336#ifdef PDEBUG
4337        Print("gen[%d] -> gen(%d)\n", c, MIN + (*V)[ c - MIN - 1 ]);
4338#endif
4339
4340        p_SetComp( p, MIN + (*V)[ c - MIN - 1 ], r );
4341      }
4342    }
4343#ifdef PDEBUG
4344    Print("new F[%d]:", j);
4345    p_Test(F->m[j], r);
4346    p_DebugPrint(F->m[j], r, r, 0);
4347#endif
4348  }
4349
4350}
4351
4352
4353
4354
4355/*2
4356* asssume that rComplete was called with r
4357* assume that the first block ist ringorder_S
4358* change the block to reflect the sequence given by appending v
4359*/
4360
4361#ifdef PDEBUG
4362void rDBChangeSComps(int* currComponents,
4363                     long* currShiftedComponents,
4364                     int length,
4365                     ring r)
4366{
4367  assume(r->typ[1].ord_typ == ro_syzcomp);
4368     
4369  r->typ[1].data.syzcomp.length = length;
4370  rNChangeSComps( currComponents, currShiftedComponents, r);
4371}
4372void rDBGetSComps(int** currComponents,
4373                 long** currShiftedComponents,
4374                 int *length,
4375                 ring r)
4376{
4377  assume(r->typ[1].ord_typ == ro_syzcomp);
4378 
4379  *length = r->typ[1].data.syzcomp.length;
4380  rNGetSComps( currComponents, currShiftedComponents, r);
4381}
4382#endif
4383
4384void rNChangeSComps(int* currComponents, long* currShiftedComponents, ring r)
4385{
4386  assume(r->typ[1].ord_typ == ro_syzcomp);
4387
4388  r->typ[1].data.syzcomp.ShiftedComponents = currShiftedComponents;
4389  r->typ[1].data.syzcomp.Components = currComponents;
4390}
4391
4392void rNGetSComps(int** currComponents, long** currShiftedComponents, ring r)
4393{
4394  assume(r->typ[1].ord_typ == ro_syzcomp);
4395
4396  *currShiftedComponents = r->typ[1].data.syzcomp.ShiftedComponents;
4397  *currComponents =   r->typ[1].data.syzcomp.Components;
4398}
4399
4400/////////////////////////////////////////////////////////////////////////////
4401//
4402// The following routines all take as input a ring r, and return R
4403// where R has a certain property. R might be equal r in which case r
4404// had already this property
4405//
4406ring rAssure_SyzComp(const ring r, BOOLEAN complete)
4407{
4408  if ( (r->order[0] == ringorder_s) ) return r;
4409
4410  if ( (r->order[0] == ringorder_IS) )
4411  {
4412#ifndef NDEBUG
4413    WarnS("rAssure_SyzComp: input ring has an IS-ordering!");
4414#endif
4415//    return r;
4416  }
4417  ring res=rCopy0(r, FALSE, FALSE);
4418  int i=rBlocks(r);
4419  int j;
4420
4421  res->order=(int *)omAlloc((i+1)*sizeof(int));
4422  res->block0=(int *)omAlloc0((i+1)*sizeof(int));
4423  res->block1=(int *)omAlloc0((i+1)*sizeof(int));
4424  int ** wvhdl =(int **)omAlloc0((i+1)*sizeof(int**));
4425  for(j=i;j>0;j--)
4426  {
4427    res->order[j]=r->order[j-1];
4428    res->block0[j]=r->block0[j-1];
4429    res->block1[j]=r->block1[j-1];
4430    if (r->wvhdl[j-1] != NULL)
4431    {
4432      wvhdl[j] = (int*) omMemDup(r->wvhdl[j-1]);
4433    }
4434  }
4435  res->order[0]=ringorder_s;
4436
4437  res->wvhdl = wvhdl;
4438
4439  if (complete)
4440  {
4441    rComplete(res, 1);
4442
4443#ifdef HAVE_PLURAL
4444    if (rIsPluralRing(r))
4445    {
4446      if ( nc_rComplete(r, res, false) ) // no qideal!
4447      {
4448#ifndef NDEBUG
4449        WarnS("error in nc_rComplete");      // cleanup?//      rDelete(res);//      return r;      // just go on..
4450#endif
4451      }
4452    }
4453    assume(rIsPluralRing(r) == rIsPluralRing(res));
4454#endif
4455
4456
4457#ifdef HAVE_PLURAL
4458    ring old_ring = r;
4459#endif
4460
4461    if (r->qideal!=NULL)
4462    {
4463      res->qideal= idrCopyR_NoSort(r->qideal, r, res);
4464
4465      assume(id_RankFreeModule(res->qideal, res) == 0);
4466
4467#ifdef HAVE_PLURAL
4468      if( rIsPluralRing(res) )
4469        if( nc_SetupQuotient(res, r, true) )
4470        {
4471//          WarnS("error in nc_SetupQuotient"); // cleanup?      rDelete(res);       return r;  // just go on...?
4472        }
4473
4474#endif
4475      assume(id_RankFreeModule(res->qideal, res) == 0);
4476    }
4477
4478#ifdef HAVE_PLURAL
4479    assume((res->qideal==NULL) == (old_ring->qideal==NULL));
4480    assume(rIsPluralRing(res) == rIsPluralRing(old_ring));
4481    assume(rIsSCA(res) == rIsSCA(old_ring));
4482    assume(ncRingType(res) == ncRingType(old_ring));
4483#endif
4484  }
4485
4486  return res;
4487}
4488
4489ring rAssure_TDeg(ring r, int start_var, int end_var, int &pos)
4490{
4491  int i;
4492  if (r->typ!=NULL)
4493  {
4494    for(i=r->OrdSize-1;i>=0;i--)
4495    {
4496      if ((r->typ[i].ord_typ==ro_dp)
4497      && (r->typ[i].data.dp.start==start_var)
4498      && (r->typ[i].data.dp.end==end_var))
4499      {
4500        pos=r->typ[i].data.dp.place;
4501        //printf("no change, pos=%d\n",pos);
4502        return r;
4503      }
4504    }
4505  }
4506
4507#ifdef HAVE_PLURAL
4508  nc_struct* save=r->GetNC();
4509  r->GetNC()=NULL;
4510#endif
4511  ring res=rCopy(r);
4512
4513  i=rBlocks(r);
4514  int j;
4515
4516  res->ExpL_Size=r->ExpL_Size+1; // one word more in each monom
4517  res->PolyBin=omGetSpecBin(POLYSIZE + (res->ExpL_Size)*sizeof(long));
4518  omFree((ADDRESS)res->ordsgn);
4519  res->ordsgn=(long *)omAlloc0(res->ExpL_Size*sizeof(long));
4520  for(j=0;j<r->CmpL_Size;j++)
4521  {
4522    res->ordsgn[j] = r->ordsgn[j];
4523  }
4524  res->OrdSize=r->OrdSize+1;   // one block more for pSetm
4525  if (r->typ!=NULL)
4526    omFree((ADDRESS)res->typ);
4527  res->typ=(sro_ord*)omAlloc0(res->OrdSize*sizeof(sro_ord));
4528  if (r->typ!=NULL)
4529    memcpy(res->typ,r->typ,r->OrdSize*sizeof(sro_ord));
4530  // the additional block for pSetm: total degree at the last word
4531  // but not included in the compare part
4532  res->typ[res->OrdSize-1].ord_typ=ro_dp;
4533  res->typ[res->OrdSize-1].data.dp.start=start_var;
4534  res->typ[res->OrdSize-1].data.dp.end=end_var;
4535  res->typ[res->OrdSize-1].data.dp.place=res->ExpL_Size-1;
4536  pos=res->ExpL_Size-1;
4537  //if ((start_var==1) && (end_var==res->N)) res->pOrdIndex=pos;
4538  extern void p_Setm_General(poly p, ring r);
4539  res->p_Setm=p_Setm_General;
4540  // ----------------------------
4541  omFree((ADDRESS)res->p_Procs);
4542  res->p_Procs = (p_Procs_s*)omAlloc(sizeof(p_Procs_s));
4543
4544  p_ProcsSet(res, res->p_Procs);
4545  if (res->qideal!=NULL) id_Delete(&res->qideal,res);
4546#ifdef HAVE_PLURAL
4547  r->GetNC()=save;
4548  if (rIsPluralRing(r))
4549  {
4550    if ( nc_rComplete(r, res, false) ) // no qideal!
4551    {
4552#ifndef NDEBUG
4553      WarnS("error in nc_rComplete");
4554#endif
4555      // just go on..
4556    }
4557  }
4558#endif
4559  if (r->qideal!=NULL)
4560  {
4561     res->qideal=idrCopyR_NoSort(r->qideal,r, res);
4562#ifdef HAVE_PLURAL
4563     if (rIsPluralRing(res))
4564     {
4565//       nc_SetupQuotient(res, currRing);
4566       nc_SetupQuotient(res, r); // ?
4567     }
4568     assume((res->qideal==NULL) == (r->qideal==NULL));
4569#endif
4570  }
4571
4572#ifdef HAVE_PLURAL
4573  assume(rIsPluralRing(res) == rIsPluralRing(r));
4574  assume(rIsSCA(res) == rIsSCA(r));
4575  assume(ncRingType(res) == ncRingType(r));
4576#endif
4577
4578  return res;
4579}
4580
4581ring rAssure_HasComp(const ring r)
4582{
4583  int last_block;
4584  int i=0;
4585  do
4586  {
4587     if (r->order[i] == ringorder_c ||
4588        r->order[i] == ringorder_C) return r;
4589     if (r->order[i] == 0)
4590        break;
4591     i++;
4592  } while (1);
4593  //WarnS("re-creating ring with comps");
4594  last_block=i-1;
4595
4596  ring new_r = rCopy0(r, FALSE, FALSE);
4597  i+=2;
4598  new_r->wvhdl=(int **)omAlloc0(i * sizeof(int *));
4599  new_r->order   = (int *) omAlloc0(i * sizeof(int));
4600  new_r->block0   = (int *) omAlloc0(i * sizeof(int));
4601  new_r->block1   = (int *) omAlloc0(i * sizeof(int));
4602  memcpy(new_r->order,r->order,(i-1) * sizeof(int));
4603  memcpy(new_r->block0,r->block0,(i-1) * sizeof(int));
4604  memcpy(new_r->block1,r->block1,(i-1) * sizeof(int));
4605  for (int j=0; j<=last_block; j++)
4606  {
4607    if (r->wvhdl[j]!=NULL)
4608    {
4609      new_r->wvhdl[j] = (int*) omMemDup(r->wvhdl[j]);
4610    }
4611  }
4612  last_block++;
4613  new_r->order[last_block]=ringorder_C;
4614  //new_r->block0[last_block]=0;
4615  //new_r->block1[last_block]=0;
4616  //new_r->wvhdl[last_block]=NULL;
4617
4618  rComplete(new_r, 1);
4619
4620#ifdef HAVE_PLURAL
4621  if (rIsPluralRing(r))
4622  {
4623    if ( nc_rComplete(r, new_r, false) ) // no qideal!
4624    {
4625#ifndef NDEBUG
4626      WarnS("error in nc_rComplete");      // cleanup?//      rDelete(res);//      return r;      // just go on..
4627#endif
4628    }
4629  }
4630  assume(rIsPluralRing(r) == rIsPluralRing(new_r));
4631#endif
4632
4633  return new_r;
4634}
4635
4636static ring rAssure_CompLastBlock(ring r, BOOLEAN complete = TRUE)
4637{
4638  int last_block = rBlocks(r) - 2;
4639  if (r->order[last_block] != ringorder_c &&
4640      r->order[last_block] != ringorder_C)
4641  {
4642    int c_pos = 0;
4643    int i;
4644
4645    for (i=0; i< last_block; i++)
4646    {
4647      if (r->order[i] == ringorder_c || r->order[i] == ringorder_C)
4648      {
4649        c_pos = i;
4650        break;
4651      }
4652    }
4653    if (c_pos != -1)
4654    {
4655      ring new_r = rCopy0(r, FALSE, TRUE);
4656      for (i=c_pos+1; i<=last_block; i++)
4657      {
4658        new_r->order[i-1] = new_r->order[i];
4659        new_r->block0[i-1] = new_r->block0[i];
4660        new_r->block1[i-1] = new_r->block1[i];
4661        new_r->wvhdl[i-1] = new_r->wvhdl[i];
4662      }
4663      new_r->order[last_block] = r->order[c_pos];
4664      new_r->block0[last_block] = r->block0[c_pos];
4665      new_r->block1[last_block] = r->block1[c_pos];
4666      new_r->wvhdl[last_block] = r->wvhdl[c_pos];
4667      if (complete)
4668      {
4669        rComplete(new_r, 1);
4670
4671#ifdef HAVE_PLURAL
4672        if (rIsPluralRing(r))
4673        {
4674          if ( nc_rComplete(r, new_r, false) ) // no qideal!
4675          {
4676#ifndef NDEBUG
4677            WarnS("error in nc_rComplete");   // cleanup?//      rDelete(res);//      return r;      // just go on..
4678#endif
4679          }
4680        }
4681        assume(rIsPluralRing(r) == rIsPluralRing(new_r));
4682#endif
4683      }
4684      return new_r;
4685    }
4686  }
4687  return r;
4688}
4689
4690// Moves _c or _C ordering to the last place AND adds _s on the 1st place
4691ring rAssure_SyzComp_CompLastBlock(const ring r)
4692{
4693  ring new_r_1 = rAssure_CompLastBlock(r, FALSE); // due to this FALSE - no completion!
4694  ring new_r = rAssure_SyzComp(new_r_1, FALSE); // new_r_1 is used only here!!!
4695
4696  if (new_r != r)
4697  {
4698    ring old_r = r;
4699    if (new_r_1 != new_r && new_r_1 != old_r) rDelete(new_r_1);
4700    rComplete(new_r, 1);
4701#ifdef HAVE_PLURAL
4702    if (rIsPluralRing(old_r))
4703    {
4704      if ( nc_rComplete(old_r, new_r, false) ) // no qideal!
4705      {
4706#ifndef NDEBUG
4707        WarnS("error in nc_rComplete"); // cleanup?      rDelete(res);       return r;  // just go on...?
4708#endif
4709        }
4710    }
4711    assume(rIsPluralRing(new_r) == rIsPluralRing(old_r));
4712#endif
4713///?    rChangeCurrRing(new_r);
4714    if (old_r->qideal != NULL)
4715    {
4716      new_r->qideal = idrCopyR(old_r->qideal, old_r, new_r);
4717      //currQuotient = new_r->qideal;
4718
4719#ifdef HAVE_PLURAL
4720      if( rIsPluralRing(old_r) )
4721        if( nc_SetupQuotient(new_r, old_r, true) )
4722        {
4723#ifndef NDEBUG
4724          WarnS("error in nc_SetupQuotient"); // cleanup?      rDelete(res);       return r;  // just go on...?
4725#endif
4726        }
4727#endif
4728    }
4729
4730#ifdef HAVE_PLURAL
4731    assume((new_r->qideal==NULL) == (old_r->qideal==NULL));
4732    assume(rIsPluralRing(new_r) == rIsPluralRing(old_r));
4733    assume(rIsSCA(new_r) == rIsSCA(old_r));
4734    assume(ncRingType(new_r) == ncRingType(old_r));
4735#endif
4736
4737    rTest(new_r);
4738    rTest(old_r);
4739  }
4740  return new_r;
4741}
4742
4743// use this for global orderings consisting of two blocks
4744static ring rAssure_Global(rRingOrder_t b1, rRingOrder_t b2, const ring r)
4745{
4746  int r_blocks = rBlocks(r);
4747
4748  assume(b1 == ringorder_c || b1 == ringorder_C ||
4749         b2 == ringorder_c || b2 == ringorder_C ||
4750         b2 == ringorder_S);
4751  if ((r_blocks == 3) &&
4752      (r->order[0] == b1) &&
4753      (r->order[1] == b2) &&
4754      (r->order[2] == 0))
4755    return r;
4756  ring res = rCopy0(r, TRUE, FALSE);
4757  res->order = (int*)omAlloc0(3*sizeof(int));
4758  res->block0 = (int*)omAlloc0(3*sizeof(int));
4759  res->block1 = (int*)omAlloc0(3*sizeof(int));
4760  res->wvhdl = (int**)omAlloc0(3*sizeof(int*));
4761  res->order[0] = b1;
4762  res->order[1] = b2;
4763  if (b1 == ringorder_c || b1 == ringorder_C)
4764  {
4765    res->block0[1] = 1;
4766    res->block1[1] = r->N;
4767  }
4768  else
4769  {
4770    res->block0[0] = 1;
4771    res->block1[0] = r->N;
4772  }
4773  // HANNES: This sould be set in rComplete
4774  res->OrdSgn = 1;
4775  rComplete(res, 1);
4776#ifdef HAVE_PLURAL
4777  if (rIsPluralRing(r))
4778  {
4779    if ( nc_rComplete(r, res, false) ) // no qideal!
4780    {
4781#ifndef NDEBUG
4782      WarnS("error in nc_rComplete");
4783#endif
4784    }
4785  }
4786#endif
4787//  rChangeCurrRing(res);
4788  return res;
4789}
4790
4791ring rAssure_InducedSchreyerOrdering(const ring r, BOOLEAN complete = TRUE, int sgn = 1)
4792{ // TODO: ???? Add leading Syz-comp ordering here...????
4793
4794#if MYTEST
4795    Print("rAssure_InducedSchreyerOrdering(r, complete = %d, sgn = %d): r: \n", complete, sgn);
4796    rWrite(r);
4797#ifdef RDEBUG
4798    rDebugPrint(r);
4799#endif
4800    PrintLn();
4801#endif
4802  assume((sgn == 1) || (sgn == -1));
4803
4804  ring res=rCopy0(r, FALSE, FALSE); // No qideal & ordering copy.
4805
4806  int n = rBlocks(r); // Including trailing zero!
4807
4808  // Create 2 more blocks for prefix/suffix:
4809  res->order=(int *)omAlloc0((n+2)*sizeof(int)); // 0  ..  n+1
4810  res->block0=(int *)omAlloc0((n+2)*sizeof(int));
4811  res->block1=(int *)omAlloc0((n+2)*sizeof(int));
4812  int ** wvhdl =(int **)omAlloc0((n+2)*sizeof(int**));
4813
4814  // Encapsulate all existing blocks between induced Schreyer ordering markers: prefix and suffix!
4815  // Note that prefix and suffix have the same ringorder marker and only differ in block[] parameters!
4816
4817  // new 1st block
4818  int j = 0;
4819  res->order[j] = ringorder_IS; // Prefix
4820  res->block0[j] = res->block1[j] = 0;
4821  // wvhdl[j] = NULL;
4822  j++;
4823
4824  for(int i = 0; (i <= n) && (r->order[i] != 0); i++, j++) // i = [0 .. n-1] <- non-zero old blocks
4825  {
4826    res->order [j] = r->order [i];
4827    res->block0[j] = r->block0[i];
4828    res->block1[j] = r->block1[i];
4829
4830    if (r->wvhdl[i] != NULL)
4831    {
4832      wvhdl[j] = (int*) omMemDup(r->wvhdl[i]);
4833    } // else wvhdl[j] = NULL;
4834  }
4835
4836  // new last block
4837  res->order [j] = ringorder_IS; // Suffix
4838  res->block0[j] = res->block1[j] = sgn; // Sign of v[o]: 1 for C, -1 for c
4839  // wvhdl[j] = NULL;
4840  j++;
4841
4842  // res->order [j] = 0; // The End!
4843  res->wvhdl = wvhdl;
4844
4845  // j == the last zero block now!
4846  assume(j == (n+1));
4847  assume(res->order[0]==ringorder_IS);
4848  assume(res->order[j-1]==ringorder_IS);
4849  assume(res->order[j]==0);
4850
4851
4852  if (complete)
4853  {
4854    rComplete(res, 1);
4855
4856#ifdef HAVE_PLURAL
4857    if (rIsPluralRing(r))
4858    {
4859      if ( nc_rComplete(r, res, false) ) // no qideal!
4860      {
4861#ifndef NDEBUG
4862        WarnS("error in nc_rComplete");      // cleanup?//      rDelete(res);//      return r;      // just go on..
4863#endif
4864      }
4865    }
4866    assume(rIsPluralRing(r) == rIsPluralRing(res));
4867#endif
4868
4869
4870#ifdef HAVE_PLURAL
4871    ring old_ring = r;
4872#endif
4873
4874    if (r->qideal!=NULL)
4875    {
4876      res->qideal= idrCopyR_NoSort(r->qideal, r, res);
4877
4878      assume(id_RankFreeModule(res->qideal, res) == 0);
4879
4880#ifdef HAVE_PLURAL
4881      if( rIsPluralRing(res) )
4882        if( nc_SetupQuotient(res, r, true) )
4883        {
4884//          WarnS("error in nc_SetupQuotient"); // cleanup?      rDelete(res);       return r;  // just go on...?
4885        }
4886
4887#endif
4888      assume(id_RankFreeModule(res->qideal, res) == 0);
4889    }
4890
4891#ifdef HAVE_PLURAL
4892    assume((res->qideal==NULL) == (old_ring->qideal==NULL));
4893    assume(rIsPluralRing(res) == rIsPluralRing(old_ring));
4894    assume(rIsSCA(res) == rIsSCA(old_ring));
4895    assume(ncRingType(res) == ncRingType(old_ring));
4896#endif
4897  }
4898
4899  return res;
4900}
4901
4902ring rAssure_dp_S(const ring r)
4903{
4904  return rAssure_Global(ringorder_dp, ringorder_S,r);
4905}
4906
4907ring rAssure_dp_C(const ring r)
4908{
4909  return rAssure_Global(ringorder_dp, ringorder_C,r);
4910}
4911
4912ring rAssure_C_dp(const ring r)
4913{
4914  return rAssure_Global(ringorder_C, ringorder_dp,r);
4915}
4916
4917
4918
4919/// Finds p^th IS ordering, and returns its position in r->typ[]
4920/// returns -1 if something went wrong!
4921int rGetISPos(const int p, const ring r)
4922{
4923  // Put the reference set F into the ring -ordering -recor
4924#if MYTEST
4925  Print("rIsIS(p: %d)\nF:", p);
4926  PrintLn();
4927#endif
4928
4929  if (r->typ==NULL)
4930  {
4931    dReportError("'rIsIS:' Error: wrong ring! (typ == NULL)");
4932    return -1;
4933  }
4934
4935  int j = p; // Which IS record to use...
4936  for( int pos = 0; pos < r->OrdSize; pos++ )
4937    if( r->typ[pos].ord_typ == ro_is)
4938      if( j-- == 0 )
4939      {
4940        return pos;
4941      }
4942
4943  return -1;
4944}
4945
4946
4947
4948
4949
4950
4951/// Changes r by setting induced ordering parameters: limit and reference leading terms
4952/// F belong to r, we will DO a copy! (same to componentWeights)
4953/// We will use it AS IS!
4954/// returns true is everything was allright!
4955BOOLEAN rSetISReference(const ring r, const ideal F, const int i, const int p, const intvec * componentWeights)
4956{
4957  // Put the reference set F into the ring -ordering -recor
4958
4959  // TEST THAT THERE ARE DEGs!!!
4960  // assume( componentWeights == NULL  ); // ???
4961  if( componentWeights != NULL )
4962  {
4963//    assure that the ring r has degrees!!!
4964//    Add weights to degrees of F[i]
4965  }
4966
4967  if (r->typ==NULL)
4968  {
4969    dReportError("Error: WRONG USE of rSetISReference: wrong ring! (typ == NULL)");
4970    return FALSE;
4971  }
4972
4973
4974  int pos = rGetISPos(p, r);
4975
4976  if( pos == -1 )
4977  {
4978    dReportError("Error: WRONG USE of rSetISReference: specified ordering block was not found!!!" );
4979    return FALSE;
4980  }
4981
4982#if MYTEST
4983  if( i != r->typ[pos].data.is.limit )
4984    Print("Changing record on pos: %d\nOld limit: %d --->> New Limit: %d\n", pos, r->typ[pos].data.is.limit, i);
4985#endif
4986
4987  const ideal FF = id_Copy(F, r); // idrHeadR(F, r, r);
4988
4989
4990  if( r->typ[pos].data.is.F != NULL)
4991  {
4992#if MYTEST
4993    PrintS("Deleting old reference set F... \n");        // idShow(r->typ[pos].data.is.F, r);         PrintLn();
4994#endif
4995    id_Delete(&r->typ[pos].data.is.F, r);
4996    r->typ[pos].data.is.F = NULL;
4997  }
4998
4999  assume(r->typ[pos].data.is.F == NULL);
5000
5001  r->typ[pos].data.is.F = FF; // F is owened by ring now! TODO: delete at the end!
5002
5003  if(r->typ[pos].data.is.componentWeights != NULL)
5004  {
5005#if MYTEST
5006    PrintS("Deleting old componentWeights: "); r->typ[pos].data.is.componentWeights->show(); PrintLn();
5007#endif
5008    delete r->typ[pos].data.is.componentWeights;
5009    r->typ[pos].data.is.componentWeights = NULL;
5010  }
5011
5012
5013  assume(r->typ[pos].data.is.componentWeights == NULL);
5014
5015  if( componentWeights != NULL )
5016    componentWeights = ivCopy(componentWeights); // componentWeights is owened by ring now! TODO: delete at the end!
5017
5018  r->typ[pos].data.is.componentWeights = componentWeights;
5019
5020  r->typ[pos].data.is.limit = i; // First induced component
5021
5022#if MYTEST
5023  PrintS("New reference set FF : \n");        idShow(FF, r, r, 1);         PrintLn();
5024#endif
5025
5026  return TRUE;
5027}
5028
5029
5030void rSetSyzComp(int k, const ring r)
5031{
5032  if(k < 0)
5033  {
5034    dReportError("rSetSyzComp with negative limit!");
5035    return;
5036  }
5037
5038  assume( k >= 0 );
5039  if (TEST_OPT_PROT) Print("{%d}", k);
5040  if ((r->typ!=NULL) && (r->typ[0].ord_typ==ro_syz))
5041  {
5042    if( k == r->typ[0].data.syz.limit )
5043      return; // nothing to do
5044
5045    int i;
5046    if (r->typ[0].data.syz.limit == 0)
5047    {
5048      r->typ[0].data.syz.syz_index = (int*) omAlloc0((k+1)*sizeof(int));
5049      r->typ[0].data.syz.syz_index[0] = 0;
5050      r->typ[0].data.syz.curr_index = 1;
5051    }
5052    else
5053    {
5054      r->typ[0].data.syz.syz_index = (int*)
5055        omReallocSize(r->typ[0].data.syz.syz_index,
5056                (r->typ[0].data.syz.limit+1)*sizeof(int),
5057                (k+1)*sizeof(int));
5058    }
5059    for (i=r->typ[0].data.syz.limit + 1; i<= k; i++)
5060    {
5061      r->typ[0].data.syz.syz_index[i] =
5062        r->typ[0].data.syz.curr_index;
5063    }
5064    if(k < r->typ[0].data.syz.limit) // ?
5065    {
5066#ifndef NDEBUG
5067      Warn("rSetSyzComp called with smaller limit (%d) as before (%d)", k, r->typ[0].data.syz.limit);
5068#endif
5069      r->typ[0].data.syz.curr_index = 1 + r->typ[0].data.syz.syz_index[k];
5070    }
5071
5072
5073    r->typ[0].data.syz.limit = k;
5074    r->typ[0].data.syz.curr_index++;
5075  }
5076  else if(
5077            (r->typ!=NULL) &&
5078            (r->typ[0].ord_typ==ro_isTemp)
5079           )
5080  {
5081//      (r->typ[currRing->typ[0].data.isTemp.suffixpos].data.is.limit == k)
5082#ifndef NDEBUG
5083    Warn("rSetSyzComp(%d) in an IS ring! Be careful!", k);
5084#endif
5085  }
5086  else
5087  if ((r->order[0]!=ringorder_c) && (k!=0)) // ???
5088  {
5089    dReportError("syzcomp in incompatible ring");
5090  }
5091#ifdef PDEBUG
5092  extern int pDBsyzComp;
5093  pDBsyzComp=k;
5094#endif
5095}
5096
5097// return the max-comonent wchich has syzIndex i
5098int rGetMaxSyzComp(int i, const ring r)
5099{
5100  if ((r->typ!=NULL) && (r->typ[0].ord_typ==ro_syz) &&
5101      r->typ[0].data.syz.limit > 0 && i > 0)
5102  {
5103    assume(i <= r->typ[0].data.syz.limit);
5104    int j;
5105    for (j=0; j<r->typ[0].data.syz.limit; j++)
5106    {
5107      if (r->typ[0].data.syz.syz_index[j] == i  &&
5108          r->typ[0].data.syz.syz_index[j+1] != i)
5109      {
5110        assume(r->typ[0].data.syz.syz_index[j+1] == i+1);
5111        return j;
5112      }
5113    }
5114    return r->typ[0].data.syz.limit;
5115  }
5116  else
5117  {
5118    return 0;
5119  }
5120}
5121
5122BOOLEAN rRing_is_Homog(ring r)
5123{
5124  if (r == NULL) return FALSE;
5125  int i, j, nb = rBlocks(r);
5126  for (i=0; i<nb; i++)
5127  {
5128    if (r->wvhdl[i] != NULL)
5129    {
5130      int length = r->block1[i] - r->block0[i];
5131      int* wvhdl = r->wvhdl[i];
5132      if (r->order[i] == ringorder_M) length *= length;
5133      assume(omSizeOfAddr(wvhdl) >= length*sizeof(int));
5134
5135      for (j=0; j< length; j++)
5136      {
5137        if (wvhdl[j] != 0 && wvhdl[j] != 1) return FALSE;
5138      }
5139    }
5140  }
5141  return TRUE;
5142}
5143
5144BOOLEAN rRing_has_CompLastBlock(ring r)
5145{
5146  assume(r != NULL);
5147  int lb = rBlocks(r) - 2;
5148  return (r->order[lb] == ringorder_c || r->order[lb] == ringorder_C);
5149}
5150
5151n_coeffType rFieldType(ring r)
5152{
5153  return (r->cf->type);
5154  if (rField_is_Zp(r))     return n_Zp;
5155  if (rField_is_Q(r))      return n_Q;
5156  if (rField_is_R(r))      return n_R;
5157  if (rField_is_GF(r))     return n_GF;
5158  if (rField_is_long_R(r)) return n_long_R;
5159  if (rField_is_Zp_a(r))   return n_Zp_a;
5160  if (rField_is_Q_a(r))    return n_Q_a;
5161  if (rField_is_long_C(r)) return n_long_C;
5162  #ifdef HAVE_RINGS
5163   if (rField_is_Ring_Z(r)) return n_Z;
5164   if (rField_is_Ring_ModN(r)) return n_Zn;
5165   if (rField_is_Ring_PtoM(r)) return n_Zpn;
5166   if (rField_is_Ring_2toM(r)) return  n_Z2m;
5167  #endif
5168
5169  return n_unknown;
5170}
5171
5172int64 * rGetWeightVec(ring r)
5173{
5174  assume(r!=NULL);
5175  assume(r->OrdSize>0);
5176  int i=0;
5177  while((r->typ[i].ord_typ!=ro_wp64) && (r->typ[i].ord_typ>0)) i++;
5178  assume(r->typ[i].ord_typ==ro_wp64);
5179  return (int64*)(r->typ[i].data.wp64.weights64);
5180}
5181
5182void rSetWeightVec(ring r, int64 *wv)
5183{
5184  assume(r!=NULL);
5185  assume(r->OrdSize>0);
5186  assume(r->typ[0].ord_typ==ro_wp64);
5187  memcpy(r->typ[0].data.wp64.weights64,wv,r->N*sizeof(int64));
5188}
5189
5190#include <ctype.h>
5191
5192static int rRealloc1(ring r, int size, int pos)
5193{
5194  r->order=(int*)omReallocSize(r->order, size*sizeof(int), (size+1)*sizeof(int));
5195  r->block0=(int*)omReallocSize(r->block0, size*sizeof(int), (size+1)*sizeof(int));
5196  r->block1=(int*)omReallocSize(r->block1, size*sizeof(int), (size+1)*sizeof(int));
5197  r->wvhdl=(int **)omReallocSize(r->wvhdl,size*sizeof(int *), (size+1)*sizeof(int *));
5198  for(int k=size; k>pos; k--) r->wvhdl[k]=r->wvhdl[k-1];
5199  r->order[size]=0;
5200  size++;
5201  return size;
5202}
5203//static int rReallocM1(ring r, ring src, int size, int pos)
5204//{
5205//  r->order=(int*)omReallocSize(r->order, size*sizeof(int), (size-1)*sizeof(int));
5206//  r->block0=(int*)omReallocSize(r->block0, size*sizeof(int), (size-1)*sizeof(int));
5207//  r->block1=(int*)omReallocSize(r->block1, size*sizeof(int), (size-1)*sizeof(int));
5208//  r->wvhdl=(int **)omReallocSize(r->wvhdl,size*sizeof(int *), (size-1)*sizeof(int *));
5209//  for(int k=pos+1; k<size; k++) r->wvhdl[k]=r->wvhdl[k+1];
5210//  size--;
5211//  return size;
5212//}
5213static void rOppWeight(int *w, int l)
5214{
5215  int i2=(l+1)/2;
5216  for(int j=0; j<=i2; j++)
5217  {
5218    int t=w[j];
5219    w[j]=w[l-j];
5220    w[l-j]=t;
5221  }
5222}
5223
5224#define rOppVar(R,I) (rVar(R)+1-I)
5225
5226ring rOpposite(ring src)
5227  /* creates an opposite algebra of R */
5228  /* that is R^opp, where f (*^opp) g = g*f  */
5229  /* treats the case of qring */
5230{
5231  if (src == NULL) return(NULL);
5232
5233#ifdef RDEBUG
5234  rTest(src);
5235#endif
5236
5237  //rChangeCurrRing(src);
5238
5239#ifdef RDEBUG
5240  rTest(src);
5241//  rWrite(src);
5242//  rDebugPrint(src);
5243#endif
5244
5245
5246  ring r = rCopy0(src,FALSE); /* qideal will be deleted later on!!! */
5247
5248  // change vars v1..vN -> vN..v1
5249  int i;
5250  int i2 = (rVar(r)-1)/2;
5251  for(i=i2; i>=0; i--)
5252  {
5253    // index: 0..N-1
5254    //Print("ex var names: %d <-> %d\n",i,rOppVar(r,i));
5255    // exchange names
5256    char *p;
5257    p = r->names[rVar(r)-1-i];
5258    r->names[rVar(r)-1-i] = r->names[i];
5259    r->names[i] = p;
5260  }
5261//  i2=(rVar(r)+1)/2;
5262//  for(int i=i2; i>0; i--)
5263//  {
5264//    // index: 1..N
5265//    //Print("ex var places: %d <-> %d\n",i,rVar(r)+1-i);
5266//    // exchange VarOffset
5267//    int t;
5268//    t=r->VarOffset[i];
5269//    r->VarOffset[i]=r->VarOffset[rOppVar(r,i)];
5270//    r->VarOffset[rOppVar(r,i)]=t;
5271//  }
5272  // change names:
5273  for (i=rVar(r)-1; i>=0; i--)
5274  {
5275    char *p=r->names[i];
5276    if(isupper(*p)) *p = tolower(*p);
5277    else            *p = toupper(*p);
5278  }
5279  // change ordering: listing
5280  // change ordering: compare
5281//  for(i=0; i<r->OrdSize; i++)
5282//  {
5283//    int t,tt;
5284//    switch(r->typ[i].ord_typ)
5285//    {
5286//      case ro_dp:
5287//      //
5288//        t=r->typ[i].data.dp.start;
5289//        r->typ[i].data.dp.start=rOppVar(r,r->typ[i].data.dp.end);
5290//        r->typ[i].data.dp.end=rOppVar(r,t);
5291//        break;
5292//      case ro_wp:
5293//      case ro_wp_neg:
5294//      {
5295//        t=r->typ[i].data.wp.start;
5296//        r->typ[i].data.wp.start=rOppVar(r,r->typ[i].data.wp.end);
5297//        r->typ[i].data.wp.end=rOppVar(r,t);
5298//        // invert r->typ[i].data.wp.weights
5299//        rOppWeight(r->typ[i].data.wp.weights,
5300//                   r->typ[i].data.wp.end-r->typ[i].data.wp.start);
5301//        break;
5302//      }
5303//      //case ro_wp64:
5304//      case ro_syzcomp:
5305//      case ro_syz:
5306//         WerrorS("not implemented in rOpposite");
5307//         // should not happen
5308//         break;
5309//
5310//      case ro_cp:
5311//        t=r->typ[i].data.cp.start;
5312//        r->typ[i].data.cp.start=rOppVar(r,r->typ[i].data.cp.end);
5313//        r->typ[i].data.cp.end=rOppVar(r,t);
5314//        break;
5315//      case ro_none:
5316//      default:
5317//       Werror("unknown type in rOpposite(%d)",r->typ[i].ord_typ);
5318//       break;
5319//    }
5320//  }
5321  // Change order/block structures (needed for rPrint, rAdd etc.)
5322  int j=0;
5323  int l=rBlocks(src);
5324  for(i=0; src->order[i]!=0; i++)
5325  {
5326    switch (src->order[i])
5327    {
5328      case ringorder_c: /* c-> c */
5329      case ringorder_C: /* C-> C */
5330      case ringorder_no /*=0*/: /* end-of-block */
5331        r->order[j]=src->order[i];
5332        j++; break;
5333      case ringorder_lp: /* lp -> rp */
5334        r->order[j]=ringorder_rp;
5335        r->block0[j]=rOppVar(r, src->block1[i]);
5336        r->block1[j]=rOppVar(r, src->block0[i]);
5337        break;
5338      case ringorder_rp: /* rp -> lp */
5339        r->order[j]=ringorder_lp;
5340        r->block0[j]=rOppVar(r, src->block1[i]);
5341        r->block1[j]=rOppVar(r, src->block0[i]);
5342        break;
5343      case ringorder_dp: /* dp -> a(1..1),ls */
5344      {
5345        l=rRealloc1(r,l,j);
5346        r->order[j]=ringorder_a;
5347        r->block0[j]=rOppVar(r, src->block1[i]);
5348        r->block1[j]=rOppVar(r, src->block0[i]);
5349        r->wvhdl[j]=(int*)omAlloc((r->block1[j]-r->block0[j]+1)*sizeof(int));
5350        for(int k=r->block0[j]; k<=r->block1[j]; k++)
5351          r->wvhdl[j][k-r->block0[j]]=1;
5352        j++;
5353        r->order[j]=ringorder_ls;
5354        r->block0[j]=rOppVar(r, src->block1[i]);
5355        r->block1[j]=rOppVar(r, src->block0[i]);
5356        j++;
5357        break;
5358      }
5359      case ringorder_Dp: /* Dp -> a(1..1),rp */
5360      {
5361        l=rRealloc1(r,l,j);
5362        r->order[j]=ringorder_a;
5363        r->block0[j]=rOppVar(r, src->block1[i]);
5364        r->block1[j]=rOppVar(r, src->block0[i]);
5365        r->wvhdl[j]=(int*)omAlloc((r->block1[j]-r->block0[j]+1)*sizeof(int));
5366        for(int k=r->block0[j]; k<=r->block1[j]; k++)
5367          r->wvhdl[j][k-r->block0[j]]=1;
5368        j++;
5369        r->order[j]=ringorder_rp;
5370        r->block0[j]=rOppVar(r, src->block1[i]);
5371        r->block1[j]=rOppVar(r, src->block0[i]);
5372        j++;
5373        break;
5374      }
5375      case ringorder_wp: /* wp -> a(...),ls */
5376      {
5377        l=rRealloc1(r,l,j);
5378        r->order[j]=ringorder_a;
5379        r->block0[j]=rOppVar(r, src->block1[i]);
5380        r->block1[j]=rOppVar(r, src->block0[i]);
5381        r->wvhdl[j]=r->wvhdl[j+1]; r->wvhdl[j+1]=NULL;
5382        rOppWeight(r->wvhdl[j], r->block1[j]-r->block0[j]);
5383        j++;
5384        r->order[j]=ringorder_ls;
5385        r->block0[j]=rOppVar(r, src->block1[i]);
5386        r->block1[j]=rOppVar(r, src->block0[i]);
5387        j++;
5388        break;
5389      }
5390      case ringorder_Wp: /* Wp -> a(...),rp */
5391      {
5392        l=rRealloc1(r,l,j);
5393        r->order[j]=ringorder_a;
5394        r->block0[j]=rOppVar(r, src->block1[i]);
5395        r->block1[j]=rOppVar(r, src->block0[i]);
5396        r->wvhdl[j]=r->wvhdl[j+1]; r->wvhdl[j+1]=NULL;
5397        rOppWeight(r->wvhdl[j], r->block1[j]-r->block0[j]);
5398        j++;
5399        r->order[j]=ringorder_rp;
5400        r->block0[j]=rOppVar(r, src->block1[i]);
5401        r->block1[j]=rOppVar(r, src->block0[i]);
5402        j++;
5403        break;
5404      }
5405      case ringorder_M: /* M -> M */
5406      {
5407        r->order[j]=ringorder_M;
5408        r->block0[j]=rOppVar(r, src->block1[i]);
5409        r->block1[j]=rOppVar(r, src->block0[i]);
5410        int n=r->block1[j]-r->block0[j];
5411        /* M is a (n+1)x(n+1) matrix */
5412        for (int nn=0; nn<=n; nn++)
5413        {
5414          rOppWeight(&(r->wvhdl[j][nn*(n+1)]), n /*r->block1[j]-r->block0[j]*/);
5415        }
5416        j++;
5417        break;
5418      }
5419      case ringorder_a: /*  a(...),ls -> wp/dp */
5420      {
5421        r->block0[j]=rOppVar(r, src->block1[i]);
5422        r->block1[j]=rOppVar(r, src->block0[i]);
5423        rOppWeight(r->wvhdl[j], r->block1[j]-r->block0[j]);
5424        if (src->order[i+1]==ringorder_ls)
5425        {
5426          r->order[j]=ringorder_wp;
5427          i++;
5428          //l=rReallocM1(r,src,l,j);
5429        }
5430        else
5431        {
5432          r->order[j]=ringorder_a;
5433        }
5434        j++;
5435        break;
5436      }
5437      // not yet done:
5438      case ringorder_ls:
5439      case ringorder_rs:
5440      case ringorder_ds:
5441      case ringorder_Ds:
5442      case ringorder_ws:
5443      case ringorder_Ws:
5444      // should not occur:
5445      case ringorder_S:
5446      case ringorder_IS:
5447      case ringorder_s:
5448      case ringorder_aa:
5449      case ringorder_L:
5450      case ringorder_unspec:
5451        Werror("order %s not (yet) supported", rSimpleOrdStr(src->order[i]));
5452        break;
5453    }
5454  }
5455  rComplete(r);
5456
5457
5458#ifdef RDEBUG
5459  rTest(r);
5460#endif
5461
5462  //rChangeCurrRing(r);
5463
5464#ifdef RDEBUG
5465  rTest(r);
5466//  rWrite(r);
5467//  rDebugPrint(r);
5468#endif
5469
5470
5471#ifdef HAVE_PLURAL
5472  // now, we initialize a non-comm structure on r
5473  if (rIsPluralRing(src))
5474  {
5475//    assume( currRing == r);
5476
5477    int *perm       = (int *)omAlloc0((rVar(r)+1)*sizeof(int));
5478    int *par_perm   = NULL;
5479    nMapFunc nMap   = n_SetMap(src->cf,r->cf);
5480    int ni,nj;
5481    for(i=1; i<=r->N; i++)
5482    {
5483      perm[i] = rOppVar(r,i);
5484    }
5485
5486    matrix C = mpNew(rVar(r),rVar(r));
5487    matrix D = mpNew(rVar(r),rVar(r));
5488
5489    for (i=1; i< rVar(r); i++)
5490    {
5491      for (j=i+1; j<=rVar(r); j++)
5492      {
5493        ni = r->N +1 - i;
5494        nj = r->N +1 - j; /* i<j ==>   nj < ni */
5495
5496        assume(MATELEM(src->GetNC()->C,i,j) != NULL);
5497        MATELEM(C,nj,ni) = p_PermPoly(MATELEM(src->GetNC()->C,i,j),perm,src,r, nMap,par_perm,src->P);
5498
5499        if(MATELEM(src->GetNC()->D,i,j) != NULL)
5500          MATELEM(D,nj,ni) = p_PermPoly(MATELEM(src->GetNC()->D,i,j),perm,src,r, nMap,par_perm,src->P);
5501      }
5502    }
5503
5504    id_Test((ideal)C, r);
5505    id_Test((ideal)D, r);
5506
5507    if (nc_CallPlural(C, D, NULL, NULL, r, false, false, true, r)) // no qring setup!
5508      WarnS("Error initializing non-commutative multiplication!");
5509
5510#ifdef RDEBUG
5511    rTest(r);
5512//    rWrite(r);
5513//    rDebugPrint(r);
5514#endif
5515
5516    assume( r->GetNC()->IsSkewConstant == src->GetNC()->IsSkewConstant);
5517
5518    omFreeSize((ADDRESS)perm,(rVar(r)+1)*sizeof(int));
5519  }
5520#endif /* HAVE_PLURAL */
5521
5522  /* now oppose the qideal for qrings */
5523  if (src->qideal != NULL)
5524  {
5525    id_Delete(&(r->qideal), r);
5526
5527#ifdef HAVE_PLURAL
5528    r->qideal = idOppose(src, src->qideal, r); // into the currRing: r
5529#else
5530    r->qideal = id_Copy(src->qideal, r); // ?
5531#endif
5532
5533#ifdef HAVE_PLURAL
5534    if( rIsPluralRing(r) )
5535    {
5536      nc_SetupQuotient(r);
5537#ifdef RDEBUG
5538      rTest(r);
5539//      rWrite(r);
5540//      rDebugPrint(r);
5541#endif
5542    }
5543#endif
5544  }
5545#ifdef HAVE_PLURAL
5546  if( rIsPluralRing(r) )
5547    assume( ncRingType(r) == ncRingType(src) );
5548#endif
5549  rTest(r);
5550
5551  return r;
5552}
5553
5554ring rEnvelope(ring R)
5555  /* creates an enveloping algebra of R */
5556  /* that is R^e = R \tensor_K R^opp */
5557{
5558  ring Ropp = rOpposite(R);
5559  ring Renv = NULL;
5560  int stat = rSum(R, Ropp, Renv); /* takes care of qideals */
5561  if ( stat <=0 )
5562    WarnS("Error in rEnvelope at rSum");
5563  rTest(Renv);
5564  return Renv;
5565}
5566
5567#ifdef HAVE_PLURAL
5568BOOLEAN nc_rComplete(const ring src, ring dest, bool bSetupQuotient)
5569/* returns TRUE is there were errors */
5570/* dest is actualy equals src with the different ordering */
5571/* we map src->nc correctly to dest->src */
5572/* to be executed after rComplete, before rChangeCurrRing */
5573{
5574// NOTE: Originally used only by idElimination to transfer NC structure to dest
5575// ring created by dirty hack (without nc_CallPlural)
5576  rTest(src);
5577
5578  assume(!rIsPluralRing(dest)); // destination must be a newly constructed commutative ring
5579
5580  if (!rIsPluralRing(src))
5581  {
5582    return FALSE;
5583  }
5584
5585  const int N = dest->N;
5586
5587  assume(src->N == N);
5588
5589//  ring save = currRing;
5590
5591//  if (dest != save)
5592//    rChangeCurrRing(dest);
5593
5594  const ring srcBase = src;
5595
5596  assume( n_SetMap(srcBase->cf,dest->cf) == n_SetMap(dest->cf,dest->cf) ); // currRing is important here!
5597
5598  matrix C = mpNew(N,N); // ring independent
5599  matrix D = mpNew(N,N);
5600
5601  matrix C0 = src->GetNC()->C;
5602  matrix D0 = src->GetNC()->D;
5603
5604  // map C and D into dest
5605  for (int i = 1; i < N; i++)
5606  {
5607    for (int j = i + 1; j <= N; j++)
5608    {
5609      const number n = n_Copy(p_GetCoeff(MATELEM(C0,i,j), srcBase), srcBase->cf); // src, mapping for coeffs into currRing = dest!
5610      const poly   p = p_NSet(n, dest);
5611      MATELEM(C,i,j) = p;
5612      if (MATELEM(D0,i,j) != NULL)
5613        MATELEM(D,i,j) = prCopyR(MATELEM(D0,i,j), srcBase, dest); // ?
5614    }
5615  }
5616  /* One must test C and D _only_ in r->GetNC()->basering!!! not in r!!! */
5617
5618  id_Test((ideal)C, dest); 
5619  id_Test((ideal)D, dest);
5620
5621  if (nc_CallPlural(C, D, NULL, NULL, dest, bSetupQuotient, false, true, dest)) // also takes care about quotient ideal
5622  {
5623    //WarnS("Error transferring non-commutative structure");
5624    // error message should be in the interpreter interface
5625
5626    mp_Delete(&C, dest);
5627    mp_Delete(&D, dest);
5628
5629//    if (currRing != save)
5630//       rChangeCurrRing(save);
5631
5632    return TRUE;
5633  }
5634
5635//  mp_Delete(&C, dest); // used by nc_CallPlural!
5636//  mp_Delete(&D, dest);
5637
5638//  if (dest != save)
5639//    rChangeCurrRing(save);
5640
5641  assume(rIsPluralRing(dest));
5642  return FALSE;
5643}
5644#endif
5645
5646void rModify_a_to_A(ring r)
5647// to be called BEFORE rComplete:
5648// changes every Block with a(...) to A(...)
5649{
5650   int i=0;
5651   int j;
5652   while(r->order[i]!=0)
5653   {
5654      if (r->order[i]==ringorder_a)
5655      {
5656        r->order[i]=ringorder_a64;
5657        int *w=r->wvhdl[i];
5658        int64 *w64=(int64 *)omAlloc((r->block1[i]-r->block0[i]+1)*sizeof(int64));
5659        for(j=r->block1[i]-r->block0[i];j>=0;j--)
5660                w64[j]=(int64)w[j];
5661        r->wvhdl[i]=(int*)w64;
5662        omFreeSize(w,(r->block1[i]-r->block0[i]+1)*sizeof(int));
5663      }
5664      i++;
5665   }
5666}
Note: See TracBrowser for help on using the repository browser.