source: git/kernel/sca.cc @ 0f7420d

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