source: git/kernel/sca.cc @ b1a5c1

spielwiese
Last change on this file since b1a5c1 was b1a5c1, checked in by Hans Schönemann <hannes@…>, 15 years ago
*hannes: some fixes for mres (sheafcoh_s.tst) git-svn-id: file:///usr/local/Singular/svn/trunk@10772 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.19 2008-06-20 16:54:45 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
938    strat->enterS(h, pos, strat, -1);
939
940    if (h.lcm!=NULL) pLmFree(h.lcm);
941  }
942
943
944}
945
946#ifndef NDEBUG
947
948// set it here if needed.
949#define MYTEST 0
950
951#else
952
953#define MYTEST 0
954
955#endif
956
957#define OUTPUT 0
958
959
960ideal sca_gr_bba(const ideal F, const ideal Q, const intvec *, const intvec *, kStrategy strat)
961{
962#if MYTEST
963//   PrintS("<sca_gr_bba>\n");
964#endif
965
966  assume(rIsSCA(currRing));
967
968#ifndef NDEBUG
969  idTest(F);
970  idTest(Q);
971#endif
972
973#ifdef HAVE_PLURAL
974#if MYTEST
975  PrintS("currRing: \n");
976  rWrite(currRing);
977#ifdef RDEBUG
978  rDebugPrint(currRing);
979#endif
980
981  PrintS("F: \n");
982  idPrint(F);
983  PrintS("Q: \n");
984  idPrint(Q);
985#endif
986#endif
987
988
989  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
990  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
991
992  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
993  ideal tempQ = Q;
994
995  if(Q == currQuotient)
996    tempQ = SCAQuotient(currRing);
997
998  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
999
1000  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
1001
1002  strat->homog = strat->homog && bIdHomog;
1003
1004  assume( strat->homog == bIdHomog );
1005
1006#if MYTEST
1007  {
1008    Print("ideal tempF: \n");
1009    idPrint(tempF);
1010    Print("ideal tempQ: \n");
1011    idPrint(tempQ);
1012  }
1013#endif
1014
1015  int srmax, lrmax;
1016  int olddeg, reduc;
1017  int red_result = 1;
1018//  int hilbeledeg = 1, minimcnt = 0;
1019  int hilbcount = 0;
1020
1021  initBuchMoraCrit(strat); // set Gebauer, honey, sugarCrit
1022
1023  nc_gr_initBba(tempF,strat); // set enterS, red, initEcart, initEcartPair
1024
1025  initBuchMoraPos(strat);
1026
1027
1028  // ?? set spSpolyShort, reduce ???
1029
1030  initBuchMora(tempF, tempQ, strat); // SCAQuotient(currRing) instead of Q == squares!!!!!!!
1031
1032  strat->posInT=posInT110; // !!!
1033
1034  srmax = strat->sl;
1035  reduc = olddeg = lrmax = 0;
1036
1037
1038  /* compute------------------------------------------------------- */
1039  for(; strat->Ll >= 0;
1040#ifdef KDEBUG
1041    strat->P.lcm = NULL,
1042#endif
1043    kTest(strat)
1044    )
1045  {
1046    if (strat->Ll > lrmax) lrmax =strat->Ll;/*stat.*/
1047
1048#ifdef KDEBUG
1049    if (TEST_OPT_DEBUG) messageSets(strat);
1050#endif
1051
1052    if (strat->Ll== 0) strat->interpt=TRUE;
1053
1054    if (TEST_OPT_DEGBOUND
1055    && ((strat->honey
1056    && (strat->L[strat->Ll].ecart+pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))
1057       || ((!strat->honey) && (pFDeg(strat->L[strat->Ll].p,currRing)>Kstd1_deg))))
1058    {
1059      /*
1060      *stops computation if
1061      * 24 IN test and the degree +ecart of L[strat->Ll] is bigger then
1062      *a predefined number Kstd1_deg
1063      */
1064      while (strat->Ll >= 0) deleteInL(strat->L,&strat->Ll,strat->Ll,strat);
1065      break;
1066    }
1067
1068    /* picks the last element from the lazyset L */
1069    strat->P = strat->L[strat->Ll];
1070    strat->Ll--;
1071
1072    //kTest(strat);
1073
1074//    assume(pNext(strat->P.p) != strat->tail); // !???
1075    if(strat->P.IsNull()) continue;
1076
1077
1078    if( pNext(strat->P.p) == strat->tail )
1079    {
1080      // deletes the int spoly and computes SPoly
1081      pLmFree(strat->P.p); // ???
1082      strat->P.p = nc_CreateSpoly(strat->P.p1, strat->P.p2, currRing);
1083    }
1084
1085    if(strat->P.IsNull()) continue;
1086
1087//     poly save = NULL;
1088//
1089//     if(pNext(strat->P.p) != NULL)
1090//       save = p_Copy(strat->P.p, currRing);
1091
1092    strat->initEcart(&strat->P); // remove it?
1093
1094    if (TEST_OPT_PROT)
1095      message((strat->honey ? strat->P.ecart : 0) + strat->P.pFDeg(), &olddeg,&reduc,strat, red_result);
1096
1097    // reduction of the element chosen from L wrt S
1098    strat->red(&strat->P,strat);
1099
1100    if(strat->P.IsNull()) continue;
1101
1102    addLObject(strat->P, strat);
1103
1104    if (strat->sl > srmax) srmax = strat->sl;
1105
1106    const poly save = strat->P.p;
1107
1108#ifdef PDEBUG
1109      p_Test(save, currRing);
1110#endif
1111    assume( save != NULL );
1112
1113    // SCA Specials:
1114
1115    {
1116      const poly pNext = pNext(save);
1117
1118      if( pNext != NULL )
1119      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
1120      if( p_GetExp(save, i, currRing) != 0 )
1121      {
1122        assume(p_GetExp(save, i, currRing) == 1);
1123
1124        const poly tt = sca_pp_Mult_xi_pp(i, pNext, currRing);
1125
1126#ifdef PDEBUG
1127        p_Test(tt, currRing);
1128#endif
1129
1130        if( tt == NULL) continue;
1131
1132        LObject h(tt); // h = x_i * P
1133
1134        if (TEST_OPT_INTSTRATEGY)
1135        {
1136//           h.pCleardenom(); // also does a pContent
1137          pContent(h.p);
1138        }
1139        else
1140        {
1141          h.pNorm();
1142        }
1143
1144        strat->initEcart(&h);
1145
1146
1147//         if (pOrdSgn==-1)
1148//         {
1149//           cancelunit(&h);  // tries to cancel a unit
1150//           deleteHC(&h, strat);
1151//         }
1152
1153//         if(h.IsNull()) continue;
1154
1155//         if (TEST_OPT_PROT)
1156//           message((strat->honey ? h.ecart : 0) + h.pFDeg(), &olddeg, &reduc, strat, red_result);
1157
1158//         strat->red(&h, strat); // wrt S
1159//         if(h.IsNull()) continue;
1160
1161//         poly save = p_Copy(h.p, currRing);
1162
1163        int pos;
1164
1165        if (strat->Ll==-1)
1166          pos =0;
1167        else
1168          pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1169
1170        h.sev = pGetShortExpVector(h.p);
1171        enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
1172
1173  //       h.p = save;
1174  //       addLObject(h, strat);
1175
1176//         if (strat->sl > srmax) srmax = strat->sl;
1177      }
1178
1179      // p_Delete( &save, currRing );
1180    }
1181
1182
1183  } // for(;;)
1184
1185
1186#ifdef KDEBUG
1187  if (TEST_OPT_DEBUG) messageSets(strat);
1188#endif
1189
1190  /* release temp data-------------------------------- */
1191  exitBuchMora(strat);
1192
1193  if (TEST_OPT_WEIGHTM)
1194  {
1195    pFDeg=pFDegOld;
1196    pLDeg=pLDegOld;
1197    if (ecartWeights)
1198    {
1199      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(int));
1200      ecartWeights=NULL;
1201    }
1202  }
1203
1204  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
1205
1206  if (tempQ!=NULL) updateResult(strat->Shdl,tempQ,strat);
1207
1208  id_Delete(&tempF, currRing);
1209
1210
1211  /* complete reduction of the standard basis--------- */
1212  if (TEST_OPT_REDSB){
1213//     completeReduce(strat); // ???
1214
1215    ideal I = strat->Shdl;
1216    ideal erg = kInterRed(I,tempQ);
1217    assume(I!=erg);
1218    id_Delete(&I, currRing);
1219    strat->Shdl = erg;
1220  }
1221
1222
1223#if MYTEST
1224//   PrintS("</sca_gr_bba>\n");
1225#endif
1226
1227  return (strat->Shdl);
1228}
1229
1230
1231// should be used only inside nc_SetupQuotient!
1232// Check whether this our case:
1233//  1. rG is  a commutative polynomial ring \otimes anticommutative algebra
1234//  2. factor ideal rGR->qideal contains squares of all alternating variables.
1235//
1236// if yes, make rGR a super-commutative algebra!
1237// NOTE: Factors of SuperCommutative Algebras are supported this way!
1238//
1239//  rG == NULL means that there is no separate base G-algebra in this case take rGR == rG
1240bool sca_SetupQuotient(ring rGR, ring rG)
1241{
1242//   return false; // test Plural
1243
1244  //////////////////////////////////////////////////////////////////////////
1245  // checks...
1246  //////////////////////////////////////////////////////////////////////////
1247  if( rG == NULL )
1248    rG = rGR;
1249
1250  assume(rGR != NULL);
1251  assume(rG  != NULL);
1252  assume(rIsPluralRing(rG));
1253
1254
1255#if MYTEST
1256   PrintS("sca_SetupQuotient(rGR, rG)");
1257#endif
1258
1259  const int N = rG->N;
1260
1261  if(N < 2)
1262    return false;
1263
1264
1265//  if( (ncRingType(rG) != nc_skew) || (ncRingType(rG) != nc_comm) )
1266//    return false;
1267
1268#if OUTPUT
1269  PrintS("sca_SetupQuotient: qring?\n");
1270#endif
1271
1272  if(rGR->qideal == NULL) // there will be a factor!
1273    return false;
1274
1275#if OUTPUT
1276  PrintS("sca_SetupQuotient: qideal!!!\n");
1277#endif
1278
1279  if((rG->qideal != NULL) && (rG != rGR) ) // we cannot change from factor to factor at the time, sorry!
1280    return false;
1281
1282
1283  int iAltVarEnd = -1;
1284  int iAltVarStart   = N+1;
1285
1286  const ring rBase = rG->GetNC()->basering;
1287  const matrix C   = rG->GetNC()->C; // live in rBase!
1288
1289#if OUTPUT
1290  PrintS("sca_SetupQuotient: AltVars?!\n");
1291#endif
1292
1293  for(int i = 1; i < N; i++)
1294  {
1295    for(int j = i + 1; j <= N; j++)
1296    {
1297      assume(MATELEM(C,i,j) != NULL); // after CallPlural!
1298      number c = p_GetCoeff(MATELEM(C,i,j), rBase);
1299
1300      if( n_IsMOne(c, rBase) )
1301      {
1302        if( i < iAltVarStart)
1303          iAltVarStart = i;
1304
1305        if( j > iAltVarEnd)
1306          iAltVarEnd = j;
1307      } else
1308      {
1309        if( !n_IsOne(c, rBase) )
1310        {
1311#if OUTPUT
1312          Print("Wrong Coeff at: [%d, %d]\n", i, j);
1313#endif
1314#if MYTEST
1315           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1316#endif
1317          return false;
1318        }
1319      }
1320    }
1321  }
1322
1323#if MYTEST
1324  Print("AltVars?1: [%d, %d]\n", iAltVarStart, iAltVarEnd);
1325#endif
1326
1327
1328  if( (iAltVarEnd == -1) || (iAltVarStart == (N+1)) )
1329    return false; // either no alternating varables, or a single one => we are in commutative case!
1330
1331
1332  for(int i = 1; i < N; i++)
1333  {
1334    for(int j = i + 1; j <= N; j++)
1335    {
1336      assume(MATELEM(C,i,j) != NULL); // after CallPlural!
1337      number c = p_GetCoeff(MATELEM(C,i,j), rBase);
1338
1339      if( (iAltVarStart <= i) && (j <= iAltVarEnd) ) // S <= i < j <= E
1340      { // anticommutative part
1341        if( !n_IsMOne(c, rBase) )
1342        {
1343#ifdef PDEBUG
1344#if OUTPUT
1345           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1346#endif
1347#endif
1348
1349#if MYTEST
1350           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1351#endif
1352          return false;
1353        }
1354      } else
1355      { // should commute
1356        if( !n_IsOne(c, rBase) )
1357        {
1358#ifdef PDEBUG
1359#if OUTPUT
1360          Print("Wrong Coeff at: [%d, %d]\n", i, j);
1361#endif
1362#endif
1363#if MYTEST
1364           Print("Wrong Coeff at: [%d, %d]\n", i, j);
1365#endif
1366          return false;
1367        }
1368      }
1369    }
1370  }
1371
1372#if MYTEST
1373  Print("AltVars!?: [%d, %d]\n", iAltVarStart, iAltVarEnd);
1374#endif
1375
1376  assume( 1            <= iAltVarStart );
1377  assume( iAltVarStart < iAltVarEnd   );
1378  assume( iAltVarEnd   <= N            );
1379
1380
1381
1382  bool bWeChangeRing = false;
1383  // for sanity
1384  ring rSaveRing = currRing;
1385
1386  if(rSaveRing != rG)
1387  {
1388    rChangeCurrRing(rG);
1389    bWeChangeRing = true;
1390  }
1391
1392
1393  assume(rGR->qideal != NULL);
1394//  assume(rG->qideal == NULL); // ?
1395
1396  const ideal idQuotient = rGR->qideal;
1397
1398
1399#if MYTEST
1400  Print("Analyzing quotient ideal:\n");
1401  idPrint(idQuotient); // in rG!!!
1402#endif
1403
1404
1405  // check for
1406  // y_{iAltVarStart}^2, y_{iAltVarStart+1}^2, \ldots, y_{iAltVarEnd}^2  (iAltVarEnd > iAltVarStart)
1407  // to be within quotient ideal.
1408
1409  bool bSCA = true;
1410
1411  for ( int i = iAltVarStart; (i <= iAltVarEnd) && bSCA; i++ )
1412  {
1413    poly square = p_ISet(1, rG);
1414    p_SetExp(square, i, 2, rG); // square = var(i)^2.
1415    p_Setm(square, rG);
1416
1417    // square = NF( var(i)^2 | Q )
1418    // NOTE: rSaveRing == currRing now!
1419    // NOTE: there is no better way to check this in general!
1420    square = kNF(idQuotient, NULL, square, 0, 0); // must ran in currRing == rG!
1421
1422    if( square != NULL ) // var(i)^2 is not in Q?
1423    {
1424      p_Delete(&square, rG);
1425      bSCA = false;
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      }
1712      else
1713      {
1714
1715        strat->P.pNorm();
1716        if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1717          strat->P.p = redtailBba(&(strat->P),pos-1,strat, withT);
1718      }
1719
1720
1721#ifdef KDEBUG
1722      if (TEST_OPT_DEBUG){PrintS(" ns:");strat->P.wrp();PrintLn();}
1723#endif
1724
1725//       // min_std stuff
1726//       if ((strat->P.p1==NULL) && (strat->minim>0))
1727//       {
1728//         if (strat->minim==1)
1729//         {
1730//           strat->M->m[minimcnt]=p_Copy(strat->P.p,currRing,strat->tailRing);
1731//           p_Delete(&strat->P.p2, currRing, strat->tailRing);
1732//         }
1733//         else
1734//         {
1735//           strat->M->m[minimcnt]=strat->P.p2;
1736//           strat->P.p2=NULL;
1737//         }
1738//         if (strat->tailRing!=currRing && pNext(strat->M->m[minimcnt])!=NULL)
1739//           pNext(strat->M->m[minimcnt])
1740//             = strat->p_shallow_copy_delete(pNext(strat->M->m[minimcnt]),
1741//                                            strat->tailRing, currRing,
1742//                                            currRing->PolyBin);
1743//         minimcnt++;
1744//       }
1745
1746      // enter into S, L, and T
1747      if(withT)
1748        enterT(strat->P, strat);
1749
1750      // L
1751      enterpairs(strat->P.p,strat->sl,strat->P.ecart,pos,strat, strat->tl);
1752
1753      // posInS only depends on the leading term
1754      strat->enterS(strat->P, pos, strat, strat->tl);
1755
1756//       if (hilb!=NULL) khCheck(Q,w,hilb,hilbeledeg,hilbcount,strat);
1757
1758//      Print("[%d]",hilbeledeg);
1759      if (strat->P.lcm!=NULL) pLmFree(strat->P.lcm);
1760
1761      if (strat->sl>srmax) srmax = strat->sl;
1762
1763      // //////////////////////////////////////////////////////////
1764      // SCA:
1765      const poly pSave = strat->P.p;
1766      const poly pNext = pNext(pSave);
1767
1768//       if(0)
1769      for( unsigned int i = m_iFirstAltVar; i <= m_iLastAltVar; i++ )
1770      if( p_GetExp(pSave, i, currRing) != 0 )
1771      {
1772        assume(p_GetExp(pSave, i, currRing) == 1);
1773        if (pNext!=NULL)
1774        {
1775          const poly pNew = sca_pp_Mult_xi_pp(i, pNext, currRing);
1776
1777#ifdef PDEBUG
1778          p_Test(pNew, currRing);
1779#endif
1780
1781          if( pNew == NULL) continue;
1782
1783          LObject h(pNew); // h = x_i * strat->P
1784
1785          if (TEST_OPT_INTSTRATEGY)
1786          {
1787//            h.pCleardenom(); // also does a pContent
1788            pContent(h.p);
1789          }
1790          else
1791          {
1792            h.pNorm();
1793          }
1794
1795          strat->initEcart(&h);
1796          h.sev = pGetShortExpVector(h.p);
1797
1798          int pos;
1799          if (strat->Ll==-1)
1800            pos =0;
1801          else
1802            pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1803
1804          enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);
1805/*
1806          h.sev = pGetShortExpVector(h.p);
1807          strat->initEcart(&h);
1808
1809          h.PrepareRed(strat->use_buckets);
1810
1811          // reduction of the element choosen from L(?)
1812          red_result = strat->red(&h,strat);
1813
1814          // reduction to non-zero new poly
1815          if (red_result != 1) continue;
1816
1817
1818          int pos = posInS(strat,strat->sl,h.p,h.ecart);
1819
1820          // reduce the tail and normalize poly
1821          if (TEST_OPT_INTSTRATEGY)
1822          {
1823            h.pCleardenom();
1824            if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1825            {
1826              h.p = redtailBba(&(h),pos-1,strat, withT); // !!!
1827              h.pCleardenom();
1828            }
1829          }
1830          else
1831          {
1832            h.pNorm();
1833            if ((TEST_OPT_REDSB)||(TEST_OPT_REDTAIL))
1834              h.p = redtailBba(&(h),pos-1,strat, withT);
1835          }
1836
1837#ifdef KDEBUG
1838          if (TEST_OPT_DEBUG){PrintS(" N:");h.wrp();PrintLn();}
1839#endif
1840
1841//          h.PrepareRed(strat->use_buckets); // ???
1842
1843          h.sev = pGetShortExpVector(h.p);
1844          strat->initEcart(&h);
1845
1846          if (strat->Ll==-1)
1847            pos = 0;
1848          else
1849            pos = strat->posInL(strat->L,strat->Ll,&h,strat);
1850
1851           enterL(&strat->L,&strat->Ll,&strat->Lmax,h,pos);*/
1852
1853        }
1854      } // for all x_i \in Ann(lm(P))
1855    } // if red(P) != NULL
1856
1857//     else if (strat->P.p1 == NULL && strat->minim > 0)
1858//     {
1859//       p_Delete(&strat->P.p2, currRing, strat->tailRing);
1860//     }
1861
1862// #ifdef KDEBUG
1863    memset(&(strat->P), 0, sizeof(strat->P));
1864// #endif
1865
1866    //kTest_TS(strat); // T is not used: cannot use this test
1867
1868//     Print("\n$\n");
1869
1870  }
1871
1872#ifdef KDEBUG
1873  if (TEST_OPT_DEBUG) messageSets(strat);
1874#endif
1875
1876  /* complete reduction of the standard basis--------- */
1877  if (TEST_OPT_REDSB) completeReduce(strat);
1878
1879  /* release temp data-------------------------------- */
1880  id_Delete(&tempF, currRing);
1881
1882  exitBuchMora(strat);
1883
1884  if (TEST_OPT_WEIGHTM)
1885  {
1886    pRestoreDegProcs(pFDegOld, pLDegOld);
1887    if (ecartWeights)
1888    {
1889      omFreeSize((ADDRESS)ecartWeights,(pVariables+1)*sizeof(short));
1890      ecartWeights=NULL;
1891    }
1892  }
1893
1894  if (TEST_OPT_PROT) messageStat(srmax,lrmax,hilbcount,strat);
1895
1896//   if (Q!=NULL) updateResult(strat->Shdl,Q,strat);
1897//    PrintS("</sca>\n");
1898
1899  return (strat->Shdl);
1900}
1901
1902
1903// //////////////////////////////////////////////////////////////////////////////
1904// sca mora...
1905
1906// returns TRUE if mora should use buckets, false otherwise
1907static BOOLEAN kMoraUseBucket(kStrategy strat)
1908{
1909#ifdef MORA_USE_BUCKETS
1910  if (TEST_OPT_NOT_BUCKETS)
1911    return FALSE;
1912  if (strat->red == redFirst)
1913  {
1914#ifdef NO_LDEG
1915    if (!strat->syzComp)
1916      return TRUE;
1917#else
1918    if ((strat->homog || strat->honey) && !strat->syzComp)
1919      return TRUE;
1920#endif
1921  }
1922  else
1923  {
1924    assume(strat->red == redEcart);
1925    if (strat->honey && !strat->syzComp)
1926      return TRUE;
1927  }
1928#endif
1929  return FALSE;
1930}
1931
1932
1933#ifdef HAVE_ASSUME
1934static int sca_mora_count = 0;
1935static int sca_mora_loop_count;
1936#endif
1937
1938// ideal sca_mora (ideal F, ideal Q, intvec *w, intvec *, kStrategy strat)
1939ideal sca_mora(const ideal F, const ideal Q, const intvec *w, const intvec *, kStrategy strat)
1940{
1941  assume(rIsSCA(currRing));
1942
1943  const unsigned int m_iFirstAltVar = scaFirstAltVar(currRing);
1944  const unsigned int m_iLastAltVar  = scaLastAltVar(currRing);
1945
1946  ideal tempF = id_KillSquares(F, m_iFirstAltVar, m_iLastAltVar, currRing);
1947
1948  ideal tempQ = Q;
1949
1950  if(Q == currQuotient)
1951    tempQ = SCAQuotient(currRing);
1952
1953  bool bIdHomog = id_IsSCAHomogeneous(tempF, NULL, NULL, currRing); // wCx == wCy == NULL!
1954
1955  assume( !bIdHomog || strat->homog ); //  bIdHomog =====[implies]>>>>> strat->homog
1956
1957  strat->homog = strat->homog && bIdHomog;
1958
1959#ifdef PDEBUG
1960  assume( strat->homog == bIdHomog );
1961#endif /*PDEBUG*/
1962
1963#ifdef HAVE_ASSUME
1964  sca_mora_count++;
1965  sca_mora_loop_count = 0;
1966#endif
1967
1968#ifdef KDEBUG
1969  om_Opts.MinTrack = 5;
1970#endif
1971  int srmax;
1972  int lrmax = 0;
1973  int olddeg = 0;
1974  int reduc = 0;
1975  int red_result = 1;
1976//  int hilbeledeg=1;
1977  int hilbcount=0;
1978
1979  strat->update = TRUE;
1980  /*- setting global variables ------------------- -*/
1981  initBuchMoraCrit(strat);
1982//   initHilbCrit(F,NULL,&hilb,strat); // no Q!
1983  initMora(tempF, strat);
1984  initBuchMoraPos(strat);
1985  /*Shdl=*/initBuchMora(tempF, tempQ, strat); // temp Q, F!
1986//   if (TEST_OPT_FASTHC) missingAxis(&strat->lastAxis,strat);
1987  /*updateS in initBuchMora has Hecketest
1988  * and could have put strat->kHEdgdeFound FALSE*/
1989#if 0
1990  if (ppNoether!=NULL)
1991  {
1992    strat->kHEdgeFound = TRUE;
1993  }
1994  if (strat->kHEdgeFound && strat->update)
1995  {
1996    firstUpdate(strat);
1997    updateLHC(strat);
1998    reorderL(strat);
1999  }
2000  if (TEST_OPT_FASTHC && (strat->lastAxis) && strat->posInLOldFlag)
2001  {
2002    strat->posInLOld = strat->posInL;
2003    strat->posInLOldFlag = FALSE;
2004    strat->posInL = posInL10;
2005    updateL(strat);
2006    reorderL(strat);
2007  }
2008#endif
2009
2010  srmax = strat->sl;
2011  kTest_TS(strat);
2012
2013  strat->use_buckets = kMoraUseBucket(strat);
2014  /*- compute-------------------------------------------*/
2015
2016#undef HAVE_TAIL_RING
2017
2018#ifdef HAVE_TAIL_RING
2019//  if (strat->homog && strat->red == redFirst)
2020//     kStratInitChangeTailRing(strat);
2021#endif
2022
2023
2024  while (strat->Ll >= 0)
2025  {
2026#ifdef HAVE_ASSUME
2027    sca_mora_loop_count++;
2028#endif
2029    if (lrmax< strat->Ll) lrmax=strat->Ll; /*stat*/
2030    //test_int_std(strat->kIdeal);
2031#ifdef KDEBUG
2032    if (TEST_OPT_DEBUG) messageSets(strat);
2033#endif
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.