source: git/kernel/ring.cc @ 6e66d2

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