source: git/dyn_modules/callgfanlib/gfan.h @ 81384b

jengelh-datetimespielwiese
Last change on this file since 81384b was 81384b, checked in by Yue Ren <ren@…>, 10 years ago
new: newest version of callgfanlib package as master branch also includes minor changes in bigintmat to allow empty matrices and a new function which returns the regular print output as char* (necessary for the print routines of the convex objects)
  • Property mode set to 100644
File size: 9.4 KB
Line 
1/*
2gfan.h Interface to gfan.cc
3
4Author: monerjan
5*/
6#ifdef HAVE_FANS
7
8#ifndef GFAN_H
9#define GFAN_H
10
11#include <libpolys/misc/int64vec.h>
12
13#include <setoper.h>
14#include <cdd.h>
15#include <cddmp.h>
16#include <bbfan.h>
17#include <bbcone.h>
18extern int gfanHeuristic;
19
20#ifndef USE_ZFAN
21#define USE_ZFAN
22#endif
23#ifndef USE_ZFAN
24  lists grfan(ideal inputIdeal, int heuristic, bool singleCone);
25#else
26  #include <gfanlib/gfanlib.h>
27  gfan::ZFan *grfan(ideal inputIdeal, int h, bool singleCone);
28#endif
29// lists grcone_by_intvec(ideal inputIdeal);
30
31class facet
32{
33        private:
34                /** \brief Inner normal of the facet, describing it uniquely up to isomorphism */
35                int64vec *fNormal;
36
37                /** \brief An interior point of the facet*/
38                int64vec *interiorPoint;
39
40                /** \brief Universal Cone Number
41                 * The number of the cone the facet belongs to, Set in getConeNormals()
42                 */
43                int UCN;
44
45                /** \brief The codim of the facet
46                 */
47                short codim;
48
49                /** \brief The Groebner basis on the other side of a shared facet
50                 *
51                 * In order not to have to compute the flipped GB twice we store the basis we already get
52                 * when identifying search facets. Thus in the next step of the reverse search we can
53                 * just copy the old cone and update the facet and the gcBasis.
54                 * facet::flibGB is set via facet::setFlipGB() and printed via facet::printFlipGB
55                 */
56                ideal flipGB;           //The Groebner Basis on the other side, computed via gcone::flip
57
58        public:
59                /** \brief Boolean value to indicate whether a facet is flippable or not
60                * This is also used to mark facets that nominally are flippable but which do
61                * not intersect with the positive orthant. This check is done in gcone::getCodim2Normals
62                 */
63                bool isFlippable;       //**flippable facet? */
64                //bool isIncoming;      //Is the facet incoming or outgoing in the reverse search? No longer in use
65                facet *next;            //Pointer to next facet
66                facet *prev;            //Pointer to predecessor. Needed for the SearchList in noRevS
67                facet *codim2Ptr;       //Pointer to (codim-2)-facet. Bit of recursion here ;-)
68                int numCodim2Facets;    //#of (codim-2)-facets of this facet. Set in getCodim2Normals()
69                unsigned numRays;       //Number of spanning rays of the facet
70                ring flipRing;          //the ring on the other side of the facet
71//              int64vec **fRays;
72
73                /** The default constructor. */
74                facet();
75                /** Constructor for lower dimensional faces*/
76                facet(const int &n);
77                /**  The copy constructor */
78                facet(const facet& f);
79                /** A shallow copy of facets*/
80                facet* shallowCopy(const facet& f);
81                void shallowDelete();
82                /** The default destructor */
83                ~facet();
84                /** Comparison operator*/
85//              inline bool operator==(const facet *f,const facet *g);
86                /** \brief Comparison of facets*/
87//              inline bool areEqual(facet *f, facet *g);//Now static
88                /** Stores the facet normal \param int64vec*/
89                inline void setFacetNormal(int64vec *iv);
90                /** Returns the facet normal */
91                inline int64vec *getFacetNormal() const;
92                /** Return a reference to the facet normal*/
93                inline const int64vec *getRef2FacetNormal() const;
94                /** Method to print the facet normal*/
95                inline void printNormal() const;
96                /** Store the flipped GB*/
97                inline void setFlipGB(ideal I);
98                /** Return the flipped GB*/
99                inline ideal getFlipGB();
100                /** Print the flipped GB*/
101                inline void printFlipGB();
102                /** Set the UCN */
103                inline void setUCN(int n);
104                /** \brief Get the UCN
105                 * Returns the UCN iff this != NULL, else -1
106                 */
107                inline int getUCN();
108                /** Store an interior point of the facet */
109                inline void setInteriorPoint(int64vec *iv);
110                inline int64vec *getInteriorPoint();
111                inline const int64vec *getRef2InteriorPoint();
112                /** \brief Debugging function
113                 * prints the facet normal an all (codim-2)-facets that belong to it
114                 */
115                volatile void fDebugPrint();
116                friend class gcone;
117};
118
119
120/**
121 *\brief Implements the cone structure
122 *
123 * A cone is represented by a linked list of facet normals
124 * @see facet
125 */
126
127class gcone
128{
129        private:
130                ideal inputIdeal;       //the original
131                ring baseRing;          //the basering of the cone
132                int64vec *ivIntPt;      //an interior point of the cone
133                int UCN;                //unique number of the cone
134                int pred;               //UCN of the cone this one is derived from
135                static int counter;
136
137        public:
138                /** \brief Pointer to the first facet */
139                facet *facetPtr;        //Will hold the adress of the first facet; set by gcone::getConeNormals
140#ifdef gfanp
141                static float time_getConeNormals;
142                static float time_getCodim2Normals;
143                static float t_getExtremalRays;
144                static float t_ddPolyh;
145                static float time_flip;
146                static float time_flip2;
147                static float t_areEqual;
148                static float t_ffG;
149                static float t_markings;
150                static float t_dd;
151                static float t_kStd;
152                static float time_enqueue;
153                static float time_computeInv;
154                static float t_ddMC;
155                static float t_mI;
156                static float t_iP;
157                static float t_isParallel;
158                static unsigned parallelButNotEqual;
159                static unsigned numberOfFacetChecks;
160#endif
161                /** Matrix to contain the homogeneity/lineality space */
162                static dd_MatrixPtr dd_LinealitySpace;
163                static int lengthOfSearchList;
164                /** Maximum size of the searchlist*/
165                static int maxSize;
166                /** is the ideal homogeneous? */
167                static bool hasHomInput;
168                /** # of variables in the ring */
169                static int numVars;             //#of variables in the ring
170                /** The hilbert function - for the homogeneous case*/
171                static int64vec *hilbertFunction;
172                /** The zero vector. Needed in case of fNormal mismatch*/
173                static int64vec *ivZeroVector;
174
175                /** # of facets of the cone
176                 * This value is set by gcone::getConeNormals
177                 */
178                int numFacets;          //#of facets of the cone
179
180                /**
181                 * At least as a workaround we store the irredundant facets of a matrix here.
182                 * This is needed to compute an interior points of a cone. Note that there
183                 * will be non-flippable facets in it!
184                 */
185                dd_MatrixPtr ddFacets;  //Matrix to store irredundant facets of the cone
186
187                /** Array of intvecs representing the rays of the cone*/
188                int64vec **gcRays;
189                unsigned numRays;       //#rays of the cone
190                /** Contains the Groebner basis of the cone. Is set by gcone::getGB(ideal I)*/
191                ideal gcBasis;          //GB of the cone, set by gcone::getGB();
192                gcone *next;            //Pointer to next cone
193                gcone *prev;
194
195                gcone();
196                gcone(ring r, ideal I);
197                gcone(const gcone& gc, const facet &f);
198                ~gcone();
199                inline int getCounter();
200                inline ring getBaseRing();
201                inline ring getRef2BaseRing();
202                inline void setBaseRing(ring r);
203                inline void setIntPoint(int64vec *iv);
204                inline int64vec *getIntPoint(bool shallow=FALSE);
205                inline void showIntPoint();
206                inline void setNumFacets();
207                inline int getNumFacets();
208                inline int getUCN();
209                inline int getPredUCN();
210                volatile void showFacets(short codim=1);
211//              volatile void showSLA(facet &f);
212//              void idDebugPrint(const ideal &I);
213//              void invPrint(const ideal &I);
214//              bool isMonomial(const ideal &I);
215//              int64vec *ivNeg(const int64vec *iv);
216//              inline int dotProduct(int64vec &iva, int64vec &ivb);
217//              inline int dotProduct(const int64vec &iva, const int64vec &ivb);
218//              inline bool isParallel(const int64vec &a, const int64vec &b);
219                void noRevS(gcone &gcRoot, bool usingIntPoint=FALSE);
220//              inline int intgcd(const int &a, const int &b);
221                void writeConeToFile(const gcone &gc, bool usingIntPoints=FALSE);
222                void readConeFromFile(int gcNum, gcone *gc);
223                int64vec f2M(gcone *gc, facet *f, int n=1);
224//              inline void sortRays(gcone *gc);
225                //The real stuff
226                void getConeNormals(const ideal &I, bool compIntPoint=FALSE);
227                void getCodim2Normals(const gcone &gc);
228                void getExtremalRays(const gcone &gc);
229                void orderRays();
230                void flip(ideal gb, facet *f);
231                void flip2(const ideal &gb, facet *f);
232                void computeInv(const ideal &gb, ideal &inv, const int64vec &f);
233                //poly restOfDiv(poly const &f, ideal const &I); removed with r12286
234                inline ideal ffG(const ideal &H, const ideal &G);
235                inline void getGB(ideal const &inputIdeal);
236                void interiorPoint( dd_MatrixPtr &M, int64vec &iv);//used from flip and optionally from getConeNormals
237//              void interiorPoint2(); //removed Feb 8th, 2010, new method Feb 19th, 2010, again removed Mar 16th, 2010
238                void preprocessInequalities(dd_MatrixPtr &M);
239                ring rCopyAndAddWeight(const ring &r, int64vec *ivw);
240                ring rCopyAndAddWeight2(const ring &, const int64vec *, const int64vec *);
241//              ring rCopyAndChangeWeight(const ring &r, int64vec *ivw);        //NOTE remove
242//              void reverseSearch(gcone *gcAct); //NOTE both removed from r12286
243//              bool isSearchFacet(gcone &gcTmp, facet *testfacet); //NOTE remove
244                void makeInt(const dd_MatrixPtr &M, const int line, int64vec &n);
245//              void normalize();//NOTE REMOVE
246                facet * enqueueNewFacets(facet *f);
247                facet * enqueue2(facet *f);
248//              dd_MatrixPtr facets2Matrix(const gcone &gc);//NOTE remove
249                /** Compute the lineality space Ax=0 and return it as dd_MatrixPtr dd_LinealitySpace*/
250                dd_MatrixPtr computeLinealitySpace();
251                inline bool iv64isStrictlyPositive(const int64vec *);
252                /** Exchange 2 ordertype_a by just 1 */
253                void replaceDouble_ringorder_a_ByASingleOne();
254//              static void gcone::idPrint(ideal &I);
255//              friend class facet;
256};
257lists lprepareResult(gcone *gc, const int n);
258static int64 int64gcd(const int64 &a, const int64 &b);
259static int intgcd(const int &a, const int &b);
260static int dotProduct(const int64vec &iva, const int64vec &ivb);
261static bool isParallel(const int64vec &a, const int64vec &b);
262static int64vec *ivNeg(/*const*/ int64vec *iv);
263static void idDebugPrint(const ideal &I);
264static volatile void showSLA(facet &f);
265static bool isMonomial(const ideal &I);
266static bool ivAreEqual(const int64vec &a, const int64vec &b);
267static bool areEqual2(facet *f, facet *g);
268static bool areEqual( facet *f, facet *g);
269// bool iv64isStrictlyPositive(int64vec *);
270#endif
271#endif
Note: See TracBrowser for help on using the repository browser.