source: git/kernel/sca.cc @ 171950

spielwiese
Last change on this file since 171950 was 06879b7, checked in by Motsak Oleksandr <motsak@…>, 16 years ago
*motsak: additional controls for testing/debugging/benchmarking gring | Use NC Extensions (SCA) sca | Force Commutative git-svn-id: file:///usr/local/Singular/svn/trunk@10717 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 57.8 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/***************************************************************
5 *  File:    sca.cc
6 *  Purpose: supercommutative kernel procedures
7 *  Author:  motsak (Oleksandr Motsak)
8 *  Created: 2006/12/18
9 *  Version: $Id: sca.cc,v 1.15 2008-05-15 17:24:36 motsak Exp $
10 *******************************************************************/
11
12// #define PDEBUG 2
13#include "mod2.h"
14
15#ifdef HAVE_PLURAL
16// for
17#define PLURAL_INTERNAL_DECLARATIONS
18#include "sca.h"
19
20
21#include "febase.h"
22
23#include "p_polys.h"
24#include "kutil.h"
25#include "ideals.h"
26#include "intvec.h"
27#include "polys.h"
28
29#include "ring.h"
30#include "numbers.h"
31#include "matpol.h"
32#include "kbuckets.h"
33#include "kstd1.h"
34#include "sbuckets.h"
35#include "prCopy.h"
36#include "p_Mult_q.h"
37#include "p_MemAdd.h"
38
39#include "kutil.h"
40#include "kstd1.h"
41
42#include "weight.h"
43
44
45// poly functions defined in p_Procs :
46
47// return pPoly * pMonom; preserve pPoly and pMonom.
48poly sca_pp_Mult_mm(const poly pPoly, const poly pMonom, const ring rRing, poly &);
49
50// return pMonom * pPoly; preserve pPoly and pMonom.
51poly sca_mm_Mult_pp(const poly pMonom, const poly pPoly, const ring rRing);
52
53// return pPoly * pMonom; preserve pMonom, destroy or reuse pPoly.
54poly sca_p_Mult_mm(poly pPoly, const poly pMonom, const ring rRing);
55
56// return pMonom * pPoly; preserve pMonom, destroy or reuse pPoly.
57poly sca_mm_Mult_p(const poly pMonom, poly pPoly, const ring rRing);
58
59
60// compute the spolynomial of p1 and p2
61poly sca_SPoly(const poly p1, const poly p2, const ring r);
62poly sca_ReduceSpoly(const poly p1, poly p2, const ring r);
63
64// Modified Plural's Buchberger's algorithmus.
65ideal sca_gr_bba(const ideal F, const ideal Q, const intvec *, const intvec *, kStrategy strat);
66
67// Modified modern Sinuglar Buchberger's algorithm.
68ideal sca_bba(const ideal F, const ideal Q, const intvec *w, const intvec *, kStrategy strat);
69
70// Modified modern Sinuglar Mora's algorithm.
71ideal sca_mora(const ideal F, const ideal Q, const intvec *w, const intvec *, kStrategy strat);
72
73
74
75////////////////////////////////////////////////////////////////////////////////////////////////////
76// Super Commutative Algabra extension by Motsak
77////////////////////////////////////////////////////////////////////////////////////////////////////
78
79
80
81
82
83// returns the sign of: lm(pMonomM) * lm(pMonomMM),
84// preserves input, may return +/-1, 0
85inline int sca_Sign_mm_Mult_mm( const poly pMonomM, const poly pMonomMM, const ring rRing )
86{
87#ifdef PDEBUG
88    p_Test(pMonomM,  rRing);
89    p_Test(pMonomMM, rRing);
90#endif
91
92    const unsigned int iFirstAltVar = scaFirstAltVar(rRing);
93    const unsigned int iLastAltVar  = scaLastAltVar(rRing);
94
95    unsigned int tpower = 0;
96    unsigned int cpower = 0;
97
98    for( unsigned int j = iLastAltVar; j >= iFirstAltVar; j-- )
99    {
100      const unsigned int iExpM  = p_GetExp(pMonomM,  j, rRing);
101      const unsigned int iExpMM = p_GetExp(pMonomMM, j, rRing);
102
103      if( iExpMM != 0 )
104      {
105        if( iExpM != 0 )
106        {
107          return 0; // lm(pMonomM) * lm(pMonomMM) == 0
108        }
109        tpower += cpower; // compute degree of (-1).
110      }
111
112      cpower += iExpM;
113    }
114
115    if( (tpower&1) != 0 ) // degree is odd => negate coeff.
116      return -1;
117
118    return(1);
119}
120
121
122
123
124// returns and changes pMonomM: lt(pMonomM) = lt(pMonomM) * lt(pMonomMM),
125// preserves pMonomMM. may return NULL!
126// if result != NULL => next(result) = next(pMonomM), lt(result) = lt(pMonomM) * lt(pMonomMM)
127// if result == NULL => pMonomM MUST BE deleted manually!
128inline poly sca_m_Mult_mm( poly pMonomM, const poly pMonomMM, const ring rRing )
129{
130#ifdef PDEBUG
131    p_Test(pMonomM,  rRing);
132    p_Test(pMonomMM, rRing);
133#endif
134
135    const unsigned int iFirstAltVar = scaFirstAltVar(rRing);
136    const unsigned int iLastAltVar = scaLastAltVar(rRing);
137
138    unsigned int tpower = 0;
139    unsigned int cpower = 0;
140
141    for( unsigned int j = iLastAltVar; j >= iFirstAltVar; j-- )
142    {
143      const unsigned int iExpM  = p_GetExp(pMonomM,  j, rRing);
144      const unsigned int iExpMM = p_GetExp(pMonomMM, j, rRing);
145
146      if( iExpMM != 0 )
147      {
148        if( iExpM != 0 ) // result is zero!
149        {
150          return NULL; // we do nothing with pMonomM in this case!
151        }
152
153        tpower += cpower; // compute degree of (-1).
154      }
155
156      cpower += iExpM;
157    }
158
159    p_ExpVectorAdd(pMonomM, pMonomMM, rRing); // "exponents" are additive!!!
160
161    number nCoeffM = p_GetCoeff(pMonomM, rRing); // no new copy! should be deleted!
162
163    if( (tpower&1) != 0 ) // degree is odd => negate coeff.
164      nCoeffM = n_Neg(nCoeffM, rRing); // negate nCoeff (will destroy the original number)
165
166    const number nCoeffMM = p_GetCoeff(pMonomMM, rRing); // no new copy!
167
168    number nCoeff = n_Mult(nCoeffM, nCoeffMM, rRing); // new number!
169
170    p_SetCoeff(pMonomM, nCoeff, rRing); // delete lc(pMonomM) and set lc(pMonomM) = nCoeff
171
172#ifdef PDEBUG
173    p_Test(pMonomM, rRing);
174#endif
175
176    return(pMonomM);
177}
178
179// returns and changes pMonomM: lt(pMonomM) = lt(pMonomMM) * lt(pMonomM),
180// preserves pMonomMM. may return NULL!
181// if result != NULL => next(result) = next(pMonomM), lt(result) = lt(pMonomMM) * lt(pMonomM)
182// if result == NULL => pMonomM MUST BE deleted manually!
183inline poly sca_mm_Mult_m( const poly pMonomMM, poly pMonomM, const ring rRing )
184{
185#ifdef PDEBUG
186    p_Test(pMonomM,  rRing);
187    p_Test(pMonomMM, rRing);
188#endif
189
190    const unsigned int iFirstAltVar = scaFirstAltVar(rRing);
191    const unsigned int iLastAltVar = scaLastAltVar(rRing);
192
193    unsigned int tpower = 0;
194    unsigned int cpower = 0;
195
196    for( unsigned int j = iLastAltVar; j >= iFirstAltVar; j-- )
197    {
198      const unsigned int iExpMM = p_GetExp(pMonomMM, j, rRing);
199      const unsigned int iExpM  = p_GetExp(pMonomM,  j, rRing);
200
201      if( iExpM != 0 )
202      {
203        if( iExpMM != 0 ) // result is zero!
204        {
205          return NULL; // we do nothing with pMonomM in this case!
206        }
207
208        tpower += cpower; // compute degree of (-1).
209      }
210
211      cpower += iExpMM;
212    }
213
214    p_ExpVectorAdd(pMonomM, pMonomMM, rRing); // "exponents" are additive!!!
215
216    number nCoeffM = p_GetCoeff(pMonomM, rRing); // no new copy! should be deleted!
217
218    if( (tpower&1) != 0 ) // degree is odd => negate coeff.
219      nCoeffM = n_Neg(nCoeffM, rRing); // negate nCoeff (will destroy the original number), creates new number!
220
221    const number nCoeffMM = p_GetCoeff(pMonomMM, rRing); // no new copy!
222
223    number nCoeff = n_Mult(nCoeffM, nCoeffMM, rRing); // new number!
224
225    p_SetCoeff(pMonomM, nCoeff, rRing); // delete lc(pMonomM) and set lc(pMonomM) = nCoeff
226
227#ifdef PDEBUG
228    p_Test(pMonomM, rRing);
229#endif
230
231    return(pMonomM);
232}
233
234
235
236// returns: result = lt(pMonom1) * lt(pMonom2),
237// preserves pMonom1, pMonom2. may return NULL!
238// if result != NULL => next(result) = NULL, lt(result) = lt(pMonom1) * lt(pMonom2)
239inline poly sca_mm_Mult_mm( poly pMonom1, const poly pMonom2, const ring rRing )
240{
241#ifdef PDEBUG
242    p_Test(pMonom1, rRing);
243    p_Test(pMonom2, rRing);
244#endif
245
246    const unsigned int iFirstAltVar = scaFirstAltVar(rRing);
247    const unsigned int iLastAltVar = scaLastAltVar(rRing);
248
249    unsigned int tpower = 0;
250    unsigned int cpower = 0;
251
252    for( unsigned int j = iLastAltVar; j >= iFirstAltVar; j-- )
253    {
254      const unsigned int iExp1 = p_GetExp(pMonom1, j, rRing);
255      const unsigned int iExp2 = p_GetExp(pMonom2, j, rRing);
256
257      if( iExp2 != 0 )
258      {
259        if( iExp1 != 0 ) // result is zero!
260        {
261          return NULL;
262        }
263
264        tpower += cpower; // compute degree of (-1).
265      }
266
267      cpower += iExp1;
268    }
269
270    poly pResult;
271    omTypeAllocBin(poly, pResult, rRing->PolyBin);
272    p_SetRingOfLm(pResult, rRing);
273
274    p_ExpVectorSum(pResult, pMonom1, pMonom2, rRing); // "exponents" are additive!!!
275
276    pNext(pResult) = NULL;
277
278    const number nCoeff1 = p_GetCoeff(pMonom1, rRing); // no new copy!
279    const number nCoeff2 = p_GetCoeff(pMonom2, rRing); // no new copy!
280
281    number nCoeff = n_Mult(nCoeff1, nCoeff2, rRing); // new number!
282
283    if( (tpower&1) != 0 ) // degree is odd => negate coeff.
284      nCoeff = n_Neg(nCoeff, rRing); // negate nCoeff (will destroy the original number)
285
286    p_SetCoeff0(pResult, nCoeff, rRing); // set lc(pResult) = nCoeff, no destruction!
287
288#ifdef PDEBUG
289    p_Test(pResult, rRing);
290#endif
291
292    return(pResult);
293}
294
295// returns: result =  x_i * lt(pMonom),
296// preserves pMonom. may return NULL!
297// if result != NULL => next(result) = NULL, lt(result) = x_i * lt(pMonom)
298inline poly sca_xi_Mult_mm(unsigned int i, const poly pMonom, const ring rRing)
299{
300#ifdef PDEBUG
301    p_Test(pMonom, rRing);
302#endif
303
304    assume( i <= scaLastAltVar(rRing));
305    assume( scaFirstAltVar(rRing) <= i );
306
307    if( p_GetExp(pMonom, i, rRing) != 0 ) // => result is zero!
308      return NULL;
309
310    const unsigned int iFirstAltVar = scaFirstAltVar(rRing);
311
312    unsigned int cpower = 0;
313
314    for( unsigned int j = iFirstAltVar; j < i ; j++ )
315      cpower += p_GetExp(pMonom, j, rRing);
316
317    poly pResult = p_LmInit(pMonom, rRing);
318
319    p_SetExp(pResult, i, 1, rRing); // pResult*=X_i &&
320    p_Setm(pResult, rRing);         // addjust degree after previous step!
321
322    number nCoeff = n_Copy(p_GetCoeff(pMonom, rRing), rRing); // new number!
323
324    if( (cpower&1) != 0 ) // degree is odd => negate coeff.
325      nCoeff = n_Neg(nCoeff, rRing); // negate nCoeff (will destroy the original number)
326
327    p_SetCoeff0(pResult, nCoeff, rRing); // set lc(pResult) = nCoeff, no destruction!
328
329#ifdef PDEBUG
330    p_Test(pResult, rRing);
331#endif
332
333    return(pResult);
334}
335
336//-----------------------------------------------------------------------------------//
337
338// return poly = pPoly * pMonom; preserve pMonom, destroy or reuse pPoly.
339poly sca_p_Mult_mm(poly pPoly, const poly pMonom, const ring rRing)
340{
341  assume( rIsSCA(rRing) );
342
343#ifdef PDEBUG
344//  Print("sca_p_Mult_mm\n"); // !
345
346  p_Test(pPoly, rRing);
347  p_Test(pMonom, rRing);
348#endif
349
350  if( pPoly == NULL )
351    return NULL;
352
353  if( pMonom == NULL )
354  {
355    // pPoly != NULL =>
356    p_Delete( &pPoly, rRing );
357
358    return NULL;
359  }
360
361  const int iComponentMonomM = p_GetComp(pMonom, rRing);
362
363  poly p = pPoly; poly* ppPrev = &pPoly;
364
365  loop
366  {
367#ifdef PDEBUG
368    p_Test(p, rRing);
369#endif
370    const int iComponent = p_GetComp(p, rRing);
371
372    if( iComponent!=0 )
373    {
374      if( iComponentMonomM!=0 ) // TODO: make global if on iComponentMonomM =?= 0
375      {
376        // REPORT_ERROR
377        Werror("sca_p_Mult_mm: exponent mismatch %d and %d\n", iComponent, iComponentMonomM);
378        // what should we do further?!?
379
380        p_Delete( &pPoly, rRing); // delete the result AND rest
381        return NULL;
382      }
383#ifdef PDEBUG
384      if(iComponentMonomM==0 )
385      {
386        WarnS("Multiplication in the left module from the right");
387      }
388#endif
389    }
390
391    // terms will be in the same order because of quasi-ordering!
392    poly v = sca_m_Mult_mm(p, pMonom, rRing);
393
394    if( v != NULL )
395    {
396      ppPrev = &pNext(p); // fixed!
397
398    // *p is changed if v != NULL ( p == v )
399      pIter(p);
400
401      if( p == NULL )
402        break;
403    }
404    else
405    { // Upps! Zero!!! we must kill this term!
406
407      //
408      p = p_LmDeleteAndNext(p, rRing);
409
410      *ppPrev = p;
411
412      if( p == NULL )
413        break;
414    }
415
416
417  }
418
419#ifdef PDEBUG
420  p_Test(pPoly,rRing);
421#endif
422
423  return(pPoly);
424}
425
426// return new poly = pPoly * pMonom; preserve pPoly and pMonom.
427poly sca_pp_Mult_mm(const poly pPoly, const poly pMonom, const ring rRing, poly &)
428{
429  assume( rIsSCA(rRing) );
430
431#ifdef PDEBUG
432//  Print("sca_pp_Mult_mm\n"); // !
433
434  p_Test(pPoly, rRing);
435  p_Test(pMonom, rRing);
436#endif
437
438  if( ( pPoly == NULL ) || ( pMonom == NULL ) )
439    return NULL;
440
441  const int iComponentMonomM = p_GetComp(pMonom, rRing);
442
443  poly pResult = NULL;
444  poly* ppPrev = &pResult;
445
446  for( poly p = pPoly; p!= NULL; pIter(p) )
447  {
448#ifdef PDEBUG
449    p_Test(p, rRing);
450#endif
451    const int iComponent = p_GetComp(p, rRing);
452
453    if( iComponent!=0 )
454    {
455      if( iComponentMonomM!=0 ) // TODO: make global if on iComponentMonomM =?= 0
456      {
457        // REPORT_ERROR
458        Werror("sca_pp_Mult_mm: exponent mismatch %d and %d\n", iComponent, iComponentMonomM);
459        // what should we do further?!?
460
461        p_Delete( &pResult, rRing); // delete the result
462        return NULL;
463      }
464
465#ifdef PDEBUG
466      if(iComponentMonomM==0 )
467      {
468        WarnS("Multiplication in the left module from the right");
469      }
470#endif
471    }
472
473    // terms will be in the same order because of quasi-ordering!
474    poly v = sca_mm_Mult_mm(p, pMonom, rRing);
475
476    if( v != NULL )
477    {
478      *ppPrev = v;
479      ppPrev = &pNext(v);
480    }
481
482  }
483
484#ifdef PDEBUG
485  p_Test(pResult,rRing);
486#endif
487
488  return(pResult);
489}
490
491//-----------------------------------------------------------------------------------//
492
493// return x_i * pPoly; preserve pPoly.
494inline poly sca_xi_Mult_pp(unsigned int i, const poly pPoly, const ring rRing)
495{
496  assume( rIsSCA(rRing) );
497
498#ifdef PDEBUG
499  p_Test(pPoly, rRing);
500#endif
501
502  assume(i <= scaLastAltVar(rRing));
503  assume(scaFirstAltVar(rRing) <= i);
504
505  if( pPoly == NULL )
506    return NULL;
507
508  poly pResult = NULL;
509  poly* ppPrev = &pResult;
510
511  for( poly p = pPoly; p!= NULL; pIter(p) )
512  {
513
514    // terms will be in the same order because of quasi-ordering!
515    poly v = sca_xi_Mult_mm(i, p, rRing);
516
517#ifdef PDEBUG
518    p_Test(v, rRing);
519#endif
520
521    if( v != NULL )
522    {
523      *ppPrev = v;
524      ppPrev = &pNext(*ppPrev);
525    }
526  }
527
528
529#ifdef PDEBUG
530  p_Test(pResult, rRing);
531#endif
532
533  return(pResult);
534}
535
536
537// return new poly = pMonom * pPoly; preserve pPoly and pMonom.
538poly sca_mm_Mult_pp(const poly pMonom, const poly pPoly, const ring rRing)
539{
540  assume( rIsSCA(rRing) );
541
542#ifdef PDEBUG
543//  Print("sca_mm_Mult_pp\n"); // !
544
545  p_Test(pPoly, rRing);
546  p_Test(pMonom, rRing);
547#endif
548
549  if( ( pPoly == NULL ) || ( pMonom == NULL ) )
550    return NULL;
551
552  const int iComponentMonomM = p_GetComp(pMonom, rRing);
553
554  poly pResult = NULL;
555  poly* ppPrev = &pResult;
556
557  for( poly p = pPoly; p!= NULL; pIter(p) )
558  {
559#ifdef PDEBUG
560    p_Test(p, rRing);
561#endif
562    const int iComponent = p_GetComp(p, rRing);
563
564    if( iComponentMonomM!=0 )
565    {
566      if( iComponent!=0 ) // TODO: make global if on iComponentMonomM =?= 0
567      {
568        // REPORT_ERROR
569        Werror("sca_mm_Mult_pp: exponent mismatch %d and %d\n", iComponent, iComponentMonomM);
570        // what should we do further?!?
571
572        p_Delete( &pResult, rRing); // delete the result
573        return NULL;
574      }
575#ifdef PDEBUG
576      if(iComponent==0 )
577      {
578        WarnS("Multiplication in the left module from the right");
579      }
580#endif
581    }
582
583    // terms will be in the same order because of quasi-ordering!
584    poly v = sca_mm_Mult_mm(pMonom, p, rRing);
585
586    if( v != NULL )
587    {
588      *ppPrev = v;
589      ppPrev = &pNext(*ppPrev); // nice line ;-)
590    }
591  }
592
593#ifdef PDEBUG
594  p_Test(pResult,rRing);
595#endif
596
597  return(pResult);
598}
599
600
601// return poly = pMonom * pPoly; preserve pMonom, destroy or reuse pPoly.
602poly sca_mm_Mult_p(const poly pMonom, poly pPoly, const ring rRing) // !!!!! the MOST used procedure !!!!!
603{
604  assume( rIsSCA(rRing) );
605
606#ifdef PDEBUG
607  p_Test(pPoly, rRing);
608  p_Test(pMonom, rRing);
609#endif
610
611  if( pPoly == NULL )
612    return NULL;
613
614  if( pMonom == NULL )
615  {
616    // pPoly != NULL =>
617    p_Delete( &pPoly, rRing );
618    return NULL;
619  }
620
621  const int iComponentMonomM = p_GetComp(pMonom, rRing);
622
623  poly p = pPoly; poly* ppPrev = &pPoly;
624
625  loop
626  {
627#ifdef PDEBUG
628    if( !p_Test(p, rRing) )
629    {
630      Print("p is wrong!");
631      p_Write(p,rRing);
632    }
633#endif
634
635    const int iComponent = p_GetComp(p, rRing);
636
637    if( iComponentMonomM!=0 )
638    {
639      if( iComponent!=0 ) // TODO: make global if on iComponentMonomM =?= 0
640      {
641        // REPORT_ERROR
642        Werror("sca_mm_Mult_p: exponent mismatch %d and %d\n", iComponent, iComponentMonomM);
643        // what should we do further?!?
644
645        p_Delete( &pPoly, rRing); // delete the result
646        return NULL;
647      }
648#ifdef PDEBUG
649      if(iComponent==0 )
650      {
651        WarnS("Multiplication in the left module from the right");
652      }
653#endif
654    }
655
656    // terms will be in the same order because of quasi-ordering!
657    poly v = sca_mm_Mult_m(pMonom, p, rRing);
658
659    if( v != NULL )
660    {
661      ppPrev = &pNext(p);
662
663    // *p is changed if v != NULL ( p == v )
664      pIter(p);
665
666      if( p == NULL )
667        break;
668    }
669    else
670    { // Upps! Zero!!! we must kill this term!
671      p = p_LmDeleteAndNext(p, rRing);
672
673      *ppPrev = p;
674
675      if( p == NULL )
676        break;
677    }
678  }
679
680#ifdef PDEBUG
681    if( !p_Test(pPoly, rRing) )
682    {
683      Print("pPoly is wrong!");
684      p_Write(pPoly, rRing);
685    }
686#endif
687
688  return(pPoly);
689}
690
691//-----------------------------------------------------------------------------------//
692
693#ifdef PDEBUG
694#endif
695
696
697
698
699//-----------------------------------------------------------------------------------//
700
701// GB computation routines:
702
703/*2
704* returns the LCM of the head terms of a and b
705*/
706inline poly p_Lcm(const poly a, const poly b, const long lCompM, const ring r)
707{
708  poly m = p_ISet(1, r);
709
710  const int N = r->N;
711
712  for (int i = N; i>0; i--)
713  {
714    const int lExpA = p_GetExp (a, i, r);
715    const int lExpB = p_GetExp (b, i, r);
716
717    p_SetExp (m, i, si_max(lExpA, lExpB), r);
718  }
719
720  p_SetComp (m, lCompM, r);
721
722  p_Setm(m,r);
723
724#ifdef PDEBUG
725  p_Test(m,r);
726#endif
727
728  return(m);
729}
730
731/*4
732* creates the S-polynomial of p1 and p2
733* does not destroy p1 and p2
734*/
735poly sca_SPoly( const poly p1, const poly p2, const ring r )
736{
737  assume( rIsSCA(r) );
738
739  const long lCompP1 = p_GetComp(p1,r);
740  const long lCompP2 = p_GetComp(p2,r);
741
742  if ((lCompP1!=lCompP2) && (lCompP1!=0) && (lCompP2!=0))
743  {
744#ifdef PDEBUG
745    Print("sca_SPoly: different non-zero components!\n");
746#endif
747    return(NULL);
748  }
749
750  poly pL = p_Lcm(p1, p2, si_max(lCompP1, lCompP2), r);       // pL = lcm( lm(p1), lm(p2) )
751
752  poly m1 = p_ISet(1, r);
753  p_ExpVectorDiff(m1, pL, p1, r);                  // m1 = pL / lm(p1)
754  //p_SetComp(m1,0,r);
755  //p_Setm(m1,r);
756#ifdef PDEBUG
757  p_Test(m1,r);
758#endif
759
760
761  poly m2 = p_ISet(1, r);
762  p_ExpVectorDiff (m2, pL, p2, r);                  // m2 = pL / lm(p2)
763
764  //p_SetComp(m2,0,r);
765  //p_Setm(m2,r);
766#ifdef PDEBUG
767  p_Test(m2,r);
768#endif
769
770  p_Delete(&pL,r);
771
772  number C1  = n_Copy(p_GetCoeff(p1,r),r);      // C1 = lc(p1)
773  number C2  = n_Copy(p_GetCoeff(p2,r),r);      // C2 = lc(p2)
774
775  number C = nGcd(C1,C2,r);                     // C = gcd(C1, C2)
776
777  if (!n_IsOne(C, r))                              // if C != 1
778  {
779    C1=n_Div(C1, C, r);                              // C1 = C1 / C
780    C2=n_Div(C2, C, r);                              // C2 = C2 / C
781  }
782
783  n_Delete(&C,r); // destroy the number C
784
785  const int iSignSum = sca_Sign_mm_Mult_mm (m1, p1, r) + sca_Sign_mm_Mult_mm (m2, p2, r);
786  // zero if different signs
787
788  assume( (iSignSum*iSignSum == 0) || (iSignSum*iSignSum == 4) );
789
790  if( iSignSum != 0 ) // the same sign!
791    C2=n_Neg (C2, r);
792
793  p_SetCoeff(m1, C2, r);                           // lc(m1) = C2!!!
794  p_SetCoeff(m2, C1, r);                           // lc(m2) = C1!!!
795
796  poly tmp1 = nc_mm_Mult_pp (m1, pNext(p1), r);    // tmp1 = m1 * tail(p1),
797  p_Delete(&m1,r);  //  => n_Delete(&C1,r);
798
799  poly tmp2 = nc_mm_Mult_pp (m2, pNext(p2), r);    // tmp1 = m2 * tail(p2),
800  p_Delete(&m2,r);  //  => n_Delete(&C1,r);
801
802  poly spoly = p_Add_q (tmp1, tmp2, r); // spoly = spoly(lt(p1), lt(p2)) + m1 * tail(p1), delete tmp1,2
803
804  if (spoly!=NULL) pCleardenom (spoly); // r?
805//  if (spoly!=NULL) pContent (spoly); // r?
806
807#ifdef PDEBUG
808  p_Test (spoly, r);
809#endif
810
811  return(spoly);
812}
813
814
815
816
817/*2
818* reduction of p2 with p1
819* do not destroy p1, but p2
820* p1 divides p2 -> for use in NF algorithm
821*/
822poly sca_ReduceSpoly(const poly p1, poly p2, const ring r)
823{
824  assume( rIsSCA(r) );
825
826  assume( p1 != NULL );
827
828  const long lCompP1 = p_GetComp (p1, r);
829  const long lCompP2 = p_GetComp (p2, r);
830
831  if ((lCompP1!=lCompP2) && (lCompP1!=0) && (lCompP2!=0))
832  {
833#ifdef PDEBUG
834    Print("sca_ReduceSpoly: different non-zero components!");
835#endif
836    return(NULL);
837  }
838
839  poly m = p_ISet (1, r);
840  p_ExpVectorDiff (m, p2, p1, r);                      // m = lm(p2) / lm(p1)
841  //p_Setm(m,r);
842#ifdef PDEBUG
843  p_Test (m,r);
844#endif
845
846  number C1 = n_Copy( p_GetCoeff(p1, r), r);
847  number C2 = n_Copy( p_GetCoeff(p2, r), r);
848
849  /* GCD stuff */
850  number C = nGcd(C1, C2, r);
851
852  if (!n_IsOne(C, r))
853  {
854    C1 = n_Div(C1, C, r);
855    C2 = n_Div(C2, C, r);
856  }
857  n_Delete(&C,r);
858
859  const int iSign = sca_Sign_mm_Mult_mm( m, p1, r );
860
861  if(iSign == 1)
862    C2 = n_Neg(C2,r);
863
864  p_SetCoeff(m, C2, r);
865
866#ifdef PDEBUG
867  p_Test(m,r);
868#endif
869
870  p2 = p_LmDeleteAndNext( p2, r );
871
872  p2 = p_Mult_nn(p2, C1, r); // p2 !!!
873  n_Delete(&C1,r);
874
875  poly T = nc_mm_Mult_pp(m, pNext(p1), r);
876  p_Delete(&m, r);
877
878  p2 = p_Add_q(p2, T, r);
879
880  if ( p2!=NULL ) pContent(p2); // r?
881
882#ifdef PDEBUG
883  p_Test(p2,r);
884#endif
885
886  return(p2);
887}
888
889
890void addLObject(LObject& h, kStrategy& strat)
891{
892  if(h.IsNull()) return;
893
894  strat->initEcart(&h);
895  h.sev=0; // pGetShortExpVector(h.p);
896
897  // add h into S and L
898  int pos=posInS(strat, strat->sl, h.p, h.ecart);
899
900  if ( (pos <= strat->sl) && (pComparePolys(h.p, strat->S[pos])) )
901  {
902    if (TEST_OPT_PROT)
903      PrintS("d\n");
904  }
905  else
906  {
907    if (TEST_OPT_INTSTRATEGY)
908    {
909      pCleardenom(h.p);
910    }
911    else
912    {
913      pNorm(h.p);
914      pContent(h.p);
915    }
916
917    if ((strat->syzComp==0)||(!strat->homog))
918    {
919      h.p = redtailBba(h.p,pos-1,strat);
920
921      if (TEST_OPT_INTSTRATEGY)
922      {
923//        pCleardenom(h.p);
924        pContent(h.p);
925      }
926      else
927      {
928        pNorm(h.p);
929      }
930    }
931
932    if(h.IsNull()) return;
933
934    /* statistic */
935    if (TEST_OPT_PROT)
936    {
937      PrintS("s\n");
938    }
939
940    if (TEST_OPT_DEBUG)
941    {
942      PrintS("new s:");
943      wrp(h.p);
944      PrintLn();
945    }
946
947    enterpairs(h.p, strat->sl, h.ecart, 0, strat);
948
949    pos=0;
950
951    if (strat->sl!=-1) pos = posInS(strat, strat->sl, h.p, h.ecart);
952
953    strat->enterS(h, pos, strat, -1);
954
955    if (h.lcm!=NULL) pLmFree(h.lcm);
956  }
957
958
959}
960
961#ifndef NDEBUG
962
963// set it here if needed.
964#define MYTEST 0
965
966#else
967
968#define MYTEST 0
969
970#endif
971
972
973
974ideal sca_gr_bba(const ideal F, const ideal Q, const intvec *, const intvec *, kStrategy strat)
975{
976#if MYTEST
977   PrintS("<sca_gr_bba>\n");
978#endif
979
980  assume(rIsSCA(currRing));
981
982#ifndef NDEBUG
983  idTest(F);
984  idTest(Q);
985#endif
986
987  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
988  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
989
990  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
991  ideal tempQ = Q;
992
993  if(Q == currQuotient)
994    tempQ = currRing->nc->SCAQuotient();
995
996  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
997
998  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
999
1000  strat->homog = strat->homog && bIdHomog;
1001
1002
1003  assume( strat->homog == bIdHomog );
1004
1005
1006#if MYTEST
1007/*
1008  {
1009    Print("ideal F: \n");
1010    idPrint(F);
1011    Print("ideal tempF: \n");
1012    idPrint(F);
1013  }
1014*/
1015#endif
1016
1017
1018
1019
1020
1021  int srmax, lrmax;
1022  int olddeg, reduc;
1023  int red_result = 1;
1024//  int hilbeledeg = 1, minimcnt = 0;
1025  int hilbcount = 0;
1026
1027  initBuchMoraCrit(strat); // set Gebauer, honey, sugarCrit
1028
1029  nc_gr_initBba(tempF,strat); // set enterS, red, initEcart, initEcartPair
1030
1031  initBuchMoraPos(strat);
1032
1033
1034  // ?? set spSpolyShort, reduce ???
1035
1036  initBuchMora(tempF, tempQ, strat); // currRing->nc->SCAQuotient() instead of Q == squares!!!!!!!
1037
1038  strat->posInT=posInT110; // !!!
1039
1040  srmax = strat->sl;
1041  reduc = olddeg = lrmax = 0;
1042
1043
1044  /* compute------------------------------------------------------- */
1045  for(; strat->Ll >= 0;
1046#ifdef KDEBUG
1047    strat->P.lcm = NULL,
1048#endif
1049    kTest(strat)
1050    )
1051  {
1052    if (strat->Ll > lrmax) lrmax =strat->Ll;/*stat.*/
1053
1054    if (TEST_OPT_DEBUG) messageSets(strat);
1055
1056    if (strat->Ll== 0) strat->interpt=TRUE;
1057
1058    if (TEST_OPT_DEGBOUND
1059    && ((strat->honey
1060    && (strat->L[strat->Ll].ecart+pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
1061       || ((!strat->honey) && (pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))))
1062    {
1063      /*
1064      *stops computation if
1065      * 24 IN test and the degree +ecart of L[strat->Ll] is bigger then
1066      *a predefined number Kstd1_deg
1067      */
1068      while (strat->Ll >= 0) deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
1069      break;
1070    }
1071
1072    /* picks the last element from the lazyset L */
1073    strat->P = strat->L[strat->Ll];
1074    strat->Ll--;
1075
1076    //kTest(strat);
1077
1078    assume(pNext(strat->P.p) != strat->tail);
1079
1080//     if( pNext(strat->P.p) == strat->tail )
1081//     {
1082//       // deletes the int spoly and computes SPoly
1083//       pLmFree(strat->P.p); // ???
1084//       strat->P.p = sca_SPoly(strat->P.p1, strat->P.p2, currRing);
1085//     }
1086
1087    if(strat->P.IsNull()) continue;
1088
1089//     poly save = NULL;
1090//
1091//     if(pNext(strat->P.p) != NULL)
1092//       save = p_Copy(strat->P.p, currRing);
1093
1094    strat->initEcart(&strat->P); // remove it?
1095
1096    if (TEST_OPT_PROT)
1097      message((strat->honey ? strat->P.ecart : 0) + strat->P.pFDeg(), &olddeg,&reduc,strat, red_result);
1098
1099    // reduction of the element chosen from L wrt S
1100    strat->red(&strat->P,strat);
1101
1102    if(strat->P.IsNull()) continue;
1103
1104    addLObject(strat->P, strat);
1105
1106    if (strat->sl > srmax) srmax = strat->sl;
1107
1108    const poly save = strat->P.p;
1109
1110#ifdef PDEBUG
1111      p_Test(save, currRing);
1112#endif
1113    assume( save != NULL );
1114
1115    // SCA Specials:
1116
1117    {
1118      const poly pNext = pNext(save);
1119
1120      if( pNext != NULL )
1121      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
1122      if( p_GetExp(save, i, currRing) != 0 )
1123      {
1124        assume(p_GetExp(save, i, currRing) == 1);
1125
1126        const poly tt = sca_pp_Mult_xi_pp(i, pNext, currRing);
1127
1128#ifdef PDEBUG
1129        p_Test(tt, currRing);
1130#endif
1131
1132        if( tt == NULL) continue;
1133
1134        LObject h(tt); // h = x_i * P
1135
1136        if (TEST_OPT_INTSTRATEGY)
1137        {
1138//           h.pCleardenom(); // also does a pContent
1139          pContent(h.p);
1140        }
1141        else
1142        {
1143          h.pNorm();
1144        }
1145
1146        strat->initEcart(&h);
1147
1148
1149//         if (pOrdSgn==-1)
1150//         {
1151//           cancelunit(&h);  // tries to cancel a unit
1152//           deleteHC(&h, strat);
1153//         }
1154
1155//         if(h.IsNull()) continue;
1156
1157//         if (TEST_OPT_PROT)
1158//           message((strat->honey ? h.ecart : 0) + h.pFDeg(), &olddeg, &reduc, strat, red_result);
1159
1160//         strat->red(&h, strat); // wrt S
1161//         if(h.IsNull()) continue;
1162
1163//         poly save = p_Copy(h.p, currRing);
1164
1165        int pos;
1166
1167        if (strat->Ll==-1)
1168          pos =0;
1169        else
1170          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1171
1172        h.sev = pGetShortExpVector(h.p);
1173        enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
1174
1175  //       h.p = save;
1176  //       addLObject(h, strat);
1177
1178//         if (strat->sl > srmax) srmax = strat->sl;
1179      }
1180
1181      // p_Delete( &save, currRing );
1182    }
1183
1184
1185  } // for(;;)
1186
1187
1188  if (TEST_OPT_DEBUG) messageSets(strat);
1189
1190  /* release temp data-------------------------------- */
1191  exitBuchMora(strat);
1192
1193  if (TEST_OPT_WEIGHTM)
1194  {
1195    pFDeg=pFDegOld;
1196    pLDeg=pLDegOld;
1197    if (ecartWeights)
1198    {
1199      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(int));
1200      ecartWeights=NULL;
1201    }
1202  }
1203
1204  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
1205
1206  if (tempQ!=NULL) updateResult(strat->Shdl,tempQ,strat);
1207
1208  id_Delete(&tempF, currRing);
1209
1210
1211  /* complete reduction of the standard basis--------- */
1212  if (TEST_OPT_REDSB){
1213//     completeReduce(strat); // ???
1214
1215    ideal I = strat->Shdl;
1216    ideal erg = kInterRed(I,tempQ);
1217    assume(I!=erg);
1218    id_Delete(&I, currRing);
1219    strat->Shdl = erg;
1220  }
1221
1222
1223#if MYTEST
1224   PrintS("</sca_gr_bba>\n");
1225#endif
1226
1227  return (strat->Shdl);
1228}
1229
1230
1231// should be used only inside nc_SetupQuotient!
1232// Check whether this our case:
1233//  1. rG is  a commutative polynomial ring \otimes anticommutative algebra
1234//  2. factor ideal rGR->qideal contains squares of all alternating variables.
1235//
1236// if yes, make rGR a super-commutative algebra!
1237// NOTE: Factors of SuperCommutative Algebras are supported this way!
1238bool sca_SetupQuotient(ring rGR, const ring rG)
1239{
1240//   return false; // test Plural
1241
1242  //////////////////////////////////////////////////////////////////////////
1243  // checks...
1244  //////////////////////////////////////////////////////////////////////////
1245  assume(rGR != NULL);
1246  assume(rG  != NULL);
1247
1248  assume(rIsPluralRing(rG));
1249
1250  const int N = rG->N;
1251
1252  if(N < 2)
1253    return false;
1254
1255#if MYTEST
1256   PrintS("sca_SetupQuotient(rGR, rG)");
1257#endif
1258
1259
1260//  if( (ncRingType(rG) != nc_skew) || (ncRingType(rG) != nc_comm) )
1261//    return false;
1262
1263  if(rGR->qideal == NULL) // there will be a factor!
1264    return false;
1265
1266  if(rG->qideal != NULL) // we cannot change from factor to factor!
1267    return false;
1268
1269
1270  int iAltVarEnd = -1;
1271  int iAltVarStart   = N+1;
1272
1273  const ring rBase = rG->nc->basering;
1274  const matrix C   = rG->nc->C; // live in rBase!
1275
1276  for(int i = 1; i < N; i++)
1277  {
1278    for(int j = i + 1; j <= N; j++)
1279    {
1280      assume(MATELEM(C,i,j) != NULL); // after CallPlural!
1281      number c = p_GetCoeff(MATELEM(C,i,j), rBase);
1282
1283      if( n_IsMOne(c, rBase) )
1284      {
1285        if( i < iAltVarStart)
1286          iAltVarStart = i;
1287
1288        if( j > iAltVarEnd)
1289          iAltVarEnd = j;
1290      } else
1291      {
1292        if( !n_IsOne(c, rBase) )
1293        {
1294#if MYTEST
1295           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1296#endif
1297          return false;
1298        }
1299      }
1300    }
1301  }
1302
1303  if( (iAltVarEnd == -1) || (iAltVarStart == (N+1)) )
1304    return false; // either no alternating varables, or a single one => we are in commutative case!
1305
1306  for(int i = 1; i < N; i++)
1307  {
1308    for(int j = i + 1; j <= N; j++)
1309    {
1310      assume(MATELEM(C,i,j) != NULL); // after CallPlural!
1311      number c = p_GetCoeff(MATELEM(C,i,j), rBase);
1312
1313      if( (iAltVarStart <= i) && (j <= iAltVarEnd) ) // S <= i < j <= E
1314      { // anticommutative part
1315        if( !n_IsMOne(c, rBase) )
1316        {
1317#if MYTEST
1318           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1319#endif
1320          return false;
1321        }
1322      } else
1323      { // should commute
1324        if( !n_IsOne(c, rBase) )
1325        {
1326#if MYTEST
1327           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1328#endif
1329          return false;
1330        }
1331      }
1332    }
1333  }
1334
1335  assume( 1            <= iAltVarStart );
1336  assume( iAltVarStart < iAltVarEnd   );
1337  assume( iAltVarEnd   <= N            );
1338
1339  bool bWeChangeRing = false;
1340  // for sanity
1341  ring rSaveRing = currRing;
1342
1343  if(rSaveRing != rG)
1344  {
1345    rChangeCurrRing(rG);
1346    bWeChangeRing = true;
1347  }
1348
1349  assume(rGR->qideal != NULL);
1350  assume(rG->qideal == NULL);
1351
1352  const ideal idQuotient = rGR->qideal;
1353
1354  // check for
1355  // y_{iAltVarStart}^2, y_{iAltVarStart+1}^2, \ldots, y_{iAltVarEnd}^2  (iAltVarEnd > iAltVarStart)
1356  // to be within quotient ideal.
1357
1358  bool bSCA = true;
1359
1360  for ( int i = iAltVarStart; (i <= iAltVarEnd) && bSCA; i++ )
1361  {
1362    poly square = p_ISet(1, rSaveRing);
1363    p_SetExp(square, i, 2, rSaveRing); // square = var(i)^2.
1364    p_Setm(square, rSaveRing);
1365
1366    // square = NF( var(i)^2 | Q )
1367    // NOTE: rSaveRing == currRing now!
1368    // NOTE: there is no better way to check this in general!
1369    square = kNF(idQuotient, NULL, square, 0, 0);
1370
1371    if( square != NULL ) // var(i)^2 is not in Q?
1372    {
1373      p_Delete(&square, rSaveRing);
1374      bSCA = false;
1375    }
1376  }
1377
1378
1379  if (bWeChangeRing)
1380  {
1381    rChangeCurrRing(rSaveRing);
1382  }
1383
1384  if(!bSCA) return false;
1385
1386
1387
1388#if MYTEST
1389  Print("AltVars: [%d, %d]\n", iAltVarStart, iAltVarEnd);
1390#endif
1391
1392
1393  //////////////////////////////////////////////////////////////////////////
1394  // ok... here we go. let's setup it!!!
1395  //////////////////////////////////////////////////////////////////////////
1396  ideal tempQ = id_KillSquares(idQuotient, iAltVarStart, iAltVarEnd, rG); // in rG!!!
1397
1398  idSkipZeroes( tempQ );
1399
1400  if( idIs0(tempQ) )
1401    rGR->nc->SCAQuotient() = NULL;
1402  else
1403    rGR->nc->SCAQuotient() = idrMoveR(tempQ, rG, rGR); // deletes tempQ!
1404
1405  ncRingType( rGR, nc_exterior );
1406
1407  scaFirstAltVar( rGR, iAltVarStart );
1408  scaLastAltVar( rGR, iAltVarEnd );
1409
1410
1411  sca_p_ProcsSet(rGR, rGR->p_Procs);
1412
1413
1414  return true;
1415}
1416
1417
1418bool sca_ForceCommutative(ring rGR, int b, int e)
1419{
1420  assume(rGR != NULL);
1421  assume(rIsPluralRing(rGR));
1422  assume(!rIsSCA(rGR));
1423 
1424  const int N = rGR->N;
1425
1426  ring rSaveRing = currRing;
1427
1428  if(rSaveRing != rGR)
1429    rChangeCurrRing(rGR);
1430
1431  const ideal idQuotient = rGR->qideal;
1432
1433 
1434  ideal tempQ = idQuotient;
1435
1436  if( b <= N && e >= 1 )
1437    tempQ = id_KillSquares(idQuotient, b, e, rGR); 
1438
1439  idSkipZeroes( tempQ );
1440
1441  if( idIs0(tempQ) )
1442    rGR->nc->SCAQuotient() = NULL;
1443  else
1444    rGR->nc->SCAQuotient() = tempQ;
1445 
1446  ncRingType( rGR, nc_exterior );
1447
1448  scaFirstAltVar( rGR, b );
1449  scaLastAltVar( rGR, e );
1450
1451
1452  sca_p_ProcsSet(rGR, rGR->p_Procs);
1453
1454  if(rSaveRing != rGR)
1455    rChangeCurrRing(rSaveRing);
1456 
1457  return true;
1458 
1459}
1460
1461// return x_i * pPoly; preserve pPoly.
1462poly sca_pp_Mult_xi_pp(unsigned int i, const poly pPoly, const ring rRing)
1463{
1464  assume(1 <= i);
1465  assume(i <= (unsigned int)rRing->N);
1466
1467  if(rIsSCA(rRing))
1468    return sca_xi_Mult_pp(i, pPoly, rRing);
1469
1470
1471
1472  poly xi =  p_ISet(1, rRing);
1473  p_SetExp(xi, i, 1, rRing);
1474  p_Setm(xi, rRing);
1475
1476  poly pResult = pp_Mult_qq(xi, pPoly, rRing);
1477
1478  p_Delete( &xi, rRing);
1479
1480  return pResult;
1481
1482}
1483
1484
1485///////////////////////////////////////////////////////////////////////////////////////
1486//************* SCA BBA *************************************************************//
1487///////////////////////////////////////////////////////////////////////////////////////
1488
1489// Under development!!!
1490ideal sca_bba (const ideal F, const ideal Q, const intvec *w, const intvec * /*hilb*/, kStrategy strat)
1491{
1492  assume(rIsSCA(currRing));
1493
1494  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
1495  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
1496
1497  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
1498
1499  ideal tempQ = Q;
1500
1501  if(Q == currQuotient)
1502    tempQ = currRing->nc->SCAQuotient();
1503
1504  // Q or tempQ will not be used below :(((
1505
1506  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
1507
1508  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
1509
1510  strat->homog = strat->homog && bIdHomog;
1511
1512#ifdef PDEBUG
1513  assume( strat->homog == bIdHomog );
1514#endif /*PDEBUG*/
1515
1516
1517//    PrintS("<sca>\n");
1518
1519  om_Opts.MinTrack = 5; // ???
1520
1521  int   srmax, lrmax, red_result = 1;
1522  int   olddeg, reduc;
1523
1524//  int hilbeledeg = 1, minimcnt = 0;
1525  int hilbcount = 0;
1526
1527  BOOLEAN withT = FALSE;
1528
1529  initBuchMoraCrit(strat); // sets Gebauer, honey, sugarCrit // sca - ok???
1530  initBuchMoraPos(strat); // sets strat->posInL, strat->posInT // check!! (Plural's: )
1531
1532//   initHilbCrit(F, Q, &hilb, strat);
1533
1534//  nc_gr_initBba(F,strat);
1535  initBba(tempF, strat); // set enterS, red, initEcart, initEcartPair
1536
1537  /*set enterS, spSpolyShort, reduce, red, initEcart, initEcartPair*/
1538  // ?? set spSpolyShort, reduce ???
1539  initBuchMora(tempF, NULL, strat); // Q = squares!!!
1540
1541//   if (strat->minim>0) strat->M = idInit(IDELEMS(F),F->rank);
1542
1543  srmax = strat->sl;
1544  reduc = olddeg = lrmax = 0;
1545
1546#define NO_BUCKETS
1547
1548#ifndef NO_BUCKETS
1549  if (!TEST_OPT_NOT_BUCKETS)
1550    strat->use_buckets = 1;
1551#endif
1552
1553  // redtailBBa against T for inhomogenous input
1554  if (!K_TEST_OPT_OLDSTD)
1555    withT = ! strat->homog;
1556
1557  // strat->posInT = posInT_pLength;
1558  kTest_TS(strat);
1559
1560#undef HAVE_TAIL_RING
1561
1562#ifdef HAVE_TAIL_RING
1563  kStratInitChangeTailRing(strat);
1564#endif
1565
1566  ///////////////////////////////////////////////////////////////
1567  // SCA:
1568
1569  /* compute------------------------------------------------------- */
1570  while (strat->Ll >= 0)
1571  {
1572    if (strat->Ll > lrmax) lrmax =strat->Ll;/*stat.*/
1573
1574#ifdef KDEBUG
1575//     loop_count++;
1576    if (TEST_OPT_DEBUG) messageSets(strat);
1577#endif
1578
1579    if (strat->Ll== 0) strat->interpt=TRUE;
1580
1581    if (TEST_OPT_DEGBOUND
1582        && ((strat->honey && (strat->L[strat->Ll].ecart+pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
1583            || ((!strat->honey) && (pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))))
1584    {
1585      /*
1586       *stops computation if
1587       * 24 IN test and the degree +ecart of L[strat->Ll] is bigger then
1588       *a predefined number Kstd1_deg
1589       */
1590      while ((strat->Ll >= 0)
1591  && (strat->L[strat->Ll].p1!=NULL) && (strat->L[strat->Ll].p2!=NULL)
1592        && ((strat->honey && (strat->L[strat->Ll].ecart+pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
1593            || ((!strat->honey) && (pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg)))
1594  )
1595        deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
1596      if (strat->Ll<0) break;
1597      else strat->noClearS=TRUE;
1598    }
1599
1600    /* picks the last element from the lazyset L */
1601    strat->P = strat->L[strat->Ll];
1602    strat->Ll--;
1603
1604
1605    assume(pNext(strat->P.p) != strat->tail);
1606
1607/*    if (pNext(strat->P.p) == strat->tail)
1608    {
1609      // deletes the short spoly
1610      pLmFree(strat->P.p);
1611      strat->P.p = NULL;
1612      poly m1 = NULL, m2 = NULL;
1613
1614      // check that spoly creation is ok
1615      while (strat->tailRing != currRing &&
1616             !kCheckSpolyCreation(&(strat->P), strat, m1, m2))
1617      {
1618        assume(m1 == NULL && m2 == NULL);
1619        // if not, change to a ring where exponents are at least
1620        // large enough
1621        kStratChangeTailRing(strat);
1622      }
1623
1624#ifdef PDEBUG
1625      Print("ksCreateSpoly!#?");
1626#endif
1627
1628      // create the real one
1629      ksCreateSpoly(&(strat->P), NULL, strat->use_buckets,
1630                    strat->tailRing, m1, m2, strat->R); //?????????
1631    }
1632    else */
1633
1634    if (strat->P.p1 == NULL)
1635    {
1636//       if (strat->minim > 0)
1637//         strat->P.p2=p_Copy(strat->P.p, currRing, strat->tailRing);
1638
1639
1640      // for input polys, prepare reduction
1641      strat->P.PrepareRed(strat->use_buckets);
1642    }
1643
1644    if (TEST_OPT_PROT)
1645      message((strat->honey ? strat->P.ecart : 0) + strat->P.pFDeg(),
1646              &olddeg,&reduc,strat, red_result);
1647
1648    /* reduction of the element choosen from L */
1649    red_result = strat->red(&strat->P,strat);
1650
1651
1652    // reduction to non-zero new poly
1653    if (red_result == 1)
1654    {
1655      /* statistic */
1656      if (TEST_OPT_PROT) PrintS("s");
1657
1658      // get the polynomial (canonicalize bucket, make sure P.p is set)
1659      strat->P.GetP(strat->lmBin);
1660
1661      int pos = posInS(strat,strat->sl,strat->P.p,strat->P.ecart);
1662
1663      // reduce the tail and normalize poly
1664      if (TEST_OPT_INTSTRATEGY)
1665      {
1666        strat->P.pCleardenom();
1667        if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1668        {
1669          strat->P.p = redtailBba(&(strat->P),pos-1,strat, withT); // !!!
1670          strat->P.pCleardenom();
1671
1672        }
1673      }
1674      else
1675      {
1676
1677        strat->P.pNorm();
1678        if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1679          strat->P.p = redtailBba(&(strat->P),pos-1,strat, withT);
1680      }
1681
1682
1683#ifdef KDEBUG
1684      if (TEST_OPT_DEBUG){PrintS(" ns:");strat->P.wrp();PrintLn();}
1685#endif
1686
1687//       // min_std stuff
1688//       if ((strat->P.p1==NULL) && (strat->minim>0))
1689//       {
1690//         if (strat->minim==1)
1691//         {
1692//           strat->M->m[minimcnt]=p_Copy(strat->P.p,currRing,strat->tailRing);
1693//           p_Delete(&strat->P.p2, currRing, strat->tailRing);
1694//         }
1695//         else
1696//         {
1697//           strat->M->m[minimcnt]=strat->P.p2;
1698//           strat->P.p2=NULL;
1699//         }
1700//         if (strat->tailRing!=currRing && pNext(strat->M->m[minimcnt])!=NULL)
1701//           pNext(strat->M->m[minimcnt])
1702//             = strat->p_shallow_copy_delete(pNext(strat->M->m[minimcnt]),
1703//                                            strat->tailRing, currRing,
1704//                                            currRing->PolyBin);
1705//         minimcnt++;
1706//       }
1707
1708      // enter into S, L, and T
1709      if(withT)
1710        enterT(strat->P, strat);
1711
1712      // L
1713      enterpairs(strat->P.p,strat->sl,strat->P.ecart,pos,strat, strat->tl);
1714
1715      // posInS only depends on the leading term
1716      strat->enterS(strat->P, pos, strat, strat->tl);
1717
1718//       if (hilb!=NULL) khCheck(Q,w,hilb,hilbeledeg,hilbcount,strat);
1719
1720//      Print("[%d]",hilbeledeg);
1721      if (strat->P.lcm!=NULL) pLmFree(strat->P.lcm);
1722
1723      if (strat->sl>srmax) srmax = strat->sl;
1724
1725      // //////////////////////////////////////////////////////////
1726      // SCA:
1727      const poly pSave = strat->P.p;
1728      const poly pNext = pNext(pSave);
1729
1730//       if(0)
1731      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
1732      if( p_GetExp(pSave, i, currRing) != 0 )
1733      {
1734        assume(p_GetExp(pSave, i, currRing) == 1);
1735
1736        const poly pNew = sca_pp_Mult_xi_pp(i, pNext, currRing);
1737
1738#ifdef PDEBUG
1739        p_Test(pNew, currRing);
1740#endif
1741
1742        if( pNew == NULL) continue;
1743
1744        LObject h(pNew); // h = x_i * strat->P
1745
1746        if (TEST_OPT_INTSTRATEGY)
1747        {
1748//           h.pCleardenom(); // also does a pContent
1749          pContent(h.p);
1750        }
1751        else
1752        {
1753          h.pNorm();
1754        }
1755
1756        strat->initEcart(&h);
1757        h.sev = pGetShortExpVector(h.p);
1758
1759        int pos;
1760
1761        if (strat->Ll==-1)
1762          pos =0;
1763        else
1764          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1765
1766        enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
1767
1768
1769/*
1770        h.sev = pGetShortExpVector(h.p);
1771        strat->initEcart(&h);
1772
1773        h.PrepareRed(strat->use_buckets);
1774
1775        // reduction of the element choosen from L(?)
1776        red_result = strat->red(&h,strat);
1777
1778        // reduction to non-zero new poly
1779        if (red_result != 1) continue;
1780
1781
1782        int pos = posInS(strat,strat->sl,h.p,h.ecart);
1783
1784        // reduce the tail and normalize poly
1785        if (TEST_OPT_INTSTRATEGY)
1786        {
1787          h.pCleardenom();
1788          if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1789          {
1790            h.p = redtailBba(&(h),pos-1,strat, withT); // !!!
1791            h.pCleardenom();
1792          }
1793        }
1794        else
1795        {
1796
1797          h.pNorm();
1798          if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1799            h.p = redtailBba(&(h),pos-1,strat, withT);
1800        }
1801
1802
1803  #ifdef KDEBUG
1804        if (TEST_OPT_DEBUG){PrintS(" N:");h.wrp();PrintLn();}
1805  #endif
1806
1807//         h.PrepareRed(strat->use_buckets); // ???
1808
1809        h.sev = pGetShortExpVector(h.p);
1810        strat->initEcart(&h);
1811
1812        if (strat->Ll==-1)
1813          pos = 0;
1814        else
1815          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1816
1817         enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);*/
1818
1819      } // for all x_i \in Ann(lm(P))
1820
1821
1822
1823
1824
1825    } // if red(P) != NULL
1826
1827//     else if (strat->P.p1 == NULL && strat->minim > 0)
1828//     {
1829//       p_Delete(&strat->P.p2, currRing, strat->tailRing);
1830//     }
1831
1832// #ifdef KDEBUG
1833    memset(&(strat->P), 0, sizeof(strat->P));
1834// #endif
1835
1836    kTest_TS(strat);
1837
1838//     Print("\n$\n");
1839
1840  }
1841
1842#ifdef KDEBUG
1843  if (TEST_OPT_DEBUG) messageSets(strat);
1844#endif
1845
1846  /* complete reduction of the standard basis--------- */
1847  if (TEST_OPT_REDSB) completeReduce(strat);
1848
1849  /* release temp data-------------------------------- */
1850  id_Delete(&tempF, currRing);
1851
1852  exitBuchMora(strat);
1853
1854  if (TEST_OPT_WEIGHTM)
1855  {
1856    pRestoreDegProcs(pFDegOld, pLDegOld);
1857    if (ecartWeights)
1858    {
1859      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(short));
1860      ecartWeights=NULL;
1861    }
1862  }
1863
1864  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
1865
1866//   if (Q!=NULL) updateResult(strat->Shdl,Q,strat);
1867
1868//    PrintS("</sca>\n");
1869
1870  return (strat->Shdl);
1871}
1872
1873
1874// //////////////////////////////////////////////////////////////////////////////
1875// sca mora...
1876
1877// returns TRUE if mora should use buckets, false otherwise
1878static BOOLEAN kMoraUseBucket(kStrategy strat)
1879{
1880#ifdef MORA_USE_BUCKETS
1881  if (TEST_OPT_NOT_BUCKETS)
1882    return FALSE;
1883  if (strat->red == redFirst)
1884  {
1885#ifdef NO_LDEG
1886    if (!strat->syzComp)
1887      return TRUE;
1888#else
1889    if ((strat->homog || strat->honey) && !strat->syzComp)
1890      return TRUE;
1891#endif
1892  }
1893  else
1894  {
1895    assume(strat->red == redEcart);
1896    if (strat->honey && !strat->syzComp)
1897      return TRUE;
1898  }
1899#endif
1900  return FALSE;
1901}
1902
1903
1904#ifdef HAVE_ASSUME
1905static int sca_mora_count = 0;
1906static int sca_mora_loop_count;
1907#endif
1908
1909// ideal sca_mora (ideal F, ideal Q, intvec *w, intvec *, kStrategy strat)
1910ideal sca_mora(const ideal F, const ideal Q, const intvec *w, const intvec *, kStrategy strat)
1911{
1912  assume(rIsSCA(currRing));
1913
1914  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
1915  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
1916
1917  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
1918
1919  ideal tempQ = Q;
1920
1921  if(Q == currQuotient)
1922    tempQ = currRing->nc->SCAQuotient();
1923
1924  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
1925
1926  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
1927
1928  strat->homog = strat->homog && bIdHomog;
1929
1930#ifdef PDEBUG
1931  assume( strat->homog == bIdHomog );
1932#endif /*PDEBUG*/
1933
1934#ifdef HAVE_ASSUME
1935  sca_mora_count++;
1936  sca_mora_loop_count = 0;
1937#endif
1938
1939#ifdef KDEBUG
1940  om_Opts.MinTrack = 5;
1941#endif
1942  int srmax;
1943  int lrmax = 0;
1944  int olddeg = 0;
1945  int reduc = 0;
1946  int red_result = 1;
1947//  int hilbeledeg=1;
1948  int hilbcount=0;
1949
1950  strat->update = TRUE;
1951  /*- setting global variables ------------------- -*/
1952  initBuchMoraCrit(strat);
1953//   initHilbCrit(F,NULL,&hilb,strat); // no Q!
1954  initMora(tempF, strat);
1955  initBuchMoraPos(strat);
1956  /*Shdl=*/initBuchMora(tempF, tempQ, strat); // temp Q, F!
1957//   if (TEST_OPT_FASTHC) missingAxis(&strat->lastAxis,strat);
1958  /*updateS in initBuchMora has Hecketest
1959  * and could have put strat->kHEdgdeFound FALSE*/
1960#if 0
1961  if (ppNoether!=NULL)
1962  {
1963    strat->kHEdgeFound = TRUE;
1964  }
1965  if (strat->kHEdgeFound && strat->update)
1966  {
1967    firstUpdate(strat);
1968    updateLHC(strat);
1969    reorderL(strat);
1970  }
1971  if (TEST_OPT_FASTHC && (strat->lastAxis) && strat->posInLOldFlag)
1972  {
1973    strat->posInLOld = strat->posInL;
1974    strat->posInLOldFlag = FALSE;
1975    strat->posInL = posInL10;
1976    updateL(strat);
1977    reorderL(strat);
1978  }
1979#endif
1980
1981  srmax = strat->sl;
1982  kTest_TS(strat);
1983
1984  strat->use_buckets = kMoraUseBucket(strat);
1985  /*- compute-------------------------------------------*/
1986
1987#undef HAVE_TAIL_RING
1988
1989#ifdef HAVE_TAIL_RING
1990//  if (strat->homog && strat->red == redFirst)
1991//     kStratInitChangeTailRing(strat);
1992#endif
1993
1994
1995  while (strat->Ll >= 0)
1996  {
1997#ifdef HAVE_ASSUME
1998    sca_mora_loop_count++;
1999#endif
2000    if (lrmax< strat->Ll) lrmax=strat->Ll; /*stat*/
2001    //test_int_std(strat->kIdeal);
2002    if (TEST_OPT_DEBUG) messageSets(strat);
2003    if (TEST_OPT_DEGBOUND
2004    && (strat->L[strat->Ll].ecart+strat->L[strat->Ll].GetpFDeg()> Kstd1_deg))
2005    {
2006      /*
2007      * stops computation if
2008      * - 24 (degBound)
2009      *   && upper degree is bigger than Kstd1_deg
2010      */
2011      while ((strat->Ll >= 0)
2012        && (strat->L[strat->Ll].p1!=NULL) && (strat->L[strat->Ll].p2!=NULL)
2013        && (strat->L[strat->Ll].ecart+strat->L[strat->Ll].GetpFDeg()> Kstd1_deg)
2014      )
2015      {
2016        deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
2017        //if (TEST_OPT_PROT)
2018        //{
2019        //   PrintS("D"); mflush();
2020        //}
2021      }
2022      if (strat->Ll<0) break;
2023      else strat->noClearS=TRUE;
2024    }
2025    strat->P = strat->L[strat->Ll];/*- picks the last element from the lazyset L -*/
2026    if (strat->Ll==0) strat->interpt=TRUE;
2027    strat->Ll--;
2028
2029    // create the real Spoly
2030    assume(pNext(strat->P.p) != strat->tail);
2031
2032    if (strat->P.p1 == NULL)
2033    {
2034      // for input polys, prepare reduction (buckets !)
2035      strat->P.SetLength(strat->length_pLength);
2036      strat->P.PrepareRed(strat->use_buckets);
2037    }
2038
2039    if (!strat->P.IsNull())
2040    {
2041      // might be NULL from noether !!!
2042      if (TEST_OPT_PROT)
2043        message(strat->P.ecart+strat->P.GetpFDeg(),&olddeg,&reduc,strat, red_result);
2044      // reduce
2045      red_result = strat->red(&strat->P,strat);
2046    }
2047
2048    if (! strat->P.IsNull())
2049    {
2050      strat->P.GetP();
2051      // statistics
2052      if (TEST_OPT_PROT) PrintS("s");
2053      // normalization
2054      if (!TEST_OPT_INTSTRATEGY)
2055        strat->P.pNorm();
2056      // tailreduction
2057      strat->P.p = redtail(&(strat->P),strat->sl,strat);
2058      // set ecart -- might have changed because of tail reductions
2059      if ((!strat->noTailReduction) && (!strat->honey))
2060        strat->initEcart(&strat->P);
2061      // cancel unit
2062      cancelunit(&strat->P);
2063      // for char 0, clear denominators
2064      if (TEST_OPT_INTSTRATEGY)
2065        strat->P.pCleardenom();
2066
2067      // put in T
2068      enterT(strat->P,strat);
2069      // build new pairs
2070      enterpairs(strat->P.p,strat->sl,strat->P.ecart,0,strat, strat->tl);
2071      // put in S
2072      strat->enterS(strat->P,
2073                    posInS(strat,strat->sl,strat->P.p, strat->P.ecart),
2074                    strat, strat->tl);
2075
2076
2077      // clear strat->P
2078      if (strat->P.lcm!=NULL) pLmFree(strat->P.lcm);
2079      strat->P.lcm=NULL;
2080
2081      if (strat->sl>srmax) srmax = strat->sl; /*stat.*/
2082      if (strat->Ll>lrmax) lrmax = strat->Ll;
2083
2084
2085
2086      // //////////////////////////////////////////////////////////
2087      // SCA:
2088      const poly pSave = strat->P.p;
2089      const poly pNext = pNext(pSave);
2090
2091      if(pNext != NULL)
2092      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
2093      if( p_GetExp(pSave, i, currRing) != 0 )
2094      {
2095
2096        assume(p_GetExp(pSave, i, currRing) == 1);
2097
2098        const poly pNew = sca_pp_Mult_xi_pp(i, pNext, currRing);
2099
2100#ifdef PDEBUG
2101        p_Test(pNew, currRing);
2102#endif
2103
2104        if( pNew == NULL) continue;
2105
2106        LObject h(pNew); // h = x_i * strat->P
2107
2108        if (TEST_OPT_INTSTRATEGY)
2109           h.pCleardenom(); // also does a pContent
2110        else
2111          h.pNorm();
2112
2113        strat->initEcart(&h);
2114        h.sev = pGetShortExpVector(h.p);
2115
2116        int pos;
2117
2118        if (strat->Ll==-1)
2119          pos = 0;
2120        else
2121          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
2122
2123        enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
2124
2125        if (strat->Ll>lrmax) lrmax = strat->Ll;
2126      }
2127
2128#ifdef KDEBUG
2129      // make sure kTest_TS does not complain about strat->P
2130      memset(&strat->P,0,sizeof(strat->P));
2131#endif
2132    }
2133#if 0
2134    if (strat->kHEdgeFound)
2135    {
2136      if ((BTEST1(27))
2137      || ((TEST_OPT_MULTBOUND) && (scMult0Int((strat->Shdl)) < mu)))
2138      {
2139        // obachman: is this still used ???
2140        /*
2141        * stops computation if strat->kHEdgeFound and
2142        * - 27 (finiteDeterminacyTest)
2143        * or
2144        * - 23
2145        *   (multBound)
2146        *   && multiplicity of the ideal is smaller then a predefined number mu
2147        */
2148        while (strat->Ll >= 0) deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
2149      }
2150    }
2151#endif
2152    kTest_TS(strat);
2153  }
2154  /*- complete reduction of the standard basis------------------------ -*/
2155  if (TEST_OPT_REDSB) completeReduce(strat);
2156  /*- release temp data------------------------------- -*/
2157  exitBuchMora(strat);
2158  /*- polynomials used for HECKE: HC, noether -*/
2159  if (BTEST1(27))
2160  {
2161    if (strat->kHEdge!=NULL)
2162      Kstd1_mu=pFDeg(strat->kHEdge,currRing);
2163    else
2164      Kstd1_mu=-1;
2165  }
2166  pDelete(&strat->kHEdge);
2167  strat->update = TRUE; //???
2168  strat->lastAxis = 0; //???
2169  pDelete(&strat->kNoether);
2170  omFreeSize((ADDRESS)strat->NotUsedAxis,(pVariables+1)*sizeof(BOOLEAN));
2171  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
2172  if (TEST_OPT_WEIGHTM)
2173  {
2174    pRestoreDegProcs(pFDegOld, pLDegOld);
2175    if (ecartWeights)
2176    {
2177      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(short));
2178      ecartWeights=NULL;
2179    }
2180  }
2181  if (tempQ!=NULL) updateResult(strat->Shdl,tempQ,strat);
2182  idTest(strat->Shdl);
2183
2184  id_Delete( &tempF, currRing);
2185
2186  return (strat->Shdl);
2187}
2188
2189
2190
2191
2192
2193
2194void sca_p_ProcsSet(ring rGR, p_Procs_s* p_Procs)
2195{
2196
2197  // "commutative" procedures:
2198  rGR->p_Procs->p_Mult_mm     = sca_p_Mult_mm;
2199  rGR->p_Procs->pp_Mult_mm    = sca_pp_Mult_mm;
2200
2201  p_Procs->p_Mult_mm          = sca_p_Mult_mm;
2202  p_Procs->pp_Mult_mm         = sca_pp_Mult_mm;
2203
2204  // non-commutaitve
2205  rGR->nc->p_Procs.mm_Mult_p   = sca_mm_Mult_p;
2206  rGR->nc->p_Procs.mm_Mult_pp  = sca_mm_Mult_pp;
2207
2208
2209  if (rGR->OrdSgn==-1)
2210  {
2211#ifdef PDEBUG
2212//           Print("Local case => GB == mora!\n");
2213#endif
2214    rGR->nc->p_Procs.GB          = sca_mora; // local ordering => Mora, otherwise - Buchberger!
2215  }
2216  else
2217  {
2218#ifdef PDEBUG
2219//           Print("Global case => GB == bba!\n");
2220#endif
2221    rGR->nc->p_Procs.GB          = sca_gr_bba; // sca_bba?
2222  }
2223
2224
2225//   rGR->nc->p_Procs.GlobalGB    = sca_gr_bba;
2226//   rGR->nc->p_Procs.LocalGB     = sca_mora;
2227
2228
2229//   rGR->nc->p_Procs.SPoly         = sca_SPoly;
2230//   rGR->nc->p_Procs.ReduceSPoly   = sca_ReduceSpoly;
2231
2232#if 0
2233
2234        // Multiplication procedures:
2235
2236        p_Procs->p_Mult_mm   = sca_p_Mult_mm;
2237        _p_procs->p_Mult_mm  = sca_p_Mult_mm;
2238
2239        p_Procs->pp_Mult_mm  = sca_pp_Mult_mm;
2240        _p_procs->pp_Mult_mm = sca_pp_Mult_mm;
2241
2242        r->nc->mmMultP()     = sca_mm_Mult_p;
2243        r->nc->mmMultPP()    = sca_mm_Mult_pp;
2244
2245        r->nc->GB()            = sca_gr_bba;
2246/*
2247        // ??????????????????????????????????????? coefficients swell...
2248        r->nc->SPoly()         = sca_SPoly;
2249        r->nc->ReduceSPoly()   = sca_ReduceSpoly;
2250*/
2251//         r->nc->BucketPolyRed() = gnc_kBucketPolyRed;
2252//         r->nc->BucketPolyRed_Z()= gnc_kBucketPolyRed_Z;
2253
2254#endif
2255}
2256
2257
2258// bi-Degree (x, y) of monomial "m"
2259// x-es and y-s are weighted by wx and wy resp.
2260// [optional] components have weights by wCx and wCy.
2261inline void m_GetBiDegree(const poly m,
2262  const intvec *wx, const intvec *wy,
2263  const intvec *wCx, const intvec *wCy,
2264  int& dx, int& dy, const ring r)
2265{
2266  const unsigned int N  = r->N;
2267
2268  p_Test(m, r);
2269
2270  assume( wx != NULL );
2271  assume( wy != NULL );
2272
2273  assume( wx->cols() == 1 );
2274  assume( wy->cols() == 1 );
2275
2276  assume( (unsigned int)wx->rows() >= N );
2277  assume( (unsigned int)wy->rows() >= N );
2278
2279  int x = 0;
2280  int y = 0;
2281
2282  for(int i = N; i > 0; i--)
2283  {
2284    const int d = p_GetExp(m, i, r);
2285    x += d * (*wx)[i-1];
2286    y += d * (*wy)[i-1];
2287  }
2288
2289  if( (wCx != NULL) && (wCy != NULL) )
2290  {
2291    const int c = p_GetComp(m, r);
2292
2293    if( wCx->range(c) )
2294      x += (*wCx)[c];
2295
2296    if( wCy->range(c) )
2297      x += (*wCy)[c];
2298  }
2299
2300  dx = x;
2301  dy = y;
2302}
2303
2304// returns true if polynom p is bi-homogenous with respect to the given weights
2305// simultaneously sets bi-Degree
2306bool p_IsBiHomogeneous(const poly p,
2307  const intvec *wx, const intvec *wy,
2308  const intvec *wCx, const intvec *wCy,
2309  int &dx, int &dy,
2310  const ring r)
2311{
2312  if( p == NULL )
2313  {
2314    dx = 0;
2315    dy = 0;
2316    return true;
2317  }
2318
2319  poly q = p;
2320
2321
2322  int ddx, ddy;
2323
2324  m_GetBiDegree( q, wx, wy, wCx, wCy, ddx, ddy, r); // get bi degree of lm(p)
2325
2326  pIter(q);
2327
2328  for(; q != NULL; pIter(q) )
2329  {
2330    int x, y;
2331
2332    m_GetBiDegree( q, wx, wy, wCx, wCy, x, y, r); // get bi degree of q
2333
2334    if ( (x != ddx) || (y != ddy) ) return false;
2335  }
2336
2337  dx = ddx;
2338  dy = ddy;
2339
2340  return true;
2341}
2342
2343
2344// returns true if id is bi-homogenous without respect to the given weights
2345bool id_IsBiHomogeneous(const ideal id,
2346  const intvec *wx, const intvec *wy,
2347  const intvec *wCx, const intvec *wCy,
2348  const ring r)
2349{
2350  if (id == NULL) return true; // zero ideal
2351
2352  const int iSize = IDELEMS(id);
2353
2354  if (iSize == 0) return true;
2355
2356  bool b = true;
2357  int x, y;
2358
2359  for(int i = iSize - 1; (i >= 0 ) && b; i--)
2360    b = p_IsBiHomogeneous(id->m[i], wx, wy, wCx, wCy, x, y, r);
2361
2362  return b;
2363}
2364
2365
2366// returns an intvector with [nvars(r)] integers [1/0]
2367// 1 - for commutative variables
2368// 0 - for anticommutative variables
2369intvec *ivGetSCAXVarWeights(const ring r)
2370{
2371  const unsigned int N  = r->N;
2372
2373  const int CommutativeVariable = 1;
2374  const int AntiCommutativeVariable = 0;
2375
2376  intvec* w = new intvec(N, 1, CommutativeVariable);
2377
2378  if( rIsSCA(r) )
2379  {
2380    const unsigned int m_iFirstAltVar = scaFirstAltVar(r);
2381    const unsigned int m_iLastAltVar  = scaLastAltVar(r);
2382
2383    for (unsigned int i = m_iFirstAltVar; i<= m_iLastAltVar; i++)
2384    {
2385      (*w)[i-1] = AntiCommutativeVariable;
2386    }
2387  }
2388  return w;
2389}
2390
2391
2392// returns an intvector with [nvars(r)] integers [1/0]
2393// 0 - for commutative variables
2394// 1 - for anticommutative variables
2395intvec *ivGetSCAYVarWeights(const ring r)
2396{
2397  const unsigned int N  = r->N;
2398
2399  const int CommutativeVariable = 0;
2400  const int AntiCommutativeVariable = 1;
2401
2402  intvec* w = new intvec(N, 1, CommutativeVariable);
2403
2404  if( rIsSCA(r) )
2405  {
2406    const unsigned int m_iFirstAltVar = scaFirstAltVar(r);
2407    const unsigned int m_iLastAltVar  = scaLastAltVar(r);
2408
2409    for (unsigned int i = m_iFirstAltVar; i<= m_iLastAltVar; i++)
2410    {
2411      (*w)[i-1] = AntiCommutativeVariable;
2412    }
2413  }
2414  return w;
2415}
2416
2417
2418
2419
2420// reduce term lt(m) modulo <y_i^2> , i = iFirstAltVar .. iLastAltVar:
2421// either create a copy of m or return NULL
2422inline poly m_KillSquares(const poly m,
2423  const unsigned int iFirstAltVar, const unsigned int iLastAltVar,
2424  const ring r)
2425{
2426#ifdef PDEBUG
2427  p_Test(m, r);
2428  assume( (iFirstAltVar >= 1) && (iLastAltVar <= r->N) && (iFirstAltVar <= iLastAltVar) );
2429
2430#if 0
2431  Print("m_KillSquares, m = "); // !
2432  p_Write(m, r);
2433#endif
2434#endif
2435
2436  assume( m != NULL );
2437
2438  for(unsigned int k = iFirstAltVar; k <= iLastAltVar; k++)
2439    if( p_GetExp(m, k, r) > 1 )
2440      return NULL;
2441
2442  return p_Head(m, r);
2443}
2444
2445
2446// reduce polynomial p modulo <y_i^2> , i = iFirstAltVar .. iLastAltVar
2447poly p_KillSquares(const poly p,
2448  const unsigned int iFirstAltVar, const unsigned int iLastAltVar,
2449  const ring r)
2450{
2451#ifdef PDEBUG
2452  p_Test(p, r);
2453
2454  assume( (iFirstAltVar >= 1) && (iLastAltVar <= r->N) && (iFirstAltVar <= iLastAltVar) );
2455
2456#if 0
2457  Print("p_KillSquares, p = "); // !
2458  p_Write(p, r);
2459#endif
2460#endif
2461
2462
2463  if( p == NULL )
2464    return NULL;
2465
2466  poly pResult = NULL;
2467  poly* ppPrev = &pResult;
2468
2469  for( poly q = p; q!= NULL; pIter(q) )
2470  {
2471#ifdef PDEBUG
2472    p_Test(q, r);
2473#endif
2474
2475    // terms will be in the same order because of quasi-ordering!
2476    poly v = m_KillSquares(q, iFirstAltVar, iLastAltVar, r);
2477
2478    if( v != NULL )
2479    {
2480      *ppPrev = v;
2481      ppPrev = &pNext(v);
2482    }
2483
2484  }
2485
2486#ifdef PDEBUG
2487  p_Test(pResult, r);
2488#if 0
2489  Print("p_KillSquares => "); // !
2490  p_Write(pResult, r);
2491#endif
2492#endif
2493
2494  return(pResult);
2495}
2496
2497
2498
2499
2500// reduces ideal id modulo <y_i^2> , i = iFirstAltVar .. iLastAltVar
2501// returns the reduced ideal or zero ideal.
2502ideal id_KillSquares(const ideal id,
2503  const unsigned int iFirstAltVar, const unsigned int iLastAltVar,
2504  const ring r)
2505{
2506  if (id == NULL) return id; // zero ideal
2507
2508  assume( (iFirstAltVar >= 1) && (iLastAltVar <= r->N) && (iFirstAltVar <= iLastAltVar) );
2509 
2510  const int iSize = id->idelems();
2511
2512  if (iSize == 0) return id;
2513
2514  ideal temp = idInit(iSize, id->rank);
2515
2516#if 0
2517   PrintS("<id_KillSquares>\n");
2518  {
2519    Print("ideal id: \n");
2520    for (int i = 0; i < id->idelems(); i++)
2521    {
2522      Print("; id[%d] = ", i+1);
2523      p_Write(id->m[i], r);
2524    }
2525    Print(";\n");
2526    PrintLn();
2527  }
2528#endif
2529
2530
2531  for (int j = 0; j < iSize; j++)
2532    temp->m[j] = p_KillSquares(id->m[j], iFirstAltVar, iLastAltVar, r);
2533
2534  idSkipZeroes(temp);
2535
2536#if 0
2537   PrintS("<id_KillSquares>\n");
2538  {
2539    Print("ideal temp: \n");
2540    for (int i = 0; i < temp->idelems(); i++)
2541    {
2542      Print("; temp[%d] = ", i+1);
2543      p_Write(temp->m[i], r);
2544    }
2545    Print(";\n");
2546    PrintLn();
2547  }
2548   PrintS("</id_KillSquares>\n");
2549#endif
2550
2551//  temp->rank = idRankFreeModule(temp, r);
2552
2553  return temp;
2554}
2555
2556
2557
2558
2559#endif
Note: See TracBrowser for help on using the repository browser.