source: git/kernel/ring.cc @ 22579cf

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