Context Navigation

source: git/kernel/f5gb.cc @ 66a5f8

spielwiese

Last change on this file since 66a5f8 was 66a5f8, checked in by Christian Eder, 15 years ago
trying some stuff with exponent vectors, still not working with most of the example ideals from f5ex2.lib git-svn-id: file:///usr/local/Singular/svn/trunk@11591 2c84dea3-7e68-4137-9b89-c4e89433aadc
Property mode set to `100644`
File size: 42.8 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats: