source: git/libpolys/polys/monomials/ring.cc @ 1e872a

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