Index: kernel/gfan.cc
===================================================================
--- kernel/gfan.cc	(revision f32177e9dea6c9833a7f2f8986a2cf6273f30fab)
+++ kernel/gfan.cc	(revision 0ad81f2c92a459da3968f1add8868540fe570fb8)
@@ -2,7 +2,7 @@
 Compute the Groebner fan of an ideal
 $Author: monerjan $
-$Date: 2009-04-06 14:57:18 $
-$Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.29 2009-04-06 14:57:18 monerjan Exp $
-$Id: gfan.cc,v 1.29 2009-04-06 14:57:18 monerjan Exp $
+$Date: 2009-04-07 09:44:20 $
+$Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.30 2009-04-07 09:44:20 monerjan Exp $
+$Id: gfan.cc,v 1.30 2009-04-07 09:44:20 monerjan Exp $
 */
 
@@ -267,4 +267,5 @@
 			cout << "Cols = " << ddcols << endl;
 #endif
+			
 			/*Write the normals into class facet*/
 #ifdef gfan_DEBUG
@@ -275,5 +276,5 @@
 			this->facetPtr = fRoot;		//set variable facetPtr of class gcone to first facet
 			facet *fAct; 			//instantiate pointer to active facet
-			fAct = fRoot;		//This does not seem to do the trick. fRoot and fAct have to point to the same adress!
+			fAct = fRoot;			//Seems to do the trick. fRoot and fAct have to point to the same adress!
 			std::cout << "fRoot = " << fRoot << ", fAct = " << fAct << endl;
 			for (int kk = 0; kk<ddrows; kk++)
@@ -291,13 +292,38 @@
 					if (jj<ddcols)				//Is this facet NOT the last facet? Writing while instead of if is a really bad idea :)
 					{
-						fAct->next = new facet();	//If so: instantiate new facet. Otherwise this->next=NULL due to the constructor						
+						//fAct->next = new facet();	//If so: instantiate new facet. Otherwise this->next=NULL due to the constructor						
 					}
-				}//for jj<ddcols
-				/*Now load should be full and we can call setFacetNormal*/
-				fAct->setFacetNormal(load);
-				//fAct->printNormal();
-				fAct=fAct->next;	//this should definitely not be called in the above while-loop :D
+				}//for (int jj = 1; jj <ddcols; jj++)
+				/*Quick'n'dirty hack for flippability*/	
+				bool isFlippable;			
+				for (int jj = 0; jj<this->numVars; jj++)
+				{					
+					intvec *ivCanonical = new intvec(this->numVars);
+					(*ivCanonical)[jj]=1;					
+					if (dotProduct(load,ivCanonical)>=0)
+					{
+						isFlippable=FALSE;						
+					}
+					else
+					{
+						isFlippable=TRUE;
+					}					
+					delete ivCanonical;
+				}//for (int jj = 0; jj<this->numVars; jj++)
+				if (isFlippable==FALSE)
+				{
+					cout << "Ignoring facet";
+					load->show();
+					//fAct->next=NULL;
+				}
+				else
+				{	/*Now load should be full and we can call setFacetNormal*/
+					fAct->setFacetNormal(load);
+					fAct->next = new facet();
+					//fAct->printNormal();
+					fAct=fAct->next;	//this should definitely not be called in the above while-loop :D
+				}//if (isFlippable==FALSE)
 				delete load;
-			}
+			}//for (int kk = 0; kk<ddrows; kk++)
 			/*
 			Now we should have a linked list containing the facet normals of those facets that are
@@ -316,7 +342,4 @@
 		}//method getConeNormals(ideal I)	
 		
-		/*bool isParallel(int v[], intvec iv)
-		{
-}*/
 		
 		/** \brief Compute the Groebner Basis on the other side of a shared facet 
@@ -335,8 +358,9 @@
 		void flip(ideal gb, facet *f)		//Compute "the other side"
 		{			
-			
 			intvec *fNormal = new intvec(this->numVars);	//facet normal, check for parallelity			
 			fNormal = f->getFacetNormal();	//read this->fNormal;
 #ifdef gfan_DEBUG
+			cout << "===" << endl;
+			cout << "running gcone::flip" << endl;
 			cout << "fNormal=";
 			fNormal->show();
@@ -384,6 +408,4 @@
 				cout << "---" << endl;
 				/*Now initialFormElement must be added to (ideal)initialForm */
-				//f->flibGB->m[ii]=(poly)initialFormElement[ii];
-				//(poly)initialForm->m[ii]=pAdd(initialForm[ii],initialFormElement[ii]);
 				initialForm->m[ii]=initialFormElement[ii];
 			}//for
@@ -391,5 +413,5 @@
 #ifdef gfan_DEBUG
 			cout << "Initial ideal is: " << endl;
-			//idShow(initialForm);
+			idShow(initialForm);
 			f->printFlipGB();
 			cout << "===" << endl;
@@ -397,27 +419,45 @@
 			delete check;
 			
-			/*2nd step: lift initial ideal to a GB of the neighbouring cone*/
-			ring tmpring=rCopy0(currRing);
-			tmpring->order[0]=ringorder_a;
-			tmpring->order[1]=ringorder_dp;
-			tmpring->order[2]=ringorder_C;
-			//rWrite(tmpring);
-			tmpring->wvhdl[0] =( int *)omAlloc((fNormal->length())*sizeof(int)); //found in Singular/ipshell.cc
+			/*2nd step: lift initial ideal to a GB of the neighbouring cone using minus alpha as weight*/
+			/*Substep 2.1
+			compute $G_{-\alpha}(in_v(I)) 
+			see journal p. 66
+			*/
+			ring srcRing=currRing;
+			ring dstRing=rCopy0(srcRing);
+			dstRing->order[0]=ringorder_a;
+			//tmpring->order[1]=ringorder_dp;
+			//tmpring->order[2]=ringorder_C;
+			dstRing->wvhdl[0] =( int *)omAlloc((fNormal->length())*sizeof(int)); //found in Singular/ipshell.cc
 			
 			for (int ii=0;ii<this->numVars;ii++)
-			{
-				cout << "ping" << endl;
-				tmpring->wvhdl[0][ii]=(*fNormal)[ii];	//What exactly am I doing here?
-				cout << "pong" << endl;
-				cout << *tmpring->wvhdl[ii] << endl;
+			{				
+				dstRing->wvhdl[0][ii]=-(*fNormal)[ii];	//What exactly am I doing here?
+				//cout << tmpring->wvhdl[0][ii] << endl;
 			}
-			//tmpring->wvhdl=(int**)(fNormal);
-			rComplete(tmpring);
-			rWrite(tmpring);
-			/*setring(tmpring);
-			ideal ina=initialForm;
+			rComplete(dstRing);
+			rChangeCurrRing(dstRing);
+			map theMap=(map)idMaxIdeal(1);
+			rWrite(currRing); cout << endl;
+			ideal ina;
+			ina=fast_map(initialForm,srcRing,(ideal)theMap,dstRing);			
+			cout << "ina=";
+			idShow(ina); cout << endl;
 			ideal H;
-			H=kstd(ina,NULL,testHomog,NULL);
-			idShow(H);*/
+			H=kStd(ina,NULL,isHomog,NULL);	//we know it is homogeneous
+			idSkipZeroes(H);
+			cout << "H="; idShow(H); cout << endl;
+			/*Substep 2.2
+			do the lifting
+			*/
+			rChangeCurrRing(srcRing);
+			ideal srcRing_H;
+			ideal srcRing_HH;
+			//map theMap = (map)idMaxIdeal(1);
+			srcRing_H=fast_map(H,dstRing,(ideal)theMap,srcRing);
+			srcRing_HH=srcRing_H-stdred(srcRing_H,this->gcBasis);
+			/*Substep 2.3
+			turn the minimal basis into a reduced one
+			*/
 			
 		}//void flip(ideal gb, facet *f)
@@ -429,5 +469,5 @@
 		*\return void
 		*/
-		void getGB(const ideal &inputIdeal)
+		void getGB(ideal const &inputIdeal)
 		{
 			ideal gb;
@@ -435,5 +475,5 @@
 			idSkipZeroes(gb);
 			this->gcBasis=gb;	//write the GB into gcBasis
-		}
+		}//void getGB
 		
 		ideal GenGrbWlk(ideal, ideal);	//Implementation of the Generic Groebner Walk. Needed for a) Computing the sink and b) finding search facets
@@ -442,10 +482,9 @@
 		/** \brief Compute the dot product of two intvecs
 		*
-		* THIS IS WEIRD - BUT WORKS
 		*/
-		int dotProduct(intvec **a, intvec **b)				
-		{
-			intvec iva=*a;
-			intvec ivb=*b;
+		int dotProduct(intvec const &iva, intvec const &ivb)				
+		{
+			//intvec iva=a;
+			//intvec ivb=b;
 			int res=0;
 			for (int i=0;i<this->numVars;i++)
@@ -454,5 +493,5 @@
 			}
 			return res;
-		}
+		}//int dotProduct
 
 		/** \brief Check whether two intvecs are parallel
@@ -460,9 +499,9 @@
 		* \f$ \alpha\parallel\beta\Leftrightarrow\langle\alpha,\beta\rangle^2=\langle\alpha,\alpha\rangle\langle\beta,\beta\rangle \f$
 		*/
-		bool isParallel(intvec *a, intvec *b)
+		bool isParallel(intvec const &a, intvec const &b)
 		{			
 			int lhs,rhs;
-			lhs=dotProduct(&a,&b)*dotProduct(&a,&b);
-			rhs=dotProduct(&a,&a)*dotProduct(&b,&b);
+			lhs=dotProduct(a,b)*dotProduct(a,b);
+			rhs=dotProduct(a,a)*dotProduct(b,b);
 			cout << "LHS="<<lhs<<", RHS="<<rhs<<endl;
 			if (lhs==rhs)
@@ -475,28 +514,4 @@
 			}
 		}//bool isParallel
-		
-		/** \brief Check two intvecs for equality --- PROBABLY NOT NEEDED
-		*
-		* \f$ \alpha=\beta\Leftrightarrow\langle\alpha-\beta,\alpha-\beta\rangle=0 \f$
-		*/
-		bool isEqual(intvec a, intvec b)
-		{
-			intvec *ivdiff;
-			int res;
-			
-			ivdiff=ivSub(&a,&b);
-			cout << "ivdiff=";
-			ivdiff->show();
-			cout << endl;
-			//res=dotProduct(ivdiff,ivdiff);
-			if (res==0)
-			{
-				return TRUE;
-			}
-			else
-			{
-				return FALSE;
-			}			
-		}//bool isEqual
 };//class gcone
 
@@ -566,5 +581,5 @@
 	=> Count the cones!
 	*/
-	
+	rChangeCurrRing(inputRing);
 	res=gcAct->gcBasis;
 	//cout << fRoot << endl;
