source: git/kernel/f5gb.cc @ a05c71

spielwiese
Last change on this file since a05c71 was a05c71, checked in by Christian Eder, 15 years ago
termination tests included in procedure findReducers() git-svn-id: file:///usr/local/Singular/svn/trunk@12090 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 70.7 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/*
5* ABSTRACT: f5gb interface
6*/
7
8#include "mod2.h"
9#ifdef HAVE_F5
10#include <unistd.h>
11#include "structs.h"
12#include "kutil.h"
13#include "omalloc.h"
14#include "polys.h"
15#include "p_polys.h"
16#include "p_Procs.h"
17#include "ideals.h"
18#include "febase.h"
19#include "kstd1.h"
20#include "khstd.h"
21#include "kbuckets.h"
22#include "weight.h"
23#include "intvec.h"
24#include "pInline1.h"
25#include "f5gb.h"
26#include "f5data.h"
27#include "f5lists.h"
28#include "timer.h"
29int notInG              =   0;
30int numberOfRules       =   0;
31int reductionsToZero    =   0;
32int reductionTime       =   0;
33int spolsTime           =   0;
34int highestDegree       =   0;
35int degreeBound         =   0;
36/*
37====================================================================
38sorting ideals by decreasing total degree "left" and "right" are the
39pointer of the first and last polynomial in the considered ideal
40====================================================================
41*/
42void qsortDegree(poly* left, poly* right) {
43    poly* ptr1 = left;
44    poly* ptr2 = right;
45    poly p1,p2;
46    p2 = *(left + (right - left >> 1));
47    do {
48            while(pTotaldegree(*ptr1, currRing) < pTotaldegree(p2, currRing)) {
49                    ptr1++;
50            } 
51            while(pTotaldegree(*ptr2, currRing) > pTotaldegree(p2,currRing)) {
52                    ptr2--;
53            }
54            if(ptr1 > ptr2) {
55                    break;
56            }
57            p1    = *ptr1;
58            *ptr1 = *ptr2;
59            *ptr2 = p1;
60    } while(++ptr1 <= --ptr2);
61
62    if(left < ptr2) {
63            qsortDegree(left,ptr2);
64    }
65    if(ptr1 < right) {
66            qsortDegree(ptr1,right);
67    }
68}
69
70/*!
71 * ======================================================================
72 * builds the sum of the entries of the exponent vectors, i.e. the degree
73 * of the corresponding monomial
74 * ======================================================================
75*/
76long sumVector(int* v, int k) {
77    int i;
78    long sum =  0;
79    for(i=1; i<=k; i++) {
80        Print("%d\n",v[i]);
81        Print("%ld\n",sum);
82        sum =   sum + v[i];
83    }
84    return sum;
85}
86
87/*!
88==========================================================================
89compares monomials, i.e. divisibility tests for criterion 1 and criterion 2
90==========================================================================
91*/
92bool compareMonomials(int* m1, int** m2, int numberOfRules, int k) {
93    int i,j;
94    long sumM1  =   sumVector(m1,k);
95    //int k   =   sizeof(m1) / sizeof(int);
96    for(i=0; i<numberOfRules; i++) {
97        if(sumVector(m2[i],k) <= sumM1) {
98            for(j=1; j<=k; j++) {
99                if(m1[j]>m2[i][j]) {
100                    return true;
101                }
102            }   
103        }
104    }
105    return false;
106}
107
108/*
109==================================================
110computes incrementally gbs of subsets of the input
111gb{f_m} -> gb{f_m,f_(m-1)} -> gb{f_m,...,f_1}
112==================================================
113*/
114LList* F5inc(int i, poly f_i, LList* gPrev, LList* reducers, ideal gbPrev, poly ONE, LTagList* lTag, RList* rules, RTagList* rTag, int termination) {
115    //Print("in f5inc\n");           
116    //pWrite(rules->getFirst()->getRuleTerm());
117    int j;
118    //Print("%p\n",gPrev->getFirst());
119    //pWrite(gPrev->getFirst()->getPoly());
120    poly tempNF =   kNF(gbPrev,currQuotient,f_i);
121    f_i         =   tempNF;
122    //gPrev->insert(ONE,i,f_i);
123    gPrev->insert(ONE,gPrev->getLength()+1,f_i);
124    // tag the first element in gPrev of the current index for findReductor()
125    lTag->setFirstCurrentIdx(gPrev->getLast());
126    //Print("1st gPrev: ");
127    //pWrite(gPrev->getFirst()->getPoly());
128    //Print("2nd gPrev: ");
129    //pWrite(gPrev->getFirst()->getNext()->getPoly());
130    //pWrite(gPrev->getFirst()->getNext()->getPoly());   
131    CListOld* critPairs        =   new CListOld();
132    CNode* critPairsMinDeg  =   new CNode();   
133    // computation of critical pairs with checking of criterion 1 and criterion 2 and saving them
134    // in the list critPairs
135    criticalPair(gPrev, critPairs, lTag, rTag, rules);
136    static LList* sPolyList        =   new LList();
137    //sPolyList->print();
138    // labeled polynomials which have passed reduction() and have to be added to list gPrev
139    static LList* completed        =   new LList();
140    // the reduced labeled polynomials which are returned from subalgorithm reduction()
141    static LList* reducedLPolys     =   new LList();
142    // while there are critical pairs to be further checked and deleted/computed
143    while(NULL != critPairs->getFirst()) { 
144        // critPairs->getMinDeg() deletes the first elements of minimal degree from
145        // critPairs, thus the while loop is not infinite.
146        critPairsMinDeg =   critPairs->getMinDeg();
147        // adds all to be reduced S-polynomials in the list sPolyList and adds
148        // the corresponding rules to the list rules
149        // NOTE: inside there is a second check of criterion 2 if new rules are
150        // added
151        //int timer4  =   initTimer();
152        //startTimer();
153        //critPairs->print();
154        computeSPols(critPairsMinDeg,rTag,rules,sPolyList);
155        //timer4  =   getTimer();
156        //Print("SPOLS TIMER: %d\n",timer4);
157        //spolsTime  =   spolsTime  +   timer4;
158        // DEBUG STUFF FOR SPOLYLIST
159        LNode* temp     =   sPolyList->getFirst();
160        //while(NULL != temp && NULL != temp->getLPoly()) {
161            //Print("Spolylist element: ");
162            //pWrite(temp->getPoly());
163            //temp    =   temp->getNext();
164        //}
165        // reduction process of new S-polynomials and also adds new critical pairs to critPairs
166        //int timer3  =   initTimer();
167        //startTimer();
168        //sPolyList->print();
169        //reduction(sPolyList, critPairs, gPrev, rules, lTag, rTag, gbPrev);
170        newReduction(sPolyList, critPairs, gPrev, reducers, rules, lTag, rTag, gbPrev, termination);
171        //timer3      =  getTimer();
172        //reductionTime = reductionTime + timer3;
173        //Print("REDUCTION TIMER: %d\n",timer3);
174        // DEBUG STUFF FOR GPREV
175        //temp    =   gPrev->getFirst();
176        //int number  =   1;
177        //Print("\n\n");
178        //while(NULL != temp) {
179        //    Print("%d.  ",number);
180        //    pWrite(temp->getPoly());
181        //    temp    =   temp->getNext();
182        //    number++;
183        //    Print("\n");
184        //}
185        //sleep(5);
186   
187    }
188    //Print("REDUCTION DONE\n");
189    //Print("%p\n",rules->getFirst());
190    //Print("%p\n",rTag->getFirst());
191    //if(rules->getFirst() != rTag->getFirst()) {
192        //Print("+++++++++++++++++++++++++++++++++++++NEW RULES+++++++++++++++++++++++++++++++++++++\n");
193        //rTag->insert(rules->getFirst());
194    //}
195    //else {
196        //Print("+++++++++++++++++++++++++++++++++++NO NEW RULES++++++++++++++++++++++++++++++++++++\n");
197    //}
198    lTag->insert(lTag->getFirstCurrentIdx());
199    //Print("INDEX: %d\n",tempTag->getIndex());
200    //pWrite(tempTag->getPoly());
201    //Print("COMPLETED FIRST IN F5INC: \n");
202    //Print("1st gPrev: ");
203    //pWrite(gPrev->getFirst()->getPoly());
204    //Print("2nd gPrev: ");
205    //pWrite(gPrev->getFirst()->getNext()->getPoly());
206    //Print("3rd gPrev: ");
207    //pWrite(gPrev->getFirst()->getNext()->getNext()->getPoly());
208    //delete sPolyList;
209    //critPairs->print();
210    delete critPairs;
211    //Print("IN F5INC\n");
212    /*
213    Print("\n\n\nRULES: \n");
214        RNode* tempR    =   rules->getFirst();
215        Print("%p\n",tempR);
216        int t   = 1;
217        while(NULL != tempR) {
218            Print("ADDRESS OF %d RNODE: %p\n",t,tempR);
219            Print("ADDRESS OF RULE: %p\n",tempR->getRule());
220            pWrite(tempR->getRuleTerm());
221            Print("ADDRESS OF TERM: %p\n",tempR->getRuleTerm());
222            Print("%d\n\n",tempR->getRuleIndex());
223            tempR   =   tempR->getNext();
224            t++;
225        }
226    */
227    //gPrev->print();
228    //Print("COMPLETE REDUCTION TIME UNTIL NOW: %d\n",reductionTime);
229    //Print("COMPLETE SPOLS TIME UNTIL NOW:     %d\n",spolsTime);
230    degreeBound = 0;
231    return gPrev;
232}
233
234
235
236/*
237================================================================
238computes a list of critical pairs for the next reduction process
239first element in gPrev is always the newest element which must
240build critical pairs with all other elements in gPrev
241================================================================
242*/
243inline void criticalPair(LList* gPrev, CListOld* critPairs, LTagList* lTag, RTagList* rTag, RList* rules) {
244    //Print("IN CRITPAIRS\n");
245    // initialization for usage in pLcm()
246    number nOne         =   nInit(1);
247    LNode* newElement   =   gPrev->getLast();
248    LNode* temp         =   gPrev->getFirst();
249    poly u1             =   pOne();
250    poly u2             =   pOne();
251    poly lcm            =   pOne();
252    poly t              =   pHead(newElement->getPoly());
253    RuleOld* testedRuleOld    =   NULL;
254    if(NULL != rules->getFirst()) {
255        testedRuleOld  =   rules->getFirst()->getRuleOld();
256    }
257    // computation of critical pairs
258    while( gPrev->getLast() != temp) {
259        pLcm(newElement->getPoly(), temp->getPoly(), lcm);
260        pSetCoeff(lcm,nOne);
261        // computing factors u2 for new labels
262        u1 = pDivide(lcm,t);
263        if(NULL == u1) {
264            break;
265        }
266        pSetCoeff(u1,nOne);
267        u2 = pDivide(lcm,pHead(temp->getPoly()));
268        pSetCoeff(u2,nOne);
269        // testing both new labels by the F5 Criterion
270        if(!criterion2(gPrev->getFirst()->getIndex(), u2, temp, rules, rTag)
271           && !criterion2(gPrev->getFirst()->getIndex(), u1, newElement, rules, rTag) 
272           && !criterion1(gPrev,u1,newElement,lTag) && !criterion1(gPrev,u2,temp,lTag)) {
273            // if they pass the test, add them to CListOld critPairs, having the LPolyOld with greater
274            // label as first element in the CPairOld
275            if(newElement->getIndex() == temp->getIndex() && 
276            -1 == pLmCmp(ppMult_qq(u1, newElement->getTerm()),ppMult_qq(u2, temp->getTerm()))) {
277                CPairOld* cp   =   new CPairOld(pDeg(ppMult_qq(u2,pHead(temp->getPoly()))), u2, 
278                                temp->getLPolyOld(), u1, newElement->getLPolyOld(), testedRuleOld);                   
279                critPairs->insert(cp);
280            }
281            else {
282                CPairOld* cp   =   new CPairOld(pDeg(ppMult_qq(u2,pHead(temp->getPoly()))), u1, 
283                                newElement->getLPolyOld(), u2, temp->getLPolyOld(), testedRuleOld);                   
284                critPairs->insert(cp);
285            }
286        }
287        else {
288        }
289        temp    =   temp->getNext();
290    }
291}
292
293
294
295
296
297
298
299/*
300========================================
301Criterion 1, i.e. Faugere's F5 Criterion
302========================================
303*/
304inline bool criterion1(LList* gPrev, poly t, LNode* l, LTagList* lTag) {
305    // starts at the first element in gPrev with index = (index of l)-1, these tags are saved in lTag
306        int idx =   l->getIndex();
307    int i;
308    if(idx == 1) {
309        //Print("HIER\n");
310        return false;
311    }
312    else {
313        LNode* testNode =   gPrev->getFirst();
314        // save the monom t1*label_term(l) as it is tested various times in the following
315        poly u1 = ppMult_qq(t,l->getTerm());
316        //Print("------------------------------IN CRITERION 1-----------------------------------------\n");
317        //Print("TESTED ELEMENT: ");
318        //pWrite(l->getPoly());
319        //pWrite(l->getTerm());
320        //pWrite(ppMult_qq(t,l->getTerm()));
321        //Print("%d\n\n",l->getIndex());
322        while(testNode->getIndex() < idx) { // && NULL != testNode->getLPoly()) {
323            //pWrite(testNode->getPoly());
324            //Print("%d\n",testNode->getIndex());
325            if(pLmDivisibleByNoComp(testNode->getPoly(),u1)) {
326                //Print("Criterion 1 NOT passed!\n");
327                if(idx != gPrev->getLast()->getIndex()) {
328                    //Print("DOCH!\n");
329                }
330                return true;
331            }
332            //pWrite(testNode->getNext()->getPoly());
333            testNode    =   testNode->getNext();
334        }
335        /*
336        ideal testId    =   idInit(idx-1,1);
337        for(i=0;i<idx-1;i++) {
338            testId->m[i]  =   pHead(testNode->getPoly());
339            testNode        =   testNode->getNext();
340        }
341        poly temp   =   kNF(testId,currQuotient,u1);
342        //pWrite(temp);
343        for(i=0;i<IDELEMS(testId);i++) {
344            testId->m[i]    =   NULL;
345        }
346        idDelete(&testId);
347        if(NULL == temp) {
348            //if(l->getIndex() != gPrev->getFirst()->getIndex()) {
349            //    Print("----------------------------Criterion1 not passed----------------------------------\n");
350            //}
351            return true;
352        }
353        */
354        return false;
355    }
356}
357
358
359
360/*
361=====================================
362Criterion 2, i.e. Rewritten Criterion
363=====================================
364*/
365inline bool criterion2(int idx, poly t, LNode* l, RList* rules, RTagList* rTag) {
366    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n");
367    /* 
368    Print("RULES: \n");
369        RNode* tempR    =   rules->getFirst();
370        Print("%p\n",tempR);
371        int i   = 1;
372        while(NULL != tempR) {
373            Print("ADDRESS OF %d RNODE: %p\n",i,tempR);
374            Print("ADDRESS OF RULE: %p\n",tempR->getRule());
375            pWrite(tempR->getRuleTerm());
376            Print("ADDRESS OF TERM: %p\n",tempR->getRuleTerm());
377            Print("%d\n\n",tempR->getRuleIndex());
378            tempR   =   tempR->getNext();
379            i++;
380        }
381        //Print("TESTED ELEMENT: ");
382        //pWrite(l->getPoly());
383        //pWrite(l->getTerm());
384        //pWrite(ppMult_qq(t,l->getTerm()));
385        //Print("%d\n\n",l->getIndex());
386      */
387// start at the previously added element to gPrev, as all other elements will have the same index for sure
388    if(idx > l->getIndex()) {
389        return false;
390    }
391   
392    RNode* testNode; // =   new RNode();
393    testNode    =   rules->getFirst();
394    /*
395     if(NULL == rTag->getFirst()) {
396        if(NULL != rules->getFirst()) {
397            testNode    =   rules->getFirst();
398        }
399        else {
400            return false;
401        }
402    }
403   
404     else {
405
406        if(l->getIndex() > rTag->getFirst()->getRuleIndex()) {
407            testNode    =   rules->getFirst();
408        }
409        else {
410       //Print("HIER\n");
411            //Print("DEBUG\n");
412        //Print("L INDEX: %d\n",l->getIndex());
413            *-------------------------------------
414             * TODO: WHEN INTERREDUCING THE GB THE
415             *       INDEX OF THE PREVIOUS ELEMENTS
416             *       GETS HIGHER!
417             *-----------------------------------*
418            //testNode    =   rules->getFirst();
419            testNode    =   rTag->get(l->getIndex());
420            if(NULL == testNode) {
421                testNode    =   rules->getFirst();
422            }
423            //Print("TESTNODE ADDRESS: %p\n",testNode);
424        }
425    }
426    */
427    //testNode    =   rules->getFirst();
428        // save the monom t1*label_term(l) as it is tested various times in the following
429    poly u1 = ppMult_qq(t,l->getTerm());
430    // first element added to rTag was NULL, check for this
431    //Print("%p\n",testNode->getRule());
432    // NOTE: testNode is possibly NULL as rTag->get() returns NULL for elements of index <=1!
433    //Print("TESTNODE: %p\n",testNode);
434    //pWrite(testNode->getRuleTerm());
435    if(NULL != testNode ) {   
436        //pWrite(testNode->getRuleTerm());
437    }
438    if(NULL != testNode) {
439        if(testNode->getRuleOld() == l->getRuleOld()) {
440            //Print("%p\n%p\n",testNode->getRule(),l->getRule());
441            //Print("EQUAL\n");
442        }
443        else {
444            //Print("NOT EQUAL\n");
445        }
446    }
447    while(NULL != testNode && testNode->getRuleOld() != l->getRuleOld() 
448          && l->getIndex() == testNode->getRuleOldIndex()) {
449        //Print("%p\n",testNode);
450        //pWrite(testNode->getRuleTerm());
451        //pWrite(t);
452        //pWrite(l->getTerm());
453        //pWrite(u1);
454        //Print("%d\n",testNode->getRuleIndex());
455        if(pLmDivisibleByNoComp(testNode->getRuleOldTerm(),u1)) {
456            //Print("-----------------Criterion 2 NOT passed!-----------------------------------\n");
457            //Print("INDEX: %d\n",l->getIndex());
458            pDelete(&u1);
459    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n\n");
460            return true;
461        }
462                testNode    =   testNode->getNext();
463    }
464    //delete testNode;
465    pDelete(&u1);
466    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n\n");
467    return false;
468}
469
470
471
472/*
473=================================================================================================================
474Criterion 2, i.e. Rewritten Criterion, for its second call in computeSPols(), with added lastRuleTested parameter
475=================================================================================================================
476*/
477inline bool criterion2(poly t, LPolyOld* l, RList* rules, RuleOld* testedRuleOld) {
478    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n");
479    //Print("LAST RuleOld TESTED: %p",testedRuleOld);
480    /*
481    Print("RULES: \n");
482        RNode* tempR    =   rules->getFirst();
483        while(NULL != tempR) {
484            pWrite(tempR->getRuleTerm());
485            Print("%d\n\n",tempR->getRuleIndex());
486            tempR   =   tempR->getNext();
487        }
488        //Print("TESTED ELEMENT: ");
489        //pWrite(l->getPoly());
490        //pWrite(l->getTerm());
491        //pWrite(ppMult_qq(t,l->getTerm()));
492        //Print("%d\n\n",l->getIndex());
493    */
494// start at the previously added element to gPrev, as all other elements will have the same index for sure
495        RNode* testNode =   rules->getFirst();
496    // save the monom t1*label_term(l) as it is tested various times in the following
497    poly u1 = ppMult_qq(t,l->getTerm());
498        // first element added to rTag was NULL, check for this
499        while(NULL != testNode && testNode->getRuleOld() != testedRuleOld) {
500        //pWrite(testNode->getRuleTerm());
501        if(pLmDivisibleByNoComp(testNode->getRuleOldTerm(),u1)) {
502            //Print("--------------------------Criterion 2 NOT passed!------------------------------\n");
503            //Print("INDEX: %d\n",l->getIndex());
504            pDelete(&u1);
505    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n\n");
506            return true;
507        }
508                testNode    =   testNode->getNext();
509    }
510    pDelete(&u1);
511    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n\n");
512    return false;
513}
514
515
516
517/*
518==================================
519Computation of S-Polynomials in F5
520==================================
521*/
522void computeSPols(CNode* first, RTagList* rTag, RList* rules, LList* sPolyList) { 
523    CNode* temp         =   first;
524    poly sp     =   pInit();
525    number sign =   nInit(-1);   
526    //Print("###############################IN SPOLS##############################\n");
527    //first->print();
528
529    while(NULL != temp) {
530        //Print("JA\n");
531        // only if a new RuleOld was added since the last test in subalgorithm criticalPair()
532        //if(rules->getFirst() != rTag->getFirst()) {
533        if(!criterion2(temp->getT1(),temp->getAdLp1(),rules,temp->getTestedRuleOld())) {
534                // the second component is tested only when it has the actual index, otherwise there is
535                // no new RuleOld to test since the last test in subalgorithm criticalPair()
536                if(highestDegree < pDeg(ppMult_qq(temp->getT1(),temp->getLp1Poly()))) { 
537                    highestDegree   = pDeg(ppMult_qq(temp->getT1(),temp->getLp1Poly()));
538                    //pWrite(pHead(ppMult_qq(temp->getT1(),temp->getLp1Poly())));
539                }   
540                if(temp->getLp2Index() == temp->getLp1Index()) {
541                    if(!criterion2(temp->getT2(),temp->getAdLp2(),rules,temp->getTestedRuleOld())) {
542                        // computation of S-polynomial
543                        //poly p1 =   temp->getLp1Poly();
544                        //poly p2 =   temp->getLp2Poly();
545                        //pIter(p1);
546                        //pIter(p2);
547                        //sp  =   pAdd(ppMult_qq(temp->getT1(),p1),pMult_nn(ppMult_qq(temp->getT2(),p2),sign)); 
548                        sp      =   ksOldSpolyRedNew(ppMult_qq(temp->getT1(),temp->getLp1Poly()),
549                                         ppMult_qq(temp->getT2(),temp->getLp2Poly()));
550                        //Print("BEGIN SPOLY1\n====================\n");
551                        pNorm(sp);
552                        //pWrite(sp);
553                        //Print("END SPOLY1\n====================\n");
554                        if(NULL == sp) {
555                            // as rules consist only of pointers we need to save the labeled
556                            // S-polynomial also of a zero S-polynomial
557                            //zeroList->insert(temp->getAdT1(),temp->getLp1Index(),&sp);
558                            // origin of RuleOld can be set NULL as the labeled polynomial
559                            // will never be used again as it is zero => no problems with
560                            // further criterion2() tests and termination conditions
561                            //Print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ZERO REDUCTION~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n");
562                                                        reductionsToZero++;
563                        //Print("IN SPOLS 1\n");
564                            //Rule* rNew  =  new Rule(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
565                            //rules->insertOrdered(rNew);
566                            rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
567                            numberOfRules++;
568                            //Print("RuleOld ADDED: \n");
569                            //pWrite(rules->getFirst()->getRuleTerm());
570                            //rules->print();
571                            // as sp = NULL, delete it
572                            pDelete(&sp);
573                            //Print("HIER\n");
574                        }
575                        else {
576                        //Print("IN SPOLS 2\n");
577                            //Rule* rNew  =  new Rule(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
578                            //rules->insertOrdered(rNew);
579                            rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
580                            numberOfRules++;
581                            //Print("RuleOld ADDED: \n");
582                            //pWrite(rules->getFirst()->getRuleTerm()); 
583                            //rules->print();
584                            sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rules->getFirst()->getRuleOld());
585                            //sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rNew);
586                        }
587                        // data is saved in sPolyList or zero => delete sp
588                    }
589                }
590                else { // temp->getLp2Index() < temp->getLp1Index()
591                    // computation of S-polynomial
592                        //poly p1 =   temp->getLp1Poly();
593                        //poly p2 =   temp->getLp2Poly();
594                        //pIter(p1);
595                        //pIter(p2);
596                        //sp  =   pAdd(ppMult_qq(temp->getT1(),p1),pMult_nn(ppMult_qq(temp->getT2(),p2),sign)); 
597                    sp      =   ksOldSpolyRedNew(ppMult_qq(temp->getT1(),temp->getLp1Poly()),
598                                     ppMult_qq(temp->getT2(),temp->getLp2Poly()));
599                    //Print("BEGIN SPOLY2\n====================\n");
600                    pNorm(sp);
601                    //pWrite(sp);
602                    //Print("END SPOLY2\n====================\n");
603                    if(NULL == sp) {
604                        // zeroList->insert(temp->getAdT1(),temp->getLp1Index(),&sp);
605                        // origin of RuleOld can be set NULL as the labeled polynomial
606                        // will never be used again as it is zero => no problems with
607                        // further criterion2() tests and termination conditions
608                            //Print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ZERO REDUCTION~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n");
609                        reductionsToZero++;
610                        //Print("IN SPOLS 3\n");
611                        rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
612                        numberOfRules++;
613                        //Print("RuleOld ADDED: \n");
614                        //pWrite(rules->getFirst()->getRuleTerm());
615                        //rules->print();
616                        // as sp = NULL, delete it
617                        pDelete(&sp);
618                    }
619                    else {
620                        //Print("IN SPOLS 4\n");
621                       
622                        //////////////////////////////////////////////////////////
623                        //////////////////////////////////////////////////////////
624                        // TODO: Rules inserted ordered by their label monomial!//
625                        //////////////////////////////////////////////////////////
626                        //////////////////////////////////////////////////////////
627
628                        //Rule* rNew      =   new Rule(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
629                        //RNode* rNodeNew =   new RNode(rNew);
630                        //rules->insertOrdered(rNew);
631                        rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
632                        numberOfRules++;
633                        //Print("RuleOld ADDED: \n");
634                        //pWrite(rules->getFirst()->getRuleTerm());
635                        //rules->print();
636                        //Print("%d\n",rules->getFirst()->getRuleIndex());
637                        //Print("%p\n",sPolyList->getFirst());
638                        sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rules->getFirst()->getRuleOld());
639                        //sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rNew);
640                    }
641                }
642            }
643        //}
644        //Print("%p\n",temp);
645        temp    =   temp->getNext();
646        //Print("%p\n",temp);
647        //Print("%p\n",temp->getData());
648        //pWrite(temp->getLp1Poly());
649    }
650    // these critical pairs can be deleted now as they are either useless for further computations or
651    // already saved as an S-polynomial to be reduced in the following
652    delete first;   
653}
654
655
656
657/*
658========================================================================
659reduction including subalgorithm topReduction() using Faugere's criteria
660========================================================================
661*/
662void reduction(LList* sPolyList, CListOld* critPairs, LList* gPrev, RList* rules, LTagList* lTag, RTagList* rTag, 
663                 ideal gbPrev) { 
664    //Print("##########################################In REDUCTION!########################################\n");
665    // check if sPolyList has any elements
666    // NOTE: due to initialization sPolyList always has a default NULL element
667    LNode* temp = sPolyList->getFirst();
668    while(NULL != temp) {
669        // temp is the first element in the sPolyList which should be reduced
670        // due to earlier sorting this is the element of minimal degree AND
671        // minimal label
672        // delete the above first element from sPolyList, temp will be either reduced to
673        // zero or added to gPrev, but never come back to sPolyList
674        //pWrite(sPolyList->getFirst()->getPoly());
675        //Print("LIST OF SPOLYNOMIALS!\n");
676        //sPolyList->print();
677        sPolyList->setFirst(temp->getNext());
678        poly tempNF = kNF(gbPrev,currQuotient,temp->getPoly());
679        if(NULL != tempNF) {
680            pNorm(tempNF);
681            temp->setPoly(tempNF);
682            // try further reductions of temp with polynomials in gPrev
683            // with label index = current label index: this is done such that there
684            // is no label corruption during the reduction process
685            //Print("lower label reduction:  ");
686            //pWrite(tempNF);
687            topReduction(temp,sPolyList,gPrev,critPairs,rules,lTag,rTag,gbPrev);
688       
689        }
690        else {
691            reductionsToZero++;
692            delete temp;
693        }
694        //if(NULL != temp->getPoly()) {
695        //    criticalPair(gPrev,critPairs,lTag,rTag,rules);
696        //}
697        temp =   sPolyList->getFirst();
698    }
699    //sPolyList->print();
700    //delete sPolyList;
701}   
702
703/*
704========================================================================
705reduction including subalgorithm topReduction() using Faugere's criteria
706========================================================================
707*/
708void newReduction(LList* sPolyList, CListOld* critPairs, LList* gPrev, LList* reducers, RList* rules, LTagList* lTag, RTagList* rTag, 
709                 ideal gbPrev, int termination) { 
710    //Print("##########################################In REDUCTION!########################################\n");
711    // check if sPolyList has any elements
712    // NOTE: due to initialization sPolyList always has a default NULL element
713    //Print("--1--\n");
714    LNode* temp = sPolyList->getFirst();
715    while(NULL != temp) {
716        // temp is the first element in the sPolyList which should be reduced
717        // due to earlier sorting this is the element of minimal degree AND
718        // minimal label
719        // delete the above first element from sPolyList, temp will be either reduced to
720        // zero or added to gPrev, but never come back to sPolyList
721        //Print("LIST OF SPOLYNOMIALS!\n");
722        //sPolyList->print();
723        //pWrite(sPolyList->getFirst()->getPoly());
724        sPolyList->setFirst(temp->getNext());
725        //pWrite(temp->getPoly());
726        //poly tempNF = kNF(gbPrev,currQuotient,temp->getPoly());
727        //Print("!!!\n");
728        //if(NULL != tempNF) {
729            //pNorm(tempNF);
730            //temp->setPoly(tempNF);
731            //Print("lower label reduction:  ");
732            //pWrite(tempNF);
733            // try further reductions of temp with polynomials in gPrev
734            // with label index = current label index: this is done such that there
735            // is no label corruption during the reduction process
736            findReducers(temp,sPolyList,gbPrev,gPrev,reducers,critPairs,rules,lTag,rTag, termination); 
737        //}
738        //else {
739        //    reductionsToZero++;
740        //    delete temp;
741        //}
742        //if(NULL != temp->getPoly()) {
743        //    criticalPair(gPrev,critPairs,lTag,rTag,rules);
744        //}
745        //Print("HIER AUCH ?\n");
746        //Print("--2--\n");
747        //sPolyList->print();
748        //critPairs->print();
749        temp =   sPolyList->getFirst();
750        //Print("%p\n",temp);
751    }
752    //sPolyList->print();
753    //delete sPolyList;
754    //Print("REDUCTION FERTIG\n");
755}     
756
757
758/*!
759 * ================================================================================
760 * searches for reducers of temp similar to the symbolic preprocessing of F4  and
761 * divides them into a "good" and "bad" part:
762 *
763 * the "good" ones are the reducers which do not corrupt the label of temp, with
764 * these the normal form of temp is computed
765 *
766 * the "bad" ones are the reducers which corrupt the label of temp, they are tested
767 * later on for possible new rules and S-polynomials to be added to the algorithm
768 * ================================================================================
769*/
770void findReducers(LNode* l, LList* sPolyList, ideal gbPrev, LList* gPrev, LList* reducers, CListOld* critPairs, RList* rules, LTagList* lTag, RTagList* rTag, int termination) {
771    int canonicalize    =   0;
772    //int timerRed        =   0;
773    bool addToG         =   1;
774    number sign         =   nInit(-1);
775    LList* good         =   new LList();
776    LList* bad          =   new LList();
777    LList* monomials    =   new LList(l->getLPolyOld());
778    poly u              =   pOne();
779    number nOne         =   nInit(1);
780    LNode* tempRed      =   lTag->getFirstCurrentIdx();
781    LNode* tempMon      =   monomials->getFirst();
782    poly tempPoly       =   pInit();
783    poly redPoly        =   NULL;
784    int idx             =   l->getIndex();
785    //Print("IN FIND REDUCERS:  ");
786    //pWrite(l->getPoly());
787    tempPoly    =   pCopy(l->getPoly());
788    //tempMon->setPoly(tempPoly);
789    //while(NULL != tempMon) {
790        // iteration over all monomials in tempMon
791        kBucket* bucket  =   kBucketCreate();
792        kBucketInit(bucket,tempPoly,0);
793        tempPoly    =   kBucketGetLm(bucket);
794        //Print("\n\n\nTO BE REDUCED:  ");
795        //pWrite(l->getPoly());
796        //pWrite(l->getTerm());
797        //pWrite(tempPoly);
798        while(NULL != tempPoly) {
799            // iteration of all elements in gPrev of the current index
800            tempRed =   gPrev->getFirst();
801            while(NULL != tempRed) {
802                //Print("TEMPREDPOLY:  ");
803                //pWrite(tempRed->getPoly());
804                if(pLmDivisibleByNoComp(tempRed->getPoly(),tempPoly)) {
805                    u   =   pDivideM(pHead(tempPoly),pHead(tempRed->getPoly()));
806                    //Print("U:  ");
807                    //pWrite(u);
808                    if(tempRed->getIndex() != idx) {
809                            // passing criterion1 ?
810                            if(!criterion1(gPrev,u,tempRed,lTag)) {
811                                    poly tempRedPoly    =   tempRed->getPoly();
812                                    //Print("REDUCER: ");
813                                    //pWrite(ppMult_qq(u,tempRedPoly));
814                                    pIter(tempRedPoly);
815                                    int lTempRedPoly    =   pLength(tempRedPoly);
816                                    kBucketExtractLm(bucket);
817                                    kBucket_Minus_m_Mult_p(bucket,u,tempRedPoly,&lTempRedPoly);
818                                    canonicalize++;
819                                    //Print("Reduction\n");
820                                    if(!(canonicalize % 50)) {
821                                        kBucketCanonicalize(bucket);
822                                    }
823                                    tempPoly    =   kBucketGetLm(bucket);
824                                    //Print("TEMPPOLY:  ");
825                                    //pWrite(tempPoly);
826                                    if(NULL != tempPoly) {
827                                        tempRed     =   gPrev->getFirst();
828                                        continue;
829                                    }
830                                    else {
831                                        break;
832                                    }
833                             }   
834                   
835                    }
836                    else {
837                        if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 0) {
838                         //Print("NOT ALLOWED REDUCER, EQUAL LABELS:\n");
839                                      //pWrite(u);
840                                      //pWrite(tempRed->getTerm());
841                                      //pWrite(pHead(tempRed->getPoly()));
842                                      //addToG  = 0;
843                        }
844                        if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) != 0) {
845                            // passing criterion2 ?
846                            if(!criterion2(gPrev->getFirst()->getIndex(), u,tempRed,rules,rTag)) {
847                                // passing criterion1 ?
848                                if(!criterion1(gPrev,u,tempRed,lTag)) {
849                                    if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 1) {
850                                        if(NULL == redPoly) {
851                                            bad->insert(tempRed->getLPolyOld());
852                                            addToG  = 1;
853                                            //poly tempRedPoly    =   tempRed->getPoly();
854                                            //break;
855                                        }
856                                    }
857                                    else {
858                                        poly tempRedPoly    =   tempRed->getPoly();
859                                        //Print("REDUCER: ");
860                                        //pWrite(ppMult_qq(u,tempRedPoly));
861                                        pIter(tempRedPoly);
862                                        int lTempRedPoly    =   pLength(tempRedPoly);
863                                        //Print("HEAD MONOMIAL KBUCKET: ");
864                                        //pWrite(kBucketGetLm(bucket));
865                                        kBucketExtractLm(bucket);
866                                        kBucket_Minus_m_Mult_p(bucket,u,tempRedPoly,&lTempRedPoly);
867                                        canonicalize++;
868                                        //Print("REDUCTION\n");
869                                        addToG  = 1;
870                                        if(!(canonicalize % 50)) {
871                                            kBucketCanonicalize(bucket);
872                                        }
873                                        //Print("HEAD MONOMIAL KBUCKET: ");
874                                        //pWrite(kBucketGetLm(bucket));
875                                        tempPoly    =   kBucketGetLm(bucket);
876                                        //Print("TEMPPOLY:  ");
877                                        //pWrite(tempPoly);
878                                        if(NULL != tempPoly) {
879                                            tempRed     =   gPrev->getFirst();
880                                            continue;
881                                        }
882                                        else {
883                                            break;
884                                        }
885                                    }
886                                }   
887                            }
888                            else {
889                              //Print("CRIT 1  ");
890                                   
891                                      if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 1 ) {
892                                      //Print("NOT ALLOWED REDUCER, CRIT 1:\n");
893                                      //pWrite(u);
894                                      //pWrite(tempRed->getTerm());
895                                      //pWrite(tempRed->getPoly());
896                                      if(pLmDivisibleByNoComp(tempRed->getTerm(),l->getTerm())) {
897                                        //Print("REDUCER LABEL:  ");
898                                        //pWrite(tempRed->getTerm());
899addToG  = 0;
900                                      }
901                                    }
902                            }
903                        }
904                        else {
905                          //Print("CRIT 2  ");
906                                    if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 1) {
907                                      //Print("NOT ALLOWED REDUCER, CRIT 2:\n");
908                                      //pWrite(u);
909                                      //pWrite(tempRed->getTerm());
910                                      //pWrite(tempRed->getPoly());
911                                      if(pLmDivisibleByNoComp(tempRed->getTerm(),l->getTerm())) {
912                                        //Print("REDUCER LABEL:  ");
913                                        //pWrite(tempRed->getTerm());
914addToG  = 0;
915                                      }
916                                    }
917                        }
918                    }
919                   
920                }
921                tempRed =   tempRed->getNext();
922            }
923            // reduction process with elements in LList* reducers
924          if(NULL != tempPoly) {
925            //Print("HERE\n");
926            tempRed =   reducers->getFirst();
927            while(NULL != tempRed) {
928                //Print("TEMPREDPOLY:  ");
929                //pWrite(tempRed->getPoly());
930              //pWrite(tempPoly); 
931              if(pLmDivisibleByNoComp(tempRed->getPoly(),tempPoly)) {
932                    //Print("A\n");
933                    u   =   pDivideM(pHead(tempPoly),pHead(tempRed->getPoly()));
934                    //Print("U:  ");
935                    //pWrite(u);
936                    if(tempRed->getIndex() != idx) {
937                            // passing criterion1 ?
938                            if(!criterion1(gPrev,u,tempRed,lTag)) {
939                                    poly tempRedPoly    =   tempRed->getPoly();
940                                    //Print("REDUCER: ");
941                                    //pWrite(ppMult_qq(u,tempRedPoly));
942                                    pIter(tempRedPoly);
943                                    int lTempRedPoly    =   pLength(tempRedPoly);
944                                    kBucketExtractLm(bucket);
945                                    kBucket_Minus_m_Mult_p(bucket,u,tempRedPoly,&lTempRedPoly);
946                                    canonicalize++;
947                                    //Print("Reduction\n");
948                                    if(!(canonicalize % 50)) {
949                                        kBucketCanonicalize(bucket);
950                                    }
951                                    tempPoly    =   kBucketGetLm(bucket);
952                                    //Print("TEMPPOLY:  ");
953                                    //pWrite(tempPoly);
954                                    if(NULL != tempPoly) {
955                                        tempRed     =   gPrev->getFirst();
956                                        continue;
957                                    }
958                                    else {
959                                        break;
960                                    }
961                             }   
962                   
963                    }
964                    else {
965                        if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 0) {
966                         //Print("NOT ALLOWED REDUCER, EQUAL LABELS:\n");
967                                      //pWrite(u);
968                                      //pWrite(tempRed->getTerm());
969                                      //pWrite(pHead(tempRed->getPoly()));
970                                      //addToG  = 0;
971                        }
972                        if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) != 0) {
973                            // passing criterion2 ?
974                            if(!criterion2(gPrev->getFirst()->getIndex(), u,tempRed,rules,rTag)) {
975                                // passing criterion1 ?
976                                if(!criterion1(gPrev,u,tempRed,lTag)) {
977                                    if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 1) {
978                                        if(NULL == redPoly) {
979                                            bad->insert(tempRed->getLPolyOld());
980                                            addToG  = 1;
981                                            //poly tempRedPoly    =   tempRed->getPoly();
982                                            //break;
983                                        }
984                                    }
985                                    else {
986                                        poly tempRedPoly    =   tempRed->getPoly();
987                                        //Print("REDUCER: ");
988                                        //pWrite(ppMult_qq(u,tempRedPoly));
989                                        pIter(tempRedPoly);
990                                        int lTempRedPoly    =   pLength(tempRedPoly);
991                                        //Print("HEAD MONOMIAL KBUCKET: ");
992                                        //pWrite(kBucketGetLm(bucket));
993                                        kBucketExtractLm(bucket);
994                                        kBucket_Minus_m_Mult_p(bucket,u,tempRedPoly,&lTempRedPoly);
995                                        canonicalize++;
996                                        //Print("REDUCTION\n");
997                                        addToG  = 1;
998                                        if(!(canonicalize % 50)) {
999                                            kBucketCanonicalize(bucket);
1000                                        }
1001                                        //Print("HEAD MONOMIAL KBUCKET: ");
1002                                        //pWrite(kBucketGetLm(bucket));
1003                                        tempPoly    =   kBucketGetLm(bucket);
1004                                        //Print("TEMPPOLY:  ");
1005                                        //pWrite(tempPoly);
1006                                        if(NULL != tempPoly) {
1007                                            tempRed     =   gPrev->getFirst();
1008                                            continue;
1009                                        }
1010                                        else {
1011                                            break;
1012                                        }
1013                                    }
1014                                }   
1015                            }
1016                            else {
1017                              //Print("CRIT 1  ");
1018                                   
1019                                      if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 1 ) {
1020                                      //Print("NOT ALLOWED REDUCER, CRIT 1:\n");
1021                                      //pWrite(u);
1022                                      //pWrite(tempRed->getTerm());
1023                                      //pWrite(tempRed->getPoly());
1024                                      if(pLmDivisibleByNoComp(tempRed->getTerm(),l->getTerm())) {
1025                                        //Print("REDUCER LABEL:  ");
1026                                        //pWrite(tempRed->getTerm());
1027                                        addToG  = 0;
1028                                      }
1029                                    }
1030                            }
1031                        }
1032                        else {
1033                          //Print("CRIT 2  ");
1034                                    if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 1) {
1035                                      //Print("NOT ALLOWED REDUCER, CRIT 2:\n");
1036                                      //pWrite(u);
1037                                      //pWrite(tempRed->getTerm());
1038                                      //pWrite(tempRed->getPoly());
1039                                      if(pLmDivisibleByNoComp(tempRed->getTerm(),l->getTerm())) {
1040                                        //Print("REDUCER LABEL:  ");
1041                                        //pWrite(tempRed->getTerm());
1042addToG  = 0;
1043                                      }
1044                                    }
1045                        }
1046                    }
1047                   
1048                }
1049                tempRed =   tempRed->getNext();
1050            }
1051          }
1052            if(NULL != tempPoly) {
1053                if(NULL == redPoly) {
1054                    redPoly =   kBucketExtractLm(bucket);
1055                }
1056                else {
1057                    redPoly     =   p_Merge_q(redPoly,kBucketExtractLm(bucket),currRing);
1058                }
1059               // for top-reduction only
1060               redPoly  = p_Merge_q(redPoly,kBucketClear(bucket),currRing);
1061               break;
1062               // end for top-reduction only
1063               tempPoly    =   kBucketGetLm(bucket);
1064               
1065            }
1066        }
1067        if(NULL == redPoly) {
1068            reductionsToZero++;
1069            //pDelete(&redPoly);
1070        }
1071        else {
1072            //Print("\nELEMENT ADDED TO GPREV: ");
1073            pNorm(redPoly);
1074            //pWrite(pHead(redPoly));
1075            //pWrite(l->getTerm());
1076            //Print("%d\n",canonicalize);
1077            l->setPoly(redPoly);
1078            if(addToG == 0 && termination == 1) {
1079              reducers->insert(l->getLPolyOld());
1080            }
1081            else {
1082              gPrev->insert(l->getLPolyOld());
1083            }
1084            //Print("%d\n\n",termination);
1085            if(termination == 1) {
1086            if(addToG) {
1087              //Print("----------------HERE?-----------------\n");
1088              criticalPair(gPrev,critPairs,lTag,rTag,rules);
1089            }
1090            else {
1091              notInG++;
1092              //Print("\nNONONO");
1093              //pWrite(pHead(l->getPoly()));
1094              //pWrite(l->getTerm());
1095            }
1096            }
1097            else {
1098              criticalPair(gPrev,critPairs,lTag,rTag,rules);
1099            }
1100        }
1101   
1102    // if there are "bad" reducers than try to compute new S-polynomials and rules
1103   
1104    if(NULL != bad->getFirst()) {
1105        //Print("BAD STUFF LIST:\n");
1106        //bad->print();
1107        LNode* tempBad  =   bad->getFirst();
1108        //pWrite(l->getPoly());
1109        while(NULL != tempBad) {
1110            if(pDivisibleBy(tempBad->getPoly(),l->getPoly())) {
1111                //Print("BAD STUFF\n");
1112                //pWrite(l->getPoly());
1113                //pWrite(tempBad->getPoly());
1114                u   =   pDivide(pHead(l->getPoly()),pHead(tempBad->getPoly()));
1115                //Print("MULTIPLIER:  ");
1116                //pWrite(u);
1117                pSetCoeff(u,nOne);
1118                if(pLmCmp(ppMult_qq(u,tempBad->getTerm()),l->getTerm()) != 0) {
1119                    // passing criterion2 ?
1120                    if(!criterion2(gPrev->getFirst()->getIndex(), u,tempBad,rules,rTag)) {
1121                        // passing criterion1 ?
1122                        if(!criterion1(gPrev,u,tempBad,lTag)) {
1123                            //Print("HIERHIERHIERHIERHIERHIER\n");
1124                            rules->insert(tempBad->getIndex(),ppMult_qq(u,tempBad->getTerm()));
1125                            numberOfRules++;
1126                            //gPrev->print();
1127                            //pWrite(l->getPoly());
1128                            poly temp   =   ksOldSpolyRedNew(ppMult_qq(u,tempBad->getPoly()),l->getPoly());
1129                            //pWrite(l->getPoly());
1130                            //Print("%p\n",temp);
1131                            //gPrev->print();
1132                            if(NULL != temp) {
1133                                pNorm(temp);
1134                                LNode* tempBadNew   =   new LNode(ppMult_qq(u,tempBad->getTerm()),tempBad->getIndex(),temp,rules->getFirst()->getRuleOld());
1135                                //pWrite(temp);
1136                                //pWrite(tempBadNew->getPoly());
1137                                //pWrite(tempBadNew->getTerm());
1138                                //pWrite(pHead(tempBadNew->getPoly()));
1139                                //Print("%p\n",tempBadNew->getPoly());
1140                                //tempRed->getLPoly()->setRule(rules->getFirst()->getRule());
1141                                tempBadNew->setDel(1);
1142                           
1143                                sPolyList->insertByLabel(tempBadNew);
1144                                //Print("BAD SPOLYLIST: \n");
1145                                //sPolyList->print();
1146                            }
1147                        }
1148                    }
1149                }
1150            }
1151        //Print("HIER\n");
1152            tempBad =   tempBad->getNext();
1153            //Print("%p\n",tempBad);
1154        }
1155       // Print("-------------------BAD STUFF LIST-----------------------------\n");
1156    }
1157        //Print("HIER AUCH\n");
1158        //Print("SPOLYLIST IN BAD: \n");
1159        //sPolyList->print();
1160    //Print("END FIND REDUCERS\n");
1161}
1162
1163/*
1164=======================================================================================
1165merging 2 polynomials p & q without requiring that all monomials of p & q are different
1166if there are equal monomials in p & q only one of these monomials (always that of p!)
1167is taken into account
1168=======================================================================================
1169
1170poly p_MergeEq_q(poly p, poly q, const ring r) {
1171  assume(p != NULL && q != NULL);
1172  p_Test(p, r);
1173  p_Test(q, r);
1174#if PDEBUG > 0
1175  int l = pLength(p) + pLength(q);
1176#endif
1177
1178  spolyrec rp;
1179  poly a = &rp;
1180  DECLARE_LENGTH(const unsigned long length = r->CmpL_Size);
1181  DECLARE_ORDSGN(const long* ordsgn = r->ordsgn);
1182
1183  Top:     // compare p and q w.r.t. monomial ordering
1184  p_MemCmp(p->exp, q->exp, length, ordsgn, goto Equal, goto Greater , goto Smaller);
1185
1186  Equal:
1187  a =   pNext(a)    =   p;
1188  pIter(p);
1189  pIter(q);
1190  if(NULL == p) {
1191      if(NULL == q) {
1192          goto Finish;
1193      }
1194      else {
1195          pNext(a)  =   q;
1196          goto Finish;
1197      }
1198  }
1199  goto Top;
1200
1201  Greater:
1202  a = pNext(a) = p;
1203  pIter(p);
1204  if (p==NULL) { pNext(a) = q; goto Finish;}
1205  goto Top;
1206
1207  Smaller:
1208  a = pNext(a) = q;
1209  pIter(q);
1210  if (q==NULL) { pNext(a) = p; goto Finish;}
1211  goto Top;
1212
1213  Finish:
1214
1215  p_Test(pNext(&rp), r);
1216#if PDEBUG > 0
1217  pAssume1(l - pLength(pNext(&rp)) == 0);
1218#endif
1219  return pNext(&rp);
1220}
1221*/
1222
1223/*
1224=====================================================================================
1225top reduction in F5, i.e. reduction of a given S-polynomial by labeled polynomials of
1226the same index whereas the labels are taken into account
1227=====================================================================================
1228*/
1229void topReduction(LNode* l, LList* sPolyList, LList* gPrev, CListOld* critPairs,  RList* rules, LTagList* lTag, RTagList* rTag, ideal gbPrev) {
1230    //Print("##########################################In topREDUCTION!########################################\n");
1231    // try l as long as there are reductors found by findReductor()
1232    LNode* gPrevRedCheck    =   lTag->getFirstCurrentIdx();
1233    LNode* tempRed;
1234    poly pOne               =   pOne();
1235    do {
1236        //int timer5  =   initTimer();
1237        //startTimer();
1238        //Print("TOP REDUCTION:  ");
1239        //pWrite(l->getPoly());
1240        tempRed  =   findReductor(l,sPolyList,gPrevRedCheck,gPrev,rules,lTag,rTag);
1241        //timer5  =   getTimer();
1242        //reductionTime   =   reductionTime   +   timer5;
1243        // if a reductor for l is found and saved in tempRed
1244        if(NULL != tempRed) {
1245            // if label of reductor is greater than the label of l we have to built a new element
1246            // and add it to sPolyList
1247           
1248            if(pLmCmp(tempRed->getTerm(),l->getTerm()) == 1) {
1249                // needed sinc pSub destroys the arguments!
1250                //poly temp_poly_l    =   pInit();
1251                //temp_poly_l         =   pCopy(l->getPoly());
1252                //Print("VORHER: ");
1253                //pWrite(tempRed->getPoly());
1254                //poly temp           =   pMinus_mm_Mult_qq(tempRed->getPoly(),pOne,l->getPoly());
1255                poly temp   =   ksOldSpolyRedNew(l->getPoly(),tempRed->getPoly());
1256                rules->insert(tempRed->getIndex(),tempRed->getTerm());
1257                //Print("NACHHER: ");
1258                //pWrite(tempRed->getPoly());
1259                //Print("TEMP: ");
1260                //pWrite(temp);
1261                if(NULL != temp) {
1262                    pNorm(temp);
1263                    //tempRed->setPoly(temp);
1264                    //tempRed->setDel(1);
1265                    //rules->insert(tempRed->getIndex(),tempRed->getTerm());
1266                    LNode* tempRedNew   =   new LNode(tempRed->getTerm(),tempRed->getIndex(),temp,rules->getFirst()->getRuleOld());
1267                    //tempRed->getLPoly()->setRule(rules->getFirst()->getRule());
1268                    tempRedNew->setDel(1);
1269                    sPolyList->insertByLabel(tempRedNew);
1270                }
1271                else {
1272                    pDelete(&temp);
1273                    reductionsToZero++;
1274                    //delete tempRed;
1275                }
1276            }
1277           
1278            // label of reductor is smaller than the label of l, subtract reductor from l and delete the
1279            // gPrevRedCheck pointer added to l during findReductor() as the head term of l changes
1280            // after subtraction
1281            else {
1282               
1283                //poly temp_poly_l    =   pInit();
1284                //temp_poly_l         =   pCopy(l->getPoly());
1285                //poly temp   =   pMinus_mm_Mult_qq(tempRed->getPoly(),pOne,l->getPoly());
1286                //Print("REDUCER: ");
1287                //pWrite(tempRed->getPoly());
1288                //pWrite(tempRed->getTerm());
1289                poly temp   =   ksOldSpolyRedNew(l->getPoly(),tempRed->getPoly());
1290                //Print("REDUCED ELEMENT:  ");
1291                if(NULL != temp) {
1292                    pNorm(temp);
1293                    //pWrite(temp);
1294                    poly tempNF =   kNF(gbPrev,currQuotient,temp); 
1295                    pNorm(tempNF);
1296                    if(NULL == tempNF) {
1297                        reductionsToZero++;
1298                        pDelete(&tempNF);
1299                        l->setPoly(NULL);
1300                        break;
1301                    }
1302                    l->setPoly(tempNF);
1303                   
1304                    gPrevRedCheck   =   lTag->getFirstCurrentIdx();
1305                }
1306                else {
1307                    reductionsToZero++;
1308                    pDelete(&temp);
1309                    l->setPoly(NULL);
1310                    break;
1311                }
1312            }   
1313        }
1314        else {
1315            if(NULL != l->getPoly()) {
1316                pNorm(l->getPoly());
1317                //Print("ELEMENT ADDED TO GPREV: ");
1318                //pWrite(l->getPoly());
1319                gPrev->insert(l->getLPolyOld());
1320                //Print("TEMP RED == 0  ");
1321                //pWrite(l->getPoly());
1322                //pWrite(l->getTerm());
1323                //rules->print();
1324                criticalPair(gPrev,critPairs,lTag,rTag,rules);
1325                //Print("LIST OF CRITICAL PAIRS:    \n");
1326                //critPairs->print();
1327            }
1328            break;
1329        }
1330    } while(1);
1331}
1332
1333
1334/*
1335=====================================================================
1336subalgorithm to find a possible reductor for the labeled polynomial l
1337=====================================================================
1338*/
1339LNode* findReductor(LNode* l, LList* sPolyList, LNode* gPrevRedCheck, LList* gPrev, RList* rules, LTagList* lTag,RTagList* rTag) {
1340    // allociation of memory for the possible reductor
1341    //Print("LPOLY:  ");
1342    //pWrite(l->getPoly());
1343    poly u      =   pOne();
1344    poly red;
1345    number nOne =   nInit(1);
1346    LNode* temp;
1347    // head term of actual element such that we do not have to call pHead at each new reductor test
1348    poly t      =   pHead(l->getPoly());
1349    // if l was already checked use the information in gPrevRedCheck such
1350    // that we can start searching for new reducers from this point and
1351    // not from the first element of gPrev with the current index
1352    temp    =   gPrevRedCheck;
1353    // search for reductors until we are at the end of gPrev resp. at the
1354    // end of the elements of the current index
1355    while(NULL != temp && temp->getIndex() == l->getIndex()) {
1356        // does the head of the element of gPrev divides the head of
1357        // the to be reduced element?
1358        if(pLmDivisibleByNoComp(pHead(temp->getPoly()),t)) {
1359            // get all the information needed for the following tests
1360            // of the criteria
1361            u   =   pDivide(t,pHead(temp->getPoly()));
1362            pSetCoeff(u,nOne);
1363            red =   ppMult_qq(u,temp->getPoly());
1364            pNorm(red);
1365            // check if both have the same label
1366            if(pLmCmp(ppMult_qq(u,temp->getTerm()),l->getTerm()) != 0) {
1367                // passing criterion2 ?
1368                if(!criterion2(gPrev->getFirst()->getIndex(), u,temp,rules,rTag)) {
1369                    // passing criterion1 ?
1370                    if(!criterion1(gPrev,u,temp,lTag)) {
1371                            gPrevRedCheck   =   temp->getNext();
1372                            LNode* redNode  =   new LNode(ppMult_qq(u,temp->getTerm()),temp->getIndex(),red,NULL,NULL);
1373                            return redNode;
1374                    }
1375                }
1376            }
1377            if(pLmCmp(ppMult_qq(u,temp->getTerm()),l->getTerm()) == 1) {
1378                // passing criterion2 ?
1379                if(!criterion2(gPrev->getFirst()->getIndex(), u,temp,rules,rTag)) {
1380                    // passing criterion1 ?
1381                    if(!criterion1(gPrev,u,temp,lTag)) {
1382                        poly tempSpoly  =   ksOldSpolyRedNew(red,l->getPoly());
1383                        rules->insert(temp->getIndex(),ppMult_qq(u,temp->getTerm()));
1384                        gPrevRedCheck   =   temp->getNext();
1385                        if(NULL != tempSpoly) {
1386                            pNorm(tempSpoly);
1387                            sPolyList->insertByLabel(ppMult_qq(u,temp->getTerm()),temp->getIndex(),tempSpoly,rules->getFirst()->getRuleOld());
1388                    //Print("NEW ONE: ");
1389                    //pWrite(tempSpoly);
1390                    //Print("HIER\n");
1391                            //sPolyList->print();
1392                            //sleep(5);
1393                        }
1394                    }
1395                }
1396            }
1397        }
1398        //Print("AUCH HIER\n");
1399        temp = temp->getNext();
1400    }
1401   
1402//    delete temp;
1403 return NULL;
1404}
1405
1406
1407
1408/*
1409==========================================================================
1410MAIN:computes a gb of the ideal i in the ring r with our F5 implementation
1411OPTIONS: INTEGER "opt" is to be set "0" for F5, "1" for F5R, "2" for F5C
1412==========================================================================
1413*/
1414ideal F5main(ideal id, ring r, int opt, int termination) {
1415  switch(opt) { 
1416    case 0:
1417      Print("\nComputations are done by the standard F5 Algorithm");
1418      break;
1419    case 1:
1420      Print("\nComputations are done by the variant F5R of the F5 Algorithm");
1421      break;
1422    case 2:
1423      Print("\nComputations are done by the variant F5C of the F5 Algorithm");
1424      break;
1425    default:
1426      WerrorS("\nThe option can only be set to \"0\", \"1\" or \"2\":\n\"0\": standard F5 Algorithm\n\"1\": variant F5R\n\"2\": variant F5C\nComputations are aborted!\n");
1427      return id;
1428  }
1429    int timer   =   initTimer();
1430    startTimer();
1431    int i,j,k,l;
1432    int gbLength;
1433    // 1 polynomial for defining initial labels & further tests
1434    poly ONE = pOne();
1435    poly pOne   =   pOne();
1436    number nOne =   nInit(1);
1437    // tag the first element of index i-1 for criterion 1
1438    //Print("LTAG BEGINNING: %p\n",lTag);
1439   
1440    // DEBUGGING STUFF START
1441    //Print("NUMBER: %d\n",r->N);
1442    /*
1443    int* ev = new int[r->N +1];
1444    for(i=0;i<IDELEMS(id);i++) {
1445        pGetExpV(id->m[i],ev);
1446        //ev2  =   pGetExp(id->m[i],1);
1447        pWrite(id->m[i]);
1448        Print("EXP1: %d\n",ev[1]);
1449        Print("EXP2: %d\n",ev[2]);
1450        Print("EXP3: %d\n\n",ev[3]);
1451        Print("SUM: %ld\n\n\n",sumVector(ev,r->N));
1452    }
1453    delete ev;
1454    */
1455    /*DEBUGGING STUFF END */
1456   
1457    // first element in rTag is first element of rules which is NULL RNode,
1458    // this must be done due to possible later improvements
1459    RList* rules    =   new RList();
1460    //Print("RULES FIRST: %p\n",rules->getFirst());
1461    //Print("RULES FIRST DATA: %p\n",rules->getFirst()->getRule());
1462    //RTagList* rTag  =   new RTagList(rules->getFirst());
1463    RTagList* rTag  =   NULL;
1464    i = 1;
1465    /*for(j=0; j<IDELEMS(id); j++) {
1466        if(NULL != id->m[j]) {
1467            if(pComparePolys(id->m[j],ONE)) {
1468                Print("One Polynomial in Input => Computations stopped\n");
1469                ideal idNew = idInit(1,1);
1470                idNew->m[0] = ONE;
1471                return(idNew);
1472            }   
1473        }
1474    }*/ 
1475    ideal idNew     =   kInterRed(id); 
1476    id              =   idNew;
1477    //qsortDegree(&id->m[0],&id->m[IDELEMS(id)-1]);
1478    //idShow(id);
1479    LList* gPrev    =   new LList(ONE, i, id->m[0]);
1480    LList* reducers =   new LList();
1481    //idShow(id);
1482    //Print("%p\n",id->m[0]);
1483    //pWrite(id->m[0]);
1484    //Print("%p\n",gPrev->getFirst()->getPoly());
1485    //pWrite(gPrev->getFirst()->getPoly());
1486
1487    LTagList* lTag  =   new LTagList(gPrev->getFirst());
1488    //lTag->insert(gPrev->getFirst());
1489    lTag->setFirstCurrentIdx(gPrev->getFirst());
1490    // computing the groebner basis of the elements of index < actual index
1491    gbLength    =   gPrev->getLength();
1492    //Print("Laenge der bisherigen Groebner Basis: %d\n",gPrev->getLength());
1493    ideal gbPrev    =   idInit(gbLength,1);
1494    // initializing the groebner basis of elements of index < actual index
1495    gbPrev->m[0]    =   gPrev->getFirst()->getPoly();
1496    //idShow(gbPrev);
1497    //idShow(currQuotient);
1498    for(i=2; i<=IDELEMS(id); i++) {
1499        LNode* gPrevTag =   gPrev->getLast();
1500        //Print("Last POlynomial in GPREV: ");
1501        //Print("%p\n",gPrevTag);   
1502        //pWrite(gPrevTag->getPoly());
1503        gPrev   =   F5inc(i, id->m[i-1], gPrev, reducers, gbPrev, ONE, lTag, rules, rTag, termination);
1504        //Print("%d\n",gPrev->count(gPrevTag->getNext()));
1505        //Print("%d\n",gPrev->getLength());
1506        //Print("____________________________________ITERATION STEP DONE________________________________________\n");
1507       
1508        // DEBUGGING STUFF
1509        LNode* temp    =   gPrev->getFirst();
1510   
1511
1512        /////////////////////////////////////////////////////////////////////////////////
1513        //                                                                             //
1514        // one needs to choose one of the following 3 implementations of the algorithm //
1515        // F5,F5R or F5C                                                               //
1516        //                                                                             //
1517        /////////////////////////////////////////////////////////////////////////////////                                                                           
1518       
1519       
1520        //   
1521        // standard "F5"
1522        //
1523        if(0 == opt) { 
1524          if(gPrev->getLength() > gbLength) {
1525            if(i < IDELEMS(id)) {
1526                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1527                LNode*  temp =   gPrevTag;
1528                int counter =   0;
1529                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1530                    temp        =   temp->getNext();
1531                    if(0 == temp->getDel()) {
1532                        counter++;
1533                        gbAdd->m[j] =   temp->getPoly();
1534                    }
1535                }
1536                    gbPrev          =   idAdd(gbPrev,gbAdd);
1537            }
1538            else {
1539                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1540                LNode*  temp =   gPrevTag;
1541                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1542                    temp        =   temp->getNext();
1543                    gbAdd->m[j] =   temp->getPoly();
1544                }
1545                gbPrev          =   idAdd(gbPrev,gbAdd);
1546            }
1547            if(i == IDELEMS(id)) {
1548                ideal tempId        =   kInterRed(gbPrev);
1549                gbPrev              =   tempId;
1550            }
1551        }
1552        gbLength    =   gPrev->getLength();
1553        }
1554       
1555
1556        //
1557        //  "F5R"
1558        //
1559        if(1 == opt) {
1560        if(gPrev->getLength() > gbLength) {
1561            if(i < IDELEMS(id)) {
1562                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1563                LNode*  temp =   gPrevTag;
1564                int counter =   0;
1565                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1566                    temp        =   temp->getNext();
1567                    if(0 == temp->getDel()) {
1568                        counter++;
1569                        gbAdd->m[j] =   temp->getPoly();
1570                    }
1571                }
1572                    gbPrev          =   idAdd(gbPrev,gbAdd);
1573            }
1574            else {
1575                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1576                LNode*  temp =   gPrevTag;
1577                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1578                    temp        =   temp->getNext();
1579                    gbAdd->m[j] =   temp->getPoly();
1580                }
1581                gbPrev          =   idAdd(gbPrev,gbAdd);
1582            }
1583            // interreduction stuff
1584            // comment this out if you want F5 instead of F5R
1585            //if(i<IDELEMS(id)) {
1586                ideal tempId    =   kInterRed(gbPrev);
1587                gbPrev          =   tempId;
1588            //}
1589        }
1590        gbLength    =   gPrev->getLength();
1591        }
1592       
1593
1594        //
1595        // "F5C"
1596        //
1597        if(2 == opt) { 
1598        if(gPrev->getLength() > gbLength) {
1599            if(i < IDELEMS(id)) {
1600                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1601                LNode*  temp =   gPrevTag;
1602                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1603                    temp        =   temp->getNext();
1604                        gbAdd->m[j] =   temp->getPoly();
1605                }
1606                    gbPrev          =   idAdd(gbPrev,gbAdd);
1607            }
1608            else {
1609                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1610                LNode*  temp =   gPrevTag;
1611                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1612                    temp        =   temp->getNext();
1613                    gbAdd->m[j] =   temp->getPoly();
1614                }
1615                gbPrev          =   idAdd(gbPrev,gbAdd);
1616            }
1617            if(i<IDELEMS(id)) {
1618                ideal tempId    =   kInterRed(gbPrev);
1619                gbPrev          =   tempId;
1620                delete gPrev;
1621                delete rules;
1622                gPrev    =   new LList(pOne,1,gbPrev->m[0]);
1623                gPrev->insert(pOne,1,gbPrev->m[1]);
1624                rules    =   new RList();
1625                //rTag     =   new RTagList(rules->getFirst());
1626                for(k=2; k<IDELEMS(gbPrev); k++) {
1627                    gPrev->insert(pOne,k+1,gbPrev->m[k]);
1628                    /*
1629                    for(l=0; l<k; l++) {
1630                    }
1631                    rTag->insert(rules->getFirst());
1632                    */
1633                }
1634            }
1635            gbLength    =   gPrev->getLength(); 
1636        } 
1637        }
1638
1639
1640    }
1641    timer   =   getTimer();
1642    //Print("\n\nADDING TIME IN REDUCTION: %d\n\n",reductionTime);
1643    Print("\n\nNumber of zero-reductions:           %d\n",reductionsToZero);
1644    Print("Number of rules:                     %d\n",numberOfRules); 
1645    Print("Elements not added to G:             %d\n",notInG); 
1646    Print("Highest Degree during computations:  %d\n",highestDegree);
1647    Print("Time for computations:               %d/1000 seconds\n",timer);
1648    Print("Number of elements in gb:            %d\n",IDELEMS(gbPrev));
1649    //LNode* temp    =   gPrev->getFirst();
1650    //while(NULL != temp) {
1651    //    pWrite(temp->getPoly());
1652    //    temp    =   temp->getNext();
1653    // }
1654    //gPrev->print();
1655    //delete lTag;
1656    delete rTag;
1657    delete gPrev;
1658    notInG              =   0;
1659    numberOfRules       =   0;
1660    reductionsToZero    =   0;
1661    highestDegree       =   0;
1662    reductionTime       =   0;
1663    spolsTime           =   0;
1664    return(gbPrev);
1665
1666
1667}
1668
1669#endif
Note: See TracBrowser for help on using the repository browser.