source: git/kernel/f5gb.cc @ 24ac666

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