source: git/polys/monomials/ring.cc @ 879cc3c

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