source: git/kernel/sca.cc @ 8fbdb2

fieker-DuValspielwiese
Last change on this file since 8fbdb2 was ab5a00, checked in by Motsak Oleksandr <motsak@…>, 17 years ago
*motsak: bi-weighted polynomials/vectors/modules for SCA and others... git-svn-id: file:///usr/local/Singular/svn/trunk@9866 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 56.5 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.9 2007-02-16 11:05:54 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        Print("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        Print("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        Print("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        Print("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 pVariables = r->N;
711
712  for (int i = pVariables; i; 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
962ideal sca_gr_bba(const ideal F, const ideal Q, const intvec *, const intvec *, kStrategy strat)
963{
964  assume(rIsSCA(currRing));
965
966  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
967  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
968 
969  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
970  const ideal tempQ = currRing->nc->SCAQuotient();
971
972  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
973
974  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
975
976  strat->homog = strat->homog && bIdHomog;
977
978#ifdef PDEBUG
979  assume( strat->homog == bIdHomog );
980#endif /*PDEBUG*/
981
982
983#if 0
984   PrintS("<sca::bba>\n");
985  {
986    Print("ideal F: \n");
987    for (int i = 0; i < tempF->idelems(); i++)
988    {
989      Print("; F[%d] = ", i+1);
990      p_Write(tempF->m[i], currRing);
991    }
992    Print(";\n");
993    PrintLn();
994  }
995#endif
996
997
998#ifdef PDEBUG
999//   PrintS("sca_gr_bba\n");
1000#endif /*PDEBUG*/
1001
1002
1003
1004  int srmax, lrmax;
1005  int olddeg, reduc;
1006  int red_result = 1;
1007  int hilbeledeg = 1, hilbcount = 0, minimcnt = 0;
1008
1009  initBuchMoraCrit(strat); // set Gebauer, honey, sugarCrit
1010
1011  nc_gr_initBba(tempF,strat); // set enterS, red, initEcart, initEcartPair
1012
1013  initBuchMoraPos(strat);
1014
1015
1016  // ?? set spSpolyShort, reduce ???
1017
1018  initBuchMora(tempF, tempQ, strat); // currRing->nc->SCAQuotient() instead of Q == squares!!!!!!!
1019
1020  strat->posInT=posInT110; // !!!
1021
1022  srmax = strat->sl;
1023  reduc = olddeg = lrmax = 0;
1024
1025
1026  /* compute------------------------------------------------------- */
1027  for(; strat->Ll >= 0;
1028#ifdef KDEBUG
1029    strat->P.lcm = NULL,
1030#endif
1031    kTest(strat)
1032    )
1033  {
1034    if (strat->Ll > lrmax) lrmax =strat->Ll;/*stat.*/
1035
1036    if (TEST_OPT_DEBUG) messageSets(strat);
1037
1038    if (strat->Ll== 0) strat->interpt=TRUE;
1039
1040    if (TEST_OPT_DEGBOUND
1041    && ((strat->honey
1042    && (strat->L[strat->Ll].ecart+pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
1043       || ((!strat->honey) && (pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))))
1044    {
1045      /*
1046      *stops computation if
1047      * 24 IN test and the degree +ecart of L[strat->Ll] is bigger then
1048      *a predefined number Kstd1_deg
1049      */
1050      while (strat->Ll >= 0) deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
1051      break;
1052    }
1053
1054    /* picks the last element from the lazyset L */
1055    strat->P = strat->L[strat->Ll];
1056    strat->Ll--;
1057
1058    //kTest(strat);
1059
1060    assume(pNext(strat->P.p) != strat->tail);
1061
1062//     if( pNext(strat->P.p) == strat->tail )
1063//     {
1064//       // deletes the int spoly and computes SPoly
1065//       pLmFree(strat->P.p); // ???
1066//       strat->P.p = sca_SPoly(strat->P.p1, strat->P.p2, currRing);
1067//     }
1068
1069    if(strat->P.IsNull()) continue;
1070
1071//     poly save = NULL;
1072//
1073//     if(pNext(strat->P.p) != NULL)
1074//       save = p_Copy(strat->P.p, currRing);
1075
1076    strat->initEcart(&strat->P); // remove it?
1077
1078    if (TEST_OPT_PROT)
1079      message((strat->honey ? strat->P.ecart : 0) + strat->P.pFDeg(), &olddeg,&reduc,strat, red_result);
1080
1081    // reduction of the element chosen from L wrt S
1082    strat->red(&strat->P,strat);
1083
1084    if(strat->P.IsNull()) continue;
1085
1086    addLObject(strat->P, strat);
1087
1088    if (strat->sl > srmax) srmax = strat->sl;
1089
1090    const poly save = strat->P.p;
1091
1092#ifdef PDEBUG
1093      p_Test(save, currRing);
1094#endif
1095    assume( save != NULL );
1096   
1097    // SCA Specials:
1098
1099    {
1100      const poly pNext = pNext(save);
1101     
1102      if( pNext != NULL )
1103      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
1104      if( p_GetExp(save, i, currRing) != 0 )
1105      {
1106        assume(p_GetExp(save, i, currRing) == 1);
1107       
1108        const poly tt = sca_pp_Mult_xi_pp(i, pNext, currRing);
1109
1110#ifdef PDEBUG
1111        p_Test(tt, currRing);
1112#endif
1113
1114        if( tt == NULL) continue;
1115
1116        LObject h(tt); // h = x_i * P
1117
1118        if (TEST_OPT_INTSTRATEGY)
1119        {
1120//           h.pCleardenom(); // also does a pContent
1121          pContent(h.p);
1122        }
1123        else
1124        {
1125          h.pNorm();
1126        }
1127
1128        strat->initEcart(&h);
1129
1130
1131//         if (pOrdSgn==-1)
1132//         {
1133//           cancelunit(&h);  // tries to cancel a unit
1134//           deleteHC(&h, strat);
1135//         }
1136
1137//         if(h.IsNull()) continue;
1138
1139//         if (TEST_OPT_PROT)
1140//           message((strat->honey ? h.ecart : 0) + h.pFDeg(), &olddeg, &reduc, strat, red_result);
1141
1142//         strat->red(&h, strat); // wrt S
1143//         if(h.IsNull()) continue;
1144
1145//         poly save = p_Copy(h.p, currRing);
1146
1147        int pos;
1148
1149        if (strat->Ll==-1)
1150          pos =0;
1151        else
1152          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1153
1154        h.sev = pGetShortExpVector(h.p);
1155        enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
1156
1157  //       h.p = save;
1158  //       addLObject(h, strat);
1159
1160//         if (strat->sl > srmax) srmax = strat->sl;
1161      }
1162
1163      // p_Delete( &save, currRing );
1164    }
1165
1166
1167  } // for(;;)
1168
1169
1170  if (TEST_OPT_DEBUG) messageSets(strat);
1171
1172  /* release temp data-------------------------------- */
1173  exitBuchMora(strat);
1174
1175  if (TEST_OPT_WEIGHTM)
1176  {
1177    pFDeg=pFDegOld;
1178    pLDeg=pLDegOld;
1179    if (ecartWeights)
1180    {
1181      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(int));
1182      ecartWeights=NULL;
1183    }
1184  }
1185
1186  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
1187
1188  if (tempQ!=NULL) updateResult(strat->Shdl,tempQ,strat);
1189 
1190  id_Delete(&tempF, currRing);
1191
1192
1193  /* complete reduction of the standard basis--------- */
1194  if (TEST_OPT_REDSB){
1195//     completeReduce(strat); // ???
1196
1197    ideal I = strat->Shdl;
1198    ideal erg = kInterRed(I,tempQ);
1199    assume(I!=erg);
1200    id_Delete(&I, currRing);
1201    strat->Shdl = erg;
1202  }
1203
1204
1205#if 0
1206   PrintS("</sca::bba>\n");
1207#endif
1208
1209  return (strat->Shdl);
1210}
1211
1212
1213// should be used only inside nc_SetupQuotient!
1214// Check whether this our case:
1215//  1. rG is  a commutative polynomial ring \otimes anticommutative algebra
1216//  2. factor ideal rGR->qideal contains squares of all alternating variables.
1217//
1218// if yes, make rGR a super-commutative algebra!
1219// NOTE: Factors of SuperCommutative Algebras are supported this way!
1220bool sca_SetupQuotient(ring rGR, const ring rG)
1221{
1222//   return false; // test Plural
1223
1224  //////////////////////////////////////////////////////////////////////////
1225  // checks...
1226  //////////////////////////////////////////////////////////////////////////
1227  assume(rGR != NULL);
1228  assume(rG  != NULL);
1229
1230  assume(rIsPluralRing(rG));
1231
1232  const int N = rG->N;
1233
1234  if(N < 2)
1235    return false;
1236
1237  if(ncRingType(rG) != nc_skew)
1238    return false;
1239
1240  if(rGR->qideal == NULL) // there will be a factor!
1241    return false;
1242
1243  if(rG->qideal != NULL) // we cannot change from factor to factor!
1244    return false;
1245
1246
1247  int iAltVarEnd = -1;
1248  int iAltVarStart   = N+1;
1249
1250  const ring rBase = rG->nc->basering;
1251  const matrix C   = rG->nc->C; // live in rBase!
1252
1253  for(int i = 1; i < N; i++)
1254  {
1255    for(int j = i + 1; j <= N; j++)
1256    {
1257      assume(MATELEM(C,i,j) != NULL); // after CallPlural!
1258      number c = p_GetCoeff(MATELEM(C,i,j), rBase);
1259     
1260      if( n_IsMOne(c, rBase) )
1261      {       
1262        if( i < iAltVarStart) 
1263          iAltVarStart = i;
1264         
1265        if( j > iAltVarEnd)
1266          iAltVarEnd = j;
1267      } else
1268      {
1269        if( !n_IsOne(c, rBase) )
1270        {
1271#ifdef PDEBUG
1272//           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1273#endif
1274          return false;
1275        }
1276      }
1277    }
1278  }
1279
1280  if( (iAltVarEnd == -1) || (iAltVarStart == (N+1)) )
1281    return false; // either no alternating varables, or a single one => we are in commutative case!
1282 
1283  for(int i = 1; i < N; i++)
1284  {
1285    for(int j = i + 1; j <= N; j++)
1286    {
1287      assume(MATELEM(C,i,j) != NULL); // after CallPlural!
1288      number c = p_GetCoeff(MATELEM(C,i,j), rBase);     
1289
1290      if( (iAltVarStart <= i) && (j <= iAltVarEnd) ) // S <= i < j <= E
1291      { // anticommutative part
1292        if( !n_IsMOne(c, rBase) )
1293        {
1294#ifdef PDEBUG
1295//           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1296#endif
1297          return false;
1298        }
1299      } else
1300      { // should commute
1301        if( !n_IsOne(c, rBase) )
1302        {
1303#ifdef PDEBUG
1304//           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1305#endif
1306          return false;
1307        }
1308      }
1309    }
1310  }
1311
1312  assume( 1            <= iAltVarStart );
1313  assume( iAltVarStart < iAltVarEnd   );
1314  assume( iAltVarEnd   <= N            );
1315
1316  bool bWeChangeRing = false;
1317  // for sanity
1318  ring rSaveRing = currRing;
1319
1320  if(rSaveRing != rG)
1321  {
1322    rChangeCurrRing(rG);
1323    bWeChangeRing = true;
1324  }
1325
1326  assume(rGR->qideal != NULL);
1327  assume(rG->qideal == NULL);
1328 
1329  const ideal idQuotient = rGR->qideal;
1330
1331  // check for
1332  // y_{iAltVarStart}^2, y_{iAltVarStart+1}^2, \ldots, y_{iAltVarEnd}^2  (iAltVarEnd > iAltVarStart)
1333  // to be within quotient ideal.
1334
1335  bool bSCA = true;
1336
1337  for ( int i = iAltVarStart; (i <= iAltVarEnd) && bSCA; i++ )
1338  {
1339    poly square = p_ISet(1, rSaveRing);
1340    p_SetExp(square, i, 2, rSaveRing); // square = var(i)^2.
1341    p_Setm(square, rSaveRing);
1342
1343    // square = NF( var(i)^2 | Q )
1344    // NOTE: rSaveRing == currRing now!
1345    // NOTE: there is no better way to check this in general!
1346    square = kNF(idQuotient, NULL, square, 0, 0); 
1347   
1348    if( square != NULL ) // var(i)^2 is not in Q?
1349    {
1350      p_Delete(&square, rSaveRing);
1351      bSCA = false;     
1352    }   
1353  }
1354 
1355
1356  if (bWeChangeRing)
1357  {
1358    rChangeCurrRing(rSaveRing);
1359  }
1360
1361  if(!bSCA) return false;
1362 
1363
1364
1365#ifdef PDEBUG
1366//  Print("AltVars: [%d, %d]\n", iAltVarStart, iAltVarEnd);
1367#endif
1368
1369
1370  //////////////////////////////////////////////////////////////////////////
1371  // ok... here we go. let's setup it!!!
1372  //////////////////////////////////////////////////////////////////////////
1373  ideal tempQ = id_KillSquares(idQuotient, iAltVarStart, iAltVarEnd, rG); // in rG!!!
1374
1375  rGR->nc->SCAQuotient() = idrMoveR(tempQ, rG, rGR); // deletes tempQ!
1376
1377  ncRingType( rGR, nc_exterior );
1378
1379  scaFirstAltVar( rGR, iAltVarStart );
1380  scaLastAltVar( rGR, iAltVarEnd );
1381
1382
1383
1384  sca_p_ProcsSet(rGR, rGR->p_Procs);
1385
1386 
1387  return true;
1388}
1389
1390// return x_i * pPoly; preserve pPoly.
1391poly sca_pp_Mult_xi_pp(unsigned int i, const poly pPoly, const ring rRing)
1392{
1393  assume(1 <= i);
1394  assume(i <= rRing->N);
1395
1396  if(rIsSCA(rRing))
1397    return sca_xi_Mult_pp(i, pPoly, rRing);
1398
1399
1400
1401  poly xi =  p_ISet(1, rRing);
1402  p_SetExp(xi, i, 1, rRing);
1403  p_Setm(xi, rRing);
1404
1405  poly pResult = pp_Mult_qq(xi, pPoly, rRing);
1406
1407  p_Delete( &xi, rRing);
1408
1409  return pResult;
1410
1411}
1412
1413
1414///////////////////////////////////////////////////////////////////////////////////////
1415//************* SCA BBA *************************************************************//
1416///////////////////////////////////////////////////////////////////////////////////////
1417
1418// Under development!!!
1419ideal sca_bba (const ideal F, const ideal Q, const intvec *w, const intvec *hilb, kStrategy strat)
1420{
1421  assume(rIsSCA(currRing));
1422
1423  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
1424  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
1425
1426  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
1427
1428  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
1429
1430  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
1431
1432  strat->homog = strat->homog && bIdHomog;
1433
1434#ifdef PDEBUG
1435  assume( strat->homog == bIdHomog );
1436#endif /*PDEBUG*/
1437
1438
1439//    PrintS("<sca>\n");
1440
1441  om_Opts.MinTrack = 5; // ???
1442
1443  int   srmax, lrmax, red_result = 1;
1444  int   olddeg, reduc;
1445  int hilbeledeg=1, hilbcount=0, minimcnt=0;
1446
1447  BOOLEAN withT = FALSE;
1448
1449  initBuchMoraCrit(strat); // sets Gebauer, honey, sugarCrit // sca - ok???
1450  initBuchMoraPos(strat); // sets strat->posInL, strat->posInT // check!! (Plural's: )
1451
1452//   initHilbCrit(F, Q, &hilb, strat);
1453
1454//  nc_gr_initBba(F,strat);
1455  initBba(tempF, strat); // set enterS, red, initEcart, initEcartPair
1456
1457  /*set enterS, spSpolyShort, reduce, red, initEcart, initEcartPair*/
1458  // ?? set spSpolyShort, reduce ???
1459  initBuchMora(tempF, NULL, strat); // Q = squares!!!
1460
1461//   if (strat->minim>0) strat->M = idInit(IDELEMS(F),F->rank);
1462
1463  srmax = strat->sl;
1464  reduc = olddeg = lrmax = 0;
1465
1466#define NO_BUCKETS
1467
1468#ifndef NO_BUCKETS
1469  if (!TEST_OPT_NOT_BUCKETS)
1470    strat->use_buckets = 1;
1471#endif
1472
1473  // redtailBBa against T for inhomogenous input
1474  if (!K_TEST_OPT_OLDSTD)
1475    withT = ! strat->homog;
1476
1477  // strat->posInT = posInT_pLength;
1478  kTest_TS(strat);
1479
1480#undef HAVE_TAIL_RING
1481
1482#ifdef HAVE_TAIL_RING
1483  kStratInitChangeTailRing(strat);
1484#endif
1485
1486  ///////////////////////////////////////////////////////////////
1487  // SCA:
1488
1489  /* compute------------------------------------------------------- */
1490  while (strat->Ll >= 0)
1491  {
1492    if (strat->Ll > lrmax) lrmax =strat->Ll;/*stat.*/
1493
1494#ifdef KDEBUG
1495//     loop_count++;
1496    if (TEST_OPT_DEBUG) messageSets(strat);
1497#endif
1498
1499    if (strat->Ll== 0) strat->interpt=TRUE;
1500
1501    if (TEST_OPT_DEGBOUND
1502        && ((strat->honey && (strat->L[strat->Ll].ecart+pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
1503            || ((!strat->honey) && (pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))))
1504    {
1505      /*
1506       *stops computation if
1507       * 24 IN test and the degree +ecart of L[strat->Ll] is bigger then
1508       *a predefined number Kstd1_deg
1509       */
1510      while ((strat->Ll >= 0)
1511  && (strat->L[strat->Ll].p1!=NULL) && (strat->L[strat->Ll].p2!=NULL)
1512        && ((strat->honey && (strat->L[strat->Ll].ecart+pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
1513            || ((!strat->honey) && (pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg)))
1514  )
1515        deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
1516      if (strat->Ll<0) break;
1517      else strat->noClearS=TRUE;
1518    }
1519
1520    /* picks the last element from the lazyset L */
1521    strat->P = strat->L[strat->Ll];
1522    strat->Ll--;
1523
1524
1525    assume(pNext(strat->P.p) != strat->tail);
1526
1527/*    if (pNext(strat->P.p) == strat->tail)
1528    {
1529      // deletes the short spoly
1530      pLmFree(strat->P.p);
1531      strat->P.p = NULL;
1532      poly m1 = NULL, m2 = NULL;
1533
1534      // check that spoly creation is ok
1535      while (strat->tailRing != currRing &&
1536             !kCheckSpolyCreation(&(strat->P), strat, m1, m2))
1537      {
1538        assume(m1 == NULL && m2 == NULL);
1539        // if not, change to a ring where exponents are at least
1540        // large enough
1541        kStratChangeTailRing(strat);
1542      }
1543
1544#ifdef PDEBUG
1545      Print("ksCreateSpoly!#?");
1546#endif
1547
1548      // create the real one
1549      ksCreateSpoly(&(strat->P), NULL, strat->use_buckets,
1550                    strat->tailRing, m1, m2, strat->R); //?????????
1551    }
1552    else */
1553
1554    if (strat->P.p1 == NULL)
1555    {
1556//       if (strat->minim > 0)
1557//         strat->P.p2=p_Copy(strat->P.p, currRing, strat->tailRing);
1558
1559
1560      // for input polys, prepare reduction
1561      strat->P.PrepareRed(strat->use_buckets);
1562    }
1563
1564    if (TEST_OPT_PROT)
1565      message((strat->honey ? strat->P.ecart : 0) + strat->P.pFDeg(),
1566              &olddeg,&reduc,strat, red_result);
1567
1568    /* reduction of the element choosen from L */
1569    red_result = strat->red(&strat->P,strat);
1570
1571
1572    // reduction to non-zero new poly
1573    if (red_result == 1)
1574    {
1575      /* statistic */
1576      if (TEST_OPT_PROT) PrintS("s");
1577
1578      // get the polynomial (canonicalize bucket, make sure P.p is set)
1579      strat->P.GetP(strat->lmBin);
1580
1581      int pos = posInS(strat,strat->sl,strat->P.p,strat->P.ecart);
1582
1583      // reduce the tail and normalize poly
1584      if (TEST_OPT_INTSTRATEGY)
1585      {
1586        strat->P.pCleardenom();
1587        if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1588        {
1589          strat->P.p = redtailBba(&(strat->P),pos-1,strat, withT); // !!!
1590          strat->P.pCleardenom();
1591
1592        }
1593      }
1594      else
1595      {
1596
1597        strat->P.pNorm();
1598        if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1599          strat->P.p = redtailBba(&(strat->P),pos-1,strat, withT);
1600      }
1601
1602
1603#ifdef KDEBUG
1604      if (TEST_OPT_DEBUG){PrintS(" ns:");strat->P.wrp();PrintLn();}
1605#endif
1606
1607//       // min_std stuff
1608//       if ((strat->P.p1==NULL) && (strat->minim>0))
1609//       {
1610//         if (strat->minim==1)
1611//         {
1612//           strat->M->m[minimcnt]=p_Copy(strat->P.p,currRing,strat->tailRing);
1613//           p_Delete(&strat->P.p2, currRing, strat->tailRing);
1614//         }
1615//         else
1616//         {
1617//           strat->M->m[minimcnt]=strat->P.p2;
1618//           strat->P.p2=NULL;
1619//         }
1620//         if (strat->tailRing!=currRing && pNext(strat->M->m[minimcnt])!=NULL)
1621//           pNext(strat->M->m[minimcnt])
1622//             = strat->p_shallow_copy_delete(pNext(strat->M->m[minimcnt]),
1623//                                            strat->tailRing, currRing,
1624//                                            currRing->PolyBin);
1625//         minimcnt++;
1626//       }
1627
1628      // enter into S, L, and T
1629      if(withT)
1630        enterT(strat->P, strat);
1631
1632      // L
1633      enterpairs(strat->P.p,strat->sl,strat->P.ecart,pos,strat, strat->tl);
1634
1635      // posInS only depends on the leading term
1636      strat->enterS(strat->P, pos, strat, strat->tl);
1637
1638//       if (hilb!=NULL) khCheck(Q,w,hilb,hilbeledeg,hilbcount,strat);
1639
1640//      Print("[%d]",hilbeledeg);
1641      if (strat->P.lcm!=NULL) pLmFree(strat->P.lcm);
1642
1643      if (strat->sl>srmax) srmax = strat->sl;
1644
1645      // //////////////////////////////////////////////////////////
1646      // SCA:
1647      const poly pSave = strat->P.p;
1648      const poly pNext = pNext(pSave);
1649
1650//       if(0)
1651      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
1652      if( p_GetExp(pSave, i, currRing) != 0 )
1653      {
1654        assume(p_GetExp(pSave, i, currRing) == 1);
1655       
1656        const poly pNew = sca_pp_Mult_xi_pp(i, pNext, currRing);
1657
1658#ifdef PDEBUG
1659        p_Test(pNew, currRing);
1660#endif
1661
1662        if( pNew == NULL) continue;
1663
1664        LObject h(pNew); // h = x_i * strat->P
1665
1666        if (TEST_OPT_INTSTRATEGY)
1667        {
1668//           h.pCleardenom(); // also does a pContent
1669          pContent(h.p);
1670        }
1671        else
1672        {
1673          h.pNorm();
1674        }
1675
1676        strat->initEcart(&h);
1677        h.sev = pGetShortExpVector(h.p);
1678
1679        int pos;
1680
1681        if (strat->Ll==-1)
1682          pos =0;
1683        else
1684          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1685
1686        enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
1687
1688
1689/*
1690        h.sev = pGetShortExpVector(h.p);
1691        strat->initEcart(&h);
1692
1693        h.PrepareRed(strat->use_buckets);
1694
1695        // reduction of the element choosen from L(?)
1696        red_result = strat->red(&h,strat);
1697
1698        // reduction to non-zero new poly
1699        if (red_result != 1) continue;
1700
1701
1702        int pos = posInS(strat,strat->sl,h.p,h.ecart);
1703
1704        // reduce the tail and normalize poly
1705        if (TEST_OPT_INTSTRATEGY)
1706        {
1707          h.pCleardenom();
1708          if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1709          {
1710            h.p = redtailBba(&(h),pos-1,strat, withT); // !!!
1711            h.pCleardenom();
1712          }
1713        }
1714        else
1715        {
1716
1717          h.pNorm();
1718          if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1719            h.p = redtailBba(&(h),pos-1,strat, withT);
1720        }
1721
1722
1723  #ifdef KDEBUG
1724        if (TEST_OPT_DEBUG){PrintS(" N:");h.wrp();PrintLn();}
1725  #endif
1726
1727//         h.PrepareRed(strat->use_buckets); // ???
1728
1729        h.sev = pGetShortExpVector(h.p);
1730        strat->initEcart(&h);
1731
1732        if (strat->Ll==-1)
1733          pos = 0;
1734        else
1735          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1736
1737         enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);*/
1738
1739      } // for all x_i \in Ann(lm(P))
1740
1741
1742
1743
1744
1745    } // if red(P) != NULL
1746
1747//     else if (strat->P.p1 == NULL && strat->minim > 0)
1748//     {
1749//       p_Delete(&strat->P.p2, currRing, strat->tailRing);
1750//     }
1751
1752// #ifdef KDEBUG
1753    memset(&(strat->P), 0, sizeof(strat->P));
1754// #endif
1755
1756    kTest_TS(strat);
1757
1758//     Print("\n$\n");
1759
1760  }
1761
1762#ifdef KDEBUG
1763  if (TEST_OPT_DEBUG) messageSets(strat);
1764#endif
1765
1766  /* complete reduction of the standard basis--------- */
1767  if (TEST_OPT_REDSB) completeReduce(strat);
1768
1769  /* release temp data-------------------------------- */
1770  id_Delete(&tempF, currRing);
1771
1772  exitBuchMora(strat);
1773
1774  if (TEST_OPT_WEIGHTM)
1775  {
1776    pRestoreDegProcs(pFDegOld, pLDegOld);
1777    if (ecartWeights)
1778    {
1779      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(short));
1780      ecartWeights=NULL;
1781    }
1782  }
1783
1784  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
1785
1786//   if (Q!=NULL) updateResult(strat->Shdl,Q,strat);
1787
1788//    PrintS("</sca>\n");
1789
1790  return (strat->Shdl);
1791}
1792
1793
1794// //////////////////////////////////////////////////////////////////////////////
1795// sca mora...
1796
1797// returns TRUE if mora should use buckets, false otherwise
1798static BOOLEAN kMoraUseBucket(kStrategy strat)
1799{
1800#ifdef MORA_USE_BUCKETS
1801  if (TEST_OPT_NOT_BUCKETS)
1802    return FALSE;
1803  if (strat->red == redFirst)
1804  {
1805#ifdef NO_LDEG
1806    if (!strat->syzComp)
1807      return TRUE;
1808#else
1809    if ((strat->homog || strat->honey) && !strat->syzComp)
1810      return TRUE;
1811#endif
1812  }
1813  else
1814  {
1815    assume(strat->red == redEcart);
1816    if (strat->honey && !strat->syzComp)
1817      return TRUE;
1818  }
1819#endif
1820  return FALSE;
1821}
1822
1823
1824#ifdef HAVE_ASSUME
1825static int sca_mora_count = 0;
1826static int sca_mora_loop_count;
1827#endif
1828
1829// ideal sca_mora (ideal F, ideal Q, intvec *w, intvec *, kStrategy strat)
1830ideal sca_mora(const ideal F, const ideal Q, const intvec *w, const intvec *, kStrategy strat)
1831{
1832  assume(rIsSCA(currRing));
1833
1834  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
1835  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
1836 
1837  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
1838 
1839  const ideal tempQ = currRing->nc->SCAQuotient();
1840
1841  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
1842
1843  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
1844
1845  strat->homog = strat->homog && bIdHomog;
1846
1847#ifdef PDEBUG
1848  assume( strat->homog == bIdHomog );
1849#endif /*PDEBUG*/
1850
1851#ifdef HAVE_ASSUME
1852  sca_mora_count++;
1853  sca_mora_loop_count = 0;
1854#endif
1855
1856#ifdef KDEBUG
1857  om_Opts.MinTrack = 5;
1858#endif
1859  int srmax;
1860  int lrmax = 0;
1861  int olddeg = 0;
1862  int reduc = 0;
1863  int red_result = 1;
1864  int hilbeledeg=1,hilbcount=0;
1865
1866  strat->update = TRUE;
1867  /*- setting global variables ------------------- -*/
1868  initBuchMoraCrit(strat);
1869//   initHilbCrit(F,NULL,&hilb,strat); // no Q!
1870  initMora(tempF, strat);
1871  initBuchMoraPos(strat);
1872  /*Shdl=*/initBuchMora(tempF, tempQ, strat); // temp Q, F!
1873//   if (TEST_OPT_FASTHC) missingAxis(&strat->lastAxis,strat);
1874  /*updateS in initBuchMora has Hecketest
1875  * and could have put strat->kHEdgdeFound FALSE*/
1876#if 0
1877  if (ppNoether!=NULL)
1878  {
1879    strat->kHEdgeFound = TRUE;
1880  }
1881  if (strat->kHEdgeFound && strat->update)
1882  {
1883    firstUpdate(strat);
1884    updateLHC(strat);
1885    reorderL(strat);
1886  }
1887  if (TEST_OPT_FASTHC && (strat->lastAxis) && strat->posInLOldFlag)
1888  {
1889    strat->posInLOld = strat->posInL;
1890    strat->posInLOldFlag = FALSE;
1891    strat->posInL = posInL10;
1892    updateL(strat);
1893    reorderL(strat);
1894  }
1895#endif
1896
1897  srmax = strat->sl;
1898  kTest_TS(strat);
1899
1900  strat->use_buckets = kMoraUseBucket(strat);
1901  /*- compute-------------------------------------------*/
1902
1903#undef HAVE_TAIL_RING
1904
1905#ifdef HAVE_TAIL_RING
1906//  if (strat->homog && strat->red == redFirst)
1907//     kStratInitChangeTailRing(strat);
1908#endif
1909
1910
1911  while (strat->Ll >= 0)
1912  {
1913#ifdef HAVE_ASSUME
1914    sca_mora_loop_count++;
1915#endif
1916    if (lrmax< strat->Ll) lrmax=strat->Ll; /*stat*/
1917    //test_int_std(strat->kIdeal);
1918    if (TEST_OPT_DEBUG) messageSets(strat);
1919    if (TEST_OPT_DEGBOUND
1920    && (strat->L[strat->Ll].ecart+strat->L[strat->Ll].GetpFDeg()> Kstd1_deg))
1921    {
1922      /*
1923      * stops computation if
1924      * - 24 (degBound)
1925      *   && upper degree is bigger than Kstd1_deg
1926      */
1927      while ((strat->Ll >= 0)
1928        && (strat->L[strat->Ll].p1!=NULL) && (strat->L[strat->Ll].p2!=NULL)
1929        && (strat->L[strat->Ll].ecart+strat->L[strat->Ll].GetpFDeg()> Kstd1_deg)
1930      )
1931      {
1932        deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
1933        //if (TEST_OPT_PROT)
1934        //{
1935        //   PrintS("D"); mflush();
1936        //}
1937      }
1938      if (strat->Ll<0) break;
1939      else strat->noClearS=TRUE;
1940    }
1941    strat->P = strat->L[strat->Ll];/*- picks the last element from the lazyset L -*/
1942    if (strat->Ll==0) strat->interpt=TRUE;
1943    strat->Ll--;
1944
1945    // create the real Spoly
1946    assume(pNext(strat->P.p) != strat->tail);
1947
1948    if (strat->P.p1 == NULL)
1949    {
1950      // for input polys, prepare reduction (buckets !)
1951      strat->P.SetLength(strat->length_pLength);
1952      strat->P.PrepareRed(strat->use_buckets);
1953    }
1954
1955    if (!strat->P.IsNull())
1956    {
1957      // might be NULL from noether !!!
1958      if (TEST_OPT_PROT)
1959        message(strat->P.ecart+strat->P.GetpFDeg(),&olddeg,&reduc,strat, red_result);
1960      // reduce
1961      red_result = strat->red(&strat->P,strat);
1962    }
1963
1964    if (! strat->P.IsNull())
1965    {
1966      strat->P.GetP();
1967      // statistics
1968      if (TEST_OPT_PROT) PrintS("s");
1969      // normalization
1970      if (!TEST_OPT_INTSTRATEGY)
1971        strat->P.pNorm();
1972      // tailreduction
1973      strat->P.p = redtail(&(strat->P),strat->sl,strat);
1974      // set ecart -- might have changed because of tail reductions
1975      if ((!strat->noTailReduction) && (!strat->honey))
1976        strat->initEcart(&strat->P);
1977      // cancel unit
1978      cancelunit(&strat->P);
1979      // for char 0, clear denominators
1980      if (TEST_OPT_INTSTRATEGY)
1981        strat->P.pCleardenom();
1982
1983      // put in T
1984      enterT(strat->P,strat);
1985      // build new pairs
1986      enterpairs(strat->P.p,strat->sl,strat->P.ecart,0,strat, strat->tl);
1987      // put in S
1988      strat->enterS(strat->P,
1989                    posInS(strat,strat->sl,strat->P.p, strat->P.ecart),
1990                    strat, strat->tl);
1991
1992
1993      // clear strat->P
1994      if (strat->P.lcm!=NULL) pLmFree(strat->P.lcm);
1995      strat->P.lcm=NULL;
1996
1997      if (strat->sl>srmax) srmax = strat->sl; /*stat.*/
1998      if (strat->Ll>lrmax) lrmax = strat->Ll;
1999
2000
2001
2002      // //////////////////////////////////////////////////////////
2003      // SCA:
2004      const poly pSave = strat->P.p;
2005      const poly pNext = pNext(pSave);
2006
2007      if(pNext != NULL)
2008      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
2009      if( p_GetExp(pSave, i, currRing) != 0 )
2010      {
2011
2012        assume(p_GetExp(pSave, i, currRing) == 1);
2013
2014        const poly pNew = sca_pp_Mult_xi_pp(i, pNext, currRing);
2015
2016#ifdef PDEBUG
2017        p_Test(pNew, currRing);
2018#endif
2019
2020        if( pNew == NULL) continue;
2021
2022        LObject h(pNew); // h = x_i * strat->P
2023
2024        if (TEST_OPT_INTSTRATEGY)
2025           h.pCleardenom(); // also does a pContent
2026        else
2027          h.pNorm();
2028
2029        strat->initEcart(&h);
2030        h.sev = pGetShortExpVector(h.p);
2031
2032        int pos;
2033
2034        if (strat->Ll==-1)
2035          pos = 0;
2036        else
2037          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
2038
2039        enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
2040
2041        if (strat->Ll>lrmax) lrmax = strat->Ll;
2042      }
2043
2044#ifdef KDEBUG
2045      // make sure kTest_TS does not complain about strat->P
2046      memset(&strat->P,0,sizeof(strat->P));
2047#endif
2048    }
2049#if 0
2050    if (strat->kHEdgeFound)
2051    {
2052      if ((BTEST1(27))
2053      || ((TEST_OPT_MULTBOUND) && (scMult0Int((strat->Shdl)) < mu)))
2054      {
2055        // obachman: is this still used ???
2056        /*
2057        * stops computation if strat->kHEdgeFound and
2058        * - 27 (finiteDeterminacyTest)
2059        * or
2060        * - 23
2061        *   (multBound)
2062        *   && multiplicity of the ideal is smaller then a predefined number mu
2063        */
2064        while (strat->Ll >= 0) deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
2065      }
2066    }
2067#endif
2068    kTest_TS(strat);
2069  }
2070  /*- complete reduction of the standard basis------------------------ -*/
2071  if (TEST_OPT_REDSB) completeReduce(strat);
2072  /*- release temp data------------------------------- -*/
2073  exitBuchMora(strat);
2074  /*- polynomials used for HECKE: HC, noether -*/
2075  if (BTEST1(27))
2076  {
2077    if (strat->kHEdge!=NULL)
2078      Kstd1_mu=pFDeg(strat->kHEdge,currRing);
2079    else
2080      Kstd1_mu=-1;
2081  }
2082  pDelete(&strat->kHEdge);
2083  strat->update = TRUE; //???
2084  strat->lastAxis = 0; //???
2085  pDelete(&strat->kNoether);
2086  omFreeSize((ADDRESS)strat->NotUsedAxis,(pVariables+1)*sizeof(BOOLEAN));
2087  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
2088  if (TEST_OPT_WEIGHTM)
2089  {
2090    pRestoreDegProcs(pFDegOld, pLDegOld);
2091    if (ecartWeights)
2092    {
2093      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(short));
2094      ecartWeights=NULL;
2095    }
2096  }
2097  if (tempQ!=NULL) updateResult(strat->Shdl,tempQ,strat);
2098  idTest(strat->Shdl);
2099 
2100  id_Delete( &tempF, currRing);
2101 
2102  return (strat->Shdl);
2103}
2104
2105
2106
2107
2108
2109
2110void sca_p_ProcsSet(ring rGR, p_Procs_s* p_Procs)
2111{
2112
2113  // "commutative" procedures:
2114  rGR->p_Procs->p_Mult_mm     = sca_p_Mult_mm;
2115  rGR->p_Procs->pp_Mult_mm    = sca_pp_Mult_mm;
2116
2117  p_Procs->p_Mult_mm          = sca_p_Mult_mm;
2118  p_Procs->pp_Mult_mm         = sca_pp_Mult_mm;
2119
2120  // non-commutaitve
2121  rGR->nc->p_Procs.mm_Mult_p   = sca_mm_Mult_p;
2122  rGR->nc->p_Procs.mm_Mult_pp  = sca_mm_Mult_pp;
2123
2124
2125  if (pOrdSgn==-1)
2126  {
2127#ifdef PDEBUG
2128//           Print("Local case => GB == mora!\n");
2129#endif
2130    rGR->nc->p_Procs.GB          = sca_mora; // local ordering => Mora, otherwise - Buchberger!
2131  }
2132  else
2133  {
2134#ifdef PDEBUG
2135//           Print("Global case => GB == bba!\n");
2136#endif
2137    rGR->nc->p_Procs.GB          = sca_gr_bba; // sca_bba?
2138  }
2139
2140
2141//   rGR->nc->p_Procs.GlobalGB    = sca_gr_bba;
2142//   rGR->nc->p_Procs.LocalGB     = sca_mora;
2143
2144
2145//   rGR->nc->p_Procs.SPoly         = sca_SPoly;
2146//   rGR->nc->p_Procs.ReduceSPoly   = sca_ReduceSpoly;
2147
2148#if 0
2149
2150        // Multiplication procedures:
2151
2152        p_Procs->p_Mult_mm   = sca_p_Mult_mm;
2153        _p_procs->p_Mult_mm  = sca_p_Mult_mm;
2154
2155        p_Procs->pp_Mult_mm  = sca_pp_Mult_mm;
2156        _p_procs->pp_Mult_mm = sca_pp_Mult_mm;
2157
2158        r->nc->mmMultP()     = sca_mm_Mult_p;
2159        r->nc->mmMultPP()    = sca_mm_Mult_pp;
2160
2161        r->nc->GB()            = sca_gr_bba;
2162/*
2163        // ??????????????????????????????????????? coefficients swell...
2164        r->nc->SPoly()         = sca_SPoly;
2165        r->nc->ReduceSPoly()   = sca_ReduceSpoly;
2166*/
2167//         r->nc->BucketPolyRed() = gnc_kBucketPolyRed;
2168//         r->nc->BucketPolyRed_Z()= gnc_kBucketPolyRed_Z;
2169
2170#endif
2171}
2172
2173
2174// bi-Degree (x, y) of monomial "m"
2175// x-es and y-s are weighted by wx and wy resp.
2176// [optional] components have weights by wCx and wCy.
2177inline void m_GetBiDegree(const poly m, 
2178  const intvec *wx, const intvec *wy, 
2179  const intvec *wCx, const intvec *wCy, 
2180  int& dx, int& dy, const ring r)
2181{
2182  const unsigned int N  = r->N;
2183 
2184  p_Test(m, r);
2185 
2186  assume( wx != NULL );
2187  assume( wy != NULL );
2188 
2189  assume( wx->cols() == 1 );
2190  assume( wy->cols() == 1 );
2191
2192  assume( wx->rows() >= N );
2193  assume( wy->rows() >= N );
2194
2195  int x = 0;
2196  int y = 0;
2197 
2198  for(int i = N; i > 0; i--)
2199  {
2200    const int d = p_GetExp(m, i, r);
2201    x += d * (*wx)[i-1];
2202    y += d * (*wy)[i-1];
2203  }
2204 
2205  if( (wCx != NULL) && (wCy != NULL) )
2206  {
2207    const int c = p_GetComp(m, r);
2208
2209    if( wCx->range(c) )
2210      x += (*wCx)[c];
2211
2212    if( wCy->range(c) )
2213      x += (*wCy)[c];     
2214  }
2215 
2216  dx = x;
2217  dy = y;
2218}
2219
2220// returns true if polynom p is bi-homogenous with respect to the given weights
2221// simultaneously sets bi-Degree
2222bool p_IsBiHomogeneous(const poly p, 
2223  const intvec *wx, const intvec *wy, 
2224  const intvec *wCx, const intvec *wCy, 
2225  int &dx, int &dy,
2226  const ring r)
2227{
2228  if( p == NULL ) 
2229  {
2230    dx = 0;
2231    dy = 0;
2232    return true;
2233  }
2234
2235  poly q = p;
2236
2237
2238  int ddx, ddy;
2239
2240  m_GetBiDegree( q, wx, wy, wCx, wCy, ddx, ddy, r); // get bi degree of lm(p)
2241
2242  pIter(q);
2243
2244  for(; q != NULL; pIter(q) )
2245  {
2246    int x, y;   
2247   
2248    m_GetBiDegree( q, wx, wy, wCx, wCy, x, y, r); // get bi degree of q
2249   
2250    if ( (x != ddx) || (y != ddy) ) return false;
2251  }
2252 
2253  dx = ddx;
2254  dy = ddy;
2255
2256  return true;
2257}
2258
2259
2260// returns true if id is bi-homogenous without respect to the given weights
2261bool id_IsBiHomogeneous(const ideal id, 
2262  const intvec *wx, const intvec *wy, 
2263  const intvec *wCx, const intvec *wCy, 
2264  const ring r)
2265{
2266  if (id == NULL) return true; // zero ideal
2267
2268  const int iSize = IDELEMS(id);
2269
2270  if (iSize == 0) return true;
2271
2272  bool b = true;
2273  int x, y;
2274
2275  for(int i = iSize - 1; (i >= 0 ) && b; i--)
2276    b = p_IsBiHomogeneous(id->m[i], wx, wy, wCx, wCy, x, y, r);
2277
2278  return b;
2279}
2280
2281
2282// returns an intvector with [nvars(r)] integers [1/0]
2283// 1 - for commutative variables
2284// 0 - for anticommutative variables
2285intvec *ivGetSCAXVarWeights(const ring r)
2286{
2287  const unsigned int N  = r->N;
2288
2289  const int CommutativeVariable = 1;
2290  const int AntiCommutativeVariable = 0;
2291 
2292  intvec* w = new intvec(N, 1, CommutativeVariable);
2293 
2294  if( rIsSCA(r) )
2295  {
2296    const unsigned int m_iFirstAltVar = scaFirstAltVar(r);
2297    const unsigned int m_iLastAltVar  = scaLastAltVar(r); 
2298
2299    for (unsigned int i = m_iFirstAltVar; i<= m_iLastAltVar; i++)
2300    {
2301      (*w)[i-1] = AntiCommutativeVariable;
2302    } 
2303  }
2304  return w;
2305}
2306
2307
2308// returns an intvector with [nvars(r)] integers [1/0]
2309// 0 - for commutative variables
2310// 1 - for anticommutative variables
2311intvec *ivGetSCAYVarWeights(const ring r)
2312{
2313  const unsigned int N  = r->N;
2314
2315  const int CommutativeVariable = 0;
2316  const int AntiCommutativeVariable = 1;
2317 
2318  intvec* w = new intvec(N, 1, CommutativeVariable);
2319 
2320  if( rIsSCA(r) )
2321  {
2322    const unsigned int m_iFirstAltVar = scaFirstAltVar(r);
2323    const unsigned int m_iLastAltVar  = scaLastAltVar(r); 
2324
2325    for (unsigned int i = m_iFirstAltVar; i<= m_iLastAltVar; i++)
2326    {
2327      (*w)[i-1] = AntiCommutativeVariable;
2328    } 
2329  }
2330  return w;
2331}
2332
2333
2334
2335
2336// reduce term lt(m) modulo <y_i^2> , i = iFirstAltVar .. iLastAltVar:
2337// either create a copy of m or return NULL
2338inline poly m_KillSquares(const poly m, 
2339  const unsigned int iFirstAltVar, const unsigned int iLastAltVar, 
2340  const ring r)
2341{ 
2342#ifdef PDEBUG
2343  p_Test(m, r);
2344
2345#if 0
2346  Print("m_KillSquares, m = "); // !
2347  p_Write(m, r);
2348#endif
2349#endif
2350
2351  assume( m != NULL );
2352
2353  for(int k = iFirstAltVar; k <= iLastAltVar; k++)
2354    if( p_GetExp(m, k, r) > 1 )
2355      return NULL; 
2356 
2357  return p_Head(m, r);
2358}
2359
2360
2361// reduce polynomial p modulo <y_i^2> , i = iFirstAltVar .. iLastAltVar
2362poly p_KillSquares(const poly p, 
2363  const unsigned int iFirstAltVar, const unsigned int iLastAltVar, 
2364  const ring r)
2365{ 
2366#ifdef PDEBUG
2367  p_Test(p, r);
2368
2369#if 0
2370  Print("p_KillSquares, p = "); // !
2371  p_Write(p, r);
2372#endif
2373#endif
2374
2375
2376  if( p == NULL )
2377    return NULL;
2378
2379  poly pResult = NULL;
2380  poly* ppPrev = &pResult;
2381
2382  for( poly q = p; q!= NULL; pIter(q) )
2383  {
2384#ifdef PDEBUG
2385    p_Test(q, r);
2386#endif
2387
2388    // terms will be in the same order because of quasi-ordering!   
2389    poly v = m_KillSquares(q, iFirstAltVar, iLastAltVar, r);
2390
2391    if( v != NULL )
2392    {
2393      *ppPrev = v;
2394      ppPrev = &pNext(v);
2395    }
2396
2397  }
2398
2399#ifdef PDEBUG
2400  p_Test(pResult, r);
2401#if 0
2402  Print("p_KillSquares => "); // !
2403  p_Write(pResult, r);
2404#endif
2405#endif
2406
2407  return(pResult);
2408}
2409 
2410
2411
2412
2413// reduces ideal id modulo <y_i^2> , i = iFirstAltVar .. iLastAltVar
2414// returns the reduced ideal or zero ideal.
2415ideal id_KillSquares(const ideal id, 
2416  const unsigned int iFirstAltVar, const unsigned int iLastAltVar, 
2417  const ring r)
2418{
2419  if (id == NULL) return id; // zero ideal
2420
2421  const int iSize = id->idelems();
2422
2423  if (iSize == 0) return id;
2424 
2425  ideal temp = idInit(iSize, 1);
2426 
2427#if 0
2428   PrintS("<id_KillSquares>\n");
2429  {
2430    Print("ideal id: \n");
2431    for (int i = 0; i < id->idelems(); i++)
2432    {
2433      Print("; id[%d] = ", i+1);
2434      p_Write(id->m[i], r);
2435    }
2436    Print(";\n");
2437    PrintLn();
2438  }
2439#endif
2440 
2441
2442  for (int j = 0; j < iSize; j++)
2443    temp->m[j] = p_KillSquares(id->m[j], iFirstAltVar, iLastAltVar, r);
2444
2445  idSkipZeroes(temp);
2446
2447#if 0
2448   PrintS("<id_KillSquares>\n");
2449  {
2450    Print("ideal temp: \n");
2451    for (int i = 0; i < temp->idelems(); i++)
2452    {
2453      Print("; temp[%d] = ", i+1);
2454      p_Write(temp->m[i], r);
2455    }
2456    Print(";\n");
2457    PrintLn();
2458  }
2459   PrintS("</id_KillSquares>\n");
2460#endif
2461
2462  return temp;
2463}
2464
2465
2466
2467
2468#endif
Note: See TracBrowser for help on using the repository browser.