source: git/kernel/f5gb.cc @ 85a342

spielwiese
Last change on this file since 85a342 was 85a342, checked in by Christian Eder, 15 years ago
criterion2() for previously computed elements deleted git-svn-id: file:///usr/local/Singular/svn/trunk@11581 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 42.5 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 "ideals.h"
17#include "febase.h"
18#include "kstd1.h"
19#include "khstd.h"
20#include "kbuckets.h"
21#include "weight.h"
22#include "intvec.h"
23#include "pInline1.h"
24#include "f5gb.h"
25#include "f5data.h"
26#include "f5lists.h"
27#include "timer.h"
28int reductionsToZero   =  0;
29
30/*
31====================================================================
32sorting ideals by decreasing total degree "left" and "right" are the
33pointer of the first and last polynomial in the considered ideal
34====================================================================
35*/
36void qsortDegree(poly* left, poly* right) {
37    poly* ptr1 = left;
38    poly* ptr2 = right;
39    poly p1,p2;
40    p2 = *(left + (right - left >> 1));
41    do {
42            while(pTotaldegree(*ptr1, currRing) < pTotaldegree(p2, currRing)) {
43                    ptr1++;
44            } 
45            while(pTotaldegree(*ptr2, currRing) > pTotaldegree(p2,currRing)) {
46                    ptr2--;
47            }
48            if(ptr1 > ptr2) {
49                    break;
50            }
51            p1    = *ptr1;
52            *ptr1 = *ptr2;
53            *ptr2 = p1;
54    } while(++ptr1 <= --ptr2);
55
56    if(left < ptr2) {
57            qsortDegree(left,ptr2);
58    }
59    if(ptr1 < right) {
60            qsortDegree(ptr1,right);
61    }
62}
63
64
65
66/*
67==================================================
68computes incrementally gbs of subsets of the input
69gb{f_m} -> gb{f_m,f_(m-1)} -> gb{f_m,...,f_1}
70==================================================
71*/
72LList* F5inc(int i, poly f_i, LList* gPrev, ideal gbPrev, poly ONE, LTagList* lTag, RList* rules, RTagList* rTag) {
73    //Print("in f5inc\n");           
74    //pWrite(rules->getFirst()->getRuleTerm());
75    int j;
76    //Print("%p\n",gPrev->getFirst());
77    //pWrite(gPrev->getFirst()->getPoly());
78    poly tempNF =   kNF(gbPrev,currQuotient,f_i);
79    f_i         =   tempNF;
80    //gPrev->insert(ONE,i,f_i);
81    gPrev->insert(ONE,gPrev->getLength()+1,f_i);
82    // tag the first element in gPrev of the current index for findReductor()
83    lTag->setFirstCurrentIdx(gPrev->getLast());
84    //Print("1st gPrev: ");
85    //pWrite(gPrev->getFirst()->getPoly());
86    //Print("2nd gPrev: ");
87    //pWrite(gPrev->getFirst()->getNext()->getPoly());
88    //pWrite(gPrev->getFirst()->getNext()->getPoly());   
89    CList* critPairs        =   new CList();
90    CNode* critPairsMinDeg  =   new CNode();   
91    // computation of critical pairs with checking of criterion 1 and criterion 2 and saving them
92    // in the list critPairs
93    criticalPair(gPrev, critPairs, lTag, rTag, rules);
94    static LList* sPolyList        =   new LList();
95    //sPolyList->print();
96    // labeled polynomials which have passed reduction() and have to be added to list gPrev
97    static LList* completed        =   new LList();
98    // the reduced labeled polynomials which are returned from subalgorithm reduction()
99    static LList* reducedLPolys     =   new LList();
100    // while there are critical pairs to be further checked and deleted/computed
101    while(NULL != critPairs->getFirst()) { 
102        // critPairs->getMinDeg() deletes the first elements of minimal degree from
103        // critPairs, thus the while loop is not infinite.
104        critPairsMinDeg =   critPairs->getMinDeg();
105        // adds all to be reduced S-polynomials in the list sPolyList and adds
106        // the corresponding rules to the list rules
107        // NOTE: inside there is a second check of criterion 2 if new rules are
108        // added
109        computeSPols(critPairsMinDeg,rTag,rules,sPolyList);
110        // DEBUG STUFF FOR SPOLYLIST
111        LNode* temp     =   sPolyList->getFirst();
112        //while(NULL != temp && NULL != temp->getLPoly()) {
113            //Print("Spolylist element: ");
114            //pWrite(temp->getPoly());
115            //temp    =   temp->getNext();
116        //}
117        // reduction process of new S-polynomials and also adds new critical pairs to critPairs
118        reduction(sPolyList, critPairs, gPrev, rules, lTag, rTag, gbPrev);
119       
120        // DEBUG STUFF FOR GPREV
121        //temp    =   gPrev->getFirst();
122        //int number  =   1;
123        //Print("\n\n");
124        //while(NULL != temp) {
125        //    Print("%d.  ",number);
126        //    pWrite(temp->getPoly());
127        //    temp    =   temp->getNext();
128        //    number++;
129        //    Print("\n");
130        //}
131        //sleep(5);
132   
133    }
134    //Print("REDUCTION DONE\n");
135    //Print("%p\n",rules->getFirst());
136    //Print("%p\n",rTag->getFirst());
137    //if(rules->getFirst() != rTag->getFirst()) {
138        //Print("+++++++++++++++++++++++++++++++++++++NEW RULES+++++++++++++++++++++++++++++++++++++\n");
139        //rTag->insert(rules->getFirst());
140    //}
141    //else {
142        //Print("+++++++++++++++++++++++++++++++++++NO NEW RULES++++++++++++++++++++++++++++++++++++\n");
143    //}
144    lTag->insert(lTag->getFirstCurrentIdx());
145    //Print("INDEX: %d\n",tempTag->getIndex());
146    //pWrite(tempTag->getPoly());
147    //Print("COMPLETED FIRST IN F5INC: \n");
148    //Print("1st gPrev: ");
149    //pWrite(gPrev->getFirst()->getPoly());
150    //Print("2nd gPrev: ");
151    //pWrite(gPrev->getFirst()->getNext()->getPoly());
152    //Print("3rd gPrev: ");
153    //pWrite(gPrev->getFirst()->getNext()->getNext()->getPoly());
154    //delete sPolyList;
155    //critPairs->print();
156    delete critPairs;
157    //Print("IN F5INC\n");
158    /*
159    Print("\n\n\nRULES: \n");
160        RNode* tempR    =   rules->getFirst();
161        Print("%p\n",tempR);
162        int t   = 1;
163        while(NULL != tempR) {
164            Print("ADDRESS OF %d RNODE: %p\n",t,tempR);
165            Print("ADDRESS OF RULE: %p\n",tempR->getRule());
166            pWrite(tempR->getRuleTerm());
167            Print("ADDRESS OF TERM: %p\n",tempR->getRuleTerm());
168            Print("%d\n\n",tempR->getRuleIndex());
169            tempR   =   tempR->getNext();
170            t++;
171        }
172    */
173    //gPrev->print();
174    return gPrev;
175}
176
177
178
179/*
180================================================================
181computes a list of critical pairs for the next reduction process
182first element in gPrev is always the newest element which must
183build critical pairs with all other elements in gPrev
184================================================================
185*/
186inline void criticalPair(LList* gPrev, CList* critPairs, LTagList* lTag, RTagList* rTag, RList* rules) {
187    //Print("IN CRITPAIRS\n");
188    // initialization for usage in pLcm()
189    number nOne         =   nInit(1);
190    LNode* newElement   =   gPrev->getLast();
191    LNode* temp         =   gPrev->getFirst();
192    poly u1             =   pOne();
193    poly u2             =   pOne();
194    poly lcm            =   pOne();
195    poly t              =   pHead(newElement->getPoly());
196    Rule* testedRule    =   NULL;
197    if(NULL != rules->getFirst()) {
198        testedRule  =   rules->getFirst()->getRule();
199    }
200    // computation of critical pairs
201    //critPairs->print();
202    while( gPrev->getLast() != temp) {
203        //pWrite( *(gPrev->getFirst()->getPoly()) );
204       // pWrite( *(l->getPoly()) );
205       //pWrite(newElement->getPoly());
206        //pWrite(temp->getPoly());
207        pLcm(newElement->getPoly(), temp->getPoly(), lcm);
208        pSetCoeff(lcm,nOne);
209        // computing factors u2 for new labels
210        //pExpVectorDiff(u1,lcm,t);
211        //Print("U1: ");
212        //pWrite(u1);
213        u1 = pDivide(lcm,t);
214        //pWrite(u1);
215        pSetCoeff(u1,nOne);
216        //Print("%p\n",u1);
217        //critPairs->print();
218        //pExpVectorDiff(u2,lcm,temp->getPoly());
219        //Print("U2: ");
220        //pWrite(u2);
221        u2 = pDivide(lcm, pHead(temp->getPoly()));
222        //pWrite(u2);
223        //Print("%p\n",u2);
224        pSetCoeff(u2,nOne);
225        //if(gPrev->getLast()->getIndex()==5) {
226            //Print("IN CRITPAIRS\n");
227        //    pWrite(u1);
228        //    Print("1st ELEMENT: ");
229        //    pWrite(newElement->getPoly());
230        //    Print("2nd ELEMENT: ");
231        //    pWrite(temp->getPoly());
232        //}
233        // testing both new labels by the F5 Criterion
234        //critPairs->print();
235        if(!criterion2(gPrev->getFirst()->getIndex(), u2, temp, rules, rTag)
236           && !criterion2(gPrev->getFirst()->getIndex(), u1, newElement, rules, rTag) 
237           && !criterion1(gPrev,u1,newElement,lTag) && !criterion1(gPrev,u2,temp,lTag)) {
238            // if they pass the test, add them to CList critPairs, having the LPoly with greater
239            // label as first element in the CPair
240            if(newElement->getIndex() == temp->getIndex() && 
241            -1 == pLmCmp(ppMult_qq(u1, newElement->getTerm()),ppMult_qq(u2, temp->getTerm()))) {
242                //Print("zweites groesser\n");
243               
244                CPair* cp   =   new CPair(pDeg(ppMult_qq(u2,pHead(temp->getPoly()))), u2, 
245                                temp->getLPoly(), u1, newElement->getLPoly(), testedRule);                   
246               
247                critPairs->insert(cp);
248                //Print("LALA %p\n",critPairs->getFirst());
249                //sleep(5);
250            }
251            else {
252                CPair* cp   =   new CPair(pDeg(ppMult_qq(u2,pHead(temp->getPoly()))), u1, 
253                                newElement->getLPoly(), u2, temp->getLPoly(), testedRule);                   
254                //Print("erstes groesser\n");
255                critPairs->insert(cp);
256                //Print("LALA %p\n",critPairs->getFirst());
257                //sleep(5);
258            }
259        }
260        else {
261        }
262       
263        //Print("\n\n");
264        temp    =   temp->getNext();
265    }
266    // for debugging
267    //if(NULL != critPairs) {
268        //critPairs->print();
269        //sleep(5);
270    //}
271    //Print("END CRITPAIRS\n");
272}
273
274
275
276
277
278
279
280/*
281========================================
282Criterion 1, i.e. Faugere's F5 Criterion
283========================================
284*/
285inline bool criterion1(LList* gPrev, poly t, LNode* l, LTagList* lTag) {
286    // starts at the first element in gPrev with index = (index of l)-1, these tags are saved in lTag
287        int idx =   l->getIndex();
288    int i;
289    if(idx == 1) {
290        //Print("HIER\n");
291        return false;
292    }
293    else {
294        LNode* testNode =   gPrev->getFirst();
295        // save the monom t1*label_term(l) as it is tested various times in the following
296        poly u1 = ppMult_qq(t,l->getTerm());
297        //Print("------------------------------IN CRITERION 1-----------------------------------------\n");
298        //Print("TESTED ELEMENT: ");
299        //pWrite(l->getPoly());
300        //pWrite(l->getTerm());
301        //pWrite(ppMult_qq(t,l->getTerm()));
302        //Print("%d\n\n",l->getIndex());
303        while(testNode->getIndex() < idx) { // && NULL != testNode->getLPoly()) {
304            //pWrite(testNode->getPoly());
305            //Print("%d\n",testNode->getIndex());
306            if(pLmDivisibleByNoComp(testNode->getPoly(),u1)) {
307                //Print("Criterion 1 NOT passed!\n");
308                return true;
309            }
310            //pWrite(testNode->getNext()->getPoly());
311            testNode    =   testNode->getNext();
312        }
313        /*
314        ideal testId    =   idInit(idx-1,1);
315        for(i=0;i<idx-1;i++) {
316            testId->m[i]  =   pHead(testNode->getPoly());
317            testNode        =   testNode->getNext();
318        }
319        poly temp   =   kNF(testId,currQuotient,u1);
320        //pWrite(temp);
321        for(i=0;i<IDELEMS(testId);i++) {
322            testId->m[i]    =   NULL;
323        }
324        idDelete(&testId);
325        if(NULL == temp) {
326            //if(l->getIndex() != gPrev->getFirst()->getIndex()) {
327            //    Print("----------------------------Criterion1 not passed----------------------------------\n");
328            //}
329            return true;
330        }
331        */
332        return false;
333    }
334}
335
336
337
338/*
339=====================================
340Criterion 2, i.e. Rewritten Criterion
341=====================================
342*/
343inline bool criterion2(int idx, poly t, LNode* l, RList* rules, RTagList* rTag) {
344    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n");
345    /* 
346    Print("RULES: \n");
347        RNode* tempR    =   rules->getFirst();
348        Print("%p\n",tempR);
349        int i   = 1;
350        while(NULL != tempR) {
351            Print("ADDRESS OF %d RNODE: %p\n",i,tempR);
352            Print("ADDRESS OF RULE: %p\n",tempR->getRule());
353            pWrite(tempR->getRuleTerm());
354            Print("ADDRESS OF TERM: %p\n",tempR->getRuleTerm());
355            Print("%d\n\n",tempR->getRuleIndex());
356            tempR   =   tempR->getNext();
357            i++;
358        }
359        //Print("TESTED ELEMENT: ");
360        //pWrite(l->getPoly());
361        //pWrite(l->getTerm());
362        //pWrite(ppMult_qq(t,l->getTerm()));
363        //Print("%d\n\n",l->getIndex());
364      */
365// start at the previously added element to gPrev, as all other elements will have the same index for sure
366    if(idx > l->getIndex()) {
367        return false;
368    }
369   
370    RNode* testNode; // =   new RNode();
371   
372
373    if(NULL == rTag->getFirst()) {
374        if(NULL != rules->getFirst()) {
375            testNode    =   rules->getFirst();
376        }
377        else {
378            return false;
379        }
380    }
381    else {
382
383        if(l->getIndex() > rTag->getFirst()->getRuleIndex()) {
384            testNode    =   rules->getFirst();
385        }
386        else {
387       //Print("HIER\n");
388            //Print("DEBUG\n");
389        //Print("L INDEX: %d\n",l->getIndex());
390            /*-------------------------------------
391             * TODO: WHEN INTERREDUCING THE GB THE
392             *       INDEX OF THE PREVIOUS ELEMENTS
393             *       GETS HIGHER!
394             *-----------------------------------*/
395            //testNode    =   rules->getFirst();
396            testNode    =   rTag->get(l->getIndex());
397            if(NULL == testNode) {
398                testNode    =   rules->getFirst();
399            }
400            //Print("TESTNODE ADDRESS: %p\n",testNode);
401        }
402    }
403    //testNode    =   rules->getFirst();
404        // save the monom t1*label_term(l) as it is tested various times in the following
405    poly u1 = ppMult_qq(t,l->getTerm());
406    // first element added to rTag was NULL, check for this
407    //Print("%p\n",testNode->getRule());
408    // NOTE: testNode is possibly NULL as rTag->get() returns NULL for elements of index <=1!
409    //Print("TESTNODE: %p\n",testNode);
410    //pWrite(testNode->getRuleTerm());
411    if(NULL != testNode ) {   
412        //pWrite(testNode->getRuleTerm());
413    }
414    if(NULL != testNode) {
415        if(testNode->getRule() == l->getRule()) {
416            //Print("%p\n%p\n",testNode->getRule(),l->getRule());
417            //Print("EQUAL\n");
418        }
419        else {
420            //Print("NOT EQUAL\n");
421        }
422    }
423    while(NULL != testNode && testNode->getRule() != l->getRule() 
424          && l->getIndex() == testNode->getRuleIndex()) {
425        //Print("%p\n",testNode);
426        //pWrite(testNode->getRuleTerm());
427        //pWrite(t);
428        //pWrite(l->getTerm());
429        //pWrite(u1);
430        //Print("%d\n",testNode->getRuleIndex());
431        if(pLmDivisibleByNoComp(testNode->getRuleTerm(),u1)) {
432            //Print("-----------------Criterion 2 NOT passed!-----------------------------------\n");
433            //Print("INDEX: %d\n",l->getIndex());
434            pDelete(&u1);
435    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n\n");
436            return true;
437        }
438                testNode    =   testNode->getNext();
439    }
440    //delete testNode;
441    pDelete(&u1);
442    //Print("------------------------------IN CRITERION 2/1-----------------------------------------\n\n");
443    return false;
444}
445
446
447
448/*
449=================================================================================================================
450Criterion 2, i.e. Rewritten Criterion, for its second call in computeSPols(), with added lastRuleTested parameter
451=================================================================================================================
452*/
453inline bool criterion2(poly t, LPoly* l, RList* rules, Rule* testedRule) {
454    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n");
455    //Print("LAST RULE TESTED: %p",testedRule);
456    /*
457    Print("RULES: \n");
458        RNode* tempR    =   rules->getFirst();
459        while(NULL != tempR) {
460            pWrite(tempR->getRuleTerm());
461            Print("%d\n\n",tempR->getRuleIndex());
462            tempR   =   tempR->getNext();
463        }
464        //Print("TESTED ELEMENT: ");
465        //pWrite(l->getPoly());
466        //pWrite(l->getTerm());
467        //pWrite(ppMult_qq(t,l->getTerm()));
468        //Print("%d\n\n",l->getIndex());
469    */
470// start at the previously added element to gPrev, as all other elements will have the same index for sure
471        RNode* testNode =   rules->getFirst();
472    // save the monom t1*label_term(l) as it is tested various times in the following
473    poly u1 = ppMult_qq(t,l->getTerm());
474        // first element added to rTag was NULL, check for this
475        while(NULL != testNode && testNode->getRule() != testedRule) {
476        //pWrite(testNode->getRuleTerm());
477        if(pLmDivisibleByNoComp(testNode->getRuleTerm(),u1)) {
478            //Print("--------------------------Criterion 2 NOT passed!------------------------------\n");
479            //Print("INDEX: %d\n",l->getIndex());
480            pDelete(&u1);
481    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n\n");
482            return true;
483        }
484                testNode    =   testNode->getNext();
485    }
486    pDelete(&u1);
487    //Print("------------------------------IN CRITERION 2/2-----------------------------------------\n\n");
488    return false;
489}
490
491
492
493/*
494==================================
495Computation of S-Polynomials in F5
496==================================
497*/
498void computeSPols(CNode* first, RTagList* rTag, RList* rules, LList* sPolyList) { 
499    CNode* temp         =   first;
500    poly sp      =   pInit();
501    //Print("###############################IN SPOLS##############################\n");
502    //first->print();
503
504    while(NULL != temp) {
505        //Print("JA\n");
506        // only if a new rule was added since the last test in subalgorithm criticalPair()
507        //if(rules->getFirst() != rTag->getFirst()) {
508        if(!criterion2(temp->getT1(),temp->getAdLp1(),rules,temp->getTestedRule())) {
509                // the second component is tested only when it has the actual index, otherwise there is
510                // no new rule to test since the last test in subalgorithm criticalPair()
511                if(temp->getLp2Index() == temp->getLp1Index()) {
512                    if(!criterion2(temp->getT2(),temp->getAdLp2(),rules,temp->getTestedRule())) {
513                        // computation of S-polynomial
514                        sp      =   pSub(ppMult_qq(temp->getT1(),temp->getLp1Poly()),
515                                         ppMult_qq(temp->getT2(),temp->getLp2Poly()));
516                        //Print("BEGIN SPOLY1\n====================\n");
517                        pNorm(sp);
518                        //pWrite(sp);
519                        //Print("END SPOLY1\n====================\n");
520                        if(NULL == sp) {
521                            // as rules consist only of pointers we need to save the labeled
522                            // S-polynomial also of a zero S-polynomial
523                            //zeroList->insert(temp->getAdT1(),temp->getLp1Index(),&sp);
524                            // origin of rule can be set NULL as the labeled polynomial
525                            // will never be used again as it is zero => no problems with
526                            // further criterion2() tests and termination conditions
527                            //Print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ZERO REDUCTION~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n");
528                                                        reductionsToZero++;
529                        //Print("IN SPOLS 1\n");
530                            rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
531                            //Print("RULE ADDED: \n");
532                            //pWrite(rules->getFirst()->getRuleTerm());
533                            // as sp = NULL, delete it
534                            pDelete(&sp);
535                            //Print("HIER\n");
536                        }
537                        else {
538                        //Print("IN SPOLS 2\n");
539                            rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
540                            //Print("RULE ADDED: \n");
541                            //pWrite(rules->getFirst()->getRuleTerm()); 
542                            sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rules->getFirst()->getRule());
543                        }
544                        // data is saved in sPolyList or zero => delete sp
545                    }
546                }
547                else { // temp->getLp2Index() < temp->getLp1Index()
548                    // computation of S-polynomial
549                    sp      =   pSub(ppMult_qq(temp->getT1(),temp->getLp1Poly()),
550                                     ppMult_qq(temp->getT2(),temp->getLp2Poly()));
551                    //Print("BEGIN SPOLY2\n====================\n");
552                    pNorm(sp);
553                    //pWrite(sp);
554                    //Print("END SPOLY2\n====================\n");
555                    if(NULL == sp) {
556                        // zeroList->insert(temp->getAdT1(),temp->getLp1Index(),&sp);
557                        // origin of rule can be set NULL as the labeled polynomial
558                        // will never be used again as it is zero => no problems with
559                        // further criterion2() tests and termination conditions
560                            //Print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ZERO REDUCTION~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n");
561                        reductionsToZero++;
562                        //Print("IN SPOLS 3\n");
563                        rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
564                        //Print("RULE ADDED: \n");
565                        //pWrite(rules->getFirst()->getRuleTerm());
566                        // as sp = NULL, delete it
567                        pDelete(&sp);
568                    }
569                    else {
570                        //Print("IN SPOLS 4\n");
571                        rules->insert(temp->getLp1Index(),ppMult_qq(temp->getT1(),temp->getLp1Term()));
572                        //Print("RULE ADDED: \n");
573                        //pWrite(rules->getFirst()->getRuleTerm());
574                        //Print("%d\n",rules->getFirst()->getRuleIndex());
575                        //Print("%p\n",sPolyList->getFirst());
576                        sPolyList->insertByLabel(ppMult_qq(temp->getT1(),temp->getLp1Term()),temp->getLp1Index(),sp,rules->getFirst()->getRule());
577                    }
578                }
579            }
580        //}
581        //Print("%p\n",temp);
582        temp    =   temp->getNext();
583        //Print("%p\n",temp);
584        //Print("%p\n",temp->getData());
585        //pWrite(temp->getLp1Poly());
586    }
587    // these critical pairs can be deleted now as they are either useless for further computations or
588    // already saved as an S-polynomial to be reduced in the following
589    delete first;   
590}
591
592
593
594/*
595========================================================================
596reduction including subalgorithm topReduction() using Faugere's criteria
597========================================================================
598*/
599void reduction(LList* sPolyList, CList* critPairs, LList* gPrev, RList* rules, LTagList* lTag, RTagList* rTag, 
600                 ideal gbPrev) { 
601    //Print("##########################################In REDUCTION!########################################\n");
602    // check if sPolyList has any elements
603    // NOTE: due to initialization sPolyList always has a default NULL element
604    while(sPolyList->getLength() > 0) {
605        //Print("SPOLYLIST LENGTH: %d\n",sPolyList->getLength());
606        if(sPolyList->getLength() > 1) {
607        //Print("%p\n",sPolyList->getFirst());
608        //Print("%p\n",sPolyList->getFirst()->getLPoly());
609        //Print("%p\n",sPolyList->getFirst()->getNext());
610        //Print("%p\n",sPolyList->getFirst()->getNext()->getLPoly());
611        //Print("%p\n",sPolyList->getFirst());
612        }
613        //if(gPrev->getLast()->getIndex() == 5) {
614            //sPolyList->print();
615            //sleep(5);
616        //}
617       
618        // temp is the first element in the sPolyList which should be reduced
619        // due to earlier sorting this is the element of minimal degree AND
620        // minimal label
621        LNode* temp =   sPolyList->getFirst();
622        //pWrite(temp->getPoly());
623        // delete the above first element from sPolyList, temp will be either reduced to
624        // zero or added to gPrev, but never come back to sPolyList
625        if(NULL != temp) { 
626            sPolyList->setFirst(temp->getNext());
627        //Print("HIER\n");
628        }
629        else {
630            break;
631        }
632        //Print("HALLO %p\n",temp->getPoly());
633        //Print("%p\n",temp->getPoly());
634        //pWrite(temp->getPoly());
635        //idShow(gbPrev);
636        poly tempNF = kNF(gbPrev,currQuotient,temp->getPoly());
637        //Print("LENGTH: %d\n",sPolyList->getLength());
638        //pWrite(tempNF);
639        //pWrite(temp->getPoly());
640        if(NULL != tempNF) {
641            pNorm(tempNF);
642            // write the reduced polynomial in temp
643            temp->setPoly(tempNF);
644            // try further reductions of temp with polynomials in gPrev
645            // with label index = current label index: this is done such that there
646            // is no label corruption during the reduction process
647            topReduction(temp,sPolyList,gPrev,rules,lTag,rTag,gbPrev);
648       
649        }
650        if(NULL != temp->getPoly()) {
651            //CList* newCritPairs = new CList;
652            //Print("##################IN CRITPAIRS IN REDUCTION#####################\n");
653            criticalPair(gPrev,critPairs,lTag,rTag,rules);
654            //Print("+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++H I E R++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n");
655        }
656        else {
657            //delete temp;
658            LNode* tempLoop = gPrev->getFirst();
659            //Print("AUSGABE IN REDUCTION:\n");       
660            //while(NULL != tempLoop) {
661                //pWrite(tempLoop->getPoly());
662                //tempLoop = tempLoop->getNext();
663            //}
664            //sleep(10);
665        }
666       
667    }
668}   
669
670
671
672/*
673=====================================================================================
674top reduction in F5, i.e. reduction of a given S-polynomial by labeled polynomials of
675the same index whereas the labels are taken into account
676=====================================================================================
677*/
678void topReduction(LNode* l, LList* sPolyList, LList* gPrev, RList* rules, LTagList* lTag, RTagList* rTag, ideal gbPrev) {
679    //Print("##########################################In topREDUCTION!########################################\n");
680    // try l as long as there are reductors found by findReductor()
681    do {
682        LNode* gPrevRedCheck    =   new LNode(lTag->getFirstCurrentIdx());
683        //Print("gPrevCheckPOLY: ");
684        //pWrite(gPrevRedCheck->getPoly());
685        LNode* tempRed          =   new LNode();
686        //Print("TESTED POLYNOMIAL IN THE FOLLOWING: ");
687        //pWrite(l->getPoly());
688        //Print("HIER\n");
689        tempRed  =   findReductor(l,gPrevRedCheck,gPrev,rules,lTag,rTag);
690        //Print("TEMPRED NEXT: %p\n",tempRed->getNext());
691        //Print("--------------------------------HIER DEBUG 2----------------------------------\n");
692        // if a reductor for l is found and saved in tempRed
693        if(NULL != tempRed) {
694            // if label of reductor is greater than the label of l we have to built a new element
695            // and add it to sPolyList
696            //Print("REDUCTOR POLYNOMIAL: ");
697            //pWrite(tempRed->getPoly());
698            //Print("TEMPRED: %p\n",tempRed);
699            //Print("TERM: ");
700            //pWrite(tempRed->getTerm());
701            if(pLmCmp(tempRed->getTerm(),l->getTerm()) == 1) {
702                // needed sinc pSub destroys the arguments!
703                //Print("H----------I------------E--------------R\n");
704                poly temp_poly_l    =   pInit();
705                temp_poly_l         =   pCopy(l->getPoly());
706                //Print("POLYNOMIAL L: ");
707                //pWrite(l->getPoly());
708                //pWrite(temp_poly_l);
709                poly temp           =   pSub(tempRed->getPoly(),temp_poly_l);
710                //Print("POLYNOMIAL L: ");
711                //pWrite(l->getPoly());
712                //pWrite(temp_poly_l);
713                //Print("AFTER REDUCTION STEP: ");
714                //pWrite(temp);
715                //sleep(20);
716                //pWrite(temp);
717                if(NULL != temp) {
718                    pNorm(temp);
719                    tempRed->setPoly(temp);
720                    //pWrite(tempRed->getPoly());
721                    // for debugging
722                    //pWrite(tempRed->getPoly());
723                    //Print("RULE ADDED\n");
724                    rules->insert(tempRed->getIndex(),tempRed->getTerm());
725                    tempRed->getLPoly()->setRule(rules->getFirst()->getRule());
726                    //Print("%p\n",sPolyList->getFirst());
727                    //Print("%p\n",sPolyList->getFirst()->getLPoly());
728                    //Print("SPOLYLIST LENGTH: %d\n",sPolyList->getLength());
729                    //sPolyList->print();
730                   
731                    sPolyList->insertByLabel(tempRed);
732                    //sPolyList->print();
733                    //Print("SPOLYLIST LENGTH: %d\n",sPolyList->getLength());
734                    //Print("OH JE\n");
735                }
736                else {
737                    pDelete(&temp);
738                    reductionsToZero++;
739                    //Print("RULE ADDED\n");
740                    //rules->insert(tempRed->getIndex(),tempRed->getTerm());
741                    //pWrite(rules->getFirst()->getRuleTerm());
742                    //Print("DELETE TEMPRED\n");
743                    delete tempRed;
744                }
745            }
746            // label of reductor is smaller than the label of l, subtract reductor from l and delete the
747            // gPrevRedCheck pointer added to l during findReductor() as the head term of l changes
748            // after subtraction
749            else {
750                poly temp_poly_l    =   pInit();
751                temp_poly_l         =   pCopy(l->getPoly());
752                poly temp   =   pSub(temp_poly_l,tempRed->getPoly());
753                //Print("AFTER REDUCTION STEP: ");
754                //pWrite(temp);
755                if(NULL != temp) {
756                    pNorm(temp);
757                    poly tempNF =   kNF(gbPrev,currQuotient,temp); 
758                    pNorm(tempNF);
759                    //pWrite(tempNF);
760                    if(NULL == tempNF) {
761                        reductionsToZero++;
762                        pDelete(&tempNF);
763                        l->setPoly(NULL);
764                        break;
765                    }
766                    l->setPoly(tempNF);
767                   
768                    //pWrite(l->getPoly());
769                    gPrevRedCheck   =   lTag->getFirstCurrentIdx();
770                }
771                else {
772                    //Print("ZERO REDUCTION!\n");
773                    reductionsToZero++;
774                    pDelete(&temp);
775                    l->setPoly(NULL);
776                    //pWrite(gPrev->getLast()->getPoly());
777                    break;
778                }
779            }   
780        }
781        else {
782            if(NULL != l->getPoly()) {
783                pNorm(l->getPoly());
784                //Print("----------------------------------ADDED TO GPREV IN TOPREDUCTION:-------------------------------------- ");
785                //pWrite(l->getPoly());
786                //pWrite(l->getTerm());
787                //Print("INDEX: %d\n\n\n", l->getIndex());
788                //sleep(20);
789                gPrev->insert(l->getLPoly());
790                //Print("GPREV: \n");
791                LNode* tempLoop = gPrev->getFirst();
792                //while(NULL != tempLoop) {
793                    //Print("HERE\n");
794                    //pWrite(tempLoop->getPoly());
795                    //tempLoop = tempLoop->getNext();
796                //}
797            }
798            break;
799        }
800    } while(1);
801}
802
803
804/*
805=====================================================================
806subalgorithm to find a possible reductor for the labeled polynomial l
807=====================================================================
808*/
809LNode* findReductor(LNode* l, LNode* gPrevRedCheck, LList* gPrev, RList* rules, LTagList* lTag,RTagList* rTag) {
810    // allociation of memory for the possible reductor
811    //Print("IN FIND REDUCTOR\n");
812    poly u      =   pOne();
813    poly red    =   pOne();
814    number nOne =   nInit(1);
815    LNode* temp =   new LNode();
816    // head term of actual element such that we do not have to call pHead at each new reductor test
817    poly t      =   pHead(l->getPoly());
818    // if l was already checked use the information in gPrevRedCheck such
819    // that we can start searching for new reducers from this point and
820    // not from the first element of gPrev with the current index
821    temp    =   gPrevRedCheck;
822    //temp    =   lTag->getFirstCurrentIdx(); 
823    //Print("GPREVREDCHECK: %p\n",gPrevRedCheck);
824    //pWrite(gPrevRedCheck->getPoly());
825    // search for reductors until we are at the end of gPrev resp. at the
826    // end of the elements of the current index
827    while(NULL != temp && temp->getIndex() == l->getIndex()) {
828        //pWrite(temp->getPoly());
829        //Print("INDEX: %d\n",temp->getIndex());
830        // does the head of the element of gPrev divides the head of
831        // the to be reduced element?
832        //Print("-------------FOUND REDUCTORS----------------------\n");
833        //Print("\n");
834        //pWrite(temp->getPoly());
835        //pWrite(temp->getTerm());
836        //pWrite(t);
837        //Print("HALLO\n");
838        if(pLmDivisibleByNoComp(pHead(temp->getPoly()),t)) {
839        //Print("HALLO\n");
840            // get all the information needed for the following tests
841            // of the criteria
842            u   =   pDivide(t,pHead(temp->getPoly()));
843            pSetCoeff(u,nOne);
844            //Print("HIER FINDRED\n");
845            //pWrite(u);
846            //Print("\n");
847            red =   ppMult_qq(u,temp->getPoly());
848            pNorm(red);
849            //u   =   ppMult_qq(u,temp->getTerm());
850            //pSetCoeff(u,nOne);
851            // check if both have the same label
852        //Print("HALLO\n");
853            if(pLmCmp(u,l->getTerm()) != 0) {
854        //Print("HALLO\n");
855                // passing criterion2 ?
856                if(!criterion2(gPrev->getFirst()->getIndex(), u,temp,rules,rTag)) {
857                    // passing criterion1 ?
858                    if(!criterion1(gPrev,u,temp,lTag)) {
859                            //Print("HIER DEBUG\n");
860                            gPrevRedCheck   =   temp->getNext();
861                            LNode* redNode  =   new LNode(ppMult_qq(u,temp->getTerm()),temp->getIndex(),red,NULL,NULL);
862                            return redNode;
863                    }
864                }
865            }
866        }
867        //Print("%p\n",temp->getNext());
868        //pWrite(temp->getPoly());
869        //Print("HALLO\n");
870        temp = temp->getNext();
871    }
872   
873//    delete temp;
874   //Print("1st gPrev: ");
875    //pWrite(gPrev->getFirst()->getPoly());
876    //Print("2nd gPrev: ");
877    //pWrite(gPrev->getFirst()->getNext()->getPoly());
878 return NULL;
879}
880
881
882
883/*
884==========================================================================
885MAIN:computes a gb of the ideal i in the ring r with our F5 implementation
886==========================================================================
887*/
888ideal F5main(ideal id, ring r) {
889    int timer   =   initTimer();
890    startTimer();
891    int i,j,k,l;
892    int gbLength;
893    // 1 polynomial for defining initial labels & further tests
894    poly ONE = pOne();
895    poly pOne   =   pOne();
896    number nOne =   nInit(1);
897    // tag the first element of index i-1 for criterion 1
898    LTagList* lTag  =   new LTagList();
899    //Print("LTAG BEGINNING: %p\n",lTag);
900   
901    // first element in rTag is first element of rules which is NULL RNode,
902    // this must be done due to possible later improvements
903    RList* rules    =   new RList();
904    //Print("RULES FIRST: %p\n",rules->getFirst());
905    //Print("RULES FIRST DATA: %p\n",rules->getFirst()->getRule());
906    RTagList* rTag  =   new RTagList(rules->getFirst());
907    i = 1;
908    for(j=0; j<IDELEMS(id); j++) {
909        if(NULL != id->m[j]) { 
910            if(pComparePolys(id->m[j],ONE)) {
911                Print("One Polynomial in Input => Computations stopped\n");
912                ideal idNew = idInit(1,1);
913                idNew->m[0] = ONE;
914                return(idNew);
915            }   
916        }
917    } 
918    ideal idNew     =   kInterRed(id); 
919    id              =   idNew;
920    //qsortDegree(&id->m[0],&id->m[IDELEMS(id)-1]);
921    idShow(id);
922    LList* gPrev    =   new LList(ONE, i, id->m[0]);
923    //idShow(id);
924    //Print("%p\n",id->m[0]);
925    //pWrite(id->m[0]);
926    //Print("%p\n",gPrev->getFirst()->getPoly());
927    //pWrite(gPrev->getFirst()->getPoly());
928
929    lTag->insert(gPrev->getFirst());
930    lTag->setFirstCurrentIdx(gPrev->getFirst());
931    // computing the groebner basis of the elements of index < actual index
932    gbLength    =   gPrev->getLength();
933    //Print("Laenge der bisherigen Groebner Basis: %d\n",gPrev->getLength());
934    ideal gbPrev    =   idInit(gbLength,1);
935    // initializing the groebner basis of elements of index < actual index
936    gbPrev->m[0]    =   gPrev->getFirst()->getPoly();
937    //idShow(gbPrev);
938    //idShow(currQuotient);
939
940    for(i=2; i<=IDELEMS(id); i++) {
941        LNode* gPrevTag =   gPrev->getLast();
942        //Print("Last POlynomial in GPREV: ");
943        //Print("%p\n",gPrevTag);   
944        //pWrite(gPrevTag->getPoly());
945        gPrev   =   F5inc(i, id->m[i-1], gPrev, gbPrev, ONE, lTag, rules, rTag);
946        //Print("____________________________________ITERATION STEP DONE________________________________________\n");
947       
948        // DEBUGGING STUFF
949        LNode* temp    =   gPrev->getFirst();
950    // computing new groebner basis gbPrev
951        if(gPrev->getLength() > gbLength) {
952            ideal gbAdd =   idInit(gPrev->getLength()-gbLength,1);
953            LNode*  temp =   gPrevTag;
954        //Print("%p\n",gPrevTag);   
955        //Print("%p\n",gPrev->getLast());   
956        //pWrite(temp->getPoly());
957        //Print("LENGTH OF GPREV LIST: %d\n",gPrev->getLength());
958        //Print("%d\n",gbLength);
959        //Print("%d\n",gPrev->getLength()-gbLength-1);
960            for(j=0;j<=gPrev->getLength()-gbLength-1;j++) {
961                //Print("YES\n");
962                temp        =   temp->getNext();
963                if(temp) {
964                    //Print("%p\n",temp);
965                    //pWrite(temp->getPoly());
966                    //Print("%p\n",temp->getNext());
967                }
968                gbAdd->m[j] =   temp->getPoly();
969                //pWrite(temp->getPoly());
970            }
971            //Print("HIER AUCH\n");
972            gbPrev          =   idAdd(gbPrev,gbAdd);
973           
974           
975            // interreduction stuff
976            if(i<IDELEMS(id)) {
977                ideal tempId    =   kInterRed(gbPrev);
978                //idShow(tempId);
979                gbPrev          =   tempId;
980                //qsortDegree(&gbPrev->m[0],&gbPrev->m[IDELEMS(gbPrev)-1]);
981                delete gPrev;
982                //sleep(5);
983                //Print("RULES FIRST NOW1: %p\n",rules->getFirst());
984                //Print("HIER\n");
985                delete rules;
986                //delete rTag;
987                //Print("HIER AUCH\n");
988                //Print("%p\n",rules->getFirst());
989                gPrev    =   new LList(pOne,1,gbPrev->m[0]);
990                gPrev->insert(pOne,1,gbPrev->m[1]);
991                //poly tempPoly = pInit();
992                //pLcm(pHead(gbPrev->m[0]),pHead(gbPrev->m[1]),tempPoly);
993                //tempPoly    =   pDivide(tempPoly,pOne());
994                //pSetCoeff(tempPoly,nOne);
995                rules    =   new RList();
996                rTag     =   new RTagList(rules->getFirst());
997               
998                //Print("%p\n",rules->getFirst());
999                //pWrite(tempPoly);
1000                //rules->insert(2,tempPoly);
1001                //rTag->insert(rules->getFirst());
1002                //Print("%p\n",rules->getFirst());
1003                //Print("%p\n",rTag->getFirst());
1004                //Print("%p\n",rules->getFirst());
1005                //Print("%p\n",rules->getFirst()->getNext()->getNext());
1006                //Print("HIERLALA\n");
1007            //pWrite(rules->getFirst()->getRuleTerm());
1008           // Print("RULES FIRST NOW2: %p\n",rules->getFirst());
1009                for(k=2; k<IDELEMS(gbPrev); k++) {
1010                    gPrev->insert(pOne,k+1,gbPrev->m[k]);
1011                    for(l=0; l<k; l++) {
1012                        //pWrite(gbPrev->m[k]);
1013                        //pWrite(gbPrev->m[l]);
1014                        //poly tempPoly2  =   pOne();
1015                        //pLcm(pHead(gbPrev->m[k]),pHead(gbPrev->m[l]),tempPoly2);
1016                        //tempPoly2   =   pDivide(tempPoly2,pOne());
1017                        //pSetCoeff(tempPoly2,nOne);
1018                        //pWrite(tempPoly2);
1019                        //rules->insert(k+1,tempPoly2);
1020                    }
1021                    //rTag->insert(rules->getFirst());
1022                }
1023            }
1024           
1025            gbLength    =   gPrev->getLength(); 
1026            //Print("HIER\n");
1027        }
1028        //gPrev->print();
1029        //int anzahl  =   1;
1030        //while(NULL != temp) {
1031        //    Print("%d. Element: ",anzahl);
1032        //    pWrite(temp->getPoly());
1033        //    Print("\n");
1034        //    temp    =   temp->getNext();
1035        //    anzahl++;
1036        //}
1037        //sleep(5);
1038        //Print("GROEBNER BASIS:\n====================================================\n");
1039        //idShow(gbPrev);
1040        //Print("===================================================\n");
1041        //Print("JA\n");
1042    } 
1043        //idShow(gbPrev);
1044    Print("\n\nNumber of zero-reductions:  %d\n",reductionsToZero);
1045    timer   =   getTimer();
1046    Print("Time for computations: %d/1000 seconds\n",timer);
1047    //LNode* temp    =   gPrev->getFirst();
1048    //while(NULL != temp) {
1049    //    pWrite(temp->getPoly());
1050    //    temp    =   temp->getNext();
1051    // }
1052    //gPrev->print();
1053    //delete lTag;
1054    //delete rTag;
1055    //delete gPrev;
1056    return(gbPrev);
1057
1058
1059}
1060
1061#endif
Note: See TracBrowser for help on using the repository browser.