source: git/kernel/f5gb.cc @ bc103b

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