source: git/kernel/ring.cc @ af3d2a

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