source: git/libpolys/polys/monomials/ring.cc @ 6bec87

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