pelqhull.h

Gambit 是一个游戏库理论软件

/* *****************************************************************

   This file, and all other pel*.h and pel*.cc files in the Gambit
source code, are derived from the source code for Pelican, an
implementation of the Huber-Sturmfels algorithm for solving sparse
systems of polynomials written by Birk Huber.  That code in turn
incorporated the source code for HOMPACK and Qhull, which are included
here after being modified in the construction of both Pelican and
Gambit.  We are grateful to the authors of Pelican, Qhull, and
HOMPACK, for code that has been enormously useful in Gambit.
    More information about Qhull can be obtained from:

http://www.geom.edu/software/qhull/

    For more information concerning HOMPACK see:

http://netlib2.cs.utk.edu/hompack/

***************************************************************** */

/* qhull.h -- user-level header file for using qhull.a library
   
   see README, qhull_a.h

   copyright (c) 1993-1994, The Geometry Center

   defines qh_qh, global data structure for qhull.  
   
   NOTE: access to qh_qh is via the 'qh' macro.  This allows
   qh_qh to be either a pointer or a structure.  An example
   of using qh is "qh DROPdim" which accesses the DROPdim
   field of qh_qh.  Similarly, access to qh_qhstat is via
   the 'qhstat' macro.
   
   includes function prototypes for qhull.c, geom.c, global.c, io.c, user.c

   use mem.h for mem.c
   use set.h for set.c
   
   see unix.c for an example of using qhull.h
   
   recompile qhull if you change this file
*/


#include <stdio.h>
#include <stdlib.h>
#include <setjmp.h>
#include <string.h>
#include <math.h>
#include <float.h>    /* some compilers will not need float.h */
#include <limits.h>
#include <time.h>
#include <ctype.h>

#ifndef qhDEFqhull
#define qhDEFqhull 1

#include <setjmp.h>
#include <float.h>
#include <time.h>

/* ============ -types- ==================== */

#define pointT coordT                 /* array of hull_dim coordinates */
#define flagT unsigned int            /* Boolean flag */
typedef enum {False, True} boolT;     /* True=1, False= 0 */

/* -center type for facet->center and CENTERtype */
enum qh_CENTER {qh_none, qh_voronoi, qh_centrum};

/* -output formats for printing (qh PRINTout) */
enum qh_PRINT {qh_PRINTnone, qh_PRINTfacets, qh_PRINTgeom, qh_PRINTnormals, 
  qh_PRINTincidences, qh_PRINTmathematica, qh_PRINToff, qh_PRINTpoints,
  qh_PRINTEND};



/*---------------------------------
-realT -- select size of floating point numbers
   normally use 'float' because convex hulls are space limited in high-d
      this saves about 8% in time and 25% in space
   could use 'double' for internal calculations (realT vs. coordT)
        but this requires many type casts, and adjusted error bounds
   C compilers will often do expressions in double
   WARNING: on some machines, can not test for equality after setting REALmax.
                        Use REALmax/2 instead.
*/

#define REALfloat 1   /* use 0 for double precision */
#if REALfloat
#define realT float
#define REALmax FLT_MAX
#define REALmin FLT_MIN
#define REALepsilon FLT_EPSILON
#define qh_REALdigits 8   /* maximum number of significant digits */
#define qh_REAL_1 "%6.8g "
#define qh_REAL_2n "%6.8g %6.8g\n"
#define qh_REAL_3n "%6.8g %6.8g %6.8g\n"
#else
#define realT double
#define REALmax DBL_MAX
#define REALmin DBL_MIN
#define REALepsilon DBL_EPSILON
#define qh_REALdigits 16    /* maximum number of significant digits */
#define qh_REAL_1 "%6.16g "
#define qh_REAL_2n "%6.16g %6.16g\n"
#define qh_REAL_3n "%6.16g %6.16g %6.16g\n"
#endif
#define coordT realT  /* for stored coordinates and coefficients */

/*------------------------------------------
        constants
*/
#define qh_DEFAULTbox 0.5    /* default box size (Geomview expects 0.5) */
#define qh_INFINITE  -10.101 /* coordinates for Voronoi center at infinity */

#define qh_ALL 	    True     /* argument for printall and checkall parameters*/


/*------------------------------------
	Conditional compilation
	
-COMPUTEfurthest computing furthest saves 4% of memory and costs ca. 3% time
                         (about 40% more distance tests for partitioning)
-MAXoutside      keep maxoutside for each facet
		   this takes a realT per facet and slightly slows down qhull
		   it should speed up pre-merging with Wn or Qm
		   better outer planes for geomview output 
*/

#define qh_COMPUTEfurthest 0    /* 1 removes facet->furthestdist */
#define qh_MAXoutside 1         /* 0 removes facet->maxoutside */


/*-------------------------------------------
-ERR - qhull exit codes, for indicating errors
*/
#define qh_ERRnone  0    /* no error occurred during qhull */
#define qh_ERRinput 1    /* input inconsistency */
#define qh_ERRsingular 2 /* singular input data */
#define qh_ERRprec  3    /* precision error */
#define qh_ERRmem   4    /* insufficient memory, matches mem.h */
#define qh_ERRqhull 5    /* internal error detected, matches mem.h */

/* ============ -structures- ====================
   each of the following structures is defined by a typedef
   all realT and coordT fields occur at the beginning of a structure
        (otherwise space may be wasted due to alignment)
   define all flags together and pack into 32-bit number
*/

typedef struct vertexT vertexT;
typedef struct ridgeT ridgeT;
typedef struct facetT facetT;
#ifndef DEFsetT
#define DEFsetT 1
typedef struct setT setT;          /* defined in set.h */
#endif

/* ----------------------------------------------
-facetT- specifies a facet
*/

struct facetT {
#if !qh_COMPUTEfurthest
  coordT   furthestdist;/* distance to furthest point of outsideset */
#endif
#if qh_MAXoutside
  coordT   maxoutside;  /* max computed distance of point to facet
  			Before QHULLfinished this is an approximation
  			since maxdist not always set for mergefacet */
#endif
  coordT   offset;      /* offset of hyperplane from origin */
  coordT  *normal;      /* normal of hyperplane, hull_dim coefficients */
  coordT  *center;      /* centrum or Voronoi center, see CENTERtype */
  facetT  *previous;    /* previous facet in the facet_list */
  facetT  *next;        /* next facet in the facet_list */
  setT    *vertices;    /* vertices for this facet, inverse sorted by id */
  setT    *ridges;      /* explicit ridges for nonsimplicial facets.
  			   for simplicial facets, neighbors defines ridge */
  setT    *neighbors;   /* neighbors of the facet.
  			   if simplicial, kth neighbor skips kth vertex.
  			   if visible, 1st neighbor (if any) is newfacet 
  			         or deleted newfacet (visible) */
  setT    *outsideset;  /* set of points outside this facet
		           if non-empty, last point is furthest */
  setT    *coplanarset; /* set of points coplanar with this facet
  			   > min_vertex and <= facet->max_outside
                           a point is assigned to the furthest facet
		           if non-empty, last point is furthest away */
  unsigned visitid;     /* visit_id, for visiting all neighbors,
			   all uses must be independent */
  unsigned id:24;       /* unique identifier, =>room for 8 flags 
  			   new facet if >= qh newfacet_id */
  flagT    toporient:1; /* True if created with top orientation
			   after merging, use ridge orientation */
  flagT    simplicial:1;/* True if simplicial facet, ->ridges may be implicit */
  flagT	   good:1;      /* True if a good facet (for output) */
  flagT    tested:1;    /* True if facet changed due to merge */
  flagT    seen:1;      /* used to perform operations only once, like visitid */
  flagT	   flipped:1;   /* True (if !NEWmerges) if facet is flipped */
  flagT	   visible:1;   /* True if visible facet */
};
#define qh_isnewmerge_(facet) ((facet)->flipped)    /* if qh NEWmerges */
#define qh_newmerge_(facet) {(facet)->flipped= True; numnewmerges++;}
#define qh_clearnewmerge_(facet) (facet)->flipped= False;


/*----------------------------------------------
-ridgeT- specifies a ridge
*/

struct ridgeT {
  setT    *vertices;    /* vertices belonging to this ridge, inverse sorted by id 
                           NULL if a degen ridge (matchsame) */
  facetT  *top;         /* top facet this ridge is part of */
  facetT  *bottom;      /* bottom facet this ridge is part of */
  unsigned id:24;       /* unique identifier, =>room for 8 flags */
  flagT    seen:1;      /* used to perform operations only once */
  flagT    tested:1;    /* True when ridge is tested for convexity */
  flagT    nonconvex:1; /* True if getmergeset detected a non-convex neighbor */
  flagT    mergeridge:1; /* True if from a qh_MERGEridge (duplicates) */
};

/* ----------------------------------------------
-vertexT- specifies a vertex
*/

struct vertexT {
  vertexT *next;        /* next vertex in vertex_list */
  vertexT *previous;    /* previous vertex in vertex_list */
  pointT  *point;       /* hull_dim coordinates (coordT) */
  setT    *neighbors;   /* neighboring facets of vertex, if qh_vertexneighbors */
  unsigned visitid;     /* for use with qh vertex_visit */
  unsigned id:24;       /* unique identifier, =>room for 8 flags */
  flagT    seen:1;      /* used to perform operations only once */
  flagT    seen2:1;     /* another seen flag */
  flagT    delridge:1;  /* vertex was part of a deleted ridge */
  flagT	   deleted:1;   /* true if vertex on qh del_vertices */
  flagT    newlist:1;   /* true if vertex on qh newvertex_list */
};

/* ======= -global variables -qh ============================ 

   all global variables for qhull are in qh, qhmem, and qhstat
   
   qhmem is defined in mem.h and qhstat is defined in stat.h
   
   set qh_QHpointer 1 to enable qh_saveqhull() and qh_restoreqhull()
     this costs about 7% in time and space.
*/

typedef struct qhT qhT;    

#define qh_QHpointer 0  /* 1 for dynamic allocation, 0 for global structure */
#if qh_QHpointer
#define qh qh_qh->
extern qhT *qh_qh;     /* allocated in global.c */
#else
#define qh qh_qh.
extern qhT qh_qh;
#endif

struct qhT {

  /*-user flags */

  boolT ALLpoints;        /* true 'Qi' if search all points for initial simplex */
  boolT APPROXhull;       /* true 'W' if MINoutside set */
  realT MINoutside;       /*   min. distance for an outside point */
  boolT BESToutside;      /* true 'Qf' if partition points into best outsideset */
  boolT CHECKfrequently;  /* true 'Tc' if checking frequently */
  realT premerge_cos;     /*   'C-n'   cos_max when pre merging */
  realT postmerge_cos;    /*   'Cn'    cos_max when post merging */
  boolT DELAUNAY;         /* true 'd' if computing DELAUNAY triangulation */
  boolT DOintersections;  /* print hyperplaneintersections in 3-d and 4-d */
  int   DROPdim;          /* drops dim 'GDn' for 4-d -> 3-d output */
  boolT ERREXITcalled;    /* true during errexit (prevents duplicate calls */
  boolT FORCEoutput;      /* true 'Po' if forcing output despite degeneracies */
  int   GOODpoint;        /* 1+n, good facet if visible/not(-) from point n*/
  pointT *GOODpointp;     /*   the actual point */
  boolT GOODthreshold;    /* true if qh lower_threshold/upper_threshold defined 
  			     false if qh SPLITthreshold */
  int   GOODvertex;       /* 1+n, good facet if vertex for point n */
  pointT *GOODvertexp;     /*   the actual point */
  int   IStracing;        /* trace execution, 0=none, 1=least, 4=most, -1=events */
  boolT KEEPcoplanar;     /* true if keeping nearest facet for coplanar points */
  boolT KEEPinside;       /* true if keeping nearest facet for inside/coplanar points */