source: git/kernel/f5gb.cc @ 667a9c

spielwiese
Last change on this file since 667a9c was 667a9c, checked in by Christian Eder, 15 years ago
bad stuff updated git-svn-id: file:///usr/local/Singular/svn/trunk@11856 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 65.2 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
388    if(NULL == rTag->getFirst()) {
389        if(NULL != rules->getFirst()) {
390            testNode    =   rules->getFirst();
391        }
392        else {
393            return false;
394        }
395    }
396    else {
397
398        if(l->getIndex() > rTag->getFirst()->getRuleIndex()) {
399            testNode    =   rules->getFirst();
400        }
401        else {
402       //Print("HIER\n");
403            //Print("DEBUG\n");
404        //Print("L INDEX: %d\n",l->getIndex());
405            /*-------------------------------------
406             * TODO: WHEN INTERREDUCING THE GB THE
407             *       INDEX OF THE PREVIOUS ELEMENTS
408             *       GETS HIGHER!
409             *-----------------------------------*/
410            //testNode    =   rules->getFirst();
411            testNode    =   rTag->get(l->getIndex());
412            if(NULL == testNode) {
413                testNode    =   rules->getFirst();
414            }
415            //Print("TESTNODE ADDRESS: %p\n",testNode);
416        }
417    }
418    //testNode    =   rules->getFirst();
419        // save the monom t1*label_term(l) as it is tested various times in the following
420    poly u1 = ppMult_qq(t,l->getTerm());
421    // first element added to rTag was NULL, check for this
422    //Print("%p\n",testNode->getRule());
423    // NOTE: testNode is possibly NULL as rTag->get() returns NULL for elements of index <=1!
424    //Print("TESTNODE: %p\n",testNode);
425    //pWrite(testNode->getRuleTerm());
426    if(NULL != testNode ) {   
427        //pWrite(testNode->getRuleTerm());
428    }
429    if(NULL != testNode) {
430        if(testNode->getRule() == l->getRule()) {
431            //Print("%p\n%p\n",testNode->getRule(),l->getRule());
432            //Print("EQUAL\n");
433        }
434        else {
435            //Print("NOT EQUAL\n");
436        }
437    }
438    while(NULL != testNode && testNode->getRule() != l->getRule() 
439          && l->getIndex() == testNode->getRuleIndex()) {
440        //Print("%p\n",testNode);
441        //pWrite(testNode->getRuleTerm());
442        //pWrite(t);
443        //pWrite(l->getTerm());
444        //pWrite(u1);
445        //Print("%d\n",testNode->getRuleIndex());
446        if(pLmDivisibleByNoComp(testNode->getRuleTerm(),u1)) {
447            //Print("-----------------Criterion 2 NOT passed!-----------------------------------\n");
448            //Print("INDEX: %d\n",l->getIndex());
449            pDelete(&u1);
450    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n\n");
451            return true;
452        }
453                testNode    =   testNode->getNext();
454    }
455    //delete testNode;
456    pDelete(&u1);
457    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n\n");
458    return false;
459}
460
461
462
463/*
464=================================================================================================================
465Criterion 2, i.e. Rewritten Criterion, for its second call in computeSPols(), with added lastRuleTested parameter
466=================================================================================================================
467*/
468inline bool criterion2(poly t, LPoly* l, RList* rules, Rule* testedRule) {
469    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n");
470    //Print("LAST RULE TESTED: %p",testedRule);
471    /*
472    Print("RULES: \n");
473        RNode* tempR    =   rules->getFirst();
474        while(NULL != tempR) {
475            pWrite(tempR->getRuleTerm());
476            Print("%d\n\n",tempR->getRuleIndex());
477            tempR   =   tempR->getNext();
478        }
479        //Print("TESTED ELEMENT: ");
480        //pWrite(l->getPoly());
481        //pWrite(l->getTerm());
482        //pWrite(ppMult_qq(t,l->getTerm()));
483        //Print("%d\n\n",l->getIndex());
484    */
485// start at the previously added element to gPrev, as all other elements will have the same index for sure
486        RNode* testNode =   rules->getFirst();
487    // save the monom t1*label_term(l) as it is tested various times in the following
488    poly u1 = ppMult_qq(t,l->getTerm());
489        // first element added to rTag was NULL, check for this
490        while(NULL != testNode && testNode->getRule() != testedRule) {
491        //pWrite(testNode->getRuleTerm());
492        if(pLmDivisibleByNoComp(testNode->getRuleTerm(),u1)) {
493            //Print("--------------------------Criterion 2 NOT passed!------------------------------\n");
494            //Print("INDEX: %d\n",l->getIndex());
495            pDelete(&u1);
496    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n\n");
497            return true;
498        }
499                testNode    =   testNode->getNext();
500    }
501    pDelete(&u1);
502    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n\n");
503    return false;
504}
505
506
507
508/*
509==================================
510Computation of S-Polynomials in F5
511==================================
512*/
513void computeSPols(CNode* first, RTagList* rTag, RList* rules, LList* sPolyList) { 
514    CNode* temp         =   first;
515    poly sp     =   pInit();
516    number sign =   nInit(-1);   
517    //Print("###############################IN SPOLS##############################\n");
518    //first->print();
519
520    while(NULL != temp) {
521        //Print("JA\n");
522        // only if a new rule was added since the last test in subalgorithm criticalPair()
523        //if(rules->getFirst() != rTag->getFirst()) {
524        if(!criterion2(temp->getT1(),temp->getAdLp1(),rules,temp->getTestedRule())) {
525                // the second component is tested only when it has the actual index, otherwise there is
526                // no new rule to test since the last test in subalgorithm criticalPair()
527                if(highestDegree < pDeg(ppMult_qq(temp->getT1(),temp->getLp1Poly()))) { 
528                    highestDegree   = pDeg(ppMult_qq(temp->getT1(),temp->getLp1Poly()));
529                    //pWrite(pHead(ppMult_qq(temp->getT1(),temp->getLp1Poly())));
530                }   
531                if(temp->getLp2Index() == temp->getLp1Index()) {
532                    if(!criterion2(temp->getT2(),temp->getAdLp2(),rules,temp->getTestedRule())) {
533                        // computation of S-polynomial
534                        //poly p1 =   temp->getLp1Poly();
535                        //poly p2 =   temp->getLp2Poly();
536                        //pIter(p1);
537                        //pIter(p2);
538                        //sp  =   pAdd(ppMult_qq(temp->getT1(),p1),pMult_nn(ppMult_qq(temp->getT2(),p2),sign)); 
539                        sp      =   ksOldSpolyRedNew(ppMult_qq(temp->getT1(),temp->getLp1Poly()),
540                                         ppMult_qq(temp->getT2(),temp->getLp2Poly()));
541                        //Print("BEGIN SPOLY1\n====================\n");
542                        pNorm(sp);
543                        //pWrite(sp);
544                        //Print("END SPOLY1\n====================\n");
545                        if(NULL == sp) {
546                            // as rules consist only of pointers we need to save the labeled
547                            // S-polynomial also of a zero S-polynomial
548                            //zeroList->insert(temp->getAdT1(),temp->getLp1Index(),&sp);
549                            // origin of rule can be set NULL as the labeled polynomial
550                            // will never be used again as it is zero => no problems with
551                            // further criterion2() tests and termination conditions
552                            //Print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ZERO REDUCTION~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n");
553                                                        reductionsToZero++;
554                        //Print("IN SPOLS 1\n");
555                            //Rule* rNew  =  new Rule(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
556                            //rules->insertOrdered(rNew);
557                            rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
558                            //Print("RULE ADDED: \n");
559                            //pWrite(rules->getFirst()->getRuleTerm());
560                            //rules->print();
561                            // as sp = NULL, delete it
562                            pDelete(&sp);
563                            //Print("HIER\n");
564                        }
565                        else {
566                        //Print("IN SPOLS 2\n");
567                            //Rule* rNew  =  new Rule(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
568                            //rules->insertOrdered(rNew);
569                            rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
570                            //Print("RULE ADDED: \n");
571                            //pWrite(rules->getFirst()->getRuleTerm()); 
572                            //rules->print();
573                            sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rules->getFirst()->getRule());
574                            //sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rNew);
575                        }
576                        // data is saved in sPolyList or zero => delete sp
577                    }
578                }
579                else { // temp->getLp2Index() < temp->getLp1Index()
580                    // computation of S-polynomial
581                        //poly p1 =   temp->getLp1Poly();
582                        //poly p2 =   temp->getLp2Poly();
583                        //pIter(p1);
584                        //pIter(p2);
585                        //sp  =   pAdd(ppMult_qq(temp->getT1(),p1),pMult_nn(ppMult_qq(temp->getT2(),p2),sign)); 
586                    sp      =   ksOldSpolyRedNew(ppMult_qq(temp->getT1(),temp->getLp1Poly()),
587                                     ppMult_qq(temp->getT2(),temp->getLp2Poly()));
588                    //Print("BEGIN SPOLY2\n====================\n");
589                    pNorm(sp);
590                    //pWrite(sp);
591                    //Print("END SPOLY2\n====================\n");
592                    if(NULL == sp) {
593                        // zeroList->insert(temp->getAdT1(),temp->getLp1Index(),&sp);
594                        // origin of rule can be set NULL as the labeled polynomial
595                        // will never be used again as it is zero => no problems with
596                        // further criterion2() tests and termination conditions
597                            //Print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ZERO REDUCTION~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n");
598                        reductionsToZero++;
599                        //Print("IN SPOLS 3\n");
600                        rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
601                        //Print("RULE ADDED: \n");
602                        //pWrite(rules->getFirst()->getRuleTerm());
603                        //rules->print();
604                        // as sp = NULL, delete it
605                        pDelete(&sp);
606                    }
607                    else {
608                        //Print("IN SPOLS 4\n");
609                       
610                        //////////////////////////////////////////////////////////
611                        //////////////////////////////////////////////////////////
612                        // TODO: Rules inserted ordered by their label monomial!//
613                        //////////////////////////////////////////////////////////
614                        //////////////////////////////////////////////////////////
615
616                        //Rule* rNew      =   new Rule(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
617                        //RNode* rNodeNew =   new RNode(rNew);
618                        //rules->insertOrdered(rNew);
619                        rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
620                        //Print("RULE ADDED: \n");
621                        //pWrite(rules->getFirst()->getRuleTerm());
622                        //rules->print();
623                        //Print("%d\n",rules->getFirst()->getRuleIndex());
624                        //Print("%p\n",sPolyList->getFirst());
625                        sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rules->getFirst()->getRule());
626                        //sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rNew);
627                    }
628                }
629            }
630        //}
631        //Print("%p\n",temp);
632        temp    =   temp->getNext();
633        //Print("%p\n",temp);
634        //Print("%p\n",temp->getData());
635        //pWrite(temp->getLp1Poly());
636    }
637    // these critical pairs can be deleted now as they are either useless for further computations or
638    // already saved as an S-polynomial to be reduced in the following
639    delete first;   
640}
641
642
643
644/*
645========================================================================
646reduction including subalgorithm topReduction() using Faugere's criteria
647========================================================================
648*/
649void reduction(LList* sPolyList, CList* critPairs, LList* gPrev, RList* rules, LTagList* lTag, RTagList* rTag, 
650                 ideal gbPrev) { 
651    //Print("##########################################In REDUCTION!########################################\n");
652    // check if sPolyList has any elements
653    // NOTE: due to initialization sPolyList always has a default NULL element
654    LNode* temp = sPolyList->getFirst();
655    while(NULL != temp) {
656        // temp is the first element in the sPolyList which should be reduced
657        // due to earlier sorting this is the element of minimal degree AND
658        // minimal label
659        // delete the above first element from sPolyList, temp will be either reduced to
660        // zero or added to gPrev, but never come back to sPolyList
661        //pWrite(sPolyList->getFirst()->getPoly());
662        //Print("LIST OF SPOLYNOMIALS!\n");
663        //sPolyList->print();
664        sPolyList->setFirst(temp->getNext());
665        poly tempNF = kNF(gbPrev,currQuotient,temp->getPoly());
666        if(NULL != tempNF) {
667            pNorm(tempNF);
668            temp->setPoly(tempNF);
669            // try further reductions of temp with polynomials in gPrev
670            // with label index = current label index: this is done such that there
671            // is no label corruption during the reduction process
672            //Print("lower label reduction:  ");
673            //pWrite(tempNF);
674            topReduction(temp,sPolyList,gPrev,critPairs,rules,lTag,rTag,gbPrev);
675       
676        }
677        else {
678            reductionsToZero++;
679            delete temp;
680        }
681        //if(NULL != temp->getPoly()) {
682        //    criticalPair(gPrev,critPairs,lTag,rTag,rules);
683        //}
684        temp =   sPolyList->getFirst();
685    }
686    //sPolyList->print();
687    //delete sPolyList;
688}   
689
690/*
691========================================================================
692reduction including subalgorithm topReduction() using Faugere's criteria
693========================================================================
694*/
695void newReduction(LList* sPolyList, CList* critPairs, LList* gPrev, RList* rules, LTagList* lTag, RTagList* rTag, 
696                 ideal gbPrev) { 
697    //Print("##########################################In REDUCTION!########################################\n");
698    // check if sPolyList has any elements
699    // NOTE: due to initialization sPolyList always has a default NULL element
700    //Print("--1--\n");
701    LNode* temp = sPolyList->getFirst();
702    while(NULL != temp) {
703        // temp is the first element in the sPolyList which should be reduced
704        // due to earlier sorting this is the element of minimal degree AND
705        // minimal label
706        // delete the above first element from sPolyList, temp will be either reduced to
707        // zero or added to gPrev, but never come back to sPolyList
708        //Print("LIST OF SPOLYNOMIALS!\n");
709        //sPolyList->print();
710        //pWrite(sPolyList->getFirst()->getPoly());
711        sPolyList->setFirst(temp->getNext());
712        //pWrite(temp->getPoly());
713        //poly tempNF = kNF(gbPrev,currQuotient,temp->getPoly());
714        //Print("!!!\n");
715        //if(NULL != tempNF) {
716            //pNorm(tempNF);
717            //temp->setPoly(tempNF);
718            //Print("lower label reduction:  ");
719            //pWrite(tempNF);
720            // try further reductions of temp with polynomials in gPrev
721            // with label index = current label index: this is done such that there
722            // is no label corruption during the reduction process
723            findReducers(temp,sPolyList,gbPrev,gPrev,critPairs,rules,lTag,rTag); 
724        //}
725        //else {
726        //    reductionsToZero++;
727        //    delete temp;
728        //}
729        //if(NULL != temp->getPoly()) {
730        //    criticalPair(gPrev,critPairs,lTag,rTag,rules);
731        //}
732        //Print("HIER AUCH ?\n");
733        //Print("--2--\n");
734        //sPolyList->print();
735        //critPairs->print();
736        temp =   sPolyList->getFirst();
737        //Print("%p\n",temp);
738    }
739    //sPolyList->print();
740    //delete sPolyList;
741    //Print("REDUCTION FERTIG\n");
742}     
743
744
745/*!
746 * ================================================================================
747 * searches for reducers of temp similar to the symbolic preprocessing of F4  and
748 * divides them into a "good" and "bad" part:
749 *
750 * the "good" ones are the reducers which do not corrupt the label of temp, with
751 * these the normal form of temp is computed
752 *
753 * the "bad" ones are the reducers which corrupt the label of temp, they are tested
754 * later on for possible new rules and S-polynomials to be added to the algorithm
755 * ================================================================================
756*/
757void findReducers(LNode* l, LList* sPolyList, ideal gbPrev, LList* gPrev, CList* critPairs, RList* rules, LTagList* lTag, RTagList* rTag) {
758    int canonicalize    =   0;
759    //int timerRed        =   0;
760    number sign         =   nInit(-1);
761    LList* good         =   new LList();
762    LList* bad          =   new LList();
763    LList* monomials    =   new LList(l->getLPoly());
764    poly u              =   pOne();
765    number nOne         =   nInit(1);
766    LNode* tempRed      =   lTag->getFirstCurrentIdx();
767    LNode* tempMon      =   monomials->getFirst();
768    poly tempPoly       =   pInit();
769    poly redPoly        =   NULL;
770    int idx             =   l->getIndex();
771    //Print("IN FIND REDUCERS:  ");
772    //pWrite(l->getPoly());
773    tempPoly    =   pCopy(l->getPoly());
774    //tempMon->setPoly(tempPoly);
775    //while(NULL != tempMon) {
776        // iteration over all monomials in tempMon
777        kBucket* bucket  =   kBucketCreate();
778        kBucketInit(bucket,tempPoly,0);
779        tempPoly    =   kBucketGetLm(bucket);
780        //Print("\n\n\nTO BE REDUCED:  ");
781        //pWrite(l->getPoly());
782        //pWrite(tempPoly);
783        while(NULL != tempPoly) {
784            // iteration of all elements in gPrev of the current index
785            tempRed =   gPrev->getFirst();
786            while(NULL != tempRed) {
787                //Print("TEMPREDPOLY:  ");
788                //pWrite(tempRed->getPoly());
789                if(pLmDivisibleByNoComp(tempRed->getPoly(),tempPoly)) {
790                    u   =   pDivideM(pHead(tempPoly),pHead(tempRed->getPoly()));
791                    //Print("U:  ");
792                    //pWrite(u);
793                    if(tempRed->getIndex() != idx) {
794                            // passing criterion1 ?
795                            if(!criterion1(gPrev,u,tempRed,lTag)) {
796                                    poly tempRedPoly    =   tempRed->getPoly();
797                                    //Print("REDUCER: ");
798                                    //pWrite(ppMult_qq(u,tempRedPoly));
799                                    pIter(tempRedPoly);
800                                    int lTempRedPoly    =   pLength(tempRedPoly);
801                                    kBucketExtractLm(bucket);
802                                    kBucket_Minus_m_Mult_p(bucket,u,tempRedPoly,&lTempRedPoly);
803                                    canonicalize++;
804                                    if(!(canonicalize % 50)) {
805                                        kBucketCanonicalize(bucket);
806                                    }
807                                    tempPoly    =   kBucketGetLm(bucket);
808                                    //Print("TEMPPOLY:  ");
809                                    //pWrite(tempPoly);
810                                    if(NULL != tempPoly) {
811                                        tempRed     =   gPrev->getFirst();
812                                        continue;
813                                    }
814                                    else {
815                                        break;
816                                    }
817                                    //tempRedPoly =   ppMult_qq(u,tempRedPoly);
818                                    //pMult_nn(tempRedPoly,sign);
819                                    //pAdd(tempPoly,tempRedPoly);
820                                    //good->insert(pOne(),1,ppMult_qq(u,tempRed->getPoly()),NULL);
821                                    //break;
822                             }   
823                   
824                    }
825                    else {
826                        if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) != 0) {
827                            // passing criterion2 ?
828                            if(!criterion2(gPrev->getFirst()->getIndex(), u,tempRed,rules,rTag)) {
829                                // passing criterion1 ?
830                                if(!criterion1(gPrev,u,tempRed,lTag)) {
831                                    if(pLmCmp(ppMult_qq(u,tempRed->getTerm()),l->getTerm()) == 1) {
832                                        if(NULL == redPoly) {
833                                            bad->insert(tempRed->getLPoly());
834                                            //poly tempRedPoly    =   tempRed->getPoly();
835                                            //break;
836                                        }
837                                        /*
838                                        // DEBUGGING STUFF
839                                        //
840                                        //Print("REDUCER: ");
841                                        //pWrite(ppMult_qq(u,tempRedPoly));
842                                        rules->insert(l->getIndex(),ppMult_qq(u,tempRed->getTerm()));
843                                        l->setRule(rules->getFirst()->getRule());
844                                        poly tempRedPoly    =   tempRed->getPoly();
845                                        pIter(tempRedPoly);
846                                        int lTempRedPoly    =   pLength(tempRedPoly);
847                                        //Print("HEAD MONOMIAL KBUCKET: ");
848                                        //pWrite(kBucketGetLm(bucket));
849                                        kBucketExtractLm(bucket);
850                                        kBucket_Minus_m_Mult_p(bucket,u,tempRedPoly,&lTempRedPoly);
851                                        canonicalize++;
852                                        if(!(canonicalize % 50)) {
853                                            kBucketCanonicalize(bucket);
854                                        }
855                                        //Print("HEAD MONOMIAL KBUCKET: ");
856                                        //pWrite(kBucketGetLm(bucket));
857                                        tempPoly    =   kBucketGetLm(bucket);
858                                        //Print("TEMPPOLY:  ");
859                                        //pWrite(tempPoly);
860                                        if(NULL != tempPoly) {
861                                            tempRed     =   gPrev->getFirst();
862                                            continue;
863                                        }
864                                        else {
865                                            break;
866                                        }
867                                        */
868                                    }
869                                    else {
870                                        poly tempRedPoly    =   tempRed->getPoly();
871                                        //Print("REDUCER: ");
872                                        //pWrite(ppMult_qq(u,tempRedPoly));
873                                        pIter(tempRedPoly);
874                                        int lTempRedPoly    =   pLength(tempRedPoly);
875                                        //Print("HEAD MONOMIAL KBUCKET: ");
876                                        //pWrite(kBucketGetLm(bucket));
877                                        kBucketExtractLm(bucket);
878                                        kBucket_Minus_m_Mult_p(bucket,u,tempRedPoly,&lTempRedPoly);
879                                        canonicalize++;
880                                        if(!(canonicalize % 50)) {
881                                            kBucketCanonicalize(bucket);
882                                        }
883                                        //Print("HEAD MONOMIAL KBUCKET: ");
884                                        //pWrite(kBucketGetLm(bucket));
885                                        tempPoly    =   kBucketGetLm(bucket);
886                                        //Print("TEMPPOLY:  ");
887                                        //pWrite(tempPoly);
888                                        if(NULL != tempPoly) {
889                                            tempRed     =   gPrev->getFirst();
890                                            continue;
891                                        }
892                                        else {
893                                            break;
894                                        }
895                                    }
896                                }   
897                            }
898                        }
899                    }
900                   
901                }
902                tempRed =   tempRed->getNext();
903            }
904            //pWrite(kBucketGetLm(bucket));
905            //Print("HIER AUCH\n");
906            if(NULL == redPoly) {
907                redPoly =   kBucketExtractLm(bucket);
908                //Print("1: ");
909                //pWrite(redPoly);
910            }
911            else {
912                redPoly     =   p_Merge_q(redPoly,kBucketExtractLm(bucket),currRing);
913                //Print("REDPOLY: ");
914                //pWrite(redPoly);
915            }
916            tempPoly    =   kBucketGetLm(bucket);
917            //Print("TEMPPOLY2:  ");
918            //pWrite(tempPoly);
919        }
920        if(NULL == redPoly) {
921            reductionsToZero++;
922            //pDelete(&redPoly);
923        }
924        else {
925            //Print("\nELEMENT ADDED TO GPREV: ");
926            //pWrite(redPoly);
927            pNorm(redPoly);
928            //Print("%d\n",canonicalize);
929            l->setPoly(redPoly);
930            gPrev->insert(l->getLPoly());
931            criticalPair(gPrev,critPairs,lTag,rTag,rules);
932        }
933    //    tempMon =   tempMon->getNext();
934    //}
935   
936    // if there are reducers than reduce l
937   
938    /*   
939    if(NULL != good->getFirst()) {
940       
941        //Print("TO BE REDUCED:\n");
942        //pWrite(l->getPoly());
943        LNode* tempGood =   good->getFirst();
944        kBucket* bucket =   kBucketCreate();
945        kBucketInit(bucket,l->getPoly(),0);
946        //Print("\nREDUCERS: \n");
947        while(NULL != tempGood) {
948            //pWrite(tempGood->getPoly());
949            int lTempGood   =   pLength(tempGood->getPoly());
950            kBucket_Minus_m_Mult_p(bucket,pOne(),tempGood->getPoly(),&lTempGood);
951            //Print("KBUCKET:  ");
952            //pWrite(kBucketGetLm(bucket));
953            tempGood    =   tempGood->getNext();
954        }
955        poly temp   =  kBucketClear(bucket);
956        //pWrite(temp);
957        //Print("\n");
958        if(NULL != temp) {
959            pNorm(temp);
960            //Print("NEW REDUCTION:  ");
961            //pWrite(temp);
962            l->setPoly(temp);
963            //Print("ELEMENT ADDED TO GPREV1: ");
964            //pWrite(l->getPoly());
965            //pWrite(l->getTerm());
966            //Print("%p\n",gPrev->getLast());
967            //pWrite(gPrev->getLast()->getPoly());
968            gPrev->insert(l->getLPoly());
969            //Print("%p\n",gPrev->getLast());
970            //pWrite(gPrev->getLast()->getPoly());
971            //rules->print();
972            criticalPair(gPrev,critPairs,lTag,rTag,rules);
973            //Print("LIST OF CRITICAL PAIRS:    \n");
974            //critPairs->print();
975            //gPrev->print();
976        }
977        else {
978            reductionsToZero++;
979            //Print("ZERO REDUCTION\n");
980            pDelete(&temp);
981        }
982         
983        /*
984        //Print("HIER IN GOOD REDUCTION\n");
985        LNode* tempGood         =   good->getFirst();
986        ideal reductionId       =   idInit(good->getLength(),1);
987        int i;
988        //Print("\n\n");
989        //good->print();
990        for(i=0;i<good->getLength();i++) {
991            reductionId->m[i]   =   tempGood->getPoly();
992            //Print("REDUCERS:");
993            //pWrite(tempGood->getPoly());
994            //Print("%p\n",tempGood);
995            tempGood            =   tempGood->getNext();
996        }
997        //idShow(reductionId);
998        //reductionId =   idAdd(reductionId,gbPrev);
999        //Print("\n\nREDUCTION PROCESS DONE TWICE!");
1000        //idDelMultiples(reductionId);
1001        //idShow(reductionId);
1002        //pWrite(l->getPoly());
1003        poly temp   =   kNF(reductionId,currQuotient,l->getPoly());
1004        //pWrite(temp);
1005        //poly temp2  =   kNF(reductionId,currQuotient,temp);
1006        //pWrite(temp2);
1007        //Print("\n\n");
1008        if(NULL != temp) {
1009            pNorm(temp);
1010            //Print("NEW REDUCTION:  ");
1011            //pWrite(temp);
1012            l->setPoly(temp);
1013            //Print("ELEMENT ADDED TO GPREV1: ");
1014            //pWrite(l->getPoly());
1015            //pWrite(l->getTerm());
1016            //Print("%p\n",gPrev->getLast());
1017            //pWrite(gPrev->getLast()->getPoly());
1018            gPrev->insert(l->getLPoly());
1019            //Print("%p\n",gPrev->getLast());
1020            //pWrite(gPrev->getLast()->getPoly());
1021            //rules->print();
1022            criticalPair(gPrev,critPairs,lTag,rTag,rules);
1023            //Print("LIST OF CRITICAL PAIRS:    \n");
1024            //critPairs->print();
1025            //gPrev->print();
1026        }
1027        else {
1028            reductionsToZero++;
1029            //Print("ZERO REDUCTION\n");
1030            pDelete(&temp);
1031        }
1032        //pWrite(temp);
1033        //for(i=0;i<IDELEMS(reductionId);i++) {
1034        //}
1035        //idShow(reductionId);
1036        //Print("HIER\n");
1037        //Print("ADDRESS OF IDEAL: %p\n",&reductionId);
1038        idDelete(&reductionId);
1039        //Print("HIER\n");
1040       
1041
1042    }
1043   
1044
1045    else {
1046        //pWrite(l->getPoly());
1047        gPrev->insert(l->getLPoly());
1048        //Print("ELEMENT ADDED TO GPREV2: ");
1049        //pWrite(l->getPoly());
1050        //pWrite(l->getTerm());
1051        //Print("GENAU HIER:  ");
1052        //pWrite(l->getPoly());
1053        criticalPair(gPrev,critPairs,lTag,rTag,rules);
1054        //Print("LIST OF CRITICAL PAIRS:    \n");
1055        //critPairs->print();
1056    }
1057    */
1058    // if there are "bad" reducers than try to compute new S-polynomials and rules
1059   
1060    if(NULL != bad->getFirst()) {
1061        //Print("BAD STUFF LIST:\n");
1062        //bad->print();
1063        LNode* tempBad  =   bad->getFirst();
1064        //pWrite(l->getPoly());
1065        while(NULL != tempBad) {
1066            if(pDivisibleBy(tempBad->getPoly(),l->getPoly())) {
1067                //Print("BAD STUFF\n");
1068                //pWrite(l->getPoly());
1069                //pWrite(tempBad->getPoly());
1070                u   =   pDivide(pHead(l->getPoly()),pHead(tempBad->getPoly()));
1071                //Print("MULTIPLIER:  ");
1072                //pWrite(u);
1073                pSetCoeff(u,nOne);
1074                if(pLmCmp(ppMult_qq(u,tempBad->getTerm()),l->getTerm()) != 0) {
1075                    // passing criterion2 ?
1076                    if(!criterion2(gPrev->getFirst()->getIndex(), u,tempBad,rules,rTag)) {
1077                        // passing criterion1 ?
1078                        if(!criterion1(gPrev,u,tempBad,lTag)) {
1079                            //Print("HIERHIERHIERHIERHIERHIER\n");
1080                            rules->insert(tempBad->getIndex(),ppMult_qq(u,tempBad->getTerm()));
1081                            //gPrev->print();
1082                            //pWrite(l->getPoly());
1083                            poly temp   =   ksOldSpolyRedNew(ppMult_qq(u,tempBad->getPoly()),l->getPoly());
1084                            //pWrite(l->getPoly());
1085                            //Print("%p\n",temp);
1086                            //gPrev->print();
1087                            if(NULL != temp) {
1088                                pNorm(temp);
1089                                LNode* tempBadNew   =   new LNode(ppMult_qq(u,tempBad->getTerm()),tempBad->getIndex(),temp,rules->getFirst()->getRule());
1090                                //pWrite(temp);
1091                                //pWrite(tempBadNew->getPoly());
1092                                //pWrite(tempBadNew->getTerm());
1093                                //pWrite(pHead(tempBadNew->getPoly()));
1094                                //Print("%p\n",tempBadNew->getPoly());
1095                                //tempRed->getLPoly()->setRule(rules->getFirst()->getRule());
1096                                tempBadNew->setDel(1);
1097                           
1098                                sPolyList->insertByLabel(tempBadNew);
1099                                //Print("BAD SPOLYLIST: \n");
1100                                //sPolyList->print();
1101                            }
1102                        }
1103                    }
1104                }
1105            }
1106        //Print("HIER\n");
1107            tempBad =   tempBad->getNext();
1108            //Print("%p\n",tempBad);
1109        }
1110       // Print("-------------------BAD STUFF LIST-----------------------------\n");
1111    }
1112        //Print("HIER AUCH\n");
1113        //Print("SPOLYLIST IN BAD: \n");
1114        //sPolyList->print();
1115    //Print("END FIND REDUCERS\n");
1116}
1117
1118/*
1119=======================================================================================
1120merging 2 polynomials p & q without requiring that all monomials of p & q are different
1121if there are equal monomials in p & q only one of these monomials (always that of p!)
1122is taken into account
1123=======================================================================================
1124
1125poly p_MergeEq_q(poly p, poly q, const ring r) {
1126  assume(p != NULL && q != NULL);
1127  p_Test(p, r);
1128  p_Test(q, r);
1129#if PDEBUG > 0
1130  int l = pLength(p) + pLength(q);
1131#endif
1132
1133  spolyrec rp;
1134  poly a = &rp;
1135  DECLARE_LENGTH(const unsigned long length = r->CmpL_Size);
1136  DECLARE_ORDSGN(const long* ordsgn = r->ordsgn);
1137
1138  Top:     // compare p and q w.r.t. monomial ordering
1139  p_MemCmp(p->exp, q->exp, length, ordsgn, goto Equal, goto Greater , goto Smaller);
1140
1141  Equal:
1142  a =   pNext(a)    =   p;
1143  pIter(p);
1144  pIter(q);
1145  if(NULL == p) {
1146      if(NULL == q) {
1147          goto Finish;
1148      }
1149      else {
1150          pNext(a)  =   q;
1151          goto Finish;
1152      }
1153  }
1154  goto Top;
1155
1156  Greater:
1157  a = pNext(a) = p;
1158  pIter(p);
1159  if (p==NULL) { pNext(a) = q; goto Finish;}
1160  goto Top;
1161
1162  Smaller:
1163  a = pNext(a) = q;
1164  pIter(q);
1165  if (q==NULL) { pNext(a) = p; goto Finish;}
1166  goto Top;
1167
1168  Finish:
1169
1170  p_Test(pNext(&rp), r);
1171#if PDEBUG > 0
1172  pAssume1(l - pLength(pNext(&rp)) == 0);
1173#endif
1174  return pNext(&rp);
1175}
1176*/
1177
1178/*
1179=====================================================================================
1180top reduction in F5, i.e. reduction of a given S-polynomial by labeled polynomials of
1181the same index whereas the labels are taken into account
1182=====================================================================================
1183*/
1184void topReduction(LNode* l, LList* sPolyList, LList* gPrev, CList* critPairs,  RList* rules, LTagList* lTag, RTagList* rTag, ideal gbPrev) {
1185    //Print("##########################################In topREDUCTION!########################################\n");
1186    // try l as long as there are reductors found by findReductor()
1187    LNode* gPrevRedCheck    =   lTag->getFirstCurrentIdx();
1188    LNode* tempRed;
1189    poly pOne               =   pOne();
1190    do {
1191        //int timer5  =   initTimer();
1192        //startTimer();
1193        //Print("TOP REDUCTION:  ");
1194        //pWrite(l->getPoly());
1195        tempRed  =   findReductor(l,sPolyList,gPrevRedCheck,gPrev,rules,lTag,rTag);
1196        //timer5  =   getTimer();
1197        //reductionTime   =   reductionTime   +   timer5;
1198        // if a reductor for l is found and saved in tempRed
1199        if(NULL != tempRed) {
1200            // if label of reductor is greater than the label of l we have to built a new element
1201            // and add it to sPolyList
1202           
1203            if(pLmCmp(tempRed->getTerm(),l->getTerm()) == 1) {
1204                // needed sinc pSub destroys the arguments!
1205                //poly temp_poly_l    =   pInit();
1206                //temp_poly_l         =   pCopy(l->getPoly());
1207                //Print("VORHER: ");
1208                //pWrite(tempRed->getPoly());
1209                //poly temp           =   pMinus_mm_Mult_qq(tempRed->getPoly(),pOne,l->getPoly());
1210                poly temp   =   ksOldSpolyRedNew(l->getPoly(),tempRed->getPoly());
1211                rules->insert(tempRed->getIndex(),tempRed->getTerm());
1212                //Print("NACHHER: ");
1213                //pWrite(tempRed->getPoly());
1214                //Print("TEMP: ");
1215                //pWrite(temp);
1216                if(NULL != temp) {
1217                    pNorm(temp);
1218                    //tempRed->setPoly(temp);
1219                    //tempRed->setDel(1);
1220                    //rules->insert(tempRed->getIndex(),tempRed->getTerm());
1221                    LNode* tempRedNew   =   new LNode(tempRed->getTerm(),tempRed->getIndex(),temp,rules->getFirst()->getRule());
1222                    //tempRed->getLPoly()->setRule(rules->getFirst()->getRule());
1223                    tempRedNew->setDel(1);
1224                    sPolyList->insertByLabel(tempRedNew);
1225                }
1226                else {
1227                    pDelete(&temp);
1228                    reductionsToZero++;
1229                    //delete tempRed;
1230                }
1231            }
1232           
1233            // label of reductor is smaller than the label of l, subtract reductor from l and delete the
1234            // gPrevRedCheck pointer added to l during findReductor() as the head term of l changes
1235            // after subtraction
1236            else {
1237               
1238                //poly temp_poly_l    =   pInit();
1239                //temp_poly_l         =   pCopy(l->getPoly());
1240                //poly temp   =   pMinus_mm_Mult_qq(tempRed->getPoly(),pOne,l->getPoly());
1241                //Print("REDUCER: ");
1242                //pWrite(tempRed->getPoly());
1243                //pWrite(tempRed->getTerm());
1244                poly temp   =   ksOldSpolyRedNew(l->getPoly(),tempRed->getPoly());
1245                //Print("REDUCED ELEMENT:  ");
1246                if(NULL != temp) {
1247                    pNorm(temp);
1248                    //pWrite(temp);
1249                    poly tempNF =   kNF(gbPrev,currQuotient,temp); 
1250                    pNorm(tempNF);
1251                    if(NULL == tempNF) {
1252                        reductionsToZero++;
1253                        pDelete(&tempNF);
1254                        l->setPoly(NULL);
1255                        break;
1256                    }
1257                    l->setPoly(tempNF);
1258                   
1259                    gPrevRedCheck   =   lTag->getFirstCurrentIdx();
1260                }
1261                else {
1262                    reductionsToZero++;
1263                    pDelete(&temp);
1264                    l->setPoly(NULL);
1265                    break;
1266                }
1267            }   
1268        }
1269        else {
1270            if(NULL != l->getPoly()) {
1271                pNorm(l->getPoly());
1272                //Print("ELEMENT ADDED TO GPREV: ");
1273                //pWrite(l->getPoly());
1274                gPrev->insert(l->getLPoly());
1275                //Print("TEMP RED == 0  ");
1276                //pWrite(l->getPoly());
1277                //pWrite(l->getTerm());
1278                //rules->print();
1279                criticalPair(gPrev,critPairs,lTag,rTag,rules);
1280                //Print("LIST OF CRITICAL PAIRS:    \n");
1281                //critPairs->print();
1282            }
1283            break;
1284        }
1285    } while(1);
1286}
1287
1288
1289/*
1290=====================================================================
1291subalgorithm to find a possible reductor for the labeled polynomial l
1292=====================================================================
1293*/
1294LNode* findReductor(LNode* l, LList* sPolyList, LNode* gPrevRedCheck, LList* gPrev, RList* rules, LTagList* lTag,RTagList* rTag) {
1295    // allociation of memory for the possible reductor
1296    //Print("LPOLY:  ");
1297    //pWrite(l->getPoly());
1298    poly u      =   pOne();
1299    poly red;
1300    number nOne =   nInit(1);
1301    LNode* temp;
1302    // head term of actual element such that we do not have to call pHead at each new reductor test
1303    poly t      =   pHead(l->getPoly());
1304    // if l was already checked use the information in gPrevRedCheck such
1305    // that we can start searching for new reducers from this point and
1306    // not from the first element of gPrev with the current index
1307    temp    =   gPrevRedCheck;
1308    // search for reductors until we are at the end of gPrev resp. at the
1309    // end of the elements of the current index
1310    while(NULL != temp && temp->getIndex() == l->getIndex()) {
1311        // does the head of the element of gPrev divides the head of
1312        // the to be reduced element?
1313        if(pLmDivisibleByNoComp(pHead(temp->getPoly()),t)) {
1314            // get all the information needed for the following tests
1315            // of the criteria
1316            u   =   pDivide(t,pHead(temp->getPoly()));
1317            pSetCoeff(u,nOne);
1318            red =   ppMult_qq(u,temp->getPoly());
1319            pNorm(red);
1320            // check if both have the same label
1321            if(pLmCmp(ppMult_qq(u,temp->getTerm()),l->getTerm()) != 0) {
1322                // passing criterion2 ?
1323                if(!criterion2(gPrev->getFirst()->getIndex(), u,temp,rules,rTag)) {
1324                    // passing criterion1 ?
1325                    if(!criterion1(gPrev,u,temp,lTag)) {
1326                            gPrevRedCheck   =   temp->getNext();
1327                            LNode* redNode  =   new LNode(ppMult_qq(u,temp->getTerm()),temp->getIndex(),red,NULL,NULL);
1328                            return redNode;
1329                    }
1330                }
1331            }
1332            if(pLmCmp(ppMult_qq(u,temp->getTerm()),l->getTerm()) == 1) {
1333                // passing criterion2 ?
1334                if(!criterion2(gPrev->getFirst()->getIndex(), u,temp,rules,rTag)) {
1335                    // passing criterion1 ?
1336                    if(!criterion1(gPrev,u,temp,lTag)) {
1337                        poly tempSpoly  =   ksOldSpolyRedNew(red,l->getPoly());
1338                        rules->insert(temp->getIndex(),ppMult_qq(u,temp->getTerm()));
1339                        gPrevRedCheck   =   temp->getNext();
1340                        if(NULL != tempSpoly) {
1341                            pNorm(tempSpoly);
1342                            sPolyList->insertByLabel(ppMult_qq(u,temp->getTerm()),temp->getIndex(),tempSpoly,rules->getFirst()->getRule());
1343                    //Print("NEW ONE: ");
1344                    //pWrite(tempSpoly);
1345                    //Print("HIER\n");
1346                            //sPolyList->print();
1347                            //sleep(5);
1348                        }
1349                    }
1350                }
1351            }
1352        }
1353        //Print("AUCH HIER\n");
1354        temp = temp->getNext();
1355    }
1356   
1357//    delete temp;
1358 return NULL;
1359}
1360
1361
1362
1363/*
1364==========================================================================
1365MAIN:computes a gb of the ideal i in the ring r with our F5 implementation
1366==========================================================================
1367*/
1368ideal F5main(ideal id, ring r) {
1369    int timer   =   initTimer();
1370    startTimer();
1371    int i,j,k,l;
1372    int gbLength;
1373    // 1 polynomial for defining initial labels & further tests
1374    poly ONE = pOne();
1375    poly pOne   =   pOne();
1376    number nOne =   nInit(1);
1377    // tag the first element of index i-1 for criterion 1
1378    //Print("LTAG BEGINNING: %p\n",lTag);
1379   
1380    // DEBUGGING STUFF START
1381    //Print("NUMBER: %d\n",r->N);
1382    /*
1383    int* ev = new int[r->N +1];
1384    for(i=0;i<IDELEMS(id);i++) {
1385        pGetExpV(id->m[i],ev);
1386        //ev2  =   pGetExp(id->m[i],1);
1387        pWrite(id->m[i]);
1388        Print("EXP1: %d\n",ev[1]);
1389        Print("EXP2: %d\n",ev[2]);
1390        Print("EXP3: %d\n\n",ev[3]);
1391        Print("SUM: %ld\n\n\n",sumVector(ev,r->N));
1392    }
1393    delete ev;
1394    */
1395    /*DEBUGGING STUFF END */
1396   
1397    // first element in rTag is first element of rules which is NULL RNode,
1398    // this must be done due to possible later improvements
1399    RList* rules    =   new RList();
1400    //Print("RULES FIRST: %p\n",rules->getFirst());
1401    //Print("RULES FIRST DATA: %p\n",rules->getFirst()->getRule());
1402    RTagList* rTag  =   new RTagList(rules->getFirst());
1403    i = 1;
1404    /*for(j=0; j<IDELEMS(id); j++) {
1405        if(NULL != id->m[j]) {
1406            if(pComparePolys(id->m[j],ONE)) {
1407                Print("One Polynomial in Input => Computations stopped\n");
1408                ideal idNew = idInit(1,1);
1409                idNew->m[0] = ONE;
1410                return(idNew);
1411            }   
1412        }
1413    }*/ 
1414    ideal idNew     =   kInterRed(id); 
1415    id              =   idNew;
1416    //qsortDegree(&id->m[0],&id->m[IDELEMS(id)-1]);
1417    //idShow(id);
1418    LList* gPrev    =   new LList(ONE, i, id->m[0]);
1419    //idShow(id);
1420    //Print("%p\n",id->m[0]);
1421    //pWrite(id->m[0]);
1422    //Print("%p\n",gPrev->getFirst()->getPoly());
1423    //pWrite(gPrev->getFirst()->getPoly());
1424
1425    LTagList* lTag  =   new LTagList(gPrev->getFirst());
1426    //lTag->insert(gPrev->getFirst());
1427    lTag->setFirstCurrentIdx(gPrev->getFirst());
1428    // computing the groebner basis of the elements of index < actual index
1429    gbLength    =   gPrev->getLength();
1430    //Print("Laenge der bisherigen Groebner Basis: %d\n",gPrev->getLength());
1431    ideal gbPrev    =   idInit(gbLength,1);
1432    // initializing the groebner basis of elements of index < actual index
1433    gbPrev->m[0]    =   gPrev->getFirst()->getPoly();
1434    //idShow(gbPrev);
1435    //idShow(currQuotient);
1436    for(i=2; i<=IDELEMS(id); i++) {
1437        LNode* gPrevTag =   gPrev->getLast();
1438        //Print("Last POlynomial in GPREV: ");
1439        //Print("%p\n",gPrevTag);   
1440        //pWrite(gPrevTag->getPoly());
1441        gPrev   =   F5inc(i, id->m[i-1], gPrev, gbPrev, ONE, lTag, rules, rTag);
1442        Print("%d\n",gPrev->count(gPrevTag->getNext()));
1443        Print("%d\n",gPrev->getLength());
1444        //Print("____________________________________ITERATION STEP DONE________________________________________\n");
1445       
1446        // DEBUGGING STUFF
1447        LNode* temp    =   gPrev->getFirst();
1448   
1449
1450        /////////////////////////////////////////////////////////////////////////////////
1451        //                                                                             //
1452        // one needs to choose one of the following 3 implementations of the algorithm //
1453        // F5,F5R or F5C                                                               //
1454        //                                                                             //
1455        /////////////////////////////////////////////////////////////////////////////////                                                                           
1456       
1457       
1458        //   
1459        // remove this comment to get "F5"
1460        //
1461        /*
1462        if(gPrev->getLength() > gbLength) {
1463            if(i < IDELEMS(id)) {
1464                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1465                LNode*  temp =   gPrevTag;
1466                int counter =   0;
1467                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1468                    temp        =   temp->getNext();
1469                    if(0 == temp->getDel()) {
1470                        counter++;
1471                        gbAdd->m[j] =   temp->getPoly();
1472                    }
1473                }
1474                    gbPrev          =   idAdd(gbPrev,gbAdd);
1475            }
1476            else {
1477                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1478                LNode*  temp =   gPrevTag;
1479                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1480                    temp        =   temp->getNext();
1481                    gbAdd->m[j] =   temp->getPoly();
1482                }
1483                gbPrev          =   idAdd(gbPrev,gbAdd);
1484            }
1485        }
1486        gbLength    =   gPrev->getLength();
1487        */
1488       
1489
1490        //
1491        // remove this comment to get "F5R"
1492        //
1493        /*
1494        if(gPrev->getLength() > gbLength) {
1495            if(i < IDELEMS(id)) {
1496                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1497                LNode*  temp =   gPrevTag;
1498                int counter =   0;
1499                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1500                    temp        =   temp->getNext();
1501                    if(0 == temp->getDel()) {
1502                        counter++;
1503                        gbAdd->m[j] =   temp->getPoly();
1504                    }
1505                }
1506                    gbPrev          =   idAdd(gbPrev,gbAdd);
1507            }
1508            else {
1509                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1510                LNode*  temp =   gPrevTag;
1511                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1512                    temp        =   temp->getNext();
1513                    gbAdd->m[j] =   temp->getPoly();
1514                }
1515                gbPrev          =   idAdd(gbPrev,gbAdd);
1516            }
1517            // interreduction stuff
1518            // comment this out if you want F5 instead of F5R
1519            if(i<IDELEMS(id)) {
1520                ideal tempId    =   kInterRed(gbPrev);
1521                gbPrev          =   tempId;
1522            }
1523        }
1524        gbLength    =   gPrev->getLength();
1525        */
1526       
1527
1528        //
1529        // Remove this comment to get "F5C"
1530        // computing new groebner basis gbPrev
1531        //
1532         
1533        if(gPrev->getLength() > gbLength) {
1534            if(i < IDELEMS(id)) {
1535                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1536                LNode*  temp =   gPrevTag;
1537                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1538                    temp        =   temp->getNext();
1539                        gbAdd->m[j] =   temp->getPoly();
1540                }
1541                    gbPrev          =   idAdd(gbPrev,gbAdd);
1542            }
1543            else {
1544                ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
1545                LNode*  temp =   gPrevTag;
1546                for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
1547                    temp        =   temp->getNext();
1548                    gbAdd->m[j] =   temp->getPoly();
1549                }
1550                gbPrev          =   idAdd(gbPrev,gbAdd);
1551            }
1552            if(i<IDELEMS(id)) {
1553                ideal tempId    =   kInterRed(gbPrev);
1554                Print("HERE\n");
1555                gbPrev          =   tempId;
1556                delete gPrev;
1557                delete rules;
1558                gPrev    =   new LList(pOne,1,gbPrev->m[0]);
1559                gPrev->insert(pOne,1,gbPrev->m[1]);
1560                rules    =   new RList();
1561                rTag     =   new RTagList(rules->getFirst());
1562                for(k=2; k<IDELEMS(gbPrev); k++) {
1563                    gPrev->insert(pOne,k+1,gbPrev->m[k]);
1564                    for(l=0; l<k; l++) {
1565                    }
1566                    rTag->insert(rules->getFirst());
1567                }
1568            }
1569            gbLength    =   gPrev->getLength(); 
1570        } 
1571   
1572
1573
1574    }
1575    //Print("\n\nADDING TIME IN REDUCTION: %d\n\n",reductionTime);
1576    Print("\n\nNumber of zero-reductions:  %d\n",reductionsToZero);
1577    timer   =   getTimer();
1578    Print("Highest Degree during computations: %d\n",highestDegree);
1579    Print("Time for computations: %d/1000 seconds\n",timer);
1580    Print("Number of elements in gb: %d\n",IDELEMS(gbPrev));
1581    //LNode* temp    =   gPrev->getFirst();
1582    //while(NULL != temp) {
1583    //    pWrite(temp->getPoly());
1584    //    temp    =   temp->getNext();
1585    // }
1586    //gPrev->print();
1587    //delete lTag;
1588    delete rTag;
1589    delete gPrev;
1590    reductionsToZero    =   0;
1591    highestDegree       =   0;
1592    reductionTime       =   0;
1593    spolsTime           =   0;
1594    return(gbPrev);
1595
1596
1597}
1598
1599#endif
Note: See TracBrowser for help on using the repository browser.