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

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