source: git/libpolys/polys/monomials/ring.cc @ 8c9912

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