source: git/kernel/sca.cc @ 151000

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