Changeset 864fc2 in git for kernel/gfan.cc

Timestamp:

Apr 8, 2009, 4:02:26 PM (15 years ago)

Author:

Martin Monerjan

Branches:

(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')

Children:

e7778adedbf4b7cd363cb88f576d7b16023f2ca3

Parents:

7e96fe079ce8711823d420d4b4ba8e959e156a77

Message:

created gcone::restOfDiv(poly, ideal) to compute f^{G} and
gcone::ffG(ideal, ideal) to compute f-f^{G}. However checks and construction
of zeropolynomial don't work. This rev will hang in restOfDiv line 488

Substituted fast_map stuff by idrCopyR


git-svn-id: file:///usr/local/Singular/svn/trunk@11647 2c84dea3-7e68-4137-9b89-c4e89433aadc

File:

• kernel/gfan.cc (modified) (8 diffs)

kernel/gfan.cc

@@ -2 +2 @@
 Compute the Groebner fan of an ideal
 $Author: monerjan $
-$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 $
+$Date: 2009-04-08 14:02:26 $
+$Header: /exports/cvsroot-2/cvsroot/kernel/gfan.cc,v 1.31 2009-04-08 14:02:26 monerjan Exp $
+$Id: gfan.cc,v 1.31 2009-04-08 14:02:26 monerjan Exp $
 */
 
@@ -12 +12 @@
 
 #include "kstd1.h"
+#include "kutil.h"
 #include "intvec.h"
 #include "polys.h"
@@ -18 +19 @@
 #include "fast_maps.h"  //Mapping of ideals
 #include "maps.h"
+#include "ring.h"
+#include "prCopy.h"
 #include "iostream.h"   //deprecated
 
@@ -438 +441 @@
                         rComplete(dstRing);
                         rChangeCurrRing(dstRing);
-                        map theMap=(map)idMaxIdeal(1);
+                        //map theMap=(map)idMaxIdeal(1);
                         rWrite(currRing); cout << endl;
                         ideal ina;
-                        ina=fast_map(initialForm,srcRing,(ideal)theMap,dstRing);                        
+                        //ina=fast_map(initialForm,srcRing,(ideal)theMap,dstRing);                      
+                        ina=idrCopyR(initialForm,srcRing);
                         cout << "ina=";
                         idShow(ina); cout << endl;
@@ -455 +459 @@
                         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);
+                        //srcRing_H=fast_map(H,dstRing,(ideal)theMap,srcRing);
+                        srcRing_H=idrCopyR(H,dstRing);
+                        idShow(srcRing_H);
+                        srcRing_HH=ffG(srcRing_H,this->gcBasis);                        
+                        idShow(srcRing_HH);
                         /*Substep 2.3
                         turn the minimal basis into a reduced one
@@ -463 +470 @@
                 }//void flip(ideal gb, facet *f)
                                 
+                /** \brief Compute the remainder of a polynomial by a given ideal
+                *
+                * Compute \f$ f^{\mathcal{G}} \f$
+                * Algorithm is taken from Cox, Little, O'Shea, IVA 2nd Ed. p 62
+                * However, since we are only interested in the remainder, there is no need to
+                * compute the factors \f$ a_i \f$
+                */
+                poly restOfDiv(poly const &f, ideal const &I)
+                {
+                        cout << "Entering restOfDiv" << endl;
+                        poly p=f;
+                        pWrite(p);
+                        //poly r=kCreateZeroPoly(,currRing,currRing);   //The 0-polynomial, hopefully
+                        poly r=0;
+                        int ii;
+                        bool divOccured;
+                        cout << "TICK 1" << endl;       //Hangs after this line. Zeropoly stuff
+                        while (pHead(p)!=NULL)
+                        {
+                                ii=1;
+                                divOccured=TRUE;
+                                while( (ii<=IDELEMS(I) && (divOccured==FALSE) ))
+                                {
+                                        cout << "TICK 2" << endl;
+                                        if (pLmDivisibleBy(I->m[ii],p))
+                                        {
+                                                cout << "TICK 3" << endl;
+                                                p=pSub(p,pMult( pDivide(pHead(p),pHead(I->m[ii])), I->m[ii]));
+                                                cout << "TICK 4" << endl;
+                                                pWrite(p);
+                                                pSetm(p);
+                                                divOccured=TRUE;
+                                        }
+                                        else
+                                        {
+                                                ii += 1;
+                                        }//if (pLmDivisibleBy(I->m[ii],p,currRing))
+                                        if (divOccured==FALSE)
+                                        {
+                                                cout << "TICK 5" << endl;
+                                                r=pAdd(r,pHead(p));
+                                                cout << "TICK 6" << endl;
+                                                p=pSub(p,pHead(p));
+                                                cout << "TICK 7" << endl;
+                                                pWrite(p);
+                                        }//if (divOccured==FALSE)
+                                }//while( (ii<IDELEMS(I) && (divOccured==FALSE) ))
+                        }//while (p!=0)
+                        return r;
+                }//poly restOfDiv(poly const &f, ideal const &I)
+                
+                /** \brief Compute \f$ f-f^{\mathcal{G}} \f$
+                */
+                ideal ffG(ideal const &H, ideal const &G)
+                {
+                        cout << "Entering ffG" << endl;
+                        int size=IDELEMS(G);
+                        ideal res=idInit(size,1);
+                        for (int ii=0;ii<size;ii++)
+                        {
+                                res->m[ii]=restOfDiv(H->m[ii],G);
+                        }
+                        return res;
+                }
+                
                 /** \brief Compute a Groebner Basis
                 *
@@ -537 +609 @@
         /* Construct a new ring which will serve as our root
         Does not yet work as expected. Will work fine with order dp,Dp but otherwise hangs in getGB
+        resolved 07.04.2009 MM
         */
         rootRing=rCopy0(currRing);
@@ -545 +618 @@
         /* Fetch the inputIdeal into our rootRing */
         map theMap=(map)idMaxIdeal(1);  //evil hack!
-        //idShow(idMaxIdeal(1));
-        /*for (int ii=0;ii<pVariables;ii++)
-        {
-                theMap->m[ii]=inputIdeal->m[ii];
-        }*/
-        theMap->preimage=NULL;
+        theMap->preimage=NULL;  //neccessary?
         ideal rootIdeal;
         rootIdeal=fast_map(inputIdeal,inputRing,(ideal)theMap, currRing);