source: git/libpolys/polys/nc/sca.cc @ 32d07a5

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