source: git/kernel/f5gb.cc @ febf7e

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