source: git/kernel/ring.cc @ 959263

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