1 /*<html><pre> -<a href="qh-user.htm" 2 >-------------------------------</a><a name="TOP">-</a> 3 4 user.h 5 user redefinable constants 6 7 see qh-user.htm. see COPYING for copyright information. 8 9 before reading any code, review qhull.h for data structure definitions and 10 the "qh" macro. 11 */ 12 13 #ifndef qhDEFuser 14 #define qhDEFuser 1 15 16 /*============= data types and configuration macros ==========*/ 17 18 /*-<a href="qh-user.htm#TOC" 19 >--------------------------------</a><a name="realT">-</a> 20 21 realT 22 set the size of floating point numbers 23 24 qh_REALdigits 25 maximimum number of significant digits 26 27 qh_REAL_1, qh_REAL_2n, qh_REAL_3n 28 format strings for printf 29 30 qh_REALmax, qh_REALmin 31 maximum and minimum (near zero) values 32 33 qh_REALepsilon 34 machine roundoff. Maximum roundoff error for addition and multiplication. 35 36 notes: 37 Select whether to store floating point numbers in single precision (float) 38 or double precision (double). 39 40 Use 'float' to save about 8% in time and 25% in space. This is particularly 41 help if high-d where convex hulls are space limited. Using 'float' also 42 reduces the printed size of Qhull's output since numbers have 8 digits of 43 precision. 44 45 Use 'double' when greater arithmetic precision is needed. This is needed 46 for Delaunay triangulations and Voronoi diagrams when you are not merging 47 facets. 48 49 If 'double' gives insufficient precision, your data probably includes 50 degeneracies. If so you should use facet merging (done by default) 51 or exact arithmetic (see imprecision section of manual, qh-impre.htm). 52 You may also use option 'Po' to force output despite precision errors. 53 54 You may use 'long double', but many format statements need to be changed 55 and you may need a 'long double' square root routine. S. Grundmann 56 (sg@eeiwzb.et.tu-dresden.de) has done this. He reports that the code runs 57 much slower with little gain in precision. 58 59 WARNING: on some machines, int f(){realT a= REALmax;return (a == REALmax);} 60 returns False. Use (a > REALmax/2) instead of (a == REALmax). 61 62 REALfloat = 1 all numbers are 'float' type 63 = 0 all numbers are 'double' type 64 */ 65 #define REALfloat 0 66 67 #if (REALfloat == 1) 68 #define realT float 69 #define REALmax FLT_MAX 70 #define REALmin FLT_MIN 71 #define REALepsilon FLT_EPSILON 72 #define qh_REALdigits 8 /* maximum number of significant digits */ 73 #define qh_REAL_1 "%6.8g " 74 #define qh_REAL_2n "%6.8g %6.8g\n" 75 #define qh_REAL_3n "%6.8g %6.8g %6.8g\n" 76 77 #elif (REALfloat == 0) 78 #define realT double 79 #define REALmax DBL_MAX 80 #define REALmin DBL_MIN 81 #define REALepsilon DBL_EPSILON 82 #define qh_REALdigits 16 /* maximum number of significant digits */ 83 #define qh_REAL_1 "%6.16g " 84 #define qh_REAL_2n "%6.16g %6.16g\n" 85 #define qh_REAL_3n "%6.16g %6.16g %6.16g\n" 86 87 #else 88 #error unknown float option 89 #endif 90 91 /*-<a href="qh-user.htm#TOC" 92 >--------------------------------</a><a name="CPUclock">-</a> 93 94 qh_CPUclock 95 define the clock() function for reporting the total time spent by Qhull 96 returns CPU ticks as a 'long int' 97 qh_CPUclock is only used for reporting the total time spent by Qhull 98 99 qh_SECticks 100 the number of clock ticks per second 101 102 notes: 103 looks for CLOCKS_PER_SEC, CLOCKS_PER_SECOND, or assumes microseconds 104 to define a custom clock, set qh_CLOCKtype to 0 105 106 if your system does not use clock() to return CPU ticks, replace 107 qh_CPUclock with the corresponding function. It is converted 108 to unsigned long to prevent wrap-around during long runs. 109 110 111 Set qh_CLOCKtype to 112 113 1 for CLOCKS_PER_SEC, CLOCKS_PER_SECOND, or microsecond 114 Note: may fail if more than 1 hour elapsed time 115 116 2 use qh_clock() with POSIX times() (see global.c) 117 */ 118 #define qh_CLOCKtype 1 /* change to the desired number */ 119 120 #if (qh_CLOCKtype == 1) 121 122 #if defined (CLOCKS_PER_SECOND) 123 #define qh_CPUclock ((unsigned long)clock()) /* return CPU clock */ 124 #define qh_SECticks CLOCKS_PER_SECOND 125 126 #elif defined (CLOCKS_PER_SEC) 127 #define qh_CPUclock ((unsigned long)clock()) /* return CPU clock */ 128 #define qh_SECticks CLOCKS_PER_SEC 129 130 #elif defined (CLK_TCK) 131 #define qh_CPUclock ((unsigned long)clock()) /* return CPU clock */ 132 #define qh_SECticks CLK_TCK 133 134 #else 135 #define qh_CPUclock ((unsigned long)clock()) /* return CPU clock */ 136 #define qh_SECticks 1E6 137 #endif 138 139 #elif (qh_CLOCKtype == 2) 140 #define qh_CPUclock qh_clock() /* return CPU clock */ 141 #define qh_SECticks 100 142 143 #else /* qh_CLOCKtype == ? */ 144 #error unknown clock option 145 #endif 146 147 /*-<a href="qh-user.htm#TOC" 148 >--------------------------------</a><a name="RANDOM">-</a> 149 150 qh_RANDOMtype, qh_RANDOMmax, qh_RANDOMseed 151 define random number generator 152 153 qh_RANDOMint generates a random integer between 0 and qh_RANDOMmax. 154 qh_RANDOMseed sets the random number seed for qh_RANDOMint 155 156 Set qh_RANDOMtype (default 5) to: 157 1 for random() with 31 bits (UCB) 158 2 for rand() with RAND_MAX or 15 bits (system 5) 159 3 for rand() with 31 bits (Sun) 160 4 for lrand48() with 31 bits (Solaris) 161 5 for qh_rand() with 31 bits (included with Qhull) 162 163 notes: 164 Random numbers are used by rbox to generate point sets. Random 165 numbers are used by Qhull to rotate the input ('QRn' option), 166 simulate a randomized algorithm ('Qr' option), and to simulate 167 roundoff errors ('Rn' option). 168 169 Random number generators differ between systems. Most systems provide 170 rand() but the period varies. The period of rand() is not critical 171 since qhull does not normally use random numbers. 172 173 The default generator is Park & Miller's minimal standard random 174 number generator [CACM 31:1195 '88]. It is included with Qhull. 175 176 If qh_RANDOMmax is wrong, qhull will report a warning and Geomview 177 output will likely be invisible. 178 */ 179 #define qh_RANDOMtype 5 /* *** change to the desired number *** */ 180 181 #if (qh_RANDOMtype == 1) 182 #define qh_RANDOMmax ((realT)0x7fffffffUL) /* 31 bits, random()/MAX */ 183 #define qh_RANDOMint random() 184 #define qh_RANDOMseed_(seed) srandom(seed); 185 186 #elif (qh_RANDOMtype == 2) 187 #ifdef RAND_MAX 188 #define qh_RANDOMmax ((realT)RAND_MAX) 189 #else 190 #define qh_RANDOMmax ((realT)32767) /* 15 bits (System 5) */ 191 #endif 192 #define qh_RANDOMint rand() 193 #define qh_RANDOMseed_(seed) srand((unsigned)seed); 194 195 #elif (qh_RANDOMtype == 3) 196 #define qh_RANDOMmax ((realT)0x7fffffffUL) /* 31 bits, Sun */ 197 #define qh_RANDOMint rand() 198 #define qh_RANDOMseed_(seed) srand((unsigned)seed); 199 200 #elif (qh_RANDOMtype == 4) 201 #define qh_RANDOMmax ((realT)0x7fffffffUL) /* 31 bits, lrand38()/MAX */ 202 #define qh_RANDOMint lrand48() 203 #define qh_RANDOMseed_(seed) srand48(seed); 204 205 #elif (qh_RANDOMtype == 5) 206 #define qh_RANDOMmax ((realT)2147483646UL) /* 31 bits, qh_rand/MAX */ 207 #define qh_RANDOMint qh_rand() 208 #define qh_RANDOMseed_(seed) qh_srand(seed); 209 /* unlike rand(), never returns 0 */ 210 211 #else 212 #error: unknown random option 213 #endif 214 215 /*-<a href="qh-user.htm#TOC" 216 >--------------------------------</a><a name="ORIENTclock">-</a> 217 218 qh_ORIENTclock 219 0 for inward pointing normals by Geomview convention 220 */ 221 #define qh_ORIENTclock 0 222 223 224 /*========= performance related constants =========*/ 225 226 /*-<a href="qh-user.htm#TOC" 227 >--------------------------------</a><a name="HASHfactor">-</a> 228 229 qh_HASHfactor 230 total hash slots / used hash slots. Must be at least 1.1. 231 232 notes: 233 =2 for at worst 50% occupancy for qh hash_table and normally 25% occupancy 234 */ 235 #define qh_HASHfactor 2 236 237 /*-<a href="qh-user.htm#TOC" 238 >--------------------------------</a><a name="VERIFYdirect">-</a> 239 240 qh_VERIFYdirect 241 with 'Tv' verify all points against all facets if op count is smaller 242 243 notes: 244 if greater, calls qh_check_bestdist() instead 245 */ 246 #define qh_VERIFYdirect 1000000 247 248 /*-<a href="qh-user.htm#TOC" 249 >--------------------------------</a><a name="INITIALsearch">-</a> 250 251 qh_INITIALsearch 252 if qh_INITIALmax, search points up to this dimension 253 */ 254 #define qh_INITIALsearch 6 255 256 /*-<a href="qh-user.htm#TOC" 257 >--------------------------------</a><a name="INITIALmax">-</a> 258 259 qh_INITIALmax 260 if dim >= qh_INITIALmax, use min/max coordinate points for initial simplex 261 262 notes: 263 from points with non-zero determinants 264 use option 'Qs' to override (much slower) 265 */ 266 #define qh_INITIALmax 8 267 268 /*-<a href="qh-user.htm#TOC" 269 >--------------------------------</a><a name="JOGGLEdefault">-</a> 270 271 qh_JOGGLEdefault 272 default qh.JOGGLEmax is qh.DISTround * qh_JOGGLEdefault 273 274 notes: 275 rbox s r 100 | qhull QJ1e-15 QR0 generates 90% faults at distround 7e-16 276 rbox s r 100 | qhull QJ1e-14 QR0 generates 70% faults 277 rbox s r 100 | qhull QJ1e-13 QR0 generates 35% faults 278 rbox s r 100 | qhull QJ1e-12 QR0 generates 8% faults 279 rbox s r 100 | qhull QJ1e-11 QR0 generates 1% faults 280 rbox s r 100 | qhull QJ1e-10 QR0 generates 0% faults 281 rbox 1000 W0 | qhull QJ1e-12 QR0 generates 86% faults 282 rbox 1000 W0 | qhull QJ1e-11 QR0 generates 20% faults 283 rbox 1000 W0 | qhull QJ1e-10 QR0 generates 2% faults 284 the later have about 20 points per facet, each of which may interfere 285 286 pick a value large enough to avoid retries on most inputs 287 */ 288 #define qh_JOGGLEdefault 30000.0 289 290 /*-<a href="qh-user.htm#TOC" 291 >--------------------------------</a><a name="JOGGLEincrease">-</a> 292 293 qh_JOGGLEincrease 294 factor to increase qh.JOGGLEmax on qh_JOGGLEretry or qh_JOGGLEagain 295 */ 296 #define qh_JOGGLEincrease 10.0 297 298 /*-<a href="qh-user.htm#TOC" 299 >--------------------------------</a><a name="JOGGLEretry">-</a> 300 301 qh_JOGGLEretry 302 if ZZretry = qh_JOGGLEretry, increase qh.JOGGLEmax 303 304 notes: 305 try twice at the original value in case of bad luck the first time 306 */ 307 #define qh_JOGGLEretry 2 308 309 /*-<a href="qh-user.htm#TOC" 310 >--------------------------------</a><a name="JOGGLEagain">-</a> 311 312 qh_JOGGLEagain 313 every following qh_JOGGLEagain, increase qh.JOGGLEmax 314 315 notes: 316 1 is OK since it's already failed qh_JOGGLEretry times 317 */ 318 #define qh_JOGGLEagain 1 319 320 /*-<a href="qh-user.htm#TOC" 321 >--------------------------------</a><a name="JOGGLEmaxincrease">-</a> 322 323 qh_JOGGLEmaxincrease 324 maximum qh.JOGGLEmax due to qh_JOGGLEincrease 325 relative to qh.MAXwidth 326 327 notes: 328 qh.joggleinput will retry at this value until qh_JOGGLEmaxretry 329 */ 330 #define qh_JOGGLEmaxincrease 1e-2 331 332 /*-<a href="qh-user.htm#TOC" 333 >--------------------------------</a><a name="JOGGLEmaxretry">-</a> 334 335 qh_JOGGLEmaxretry 336 stop after qh_JOGGLEmaxretry attempts 337 */ 338 #define qh_JOGGLEmaxretry 100 339 340 /*========= memory constants =========*/ 341 342 /*-<a href="qh-user.htm#TOC" 343 >--------------------------------</a><a name="MEMalign">-</a> 344 345 qh_MEMalign 346 memory alignment for qh_meminitbuffers() in global.c 347 348 notes: 349 to avoid bus errors, memory allocation must consider alignment requirements. 350 malloc() automatically takes care of alignment. Since mem.c manages 351 its own memory, we need to explicitly specify alignment in 352 qh_meminitbuffers(). 353 354 A safe choice is sizeof(double). sizeof(float) may be used if doubles 355 do not occur in data structures and pointers are the same size. Be careful 356 of machines (e.g., DEC Alpha) with large pointers. 357 358 If using gcc, best alignment is 359 #define qh_MEMalign fmax_(__alignof__(realT),__alignof__(void *)) 360 */ 361 #define qh_MEMalign fmax_(sizeof(realT), sizeof(void *)) 362 363 /*-<a href="qh-user.htm#TOC" 364 >--------------------------------</a><a name="MEMbufsize">-</a> 365 366 qh_MEMbufsize 367 size of additional memory buffers 368 369 notes: 370 used for qh_meminitbuffers() in global.c 371 */ 372 #define qh_MEMbufsize 0x10000 /* allocate 64K memory buffers */ 373 374 /*-<a href="qh-user.htm#TOC" 375 >--------------------------------</a><a name="MEMinitbuf">-</a> 376 377 qh_MEMinitbuf 378 size of initial memory buffer 379 380 notes: 381 use for qh_meminitbuffers() in global.c 382 */ 383 #define qh_MEMinitbuf 0x20000 /* initially allocate 128K buffer */ 384 385 /*-<a href="qh-user.htm#TOC" 386 >--------------------------------</a><a name="INFINITE">-</a> 387 388 qh_INFINITE 389 on output, indicates Voronoi center at infinity 390 */ 391 #define qh_INFINITE -10.101 392 393 /*-<a href="qh-user.htm#TOC" 394 >--------------------------------</a><a name="DEFAULTbox">-</a> 395 396 qh_DEFAULTbox 397 default box size (Geomview expects 0.5) 398 */ 399 #define qh_DEFAULTbox 0.5 400 401 /*======= conditional compilation ============================*/ 402 403 /*-<a href="qh-user.htm#TOC" 404 >--------------------------------</a><a name="compiler">-</a> 405 406 __cplusplus 407 defined by C++ compilers 408 409 __MSC_VER 410 defined by Microsoft Visual C++ 411 412 __MWERKS__ && __POWERPC__ 413 defined by Metrowerks when compiling for the Power Macintosh 414 415 __STDC__ 416 defined for strict ANSI C 417 */ 418 419 /*-<a href="qh-user.htm#TOC" 420 >--------------------------------</a><a name="COMPUTEfurthest">-</a> 421 422 qh_COMPUTEfurthest 423 compute furthest distance to an outside point instead of storing it with the facet 424 =1 to compute furthest 425 426 notes: 427 computing furthest saves memory but costs time 428 about 40% more distance tests for partitioning 429 removes facet->furthestdist 430 */ 431 #define qh_COMPUTEfurthest 0 432 433 /*-<a href="qh-user.htm#TOC" 434 >--------------------------------</a><a name="KEEPstatistics">-</a> 435 436 qh_KEEPstatistics 437 =0 removes most of statistic gathering and reporting 438 439 notes: 440 if 0, code size is reduced by about 4%. 441 */ 442 #define qh_KEEPstatistics 1 443 444 /*-<a href="qh-user.htm#TOC" 445 >--------------------------------</a><a name="MAXoutside">-</a> 446 447 qh_MAXoutside 448 record outer plane for each facet 449 =1 to record facet->maxoutside 450 451 notes: 452 this takes a realT per facet and slightly slows down qhull 453 it produces better outer planes for geomview output 454 */ 455 #define qh_MAXoutside 1 456 457 /*-<a href="qh-user.htm#TOC" 458 >--------------------------------</a><a name="NOmerge">-</a> 459 460 qh_NOmerge 461 disables facet merging if defined 462 463 notes: 464 This saves about 10% space. 465 466 Unless 'Q0' 467 qh_NOmerge sets 'QJ' to avoid precision errors 468 469 #define qh_NOmerge 470 471 see: 472 <a href="mem.h#NOmem">qh_NOmem</a> in mem.c 473 474 see user.c/user_eg.c for removing io.o 475 */ 476 477 /*-<a href="qh-user.htm#TOC" 478 >--------------------------------</a><a name="NOtrace">-</a> 479 480 qh_NOtrace 481 no tracing if defined 482 483 notes: 484 This saves about 5% space. 485 486 #define qh_NOtrace 487 */ 488 489 /*-<a href="qh-user.htm#TOC" 490 >--------------------------------</a><a name="QHpointer">-</a> 491 492 qh_QHpointer 493 access global data with pointer or static structure 494 495 qh_QHpointer = 1 access globals via a pointer to allocated memory 496 enables qh_saveqhull() and qh_restoreqhull() 497 costs about 8% in time and 2% in space 498 499 = 0 qh_qh and qh_qhstat are static data structures 500 only one instance of qhull() can be active at a time 501 default value 502 503 notes: 504 all global variables for qhull are in qh, qhmem, and qhstat 505 qh is defined in qhull.h 506 qhmem is defined in mem.h 507 qhstat is defined in stat.h 508 509 see: 510 user_eg.c for an example 511 */ 512 #define qh_QHpointer 0 513 #if 0 /* sample code */ 514 qhT *oldqhA, *oldqhB; 515 516 exitcode= qh_new_qhull (dim, numpoints, points, ismalloc, 517 flags, outfile, errfile); 518 /* use results from first call to qh_new_qhull */ 519 oldqhA= qh_save_qhull(); 520 exitcode= qh_new_qhull (dimB, numpointsB, pointsB, ismalloc, 521 flags, outfile, errfile); 522 /* use results from second call to qh_new_qhull */ 523 oldqhB= qh_save_qhull(); 524 qh_restore_qhull (&oldqhA); 525 /* use results from first call to qh_new_qhull */ 526 qh_freeqhull (qh_ALL); /* frees all memory used by first call */ 527 qh_restore_qhull (&oldqhB); 528 /* use results from second call to qh_new_qhull */ 529 qh_freeqhull (!qh_ALL); /* frees long memory used by second call */ 530 qh_memfreeshort (&curlong, &totlong); /* frees short memory and memory allocator */ 531 #endif 532 533 /*-<a href="qh-user.htm#TOC" 534 >--------------------------------</a><a name="QUICKhelp">-</a> 535 536 qh_QUICKhelp 537 =1 to use abbreviated help messages, e.g., for degenerate inputs 538 */ 539 #define qh_QUICKhelp 0 540 541 /* ============ -merge constants- ==================== 542 543 These constants effect facet merging. You probably will not need 544 to modify these. They effect the performance of facet merging. 545 */ 546 547 /*-<a href="qh-user.htm#TOC" 548 >--------------------------------</a><a name="DIMmergeVertex">-</a> 549 550 qh_DIMmergeVertex 551 max dimension for vertex merging (it is not effective in high-d) 552 */ 553 #define qh_DIMmergeVertex 6 554 555 /*-<a href="qh-user.htm#TOC" 556 >--------------------------------</a><a name="DIMreduceBuild">-</a> 557 558 qh_DIMreduceBuild 559 max dimension for vertex reduction during build (slow in high-d) 560 */ 561 #define qh_DIMreduceBuild 5 562 563 /*-<a href="qh-user.htm#TOC" 564 >--------------------------------</a><a name="BESTcentrum">-</a> 565 566 qh_BESTcentrum 567 if > 2*dim+n vertices, qh_findbestneighbor() tests centrums (faster) 568 else, qh_findbestneighbor() tests all vertices (much better merges) 569 570 qh_BESTcentrum2 571 if qh_BESTcentrum2 * DIM3 + BESTcentrum < #vertices tests centrums 572 */ 573 #define qh_BESTcentrum 20 574 #define qh_BESTcentrum2 2 575 576 /*-<a href="qh-user.htm#TOC" 577 >--------------------------------</a><a name="BESTnonconvex">-</a> 578 579 qh_BESTnonconvex 580 if > dim+n neighbors, qh_findbestneighbor() tests nonconvex ridges. 581 582 notes: 583 It is needed because qh_findbestneighbor is slow for large facets 584 */ 585 #define qh_BESTnonconvex 15 586 587 /*-<a href="qh-user.htm#TOC" 588 >--------------------------------</a><a name="MAXnewmerges">-</a> 589 590 qh_MAXnewmerges 591 if >n newmerges, qh_merge_nonconvex() calls qh_reducevertices_centrums. 592 593 notes: 594 It is needed because postmerge can merge many facets at once 595 */ 596 #define qh_MAXnewmerges 2 597 598 /*-<a href="qh-user.htm#TOC" 599 >--------------------------------</a><a name="MAXnewcentrum">-</a> 600 601 qh_MAXnewcentrum 602 if <= dim+n vertices (n approximates the number of merges), 603 reset the centrum in qh_updatetested() and qh_mergecycle_facets() 604 605 notes: 606 needed to reduce cost and because centrums may move too much if 607 many vertices in high-d 608 */ 609 #define qh_MAXnewcentrum 5 610 611 /*-<a href="qh-user.htm#TOC" 612 >--------------------------------</a><a name="COPLANARratio">-</a> 613 614 qh_COPLANARratio 615 for 3-d+ merging, qh.MINvisible is n*premerge_centrum 616 617 notes: 618 for non-merging, it's DISTround 619 */ 620 #define qh_COPLANARratio 3 621 622 /*-<a href="qh-user.htm#TOC" 623 >--------------------------------</a><a name="DISToutside">-</a> 624 625 qh_DISToutside 626 When is a point clearly outside of a facet? 627 Stops search in qh_findbestnew or qh_partitionall 628 qh_findbest uses qh.MINoutside since since it is only called if no merges. 629 630 notes: 631 'Qf' always searches for best facet 632 if !qh.MERGING, same as qh.MINoutside. 633 if qh_USEfindbestnew, increase value since neighboring facets may be ill-behaved 634 [Note: Zdelvertextot occurs normally with interior points] 635 RBOX 1000 s Z1 G1e-13 t1001188774 | QHULL Tv 636 When there is a sharp edge, need to move points to a 637 clearly good facet; otherwise may be lost in another partitioning. 638 if too big then O(n^2) behavior for partitioning in cone 639 if very small then important points not processed 640 Needed in qh_partitionall for 641 RBOX 1000 s Z1 G1e-13 t1001032651 | QHULL Tv 642 Needed in qh_findbestnew for many instances of 643 RBOX 1000 s Z1 G1e-13 t | QHULL Tv 644 645 See: 646 qh_DISToutside -- when is a point clearly outside of a facet 647 qh_SEARCHdist -- when is facet coplanar with the best facet? 648 qh_USEfindbestnew -- when to use qh_findbestnew for qh_partitionpoint() 649 */ 650 #define qh_DISToutside ((qh_USEfindbestnew ? 2 : 1) * \ 651 fmax_((qh MERGING ? 2 : 1)*qh MINoutside, qh max_outside)) 652 653 /*-<a href="qh-user.htm#TOC" 654 >--------------------------------</a><a name="RATIOnearinside">-</a> 655 656 qh_RATIOnearinside 657 ratio of qh.NEARinside to qh.ONEmerge for retaining inside points for 658 qh_check_maxout(). 659 660 notes: 661 This is overkill since do not know the correct value. 662 It effects whether 'Qc' reports all coplanar points 663 Not used for 'd' since non-extreme points are coplanar 664 */ 665 #define qh_RATIOnearinside 5 666 667 /*-<a href="qh-user.htm#TOC" 668 >--------------------------------</a><a name="SEARCHdist">-</a> 669 670 qh_SEARCHdist 671 When is a facet coplanar with the best facet? 672 qh_findbesthorizon: all coplanar facets of the best facet need to be searched. 673 674 See: 675 qh_DISToutside -- when is a point clearly outside of a facet 676 qh_SEARCHdist -- when is facet coplanar with the best facet? 677 qh_USEfindbestnew -- when to use qh_findbestnew for qh_partitionpoint() 678 */ 679 #define qh_SEARCHdist ((qh_USEfindbestnew ? 2 : 1) * \ 680 (qh max_outside + 2 * qh DISTround + fmax_( qh MINvisible, qh MAXcoplanar))); 681 682 /*-<a href="qh-user.htm#TOC" 683 >--------------------------------</a><a name="USEfindbestnew">-</a> 684 685 qh_USEfindbestnew 686 Always use qh_findbestnew for qh_partitionpoint, otherwise use 687 qh_findbestnew if merged new facet or sharpnewfacets. 688 689 See: 690 qh_DISToutside -- when is a point clearly outside of a facet 691 qh_SEARCHdist -- when is facet coplanar with the best facet? 692 qh_USEfindbestnew -- when to use qh_findbestnew for qh_partitionpoint() 693 */ 694 #define qh_USEfindbestnew (zzval_(Ztotmerge) > 50) 695 696 /*-<a href="qh-user.htm#TOC" 697 >--------------------------------</a><a name="WIDEcoplanar">-</a> 698 699 qh_WIDEcoplanar 700 n*MAXcoplanar or n*MINvisible for a WIDEfacet 701 702 if vertex is further than qh.WIDEfacet from the hyperplane 703 then its ridges are not counted in computing the area, and 704 the facet's centrum is frozen. 705 706 notes: 707 qh.WIDEfacet= max(qh.MAXoutside,qh_WIDEcoplanar*qh.MAXcoplanar, 708 qh_WIDEcoplanar * qh.MINvisible); 709 */ 710 #define qh_WIDEcoplanar 6 711 712 /*-<a href="qh-user.htm#TOC" 713 >--------------------------------</a><a name="MAXnarrow">-</a> 714 715 qh_MAXnarrow 716 max. cosine in initial hull that sets qh.NARROWhull 717 718 notes: 719 If qh.NARROWhull, the initial partition does not make 720 coplanar points. If narrow, a coplanar point can be 721 coplanar to two facets of opposite orientations and 722 distant from the exact convex hull. 723 724 Conservative estimate. Don't actually see problems until it is -1.0 725 */ 726 #define qh_MAXnarrow -0.99999999 727 728 /*-<a href="qh-user.htm#TOC" 729 >--------------------------------</a><a name="WARNnarrow">-</a> 730 731 qh_WARNnarrow 732 max. cosine in initial hull to warn about qh.NARROWhull 733 734 notes: 735 this is a conservative estimate. 736 Don't actually see problems until it is -1.0. See qh-impre.htm 737 */ 738 #define qh_WARNnarrow -0.999999999999999 739 740 /*-<a href="qh-user.htm#TOC" 741 >--------------------------------</a><a name="ZEROdelaunay">-</a> 742 743 qh_ZEROdelaunay 744 a zero Delaunay facet occurs for input sites coplanar with their convex hull 745 the last normal coefficient of a zero Delaunay facet is within 746 qh_ZEROdelaunay * qh.ANGLEround of 0 747 748 notes: 749 qh_ZEROdelaunay does not allow for joggled input ('QJ'). 750 751 You can avoid zero Delaunay facets by surrounding the input with a box. 752 753 Use option 'PDk:-n' to explicitly define zero Delaunay facets 754 k= dimension of input sites (e.g., 3 for 3-d Delaunay triangulation) 755 n= the cutoff for zero Delaunay facets (e.g., 'PD3:-1e-12') 756 */ 757 #define qh_ZEROdelaunay 2 758 759 #endif /* qh_DEFuser */ 760 761 762 763