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