source: git/kernel/f5gb.cc @ ae5177

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