source: git/kernel/f5gb.cc @ 1534d9

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