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