1 /*<html><pre> -<a href="qh-merge.htm#TOC"
2 >-------------------------------</a><a name="TOP">-</a>
3
4 merge.c
5 merges non-convex facets
6
7 see qh-merge.htm and merge.h
8
9 other modules call qh_premerge() and qh_postmerge()
10
11 the user may call qh_postmerge() to perform additional merges.
12
13 To remove deleted facets and vertices (qhull() in libqhull.c):
14 qh_partitionvisible(!qh_ALL, &numoutside); // visible_list, newfacet_list
15 qh_deletevisible(); // qh.visible_list
16 qh_resetlists(False, qh_RESETvisible); // qh.visible_list newvertex_list newfacet_list
17
18 assumes qh.CENTERtype= centrum
19
20 merges occur in qh_mergefacet and in qh_mergecycle
21 vertex->neighbors not set until the first merge occurs
22
23 Copyright (c) 1993-2019 C.B. Barber.
24 $Id: //main/2019/qhull/src/libqhull/merge.c#9 $$Change: 2712 $
25 $DateTime: 2019/06/28 12:57:00 $$Author: bbarber $
26 */
27
28 #include "qhull_a.h"
29
30 #ifndef qh_NOmerge
31
32 /* MRGnone, etc. */
33 const char *mergetypes[]= {
34 "none",
35 "coplanar",
36 "anglecoplanar",
37 "concave",
38 "concavecoplanar",
39 "twisted",
40 "flip",
41 "dupridge",
42 "subridge",
43 "vertices",
44 "degen",
45 "redundant",
46 "mirror",
47 "coplanarhorizon",
48 };
49
50 /*===== functions(alphabetical after premerge and postmerge) ======*/
51
52 /*-<a href="qh-merge.htm#TOC"
53 >-------------------------------</a><a name="premerge">-</a>
54
55 qh_premerge( apexpointid, maxcentrum )
56 pre-merge nonconvex facets in qh.newfacet_list for apexpointid
57 maxcentrum defines coplanar and concave (qh_test_appendmerge)
58
59 returns:
60 deleted facets added to qh.visible_list with facet->visible set
61
62 notes:
63 only called by qh_addpoint
64 uses globals, qh.MERGEexact, qh.PREmerge
65
66 design:
67 mark dupridges in qh.newfacet_list
68 merge facet cycles in qh.newfacet_list
69 merge dupridges and concave facets in qh.newfacet_list
70 check merged facet cycles for degenerate and redundant facets
71 merge degenerate and redundant facets
72 collect coplanar and concave facets
73 merge concave, coplanar, degenerate, and redundant facets
74 */
qh_premerge(int apexpointid,realT maxcentrum,realT maxangle)75 void qh_premerge(int apexpointid, realT maxcentrum, realT maxangle /* qh.newfacet_list */) {
76 boolT othermerge= False;
77
78 if (qh ZEROcentrum && qh_checkzero(!qh_ALL))
79 return;
80 trace2((qh ferr, 2008, "qh_premerge: premerge centrum %2.2g angle %4.4g for apex p%d newfacet_list f%d\n",
81 maxcentrum, maxangle, apexpointid, getid_(qh newfacet_list)));
82 if (qh IStracing >= 4 && qh num_facets < 100)
83 qh_printlists();
84 qh centrum_radius= maxcentrum;
85 qh cos_max= maxangle;
86 if (qh hull_dim >=3) {
87 qh_mark_dupridges(qh newfacet_list, qh_ALL); /* facet_mergeset */
88 qh_mergecycle_all(qh newfacet_list, &othermerge);
89 qh_forcedmerges(&othermerge /* qh.facet_mergeset */);
90 }else /* qh.hull_dim == 2 */
91 qh_mergecycle_all(qh newfacet_list, &othermerge);
92 qh_flippedmerges(qh newfacet_list, &othermerge);
93 if (!qh MERGEexact || zzval_(Ztotmerge)) {
94 zinc_(Zpremergetot);
95 qh POSTmerging= False;
96 qh_getmergeset_initial(qh newfacet_list);
97 qh_all_merges(othermerge, False);
98 }
99 } /* premerge */
100
101 /*-<a href="qh-merge.htm#TOC"
102 >-------------------------------</a><a name="postmerge">-</a>
103
104 qh_postmerge( reason, maxcentrum, maxangle, vneighbors )
105 post-merge nonconvex facets as defined by maxcentrum and maxangle
106 'reason' is for reporting progress
107 if vneighbors ('Qv'),
108 calls qh_test_vneighbors at end of qh_all_merge from qh_postmerge
109
110 returns:
111 if first call (qh.visible_list != qh.facet_list),
112 builds qh.facet_newlist, qh.newvertex_list
113 deleted facets added to qh.visible_list with facet->visible
114 qh.visible_list == qh.facet_list
115
116 notes:
117 called by qh_qhull after qh_buildhull
118 called if a merge may be needed due to
119 qh.MERGEexact ('Qx'), qh_DIMreduceBuild, POSTmerge (e.g., 'Cn'), or TESTvneighbors ('Qv')
120 if firstmerge,
121 calls qh_reducevertices before qh_getmergeset
122
123 design:
124 if first call
125 set qh.visible_list and qh.newfacet_list to qh.facet_list
126 add all facets to qh.newfacet_list
127 mark non-simplicial facets, facet->newmerge
128 set qh.newvertext_list to qh.vertex_list
129 add all vertices to qh.newvertex_list
130 if a pre-merge occurred
131 set vertex->delridge {will retest the ridge}
132 if qh.MERGEexact
133 call qh_reducevertices()
134 if no pre-merging
135 merge flipped facets
136 determine non-convex facets
137 merge all non-convex facets
138 */
qh_postmerge(const char * reason,realT maxcentrum,realT maxangle,boolT vneighbors)139 void qh_postmerge(const char *reason, realT maxcentrum, realT maxangle,
140 boolT vneighbors) {
141 facetT *newfacet;
142 boolT othermerges= False;
143 vertexT *vertex;
144
145 if (qh REPORTfreq || qh IStracing) {
146 qh_buildtracing(NULL, NULL);
147 qh_printsummary(qh ferr);
148 if (qh PRINTstatistics)
149 qh_printallstatistics(qh ferr, "reason");
150 qh_fprintf(qh ferr, 8062, "\n%s with 'C%.2g' and 'A%.2g'\n",
151 reason, maxcentrum, maxangle);
152 }
153 trace2((qh ferr, 2009, "qh_postmerge: postmerge. test vneighbors? %d\n",
154 vneighbors));
155 qh centrum_radius= maxcentrum;
156 qh cos_max= maxangle;
157 qh POSTmerging= True;
158 if (qh visible_list != qh facet_list) { /* first call due to qh_buildhull, multiple calls if qh.POSTmerge */
159 qh NEWfacets= True;
160 qh visible_list= qh newfacet_list= qh facet_list;
161 FORALLnew_facets { /* all facets are new facets for qh_postmerge */
162 newfacet->newfacet= True;
163 if (!newfacet->simplicial)
164 newfacet->newmerge= True; /* test f.vertices for 'delridge'. 'newmerge' was cleared at end of qh_all_merges */
165 zinc_(Zpostfacets);
166 }
167 qh newvertex_list= qh vertex_list;
168 FORALLvertices
169 vertex->newfacet= True;
170 if (qh VERTEXneighbors) { /* a merge has occurred */
171 if (qh MERGEexact && qh hull_dim <= qh_DIMreduceBuild)
172 qh_reducevertices(); /* qh_all_merges did not call qh_reducevertices for v.delridge */
173 }
174 if (!qh PREmerge && !qh MERGEexact)
175 qh_flippedmerges(qh newfacet_list, &othermerges);
176 }
177 qh_getmergeset_initial(qh newfacet_list);
178 qh_all_merges(False, vneighbors); /* calls qh_reducevertices before exiting */
179 FORALLnew_facets
180 newfacet->newmerge= False; /* Was True if no vertex in f.vertices was 'delridge' */
181 } /* post_merge */
182
183 /*-<a href="qh-merge.htm#TOC"
184 >-------------------------------</a><a name="all_merges">-</a>
185
186 qh_all_merges( othermerge, vneighbors )
187 merge all non-convex facets
188
189 set othermerge if already merged facets (calls qh_reducevertices)
190 if vneighbors ('Qv' at qh.POSTmerge)
191 tests vertex neighbors for convexity at end (qh_test_vneighbors)
192 qh.facet_mergeset lists the non-convex ridges in qh_newfacet_list
193 qh.degen_mergeset is defined
194 if qh.MERGEexact && !qh.POSTmerging,
195 does not merge coplanar facets
196
197 returns:
198 deleted facets added to qh.visible_list with facet->visible
199 deleted vertices added qh.delvertex_list with vertex->delvertex
200
201 notes:
202 unless !qh.MERGEindependent,
203 merges facets in independent sets
204 uses qh.newfacet_list as implicit argument since merges call qh_removefacet()
205 [apr'19] restored qh_setdellast in place of qh_next_facetmerge. Much faster for post-merge
206
207 design:
208 while merges occur
209 for each merge in qh.facet_mergeset
210 unless one of the facets was already merged in this pass
211 merge the facets
212 test merged facets for additional merges
213 add merges to qh.facet_mergeset
214 if qh.POSTmerging
215 periodically call qh_reducevertices to reduce extra vertices and redundant vertices
216 after each pass, if qh.VERTEXneighbors
217 if qh.POSTmerging or was a merge with qh.hull_dim<=5
218 call qh_reducevertices
219 update qh.facet_mergeset if degenredundant merges
220 if 'Qv' and qh.POSTmerging
221 test vertex neighbors for convexity
222 */
qh_all_merges(boolT othermerge,boolT vneighbors)223 void qh_all_merges(boolT othermerge, boolT vneighbors) {
224 facetT *facet1, *facet2, *newfacet;
225 mergeT *merge;
226 boolT wasmerge= False, isreduce;
227 void **freelistp; /* used if !qh_NOmem by qh_memfree_() */
228 vertexT *vertex;
229 realT angle, distance;
230 mergeType mergetype;
231 int numcoplanar=0, numconcave=0, numconcavecoplanar= 0, numdegenredun= 0, numnewmerges= 0, numtwisted= 0;
232
233 trace2((qh ferr, 2010, "qh_all_merges: starting to merge %d facet and %d degenerate merges for new facets f%d, othermerge? %d\n",
234 qh_setsize(qh facet_mergeset), qh_setsize(qh degen_mergeset), getid_(qh newfacet_list), othermerge));
235
236 while (True) {
237 wasmerge= False;
238 while (qh_setsize(qh facet_mergeset) > 0 || qh_setsize(qh degen_mergeset) > 0) {
239 if (qh_setsize(qh degen_mergeset) > 0) {
240 numdegenredun += qh_merge_degenredundant();
241 wasmerge= True;
242 }
243 while ((merge= (mergeT *)qh_setdellast(qh facet_mergeset))) {
244 facet1= merge->facet1;
245 facet2= merge->facet2;
246 vertex= merge->vertex1; /* not used for qh.facet_mergeset*/
247 mergetype= merge->mergetype;
248 angle= merge->angle;
249 distance= merge->distance;
250 qh_memfree_(merge, (int)sizeof(mergeT), freelistp); /* 'merge' is invalid */
251 if (facet1->visible || facet2->visible) {
252 trace3((qh ferr, 3045, "qh_all_merges: drop merge of f%d (del? %d) into f%d (del? %d) mergetype %d, dist %4.4g, angle %4.4g. One or both facets is deleted\n",
253 facet1->id, facet1->visible, facet2->id, facet2->visible, mergetype, distance, angle));
254 continue;
255 }else if (mergetype == MRGcoplanar || mergetype == MRGanglecoplanar) {
256 if (qh MERGEindependent) {
257 if ((!facet1->tested && facet1->newfacet)
258 || (!facet2->tested && facet2->newfacet)) {
259 trace3((qh ferr, 3064, "qh_all_merges: drop merge of f%d (tested? %d) into f%d (tested? %d) mergetype %d, dist %2.2g, angle %4.4g. Merge independent sets of coplanar merges\n",
260 facet1->id, facet1->visible, facet2->id, facet2->visible, mergetype, distance, angle));
261 continue;
262 }
263 }
264 }
265 trace3((qh ferr, 3047, "qh_all_merges: merge f%d and f%d type %d dist %2.2g angle %4.4g\n",
266 facet1->id, facet2->id, mergetype, distance, angle));
267 if (mergetype == MRGtwisted)
268 qh_merge_twisted(facet1, facet2);
269 else
270 qh_merge_nonconvex(facet1, facet2, mergetype);
271 numnewmerges++;
272 numdegenredun += qh_merge_degenredundant();
273 wasmerge= True;
274 if (mergetype == MRGconcave)
275 numconcave++;
276 else if (mergetype == MRGconcavecoplanar)
277 numconcavecoplanar++;
278 else if (mergetype == MRGtwisted)
279 numtwisted++;
280 else if (mergetype == MRGcoplanar || mergetype == MRGanglecoplanar)
281 numcoplanar++;
282 else {
283 qh_fprintf(qh ferr, 6394, "qhull internal error (qh_all_merges): expecting concave, coplanar, or twisted merge. Got merge f%d f%d v%d mergetype %d\n",
284 getid_(facet1), getid_(facet2), getid_(vertex), mergetype);
285 qh_errexit2(qh_ERRqhull, facet1, facet2);
286 }
287 } /* while qh_setdellast */
288 if (qh POSTmerging && qh hull_dim <= qh_DIMreduceBuild
289 && numnewmerges > qh_MAXnewmerges) {
290 numnewmerges= 0;
291 wasmerge= othermerge= False;
292 qh_reducevertices(); /* otherwise large post merges too slow */
293 }
294 qh_getmergeset(qh newfacet_list); /* qh.facet_mergeset */
295 } /* while facet_mergeset or degen_mergeset */
296 if (qh VERTEXneighbors) { /* at least one merge */
297 isreduce= False;
298 if (qh POSTmerging && qh hull_dim >= 4) {
299 isreduce= True;
300 }else if (qh POSTmerging || !qh MERGEexact) {
301 if ((wasmerge || othermerge) && qh hull_dim > 2 && qh hull_dim <= qh_DIMreduceBuild)
302 isreduce= True;
303 }
304 if (isreduce) {
305 wasmerge= othermerge= False;
306 if (qh_reducevertices()) {
307 qh_getmergeset(qh newfacet_list); /* facet_mergeset */
308 continue;
309 }
310 }
311 }
312 if (vneighbors && qh_test_vneighbors(/* qh.newfacet_list */))
313 continue;
314 break;
315 } /* while (True) */
316 if (wasmerge || othermerge) {
317 trace3((qh ferr, 3033, "qh_all_merges: skip qh_reducevertices due to post-merging, no qh.VERTEXneighbors (%d), or hull_dim %d ==2 or >%d\n", qh VERTEXneighbors, qh hull_dim, qh_DIMreduceBuild))
318 FORALLnew_facets {
319 newfacet->newmerge= False;
320 }
321 }
322 if (qh CHECKfrequently && !qh MERGEexact) {
323 qh old_randomdist= qh RANDOMdist;
324 qh RANDOMdist= False;
325 qh_checkconvex(qh newfacet_list, qh_ALGORITHMfault);
326 /* qh_checkconnect(); [this is slow and it changes the facet order] */
327 qh RANDOMdist= qh old_randomdist;
328 }
329 trace1((qh ferr, 1009, "qh_all_merges: merged %d coplanar %d concave %d concavecoplanar %d twisted facets and %d degen or redundant facets.\n",
330 numcoplanar, numconcave, numconcavecoplanar, numtwisted, numdegenredun));
331 if (qh IStracing >= 4 && qh num_facets < 500)
332 qh_printlists();
333 } /* all_merges */
334
335 /*-<a href="qh-merge.htm#TOC"
336 >-------------------------------</a><a name="all_vertexmerges">-</a>
337
338 qh_all_vertexmerges( apexpointid, facet, &retryfacet )
339 merge vertices in qh.vertex_mergeset and subsequent merges
340
341 returns:
342 returns retryfacet for facet (if defined)
343 updates qh.facet_list, qh.num_facets, qh.vertex_list, qh.num_vertices
344 mergesets are empty
345 if merges, resets facet lists
346
347 notes:
348 called from qh_qhull, qh_addpoint, and qh_buildcone_mergepinched
349 vertex merges occur after facet merges and qh_resetlists
350
351 design:
352 while merges in vertex_mergeset (MRGvertices)
353 merge a pair of pinched vertices
354 update vertex neighbors
355 merge non-convex and degenerate facets and check for ridges with duplicate vertices
356 partition outside points of deleted, "visible" facets
357 */
qh_all_vertexmerges(int apexpointid,facetT * facet,facetT ** retryfacet)358 void qh_all_vertexmerges(int apexpointid, facetT *facet, facetT **retryfacet) {
359 int numpoints; /* ignore count of partitioned points. Used by qh_addpoint for Zpbalance */
360
361 if (retryfacet)
362 *retryfacet= facet;
363 while (qh_setsize(qh vertex_mergeset) > 0) {
364 trace1((qh ferr, 1057, "qh_all_vertexmerges: starting to merge %d vertex merges for apex p%d facet f%d\n",
365 qh_setsize(qh vertex_mergeset), apexpointid, getid_(facet)));
366 if (qh IStracing >= 4 && qh num_facets < 1000)
367 qh_printlists();
368 qh_merge_pinchedvertices(apexpointid /* qh.vertex_mergeset, visible_list, newvertex_list, newfacet_list */);
369 qh_update_vertexneighbors(); /* update neighbors of qh.newvertex_list from qh_newvertices for deleted facets on qh.visible_list */
370 /* test ridges and merge non-convex facets */
371 qh_getmergeset(qh newfacet_list);
372 qh_all_merges(True, False); /* calls qh_reducevertices */
373 if (qh CHECKfrequently)
374 qh_checkpolygon(qh facet_list);
375 qh_partitionvisible(!qh_ALL, &numpoints /* qh.visible_list qh.del_vertices*/);
376 if (retryfacet)
377 *retryfacet= qh_getreplacement(*retryfacet);
378 qh_deletevisible(/* qh.visible_list qh.del_vertices*/);
379 qh_resetlists(False, qh_RESETvisible /* qh.visible_list newvertex_list qh.newfacet_list */);
380 if (qh IStracing >= 4 && qh num_facets < 1000) {
381 qh_printlists();
382 qh_checkpolygon(qh facet_list);
383 }
384 }
385 } /* all_vertexmerges */
386
387 /*-<a href="qh-merge.htm#TOC"
388 >-------------------------------</a><a name="appendmergeset">-</a>
389
390 qh_appendmergeset( facet, vertex, neighbor, mergetype, dist, angle )
391 appends an entry to qh.facet_mergeset or qh.degen_mergeset
392 if 'dist' is unknown, set it to 0.0
393 if 'angle' is unknown, set it to 1.0 (coplanar)
394
395 returns:
396 merge appended to facet_mergeset or degen_mergeset
397 sets ->degenerate or ->redundant if degen_mergeset
398
399 notes:
400 caller collects statistics and/or caller of qh_mergefacet
401 see: qh_test_appendmerge()
402
403 design:
404 allocate merge entry
405 if regular merge
406 append to qh.facet_mergeset
407 else if degenerate merge and qh.facet_mergeset is all degenerate
408 append to qh.degen_mergeset
409 else if degenerate merge
410 prepend to qh.degen_mergeset (merged last)
411 else if redundant merge
412 append to qh.degen_mergeset
413 */
qh_appendmergeset(facetT * facet,facetT * neighbor,mergeType mergetype,coordT dist,realT angle)414 void qh_appendmergeset(facetT *facet, facetT *neighbor, mergeType mergetype, coordT dist, realT angle) {
415 mergeT *merge, *lastmerge;
416 void **freelistp; /* used if !qh_NOmem by qh_memalloc_() */
417 const char *mergename;
418
419 if ((facet->redundant && mergetype != MRGmirror) || neighbor->redundant) {
420 trace3((qh ferr, 3051, "qh_appendmergeset: f%d is already redundant (%d) or f%d is already redundant (%d). Ignore merge f%d and f%d type %d\n",
421 facet->id, facet->redundant, neighbor->id, neighbor->redundant, facet->id, neighbor->id, mergetype));
422 return;
423 }
424 if (facet->degenerate && mergetype == MRGdegen) {
425 trace3((qh ferr, 3077, "qh_appendmergeset: f%d is already degenerate. Ignore merge f%d type %d (MRGdegen)\n",
426 facet->id, facet->id, mergetype));
427 return;
428 }
429 if (!qh facet_mergeset || !qh degen_mergeset) {
430 qh_fprintf(qh ferr, 6403, "qhull internal error (qh_appendmergeset): expecting temp set defined for qh.facet_mergeset (0x%x) and qh.degen_mergeset (0x%x). Got NULL\n",
431 qh facet_mergeset, qh degen_mergeset);
432 /* otherwise qh_setappend creates a new set that is not freed by qh_freebuild() */
433 qh_errexit(qh_ERRqhull, NULL, NULL);
434 }
435 if (neighbor->flipped && !facet->flipped) {
436 if (mergetype != MRGdupridge) {
437 qh_fprintf(qh ferr, 6355, "qhull internal error (qh_appendmergeset): except for MRGdupridge, cannot merge a non-flipped facet f%d into flipped f%d, mergetype %d, dist %4.4g\n",
438 facet->id, neighbor->id, mergetype, dist);
439 qh_errexit(qh_ERRqhull, NULL, NULL);
440 }else {
441 trace2((qh ferr, 2106, "qh_appendmergeset: dupridge will merge a non-flipped facet f%d into flipped f%d, dist %4.4g\n",
442 facet->id, neighbor->id, dist));
443 }
444 }
445 qh_memalloc_((int)sizeof(mergeT), freelistp, merge, mergeT);
446 merge->angle= angle;
447 merge->distance= dist;
448 merge->facet1= facet;
449 merge->facet2= neighbor;
450 merge->vertex1= NULL;
451 merge->vertex2= NULL;
452 merge->ridge1= NULL;
453 merge->ridge2= NULL;
454 merge->mergetype= mergetype;
455 if(mergetype > 0 && mergetype <= sizeof(mergetypes))
456 mergename= mergetypes[mergetype];
457 else
458 mergename= mergetypes[MRGnone];
459 if (mergetype < MRGdegen)
460 qh_setappend(&(qh facet_mergeset), merge);
461 else if (mergetype == MRGdegen) {
462 facet->degenerate= True;
463 if (!(lastmerge= (mergeT *)qh_setlast(qh degen_mergeset))
464 || lastmerge->mergetype == MRGdegen)
465 qh_setappend(&(qh degen_mergeset), merge);
466 else
467 qh_setaddnth(&(qh degen_mergeset), 0, merge); /* merged last */
468 }else if (mergetype == MRGredundant) {
469 facet->redundant= True;
470 qh_setappend(&(qh degen_mergeset), merge);
471 }else /* mergetype == MRGmirror */ {
472 if (facet->redundant || neighbor->redundant) {
473 qh_fprintf(qh ferr, 6092, "qhull internal error (qh_appendmergeset): facet f%d or f%d is already a mirrored facet (i.e., 'redundant')\n",
474 facet->id, neighbor->id);
475 qh_errexit2(qh_ERRqhull, facet, neighbor);
476 }
477 if (!qh_setequal(facet->vertices, neighbor->vertices)) {
478 qh_fprintf(qh ferr, 6093, "qhull internal error (qh_appendmergeset): mirrored facets f%d and f%d do not have the same vertices\n",
479 facet->id, neighbor->id);
480 qh_errexit2(qh_ERRqhull, facet, neighbor);
481 }
482 facet->redundant= True;
483 neighbor->redundant= True;
484 qh_setappend(&(qh degen_mergeset), merge);
485 }
486 if (merge->mergetype >= MRGdegen) {
487 trace3((qh ferr, 3044, "qh_appendmergeset: append merge f%d and f%d type %d (%s) to qh.degen_mergeset (size %d)\n",
488 merge->facet1->id, merge->facet2->id, merge->mergetype, mergename, qh_setsize(qh degen_mergeset)));
489 }else {
490 trace3((qh ferr, 3027, "qh_appendmergeset: append merge f%d and f%d type %d (%s) dist %2.2g angle %4.4g to qh.facet_mergeset (size %d)\n",
491 merge->facet1->id, merge->facet2->id, merge->mergetype, mergename, merge->distance, merge->angle, qh_setsize(qh facet_mergeset)));
492 }
493 } /* appendmergeset */
494
495
496 /*-<a href="qh-merge.htm#TOC"
497 >-------------------------------</a><a name="appendvertexmerge">-</a>
498
499 qh_appendvertexmerge( vertex, vertex2, mergetype, distance, ridge1, ridge2 )
500 appends a vertex merge to qh.vertex_mergeset
501 MRGsubridge includes two ridges (from MRGdupridge)
502 MRGvertices includes two ridges
503
504 notes:
505 called by qh_getpinchedmerges for MRGsubridge
506 called by qh_maybe_duplicateridge and qh_maybe_duplicateridges for MRGvertices
507 only way to add a vertex merge to qh.vertex_mergeset
508 checked by qh_next_vertexmerge
509 */
qh_appendvertexmerge(vertexT * vertex,vertexT * destination,mergeType mergetype,realT distance,ridgeT * ridge1,ridgeT * ridge2)510 void qh_appendvertexmerge(vertexT *vertex, vertexT *destination, mergeType mergetype, realT distance, ridgeT *ridge1, ridgeT *ridge2) {
511 mergeT *merge;
512 void **freelistp; /* used if !qh_NOmem by qh_memalloc_() */
513 const char *mergename;
514
515 if (!qh vertex_mergeset) {
516 qh_fprintf(qh ferr, 6387, "qhull internal error (qh_appendvertexmerge): expecting temp set defined for qh.vertex_mergeset (0x%x). Got NULL\n",
517 qh vertex_mergeset);
518 /* otherwise qh_setappend creates a new set that is not freed by qh_freebuild() */
519 qh_errexit(qh_ERRqhull, NULL, NULL);
520 }
521 qh_memalloc_((int)sizeof(mergeT), freelistp, merge, mergeT);
522 merge->angle= qh_ANGLEnone;
523 merge->distance= distance;
524 merge->facet1= NULL;
525 merge->facet2= NULL;
526 merge->vertex1= vertex;
527 merge->vertex2= destination;
528 merge->ridge1= ridge1;
529 merge->ridge2= ridge2;
530 merge->mergetype= mergetype;
531 if(mergetype > 0 && mergetype <= sizeof(mergetypes))
532 mergename= mergetypes[mergetype];
533 else
534 mergename= mergetypes[MRGnone];
535 if (mergetype == MRGvertices) {
536 if (!ridge1 || !ridge2 || ridge1 == ridge2) {
537 qh_fprintf(qh ferr, 6106, "qhull internal error (qh_appendvertexmerge): expecting two distinct ridges for MRGvertices. Got r%d r%d\n",
538 getid_(ridge1), getid_(ridge2));
539 qh_errexit(qh_ERRqhull, NULL, ridge1);
540 }
541 }
542 qh_setappend(&(qh vertex_mergeset), merge);
543 trace3((qh ferr, 3034, "qh_appendvertexmerge: append merge v%d into v%d r%d r%d dist %2.2g type %d (%s)\n",
544 vertex->id, destination->id, getid_(ridge1), getid_(ridge2), distance, merge->mergetype, mergename));
545 } /* appendvertexmerge */
546
547
548 /*-<a href="qh-merge.htm#TOC"
549 >-------------------------------</a><a name="basevertices">-</a>
550
551 qh_basevertices( samecycle )
552 return temporary set of base vertices for samecycle
553 samecycle is first facet in the cycle
554 assumes apex is SETfirst_( samecycle->vertices )
555
556 returns:
557 vertices(settemp)
558 all ->seen are cleared
559
560 notes:
561 uses qh_vertex_visit;
562
563 design:
564 for each facet in samecycle
565 for each unseen vertex in facet->vertices
566 append to result
567 */
qh_basevertices(facetT * samecycle)568 setT *qh_basevertices(facetT *samecycle) {
569 facetT *same;
570 vertexT *apex, *vertex, **vertexp;
571 setT *vertices= qh_settemp(qh TEMPsize);
572
573 apex= SETfirstt_(samecycle->vertices, vertexT);
574 apex->visitid= ++qh vertex_visit;
575 FORALLsame_cycle_(samecycle) {
576 if (same->mergeridge)
577 continue;
578 FOREACHvertex_(same->vertices) {
579 if (vertex->visitid != qh vertex_visit) {
580 qh_setappend(&vertices, vertex);
581 vertex->visitid= qh vertex_visit;
582 vertex->seen= False;
583 }
584 }
585 }
586 trace4((qh ferr, 4019, "qh_basevertices: found %d vertices\n",
587 qh_setsize(vertices)));
588 return vertices;
589 } /* basevertices */
590
591 /*-<a href="qh-merge.htm#TOC"
592 >-------------------------------</a><a name="check_dupridge">-</a>
593
594 qh_check_dupridge( facet1, dist1, facet2, dist2 )
595 Check dupridge between facet1 and facet2 for wide merge
596 dist1 is the maximum distance of facet1's vertices to facet2
597 dist2 is the maximum distance of facet2's vertices to facet1
598
599 returns
600 Level 1 log of the dupridge with the minimum distance between vertices
601 Throws error if the merge will increase the maximum facet width by qh_WIDEduplicate (100x)
602
603 notes:
604 only called from qh_forcedmerges
605 */
qh_check_dupridge(facetT * facet1,realT dist1,facetT * facet2,realT dist2)606 void qh_check_dupridge(facetT *facet1, realT dist1, facetT *facet2, realT dist2) {
607 vertexT *vertex, **vertexp, *vertexA, **vertexAp;
608 realT dist, innerplane, mergedist, outerplane, prevdist, ratio, vertexratio;
609 realT minvertex= REALmax;
610
611 mergedist= fmin_(dist1, dist2);
612 qh_outerinner(NULL, &outerplane, &innerplane); /* ratio from qh_printsummary */
613 FOREACHvertex_(facet1->vertices) { /* The dupridge is between facet1 and facet2, so either facet can be tested */
614 FOREACHvertexA_(facet1->vertices) {
615 if (vertex > vertexA){ /* Test each pair once */
616 dist= qh_pointdist(vertex->point, vertexA->point, qh hull_dim);
617 minimize_(minvertex, dist);
618 /* Not quite correct. A facet may have a dupridge and another pair of nearly adjacent vertices. */
619 }
620 }
621 }
622 prevdist= fmax_(outerplane, innerplane);
623 maximize_(prevdist, qh ONEmerge + qh DISTround);
624 maximize_(prevdist, qh MINoutside + qh DISTround);
625 ratio= mergedist/prevdist;
626 vertexratio= minvertex/prevdist;
627 trace0((qh ferr, 16, "qh_check_dupridge: dupridge between f%d and f%d (vertex dist %2.2g), dist %2.2g, reverse dist %2.2g, ratio %2.2g while processing p%d\n",
628 facet1->id, facet2->id, minvertex, dist1, dist2, ratio, qh furthest_id));
629 if (ratio > qh_WIDEduplicate) {
630 qh_fprintf(qh ferr, 6271, "qhull topology error (qh_check_dupridge): wide merge (%.1fx wider) due to dupridge between f%d and f%d (vertex dist %2.2g), merge dist %2.2g, while processing p%d\n- Allow error with option 'Q12'\n",
631 ratio, facet1->id, facet2->id, minvertex, mergedist, qh furthest_id);
632 if (vertexratio < qh_WIDEpinched)
633 qh_fprintf(qh ferr, 8145, "- Experimental option merge-pinched-vertices ('Q14') may avoid this error. It merges nearly adjacent vertices.\n");
634 if (qh DELAUNAY)
635 qh_fprintf(qh ferr, 8145, "- A bounding box for the input sites may alleviate this error.\n");
636 if (!qh ALLOWwide)
637 qh_errexit2(qh_ERRwide, facet1, facet2);
638 }
639 } /* check_dupridge */
640
641 /*-<a href="qh-merge.htm#TOC"
642 >-------------------------------</a><a name="checkconnect">-</a>
643
644 qh_checkconnect( )
645 check that new facets are connected
646 new facets are on qh.newfacet_list
647
648 notes:
649 this is slow and it changes the order of the facets
650 uses qh.visit_id
651
652 design:
653 move first new facet to end of qh.facet_list
654 for all newly appended facets
655 append unvisited neighbors to end of qh.facet_list
656 for all new facets
657 report error if unvisited
658 */
qh_checkconnect(void)659 void qh_checkconnect(void /* qh.newfacet_list */) {
660 facetT *facet, *newfacet, *errfacet= NULL, *neighbor, **neighborp;
661
662 facet= qh newfacet_list;
663 qh_removefacet(facet);
664 qh_appendfacet(facet);
665 facet->visitid= ++qh visit_id;
666 FORALLfacet_(facet) {
667 FOREACHneighbor_(facet) {
668 if (neighbor->visitid != qh visit_id) {
669 qh_removefacet(neighbor);
670 qh_appendfacet(neighbor);
671 neighbor->visitid= qh visit_id;
672 }
673 }
674 }
675 FORALLnew_facets {
676 if (newfacet->visitid == qh visit_id)
677 break;
678 qh_fprintf(qh ferr, 6094, "qhull internal error (qh_checkconnect): f%d is not attached to the new facets\n",
679 newfacet->id);
680 errfacet= newfacet;
681 }
682 if (errfacet)
683 qh_errexit(qh_ERRqhull, errfacet, NULL);
684 } /* checkconnect */
685
686 /*-<a href="qh-merge.htm#TOC"
687 >-------------------------------</a><a name="checkdelfacet">-</a>
688
689 qh_checkdelfacet( facet, mergeset )
690 check that mergeset does not reference facet
691
692 */
qh_checkdelfacet(facetT * facet,setT * mergeset)693 void qh_checkdelfacet(facetT *facet, setT *mergeset) {
694 mergeT *merge, **mergep;
695
696 FOREACHmerge_(mergeset) {
697 if (merge->facet1 == facet || merge->facet2 == facet) {
698 qh_fprintf(qh ferr, 6390, "qhull internal error (qh_checkdelfacet): cannot delete f%d. It is referenced by merge f%d f%d mergetype %d\n",
699 facet->id, merge->facet1->id, getid_(merge->facet2), merge->mergetype);
700 qh_errexit2(qh_ERRqhull, merge->facet1, merge->facet2);
701 }
702 }
703 } /* checkdelfacet */
704
705 /*-<a href="qh-merge.htm#TOC"
706 >-------------------------------</a><a name="checkdelridge">-</a>
707
708 qh_checkdelridge( )
709 check that qh_delridge_merge is not needed for deleted ridges
710
711 notes:
712 called from qh_mergecycle, qh_makenewfacets, qh_attachnewfacets
713 errors if qh.vertex_mergeset is non-empty
714 errors if any visible or new facet has a ridge with r.nonconvex set
715 assumes that vertex.delfacet is not needed
716 */
qh_checkdelridge(void)717 void qh_checkdelridge(void /* qh.visible_facets, vertex_mergeset */) {
718 facetT *newfacet, *visible;
719 ridgeT *ridge, **ridgep;
720
721 if (!SETempty_(qh vertex_mergeset)) {
722 qh_fprintf(qh ferr, 6382, "qhull internal error (qh_checkdelridge): expecting empty qh.vertex_mergeset in order to avoid calling qh_delridge_merge. Got %d merges\n", qh_setsize(qh vertex_mergeset));
723 qh_errexit(qh_ERRqhull, NULL, NULL);
724 }
725
726 FORALLnew_facets {
727 FOREACHridge_(newfacet->ridges) {
728 if (ridge->nonconvex) {
729 qh_fprintf(qh ferr, 6313, "qhull internal error (qh_checkdelridge): unexpected 'nonconvex' flag for ridge r%d in newfacet f%d. Otherwise need to call qh_delridge_merge\n",
730 ridge->id, newfacet->id);
731 qh_errexit(qh_ERRqhull, newfacet, ridge);
732 }
733 }
734 }
735
736 FORALLvisible_facets {
737 FOREACHridge_(visible->ridges) {
738 if (ridge->nonconvex) {
739 qh_fprintf(qh ferr, 6385, "qhull internal error (qh_checkdelridge): unexpected 'nonconvex' flag for ridge r%d in visible facet f%d. Otherwise need to call qh_delridge_merge\n",
740 ridge->id, visible->id);
741 qh_errexit(qh_ERRqhull, visible, ridge);
742 }
743 }
744 }
745 } /* checkdelridge */
746
747
748 /*-<a href="qh-merge.htm#TOC"
749 >-------------------------------</a><a name="checkzero">-</a>
750
751 qh_checkzero( testall )
752 check that facets are clearly convex for qh.DISTround with qh.MERGEexact
753
754 if testall,
755 test all facets for qh.MERGEexact post-merging
756 else
757 test qh.newfacet_list
758
759 if qh.MERGEexact,
760 allows coplanar ridges
761 skips convexity test while qh.ZEROall_ok
762
763 returns:
764 True if all facets !flipped, !dupridge, normal
765 if all horizon facets are simplicial
766 if all vertices are clearly below neighbor
767 if all opposite vertices of horizon are below
768 clears qh.ZEROall_ok if any problems or coplanar facets
769
770 notes:
771 called by qh_premerge (qh.CHECKzero, 'C-0') and qh_qhull ('Qx')
772 uses qh.vertex_visit
773 horizon facets may define multiple new facets
774
775 design:
776 for all facets in qh.newfacet_list or qh.facet_list
777 check for flagged faults (flipped, etc.)
778 for all facets in qh.newfacet_list or qh.facet_list
779 for each neighbor of facet
780 skip horizon facets for qh.newfacet_list
781 test the opposite vertex
782 if qh.newfacet_list
783 test the other vertices in the facet's horizon facet
784 */
qh_checkzero(boolT testall)785 boolT qh_checkzero(boolT testall) {
786 facetT *facet, *neighbor;
787 facetT *horizon, *facetlist;
788 int neighbor_i, neighbor_n;
789 vertexT *vertex, **vertexp;
790 realT dist;
791
792 if (testall)
793 facetlist= qh facet_list;
794 else {
795 facetlist= qh newfacet_list;
796 FORALLfacet_(facetlist) {
797 horizon= SETfirstt_(facet->neighbors, facetT);
798 if (!horizon->simplicial)
799 goto LABELproblem;
800 if (facet->flipped || facet->dupridge || !facet->normal)
801 goto LABELproblem;
802 }
803 if (qh MERGEexact && qh ZEROall_ok) {
804 trace2((qh ferr, 2011, "qh_checkzero: skip convexity check until first pre-merge\n"));
805 return True;
806 }
807 }
808 FORALLfacet_(facetlist) {
809 qh vertex_visit++;
810 horizon= NULL;
811 FOREACHneighbor_i_(facet) {
812 if (!neighbor_i && !testall) {
813 horizon= neighbor;
814 continue; /* horizon facet tested in qh_findhorizon */
815 }
816 vertex= SETelemt_(facet->vertices, neighbor_i, vertexT);
817 vertex->visitid= qh vertex_visit;
818 zzinc_(Zdistzero);
819 qh_distplane(vertex->point, neighbor, &dist);
820 if (dist >= -2 * qh DISTround) { /* need 2x for qh_distround and 'Rn' for qh_checkconvex, same as qh.premerge_centrum */
821 qh ZEROall_ok= False;
822 if (!qh MERGEexact || testall || dist > qh DISTround)
823 goto LABELnonconvex;
824 }
825 }
826 if (!testall && horizon) {
827 FOREACHvertex_(horizon->vertices) {
828 if (vertex->visitid != qh vertex_visit) {
829 zzinc_(Zdistzero);
830 qh_distplane(vertex->point, facet, &dist);
831 if (dist >= -2 * qh DISTround) {
832 qh ZEROall_ok= False;
833 if (!qh MERGEexact || dist > qh DISTround)
834 goto LABELnonconvexhorizon;
835 }
836 break;
837 }
838 }
839 }
840 }
841 trace2((qh ferr, 2012, "qh_checkzero: testall %d, facets are %s\n", testall,
842 (qh MERGEexact && !testall) ?
843 "not concave, flipped, or dupridge" : "clearly convex"));
844 return True;
845
846 LABELproblem:
847 qh ZEROall_ok= False;
848 trace2((qh ferr, 2013, "qh_checkzero: qh_premerge is needed. New facet f%d or its horizon f%d is non-simplicial, flipped, dupridge, or mergehorizon\n",
849 facet->id, horizon->id));
850 return False;
851
852 LABELnonconvex:
853 trace2((qh ferr, 2014, "qh_checkzero: facet f%d and f%d are not clearly convex. v%d dist %.2g\n",
854 facet->id, neighbor->id, vertex->id, dist));
855 return False;
856
857 LABELnonconvexhorizon:
858 trace2((qh ferr, 2060, "qh_checkzero: facet f%d and horizon f%d are not clearly convex. v%d dist %.2g\n",
859 facet->id, horizon->id, vertex->id, dist));
860 return False;
861 } /* checkzero */
862
863 /*-<a href="qh-merge.htm#TOC"
864 >-------------------------------</a><a name="compare_anglemerge">-</a>
865
866 qh_compare_anglemerge( mergeA, mergeB )
867 used by qsort() to order qh.facet_mergeset by mergetype and angle (qh.ANGLEmerge, 'Q1')
868 lower numbered mergetypes done first (MRGcoplanar before MRGconcave)
869
870 notes:
871 qh_all_merges processes qh.facet_mergeset by qh_setdellast
872 [mar'19] evaluated various options with eg/q_benchmark and merging of pinched vertices (Q14)
873 */
qh_compare_anglemerge(const void * p1,const void * p2)874 int qh_compare_anglemerge(const void *p1, const void *p2) {
875 const mergeT *a= *((mergeT *const*)p1), *b= *((mergeT *const*)p2);
876
877 if (a->mergetype != b->mergetype)
878 return (a->mergetype < b->mergetype ? 1 : -1); /* select MRGcoplanar (1) before MRGconcave (3) */
879 else
880 return (a->angle > b->angle ? 1 : -1); /* select coplanar merge (1.0) before sharp merge (-0.5) */
881 } /* compare_anglemerge */
882
883 /*-<a href="qh-merge.htm#TOC"
884 >-------------------------------</a><a name="compare_facetmerge">-</a>
885
886 qh_compare_facetmerge( mergeA, mergeB )
887 used by qsort() to order merges by mergetype, first merge, first
888 lower numbered mergetypes done first (MRGcoplanar before MRGconcave)
889 if same merge type, flat merges are first
890
891 notes:
892 qh_all_merges processes qh.facet_mergeset by qh_setdellast
893 [mar'19] evaluated various options with eg/q_benchmark and merging of pinched vertices (Q14)
894 */
qh_compare_facetmerge(const void * p1,const void * p2)895 int qh_compare_facetmerge(const void *p1, const void *p2) {
896 const mergeT *a= *((mergeT *const*)p1), *b= *((mergeT *const*)p2);
897
898 if (a->mergetype != b->mergetype)
899 return (a->mergetype < b->mergetype ? 1 : -1); /* select MRGcoplanar (1) before MRGconcave (3) */
900 else if (a->mergetype == MRGanglecoplanar)
901 return (a->angle > b->angle ? 1 : -1); /* if MRGanglecoplanar, select coplanar merge (1.0) before sharp merge (-0.5) */
902 else
903 return (a->distance < b->distance ? 1 : -1); /* select flat (0.0) merge before wide (1e-10) merge */
904 } /* compare_facetmerge */
905
906 /*-<a href="qh-merge.htm#TOC"
907 >-------------------------------</a><a name="comparevisit">-</a>
908
909 qh_comparevisit( vertexA, vertexB )
910 used by qsort() to order vertices by their visitid
911
912 notes:
913 only called by qh_find_newvertex
914 */
qh_comparevisit(const void * p1,const void * p2)915 int qh_comparevisit(const void *p1, const void *p2) {
916 const vertexT *a= *((vertexT *const*)p1), *b= *((vertexT *const*)p2);
917
918 if (a->visitid > b->visitid)
919 return 1;
920 return -1;
921 } /* comparevisit */
922
923 /*-<a href="qh-merge.htm#TOC"
924 >-------------------------------</a><a name="copynonconvex">-</a>
925
926 qh_copynonconvex( atridge )
927 set non-convex flag on other ridges (if any) between same neighbors
928
929 notes:
930 may be faster if use smaller ridge set
931
932 design:
933 for each ridge of atridge's top facet
934 if ridge shares the same neighbor
935 set nonconvex flag
936 */
qh_copynonconvex(ridgeT * atridge)937 void qh_copynonconvex(ridgeT *atridge) {
938 facetT *facet, *otherfacet;
939 ridgeT *ridge, **ridgep;
940
941 facet= atridge->top;
942 otherfacet= atridge->bottom;
943 atridge->nonconvex= False;
944 FOREACHridge_(facet->ridges) {
945 if (otherfacet == ridge->top || otherfacet == ridge->bottom) {
946 if (ridge != atridge) {
947 ridge->nonconvex= True;
948 trace4((qh ferr, 4020, "qh_copynonconvex: moved nonconvex flag from r%d to r%d between f%d and f%d\n",
949 atridge->id, ridge->id, facet->id, otherfacet->id));
950 break;
951 }
952 }
953 }
954 } /* copynonconvex */
955
956 /*-<a href="qh-merge.htm#TOC"
957 >-------------------------------</a><a name="degen_redundant_facet">-</a>
958
959 qh_degen_redundant_facet( facet )
960 check for a degenerate (too few neighbors) or redundant (subset of vertices) facet
961
962 notes:
963 called at end of qh_mergefacet, qh_renamevertex, and qh_reducevertices
964 bumps vertex_visit
965 called if a facet was redundant but no longer is (qh_merge_degenredundant)
966 qh_appendmergeset() only appends first reference to facet (i.e., redundant)
967 see: qh_test_redundant_neighbors, qh_maydropneighbor
968
969 design:
970 test for redundant neighbor
971 test for degenerate facet
972 */
qh_degen_redundant_facet(facetT * facet)973 void qh_degen_redundant_facet(facetT *facet) {
974 vertexT *vertex, **vertexp;
975 facetT *neighbor, **neighborp;
976
977 trace3((qh ferr, 3028, "qh_degen_redundant_facet: test facet f%d for degen/redundant\n",
978 facet->id));
979 if (facet->flipped) {
980 trace2((qh ferr, 3074, "qh_degen_redundant_facet: f%d is flipped, will merge later\n", facet->id));
981 return;
982 }
983 FOREACHneighbor_(facet) {
984 if (neighbor->flipped) /* disallow merge of non-flipped into flipped, neighbor will be merged later */
985 continue;
986 if (neighbor->visible) {
987 qh_fprintf(qh ferr, 6357, "qhull internal error (qh_degen_redundant_facet): facet f%d has deleted neighbor f%d (qh.visible_list)\n",
988 facet->id, neighbor->id);
989 qh_errexit2(qh_ERRqhull, facet, neighbor);
990 }
991 qh vertex_visit++;
992 FOREACHvertex_(neighbor->vertices)
993 vertex->visitid= qh vertex_visit;
994 FOREACHvertex_(facet->vertices) {
995 if (vertex->visitid != qh vertex_visit)
996 break;
997 }
998 if (!vertex) {
999 trace2((qh ferr, 2015, "qh_degen_redundant_facet: f%d is contained in f%d. merge\n", facet->id, neighbor->id));
1000 qh_appendmergeset(facet, neighbor, MRGredundant, 0.0, 1.0);
1001 return;
1002 }
1003 }
1004 if (qh_setsize(facet->neighbors) < qh hull_dim) {
1005 qh_appendmergeset(facet, facet, MRGdegen, 0.0, 1.0);
1006 trace2((qh ferr, 2016, "qh_degen_redundant_facet: f%d is degenerate.\n", facet->id));
1007 }
1008 } /* degen_redundant_facet */
1009
1010
1011 /*-<a href="qh-merge.htm#TOC"
1012 >-------------------------------</a><a name="delridge_merge">-</a>
1013
1014 qh_delridge_merge( ridge )
1015 delete ridge due to a merge
1016
1017 notes:
1018 only called by merge.c (qh_mergeridges, qh_renameridgevertex)
1019 ridges also freed in qh_freeqhull and qh_mergecycle_ridges
1020
1021 design:
1022 if needed, moves ridge.nonconvex to another ridge
1023 sets vertex.delridge for qh_reducevertices
1024 deletes ridge from qh.vertex_mergeset
1025 deletes ridge from its neighboring facets
1026 frees up its memory
1027 */
qh_delridge_merge(ridgeT * ridge)1028 void qh_delridge_merge(ridgeT *ridge) {
1029 vertexT *vertex, **vertexp;
1030 mergeT *merge;
1031 int merge_i, merge_n;
1032
1033 trace3((qh ferr, 3036, "qh_delridge_merge: delete ridge r%d between f%d and f%d\n",
1034 ridge->id, ridge->top->id, ridge->bottom->id));
1035 if (ridge->nonconvex)
1036 qh_copynonconvex(ridge);
1037 FOREACHvertex_(ridge->vertices)
1038 vertex->delridge= True;
1039 FOREACHmerge_i_(qh vertex_mergeset) {
1040 if (merge->ridge1 == ridge || merge->ridge2 == ridge) {
1041 trace3((qh ferr, 3029, "qh_delridge_merge: drop merge of v%d into v%d (dist %2.2g r%d r%d) due to deleted, duplicated ridge r%d\n",
1042 merge->vertex1->id, merge->vertex2->id, merge->distance, merge->ridge1->id, merge->ridge2->id, ridge->id));
1043 if (merge->ridge1 == ridge)
1044 merge->ridge2->mergevertex= False;
1045 else
1046 merge->ridge1->mergevertex= False;
1047 qh_setdelnth(qh vertex_mergeset, merge_i);
1048 merge_i--; merge_n--; /* next merge after deleted */
1049 }
1050 }
1051 qh_setdel(ridge->top->ridges, ridge);
1052 qh_setdel(ridge->bottom->ridges, ridge);
1053 qh_delridge(ridge);
1054 } /* delridge_merge */
1055
1056
1057 /*-<a href="qh-merge.htm#TOC"
1058 >-------------------------------</a><a name="drop_mergevertex">-</a>
1059
1060 qh_drop_mergevertex( merge )
1061
1062 clear mergevertex flags for ridges of a vertex merge
1063 */
qh_drop_mergevertex(mergeT * merge)1064 void qh_drop_mergevertex(mergeT *merge)
1065 {
1066 if (merge->mergetype == MRGvertices) {
1067 merge->ridge1->mergevertex= False;
1068 merge->ridge1->mergevertex2= True;
1069 merge->ridge2->mergevertex= False;
1070 merge->ridge2->mergevertex2= True;
1071 trace3((qh ferr, 3032, "qh_drop_mergevertex: unset mergevertex for r%d and r%d due to dropped vertex merge v%d to v%d. Sets mergevertex2\n",
1072 merge->ridge1->id, merge->ridge2->id, merge->vertex1->id, merge->vertex2->id));
1073 }
1074 } /* drop_mergevertex */
1075
1076 /*-<a href="qh-merge.htm#TOC"
1077 >-------------------------------</a><a name="find_newvertex">-</a>
1078
1079 qh_find_newvertex( oldvertex, vertices, ridges )
1080 locate new vertex for renaming old vertex
1081 vertices is a set of possible new vertices
1082 vertices sorted by number of deleted ridges
1083
1084 returns:
1085 newvertex or NULL
1086 each ridge includes both newvertex and oldvertex
1087 vertices without oldvertex sorted by number of deleted ridges
1088 qh.vertex_visit updated
1089 sets v.seen
1090
1091 notes:
1092 called by qh_redundant_vertex due to vertex->delridge and qh_rename_sharedvertex
1093 sets vertex->visitid to 0..setsize() for vertices
1094 new vertex is in one of the ridges
1095 renaming will not cause a duplicate ridge
1096 renaming will minimize the number of deleted ridges
1097 newvertex may not be adjacent in the dual (though unlikely)
1098
1099 design:
1100 for each vertex in vertices
1101 set vertex->visitid to number of ridges
1102 remove unvisited vertices
1103 set qh.vertex_visit above all possible values
1104 sort vertices by number of ridges (minimize ridges that need renaming
1105 add each ridge to qh.hash_table
1106 for each vertex in vertices
1107 find the first vertex that would not cause a duplicate ridge after a rename
1108 */
qh_find_newvertex(vertexT * oldvertex,setT * vertices,setT * ridges)1109 vertexT *qh_find_newvertex(vertexT *oldvertex, setT *vertices, setT *ridges) {
1110 vertexT *vertex, **vertexp;
1111 setT *newridges;
1112 ridgeT *ridge, **ridgep;
1113 int size, hashsize;
1114 int hash;
1115 unsigned int maxvisit;
1116
1117 #ifndef qh_NOtrace
1118 if (qh IStracing >= 4) {
1119 qh_fprintf(qh ferr, 8063, "qh_find_newvertex: find new vertex for v%d from ",
1120 oldvertex->id);
1121 FOREACHvertex_(vertices)
1122 qh_fprintf(qh ferr, 8064, "v%d ", vertex->id);
1123 FOREACHridge_(ridges)
1124 qh_fprintf(qh ferr, 8065, "r%d ", ridge->id);
1125 qh_fprintf(qh ferr, 8066, "\n");
1126 }
1127 #endif
1128 FOREACHridge_(ridges) {
1129 FOREACHvertex_(ridge->vertices)
1130 vertex->seen= False;
1131 }
1132 FOREACHvertex_(vertices) {
1133 vertex->visitid= 0; /* v.visitid will be number of ridges */
1134 vertex->seen= True;
1135 }
1136 FOREACHridge_(ridges) {
1137 FOREACHvertex_(ridge->vertices) {
1138 if (vertex->seen)
1139 vertex->visitid++;
1140 }
1141 }
1142 FOREACHvertex_(vertices) {
1143 if (!vertex->visitid) {
1144 qh_setdelnth(vertices, SETindex_(vertices,vertex));
1145 vertexp--; /* repeat since deleted this vertex */
1146 }
1147 }
1148 maxvisit= (unsigned int)qh_setsize(ridges);
1149 maximize_(qh vertex_visit, maxvisit);
1150 if (!qh_setsize(vertices)) {
1151 trace4((qh ferr, 4023, "qh_find_newvertex: vertices not in ridges for v%d\n",
1152 oldvertex->id));
1153 return NULL;
1154 }
1155 qsort(SETaddr_(vertices, vertexT), (size_t)qh_setsize(vertices),
1156 sizeof(vertexT *), qh_comparevisit);
1157 /* can now use qh vertex_visit */
1158 if (qh PRINTstatistics) {
1159 size= qh_setsize(vertices);
1160 zinc_(Zintersect);
1161 zadd_(Zintersecttot, size);
1162 zmax_(Zintersectmax, size);
1163 }
1164 hashsize= qh_newhashtable(qh_setsize(ridges));
1165 FOREACHridge_(ridges)
1166 qh_hashridge(qh hash_table, hashsize, ridge, oldvertex);
1167 FOREACHvertex_(vertices) {
1168 newridges= qh_vertexridges(vertex, !qh_ALL);
1169 FOREACHridge_(newridges) {
1170 if (qh_hashridge_find(qh hash_table, hashsize, ridge, vertex, oldvertex, &hash)) {
1171 zinc_(Zvertexridge);
1172 break;
1173 }
1174 }
1175 qh_settempfree(&newridges);
1176 if (!ridge)
1177 break; /* found a rename */
1178 }
1179 if (vertex) {
1180 /* counted in qh_renamevertex */
1181 trace2((qh ferr, 2020, "qh_find_newvertex: found v%d for old v%d from %d vertices and %d ridges.\n",
1182 vertex->id, oldvertex->id, qh_setsize(vertices), qh_setsize(ridges)));
1183 }else {
1184 zinc_(Zfindfail);
1185 trace0((qh ferr, 14, "qh_find_newvertex: no vertex for renaming v%d (all duplicated ridges) during p%d\n",
1186 oldvertex->id, qh furthest_id));
1187 }
1188 qh_setfree(&qh hash_table);
1189 return vertex;
1190 } /* find_newvertex */
1191
1192 /*-<a href="qh-geom2.htm#TOC"
1193 >-------------------------------</a><a name="findbest_pinchedvertex">-</a>
1194
1195 qh_findbest_pinchedvertex( merge, apex, nearestp, distp )
1196 Determine the best pinched vertex to rename as its nearest neighboring vertex
1197 Renaming will remove a duplicate MRGdupridge in newfacet_list
1198
1199 returns:
1200 pinched vertex (either apex or subridge), nearest vertex (subridge or neighbor vertex), and the distance between them
1201
1202 notes:
1203 only called by qh_getpinchedmerges
1204 assumes qh.VERTEXneighbors
1205 see qh_findbest_ridgevertex
1206
1207 design:
1208 if the facets have the same vertices
1209 return the nearest vertex pair
1210 else
1211 the subridge is the intersection of the two new facets minus the apex
1212 the subridge consists of qh.hull_dim-2 horizon vertices
1213 the subridge is also a matched ridge for the new facets (its duplicate)
1214 determine the nearest vertex to the apex
1215 determine the nearest pair of subridge vertices
1216 for each vertex in the subridge
1217 determine the nearest neighbor vertex (not in the subridge)
1218 */
qh_findbest_pinchedvertex(mergeT * merge,vertexT * apex,vertexT ** nearestp,coordT * distp)1219 vertexT *qh_findbest_pinchedvertex(mergeT *merge, vertexT *apex, vertexT **nearestp, coordT *distp /* qh.newfacet_list */) {
1220 vertexT *vertex, **vertexp, *vertexA, **vertexAp;
1221 vertexT *bestvertex= NULL, *bestpinched= NULL;
1222 setT *subridge, *maybepinched;
1223 coordT dist, bestdist= REALmax;
1224 coordT pincheddist= (qh ONEmerge+qh DISTround)*qh_RATIOpinchedsubridge;
1225
1226 if (!merge->facet1->simplicial || !merge->facet2->simplicial) {
1227 qh_fprintf(qh ferr, 6351, "qhull internal error (qh_findbest_pinchedvertex): expecting merge of adjacent, simplicial new facets. f%d or f%d is not simplicial\n",
1228 merge->facet1->id, merge->facet2->id);
1229 qh_errexit2(qh_ERRqhull, merge->facet1, merge->facet2);
1230 }
1231 subridge= qh_vertexintersect_new(merge->facet1->vertices, merge->facet2->vertices); /* new setT. No error_exit() */
1232 if (qh_setsize(subridge) == qh hull_dim) { /* duplicate vertices */
1233 bestdist= qh_vertex_bestdist2(subridge, &bestvertex, &bestpinched);
1234 if(bestvertex == apex) {
1235 bestvertex= bestpinched;
1236 bestpinched= apex;
1237 }
1238 }else {
1239 qh_setdel(subridge, apex);
1240 if (qh_setsize(subridge) != qh hull_dim - 2) {
1241 qh_fprintf(qh ferr, 6409, "qhull internal error (qh_findbest_pinchedvertex): expecting subridge of qh.hull_dim-2 vertices for the intersection of new facets f%d and f%d minus their apex. Got %d vertices\n",
1242 merge->facet1->id, merge->facet2->id, qh_setsize(subridge));
1243 qh_errexit2(qh_ERRqhull, merge->facet1, merge->facet2);
1244 }
1245 FOREACHvertex_(subridge) {
1246 dist= qh_pointdist(vertex->point, apex->point, qh hull_dim);
1247 if (dist < bestdist) {
1248 bestpinched= apex;
1249 bestvertex= vertex;
1250 bestdist= dist;
1251 }
1252 }
1253 if (bestdist > pincheddist) {
1254 FOREACHvertex_(subridge) {
1255 FOREACHvertexA_(subridge) {
1256 if (vertexA->id > vertex->id) { /* once per vertex pair, do not compare addresses */
1257 dist= qh_pointdist(vertexA->point, vertex->point, qh hull_dim);
1258 if (dist < bestdist) {
1259 bestpinched= vertexA;
1260 bestvertex= vertex;
1261 bestdist= dist;
1262 }
1263 }
1264 }
1265 }
1266 }
1267 if (bestdist > pincheddist) {
1268 FOREACHvertexA_(subridge) {
1269 maybepinched= qh_neighbor_vertices(vertexA, subridge); /* subridge and apex tested above */
1270 FOREACHvertex_(maybepinched) {
1271 dist= qh_pointdist(vertex->point, vertexA->point, qh hull_dim);
1272 if (dist < bestdist) {
1273 bestvertex= vertex;
1274 bestpinched= vertexA;
1275 bestdist= dist;
1276 }
1277 }
1278 qh_settempfree(&maybepinched);
1279 }
1280 }
1281 }
1282 *distp= bestdist;
1283 qh_setfree(&subridge); /* qh_err_exit not called since allocated */
1284 if (!bestvertex) { /* should never happen if qh.hull_dim > 2 */
1285 qh_fprintf(qh ferr, 6274, "qhull internal error (qh_findbest_pinchedvertex): did not find best vertex for subridge of dupridge between f%d and f%d, while processing p%d\n", merge->facet1->id, merge->facet2->id, qh furthest_id);
1286 qh_errexit2(qh_ERRqhull, merge->facet1, merge->facet2);
1287 }
1288 *nearestp= bestvertex;
1289 trace2((qh ferr, 2061, "qh_findbest_pinchedvertex: best pinched p%d(v%d) and vertex p%d(v%d) are closest (%2.2g) for duplicate subridge between f%d and f%d\n",
1290 qh_pointid(bestpinched->point), bestpinched->id, qh_pointid(bestvertex->point), bestvertex->id, bestdist, merge->facet1->id, merge->facet2->id));
1291 return bestpinched;
1292 } /* findbest_pinchedvertex */
1293
1294 /*-<a href="qh-geom2.htm#TOC"
1295 >-------------------------------</a><a name="findbest_ridgevertex">-</a>
1296
1297 qh_findbest_ridgevertex( ridge, pinchedp, distp )
1298 Determine the best vertex/pinched-vertex to merge for ridges with the same vertices
1299
1300 returns:
1301 vertex, pinched vertex, and the distance between them
1302
1303 notes:
1304 assumes qh.hull_dim>=3
1305 see qh_findbest_pinchedvertex
1306
1307 */
qh_findbest_ridgevertex(ridgeT * ridge,vertexT ** pinchedp,coordT * distp)1308 vertexT *qh_findbest_ridgevertex(ridgeT *ridge, vertexT **pinchedp, coordT *distp) {
1309 vertexT *bestvertex;
1310
1311 *distp= qh_vertex_bestdist2(ridge->vertices, &bestvertex, pinchedp);
1312 trace4((qh ferr, 4069, "qh_findbest_ridgevertex: best pinched p%d(v%d) and vertex p%d(v%d) are closest (%2.2g) for duplicated ridge r%d (same vertices) between f%d and f%d\n",
1313 qh_pointid((*pinchedp)->point), (*pinchedp)->id, qh_pointid(bestvertex->point), bestvertex->id, *distp, ridge->id, ridge->top->id, ridge->bottom->id));
1314 return bestvertex;
1315 } /* findbest_ridgevertex */
1316
1317 /*-<a href="qh-merge.htm#TOC"
1318 >-------------------------------</a><a name="findbest_test">-</a>
1319
1320 qh_findbest_test( testcentrum, facet, neighbor, &bestfacet, &dist, &mindist, &maxdist )
1321 test neighbor of facet for qh_findbestneighbor()
1322 if testcentrum,
1323 tests centrum (assumes it is defined)
1324 else
1325 tests vertices
1326 initially *bestfacet==NULL and *dist==REALmax
1327
1328 returns:
1329 if a better facet (i.e., vertices/centrum of facet closer to neighbor)
1330 updates bestfacet, dist, mindist, and maxdist
1331
1332 notes:
1333 called by qh_findbestneighbor
1334 ignores pairs of flipped facets, unless that's all there is
1335 */
qh_findbest_test(boolT testcentrum,facetT * facet,facetT * neighbor,facetT ** bestfacet,realT * distp,realT * mindistp,realT * maxdistp)1336 void qh_findbest_test(boolT testcentrum, facetT *facet, facetT *neighbor,
1337 facetT **bestfacet, realT *distp, realT *mindistp, realT *maxdistp) {
1338 realT dist, mindist, maxdist;
1339
1340 if (facet->flipped && neighbor->flipped && *bestfacet && !(*bestfacet)->flipped)
1341 return; /* do not merge flipped into flipped facets */
1342 if (testcentrum) {
1343 zzinc_(Zbestdist);
1344 qh_distplane(facet->center, neighbor, &dist);
1345 dist *= qh hull_dim; /* estimate furthest vertex */
1346 if (dist < 0) {
1347 maxdist= 0;
1348 mindist= dist;
1349 dist= -dist;
1350 }else {
1351 mindist= 0;
1352 maxdist= dist;
1353 }
1354 }else
1355 dist= qh_getdistance(facet, neighbor, &mindist, &maxdist);
1356 if (dist < *distp) {
1357 *bestfacet= neighbor;
1358 *mindistp= mindist;
1359 *maxdistp= maxdist;
1360 *distp= dist;
1361 }
1362 } /* findbest_test */
1363
1364 /*-<a href="qh-merge.htm#TOC"
1365 >-------------------------------</a><a name="findbestneighbor">-</a>
1366
1367 qh_findbestneighbor( facet, dist, mindist, maxdist )
1368 finds best neighbor (least dist) of a facet for merging
1369
1370 returns:
1371 returns min and max distances and their max absolute value
1372
1373 notes:
1374 error if qh_ASvoronoi
1375 avoids merging old into new
1376 assumes ridge->nonconvex only set on one ridge between a pair of facets
1377 could use an early out predicate but not worth it
1378
1379 design:
1380 if a large facet
1381 will test centrum
1382 else
1383 will test vertices
1384 if a large facet
1385 test nonconvex neighbors for best merge
1386 else
1387 test all neighbors for the best merge
1388 if testing centrum
1389 get distance information
1390 */
qh_findbestneighbor(facetT * facet,realT * distp,realT * mindistp,realT * maxdistp)1391 facetT *qh_findbestneighbor(facetT *facet, realT *distp, realT *mindistp, realT *maxdistp) {
1392 facetT *neighbor, **neighborp, *bestfacet= NULL;
1393 ridgeT *ridge, **ridgep;
1394 boolT nonconvex= True, testcentrum= False;
1395 int size= qh_setsize(facet->vertices);
1396
1397 if(qh CENTERtype==qh_ASvoronoi){
1398 qh_fprintf(qh ferr, 6272, "qhull internal error: cannot call qh_findbestneighor for f%d while qh.CENTERtype is qh_ASvoronoi\n", facet->id);
1399 qh_errexit(qh_ERRqhull, facet, NULL);
1400 }
1401 *distp= REALmax;
1402 if (size > qh_BESTcentrum2 * qh hull_dim + qh_BESTcentrum) {
1403 testcentrum= True;
1404 zinc_(Zbestcentrum);
1405 if (!facet->center)
1406 facet->center= qh_getcentrum(facet);
1407 }
1408 if (size > qh hull_dim + qh_BESTnonconvex) {
1409 FOREACHridge_(facet->ridges) {
1410 if (ridge->nonconvex) {
1411 neighbor= otherfacet_(ridge, facet);
1412 qh_findbest_test(testcentrum, facet, neighbor,
1413 &bestfacet, distp, mindistp, maxdistp);
1414 }
1415 }
1416 }
1417 if (!bestfacet) {
1418 nonconvex= False;
1419 FOREACHneighbor_(facet)
1420 qh_findbest_test(testcentrum, facet, neighbor,
1421 &bestfacet, distp, mindistp, maxdistp);
1422 }
1423 if (!bestfacet) {
1424 qh_fprintf(qh ferr, 6095, "qhull internal error (qh_findbestneighbor): no neighbors for f%d\n", facet->id);
1425 qh_errexit(qh_ERRqhull, facet, NULL);
1426 }
1427 if (testcentrum)
1428 qh_getdistance(facet, bestfacet, mindistp, maxdistp);
1429 trace3((qh ferr, 3002, "qh_findbestneighbor: f%d is best neighbor for f%d testcentrum? %d nonconvex? %d dist %2.2g min %2.2g max %2.2g\n",
1430 bestfacet->id, facet->id, testcentrum, nonconvex, *distp, *mindistp, *maxdistp));
1431 return(bestfacet);
1432 } /* findbestneighbor */
1433
1434
1435 /*-<a href="qh-merge.htm#TOC"
1436 >-------------------------------</a><a name="flippedmerges">-</a>
1437
1438 qh_flippedmerges( facetlist, wasmerge )
1439 merge flipped facets into best neighbor
1440 assumes qh.facet_mergeset at top of temporary stack
1441
1442 returns:
1443 no flipped facets on facetlist
1444 sets wasmerge if merge occurred
1445 degen/redundant merges passed through
1446
1447 notes:
1448 othermerges not needed since qh.facet_mergeset is empty before & after
1449 keep it in case of change
1450
1451 design:
1452 append flipped facets to qh.facetmergeset
1453 for each flipped merge
1454 find best neighbor
1455 merge facet into neighbor
1456 merge degenerate and redundant facets
1457 remove flipped merges from qh.facet_mergeset
1458 */
qh_flippedmerges(facetT * facetlist,boolT * wasmerge)1459 void qh_flippedmerges(facetT *facetlist, boolT *wasmerge) {
1460 facetT *facet, *neighbor, *facet1;
1461 realT dist, mindist, maxdist;
1462 mergeT *merge, **mergep;
1463 setT *othermerges;
1464 int nummerge= 0, numdegen= 0;
1465
1466 trace4((qh ferr, 4024, "qh_flippedmerges: begin\n"));
1467 FORALLfacet_(facetlist) {
1468 if (facet->flipped && !facet->visible)
1469 qh_appendmergeset(facet, facet, MRGflip, 0.0, 1.0);
1470 }
1471 othermerges= qh_settemppop();
1472 if(othermerges != qh facet_mergeset) {
1473 qh_fprintf(qh ferr, 6392, "qhull internal error (qh_flippedmerges): facet_mergeset (%d merges) not at top of tempstack (%d merges)\n",
1474 qh_setsize(qh facet_mergeset), qh_setsize(othermerges));
1475 qh_errexit(qh_ERRqhull, NULL, NULL);
1476 }
1477 qh facet_mergeset= qh_settemp(qh TEMPsize);
1478 qh_settemppush(othermerges);
1479 FOREACHmerge_(othermerges) {
1480 facet1= merge->facet1;
1481 if (merge->mergetype != MRGflip || facet1->visible)
1482 continue;
1483 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
1484 qhmem.IStracing= qh IStracing= qh TRACElevel;
1485 neighbor= qh_findbestneighbor(facet1, &dist, &mindist, &maxdist);
1486 trace0((qh ferr, 15, "qh_flippedmerges: merge flipped f%d into f%d dist %2.2g during p%d\n",
1487 facet1->id, neighbor->id, dist, qh furthest_id));
1488 qh_mergefacet(facet1, neighbor, merge->mergetype, &mindist, &maxdist, !qh_MERGEapex);
1489 nummerge++;
1490 if (qh PRINTstatistics) {
1491 zinc_(Zflipped);
1492 wadd_(Wflippedtot, dist);
1493 wmax_(Wflippedmax, dist);
1494 }
1495 }
1496 FOREACHmerge_(othermerges) {
1497 if (merge->facet1->visible || merge->facet2->visible)
1498 qh_memfree(merge, (int)sizeof(mergeT)); /* invalidates merge and othermerges */
1499 else
1500 qh_setappend(&qh facet_mergeset, merge);
1501 }
1502 qh_settempfree(&othermerges);
1503 numdegen += qh_merge_degenredundant(); /* somewhat better here than after each flipped merge -- qtest.sh 10 '500 C1,2e-13 D4' 'd Qbb' */
1504 if (nummerge)
1505 *wasmerge= True;
1506 trace1((qh ferr, 1010, "qh_flippedmerges: merged %d flipped and %d degenredundant facets into a good neighbor\n",
1507 nummerge, numdegen));
1508 } /* flippedmerges */
1509
1510
1511 /*-<a href="qh-merge.htm#TOC"
1512 >-------------------------------</a><a name="forcedmerges">-</a>
1513
1514 qh_forcedmerges( wasmerge )
1515 merge dupridges
1516 calls qh_check_dupridge to report an error on wide merges
1517 assumes qh_settemppop is qh.facet_mergeset
1518
1519 returns:
1520 removes all dupridges on facet_mergeset
1521 wasmerge set if merge
1522 qh.facet_mergeset may include non-forced merges(none for now)
1523 qh.degen_mergeset includes degen/redun merges
1524
1525 notes:
1526 called by qh_premerge
1527 dupridges occur when the horizon is pinched,
1528 i.e. a subridge occurs in more than two horizon ridges.
1529 could rename vertices that pinch the horizon
1530 assumes qh_merge_degenredundant() has not be called
1531 othermerges isn't needed since facet_mergeset is empty afterwards
1532 keep it in case of change
1533
1534 design:
1535 for each dupridge
1536 find current facets by chasing f.replace links
1537 check for wide merge due to dupridge
1538 determine best direction for facet
1539 merge one facet into the other
1540 remove dupridges from qh.facet_mergeset
1541 */
qh_forcedmerges(boolT * wasmerge)1542 void qh_forcedmerges(boolT *wasmerge) {
1543 facetT *facet1, *facet2, *merging, *merged, *newfacet;
1544 mergeT *merge, **mergep;
1545 realT dist, mindist, maxdist, dist2, mindist2, maxdist2;
1546 setT *othermerges;
1547 int nummerge=0, numflip=0, numdegen= 0;
1548 boolT wasdupridge= False;
1549
1550 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
1551 qhmem.IStracing= qh IStracing= qh TRACElevel;
1552 trace3((qh ferr, 3054, "qh_forcedmerges: merge dupridges\n"));
1553 othermerges= qh_settemppop(); /* was facet_mergeset */
1554 if (qh facet_mergeset != othermerges ) {
1555 qh_fprintf(qh ferr, 6279, "qhull internal error (qh_forcedmerges): qh_settemppop (size %d) is not qh facet_mergeset (size %d)\n",
1556 qh_setsize(othermerges), qh_setsize(qh facet_mergeset));
1557 qh_errexit(qh_ERRqhull, NULL, NULL);
1558 }
1559 qh facet_mergeset= qh_settemp(qh TEMPsize);
1560 qh_settemppush(othermerges);
1561 FOREACHmerge_(othermerges) {
1562 if (merge->mergetype != MRGdupridge)
1563 continue;
1564 wasdupridge= True;
1565 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
1566 qhmem.IStracing= qh IStracing= qh TRACElevel;
1567 facet1= qh_getreplacement(merge->facet1); /* must exist, no qh_merge_degenredunant */
1568 facet2= qh_getreplacement(merge->facet2); /* previously merged facet, if any */
1569 if (facet1 == facet2)
1570 continue;
1571 if (!qh_setin(facet2->neighbors, facet1)) {
1572 qh_fprintf(qh ferr, 6096, "qhull internal error (qh_forcedmerges): f%d and f%d had a dupridge but as f%d and f%d they are no longer neighbors\n",
1573 merge->facet1->id, merge->facet2->id, facet1->id, facet2->id);
1574 qh_errexit2(qh_ERRqhull, facet1, facet2);
1575 }
1576 dist= qh_getdistance(facet1, facet2, &mindist, &maxdist);
1577 dist2= qh_getdistance(facet2, facet1, &mindist2, &maxdist2);
1578 qh_check_dupridge(facet1, dist, facet2, dist2);
1579 if (dist < dist2) {
1580 if (facet2->flipped && !facet1->flipped && dist2 < qh_WIDEdupridge*(qh ONEmerge+qh DISTround)) { /* prefer merge of flipped facet */
1581 merging= facet2;
1582 merged= facet1;
1583 dist= dist2;
1584 mindist= mindist2;
1585 maxdist= maxdist2;
1586 }else {
1587 merging= facet1;
1588 merged= facet2;
1589 }
1590 }else {
1591 if (facet1->flipped && !facet2->flipped && dist < qh_WIDEdupridge*(qh ONEmerge+qh DISTround)) { /* prefer merge of flipped facet */
1592 merging= facet1;
1593 merged= facet2;
1594 }else {
1595 merging= facet2;
1596 merged= facet1;
1597 dist= dist2;
1598 mindist= mindist2;
1599 maxdist= maxdist2;
1600 }
1601 }
1602 qh_mergefacet(merging, merged, merge->mergetype, &mindist, &maxdist, !qh_MERGEapex);
1603 numdegen += qh_merge_degenredundant(); /* better here than at end -- qtest.sh 10 '500 C1,2e-13 D4' 'd Qbb' */
1604 if (facet1->flipped) {
1605 zinc_(Zmergeflipdup);
1606 numflip++;
1607 }else
1608 nummerge++;
1609 if (qh PRINTstatistics) {
1610 zinc_(Zduplicate);
1611 wadd_(Wduplicatetot, dist);
1612 wmax_(Wduplicatemax, dist);
1613 }
1614 }
1615 FOREACHmerge_(othermerges) {
1616 if (merge->mergetype == MRGdupridge)
1617 qh_memfree(merge, (int)sizeof(mergeT)); /* invalidates merge and othermerges */
1618 else
1619 qh_setappend(&qh facet_mergeset, merge);
1620 }
1621 qh_settempfree(&othermerges);
1622 if (wasdupridge) {
1623 FORALLnew_facets {
1624 if (newfacet->dupridge) {
1625 newfacet->dupridge= False;
1626 newfacet->mergeridge= False;
1627 newfacet->mergeridge2= False;
1628 if (qh_setsize(newfacet->neighbors) < qh hull_dim) { /* not tested for MRGdupridge */
1629 qh_appendmergeset(newfacet, newfacet, MRGdegen, 0.0, 1.0);
1630 trace2((qh ferr, 2107, "qh_forcedmerges: dupridge f%d is degenerate with fewer than %d neighbors\n",
1631 newfacet->id, qh hull_dim));
1632 }
1633 }
1634 }
1635 numdegen += qh_merge_degenredundant();
1636 }
1637 if (nummerge || numflip) {
1638 *wasmerge= True;
1639 trace1((qh ferr, 1011, "qh_forcedmerges: merged %d facets, %d flipped facets, and %d degenredundant facets across dupridges\n",
1640 nummerge, numflip, numdegen));
1641 }
1642 } /* forcedmerges */
1643
1644
1645 /*-<a href="qh-merge.htm#TOC"
1646 >-------------------------------</a><a name="freemergesets">-</a>
1647
1648 qh_freemergesets( )
1649 free the merge sets
1650
1651 notes:
1652 matches qh_initmergesets
1653 */
qh_freemergesets(void)1654 void qh_freemergesets(void) {
1655
1656 if (!qh facet_mergeset || !qh degen_mergeset || !qh vertex_mergeset) {
1657 qh_fprintf(qh ferr, 6388, "qhull internal error (qh_freemergesets): expecting mergesets. Got a NULL mergeset, qh.facet_mergeset (0x%x), qh.degen_mergeset (0x%x), qh.vertex_mergeset (0x%x)\n",
1658 qh facet_mergeset, qh degen_mergeset, qh vertex_mergeset);
1659 qh_errexit(qh_ERRqhull, NULL, NULL);
1660 }
1661 if (!SETempty_(qh facet_mergeset) || !SETempty_(qh degen_mergeset) || !SETempty_(qh vertex_mergeset)) {
1662 qh_fprintf(qh ferr, 6389, "qhull internal error (qh_freemergesets): expecting empty mergesets. Got qh.facet_mergeset (%d merges), qh.degen_mergeset (%d merges), qh.vertex_mergeset (%d merges)\n",
1663 qh_setsize(qh facet_mergeset), qh_setsize(qh degen_mergeset), qh_setsize(qh vertex_mergeset));
1664 qh_errexit(qh_ERRqhull, NULL, NULL);
1665 }
1666 qh_settempfree(&qh facet_mergeset);
1667 qh_settempfree(&qh vertex_mergeset);
1668 qh_settempfree(&qh degen_mergeset);
1669 } /* freemergesets */
1670
1671 /*-<a href="qh-merge.htm#TOC"
1672 >-------------------------------</a><a name="getmergeset">-</a>
1673
1674 qh_getmergeset( facetlist )
1675 determines nonconvex facets on facetlist
1676 tests !tested ridges and nonconvex ridges of !tested facets
1677
1678 returns:
1679 returns sorted qh.facet_mergeset of facet-neighbor pairs to be merged
1680 all ridges tested
1681
1682 notes:
1683 facetlist is qh.facet_newlist, use qh_getmergeset_initial for all facets
1684 assumes no nonconvex ridges with both facets tested
1685 uses facet->tested/ridge->tested to prevent duplicate tests
1686 can not limit tests to modified ridges since the centrum changed
1687 uses qh.visit_id
1688
1689 design:
1690 for each facet on facetlist
1691 for each ridge of facet
1692 if untested ridge
1693 test ridge for convexity
1694 if non-convex
1695 append ridge to qh.facet_mergeset
1696 sort qh.facet_mergeset by mergetype and angle or distance
1697 */
qh_getmergeset(facetT * facetlist)1698 void qh_getmergeset(facetT *facetlist) {
1699 facetT *facet, *neighbor, **neighborp;
1700 ridgeT *ridge, **ridgep;
1701 int nummerges;
1702 boolT simplicial;
1703
1704 nummerges= qh_setsize(qh facet_mergeset);
1705 trace4((qh ferr, 4026, "qh_getmergeset: started.\n"));
1706 qh visit_id++;
1707 FORALLfacet_(facetlist) {
1708 if (facet->tested)
1709 continue;
1710 facet->visitid= qh visit_id;
1711 FOREACHneighbor_(facet)
1712 neighbor->seen= False;
1713 /* facet must be non-simplicial due to merge to qh.facet_newlist */
1714 FOREACHridge_(facet->ridges) {
1715 if (ridge->tested && !ridge->nonconvex)
1716 continue;
1717 /* if r.tested & r.nonconvex, need to retest and append merge */
1718 neighbor= otherfacet_(ridge, facet);
1719 if (neighbor->seen) { /* another ridge for this facet-neighbor pair was already tested in this loop */
1720 ridge->tested= True;
1721 ridge->nonconvex= False; /* only one ridge is marked nonconvex per facet-neighbor pair */
1722 }else if (neighbor->visitid != qh visit_id) {
1723 neighbor->seen= True;
1724 ridge->nonconvex= False;
1725 simplicial= False;
1726 if (ridge->simplicialbot && ridge->simplicialtop)
1727 simplicial= True;
1728 if (qh_test_appendmerge(facet, neighbor, simplicial))
1729 ridge->nonconvex= True;
1730 ridge->tested= True;
1731 }
1732 }
1733 facet->tested= True;
1734 }
1735 nummerges= qh_setsize(qh facet_mergeset);
1736 if (qh ANGLEmerge)
1737 qsort(SETaddr_(qh facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_compare_anglemerge);
1738 else
1739 qsort(SETaddr_(qh facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_compare_facetmerge);
1740 nummerges += qh_setsize(qh degen_mergeset);
1741 if (qh POSTmerging) {
1742 zadd_(Zmergesettot2, nummerges);
1743 }else {
1744 zadd_(Zmergesettot, nummerges);
1745 zmax_(Zmergesetmax, nummerges);
1746 }
1747 trace2((qh ferr, 2021, "qh_getmergeset: %d merges found\n", nummerges));
1748 } /* getmergeset */
1749
1750
1751 /*-<a href="qh-merge.htm#TOC"
1752 >-------------------------------</a><a name="getmergeset_initial">-</a>
1753
1754 qh_getmergeset_initial( facetlist )
1755 determine initial qh.facet_mergeset for facets
1756 tests all facet/neighbor pairs on facetlist
1757
1758 returns:
1759 sorted qh.facet_mergeset with nonconvex ridges
1760 sets facet->tested, ridge->tested, and ridge->nonconvex
1761
1762 notes:
1763 uses visit_id, assumes ridge->nonconvex is False
1764 see qh_getmergeset
1765
1766 design:
1767 for each facet on facetlist
1768 for each untested neighbor of facet
1769 test facet and neighbor for convexity
1770 if non-convex
1771 append merge to qh.facet_mergeset
1772 mark one of the ridges as nonconvex
1773 sort qh.facet_mergeset by mergetype and angle or distance
1774 */
qh_getmergeset_initial(facetT * facetlist)1775 void qh_getmergeset_initial(facetT *facetlist) {
1776 facetT *facet, *neighbor, **neighborp;
1777 ridgeT *ridge, **ridgep;
1778 int nummerges;
1779 boolT simplicial;
1780
1781 qh visit_id++;
1782 FORALLfacet_(facetlist) {
1783 facet->visitid= qh visit_id;
1784 FOREACHneighbor_(facet) {
1785 if (neighbor->visitid != qh visit_id) {
1786 simplicial= False; /* ignores r.simplicialtop/simplicialbot. Need to test horizon facets */
1787 if (facet->simplicial && neighbor->simplicial)
1788 simplicial= True;
1789 if (qh_test_appendmerge(facet, neighbor, simplicial)) {
1790 FOREACHridge_(neighbor->ridges) {
1791 if (facet == otherfacet_(ridge, neighbor)) {
1792 ridge->nonconvex= True;
1793 break; /* only one ridge is marked nonconvex */
1794 }
1795 }
1796 }
1797 }
1798 }
1799 facet->tested= True;
1800 FOREACHridge_(facet->ridges)
1801 ridge->tested= True;
1802 }
1803 nummerges= qh_setsize(qh facet_mergeset);
1804 if (qh ANGLEmerge)
1805 qsort(SETaddr_(qh facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_compare_anglemerge);
1806 else
1807 qsort(SETaddr_(qh facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_compare_facetmerge);
1808 nummerges += qh_setsize(qh degen_mergeset);
1809 if (qh POSTmerging) {
1810 zadd_(Zmergeinittot2, nummerges);
1811 }else {
1812 zadd_(Zmergeinittot, nummerges);
1813 zmax_(Zmergeinitmax, nummerges);
1814 }
1815 trace2((qh ferr, 2022, "qh_getmergeset_initial: %d merges found\n", nummerges));
1816 } /* getmergeset_initial */
1817
1818 /*-<a href="qh-merge.htm#TOC"
1819 >-------------------------------</a><a name="getpinchedmerges">-</a>
1820
1821 qh_getpinchedmerges( apex, maxdist, iscoplanar )
1822 get pinched merges for dupridges in qh.facet_mergeset
1823 qh.NEWtentative==True
1824 qh.newfacet_list with apex
1825 qh.horizon_list is attached to qh.visible_list instead of qh.newfacet_list
1826 maxdist for vertex-facet of a dupridge
1827 qh.facet_mergeset is empty
1828 qh.vertex_mergeset is a temporary set
1829
1830 returns:
1831 False if nearest vertex would increase facet width by more than maxdist or qh_WIDEpinched
1832 True and iscoplanar, if the pinched vertex is the apex (i.e., make the apex a coplanar point)
1833 True and !iscoplanar, if should merge a pinched vertex of a dupridge
1834 qh.vertex_mergeset contains one or more MRGsubridge with a pinched vertex and a nearby, neighboring vertex
1835 qh.facet_mergeset is empty
1836
1837 notes:
1838 called by qh_buildcone_mergepinched
1839 hull_dim >= 3
1840 a pinched vertex is in a dupridge and the horizon
1841 selects the pinched vertex that is closest to its neighbor
1842
1843 design:
1844 for each dupridge
1845 determine the best pinched vertex to be merged into a neighboring vertex
1846 if merging the pinched vertex would produce a wide merge (qh_WIDEpinched)
1847 ignore pinched vertex with a warning, and use qh_merge_degenredundant instead
1848 else
1849 append the pinched vertex to vertex_mergeset for merging
1850 */
qh_getpinchedmerges(vertexT * apex,coordT maxdupdist,boolT * iscoplanar)1851 boolT qh_getpinchedmerges(vertexT *apex, coordT maxdupdist, boolT *iscoplanar /* qh.newfacet_list, qh.vertex_mergeset */) {
1852 mergeT *merge, **mergep, *bestmerge= NULL;
1853 vertexT *nearest, *pinched, *bestvertex= NULL, *bestpinched= NULL;
1854 boolT result;
1855 coordT dist, prevdist, bestdist= REALmax/(qh_RATIOcoplanarapex+1.0); /* allow *3.0 */
1856 realT ratio;
1857
1858 trace2((qh ferr, 2062, "qh_getpinchedmerges: try to merge pinched vertices for dupridges in new facets with apex p%d(v%d) max dupdist %2.2g\n",
1859 qh_pointid(apex->point), apex->id, maxdupdist));
1860 *iscoplanar= False;
1861 prevdist= fmax_(qh ONEmerge + qh DISTround, qh MINoutside + qh DISTround);
1862 maximize_(prevdist, qh max_outside);
1863 maximize_(prevdist, -qh min_vertex);
1864 qh_mark_dupridges(qh newfacet_list, !qh_ALL); /* qh.facet_mergeset, creates ridges */
1865 /* qh_mark_dupridges is called a second time in qh_premerge */
1866 FOREACHmerge_(qh facet_mergeset) { /* read-only */
1867 if (merge->mergetype != MRGdupridge) {
1868 qh_fprintf(qh ferr, 6393, "qhull internal error (qh_getpinchedmerges): expecting MRGdupridge from qh_mark_dupridges. Got merge f%d f%d type %d\n",
1869 getid_(merge->facet1), getid_(merge->facet2), merge->mergetype);
1870 qh_errexit(qh_ERRqhull, NULL, NULL);
1871 }
1872 /* dist is distance between vertices */
1873 pinched= qh_findbest_pinchedvertex(merge, apex, &nearest, &dist /* qh.newfacet_list */);
1874 if (pinched == apex && dist < qh_RATIOcoplanarapex*bestdist) { /* prefer coplanar apex since it always works */
1875 bestdist= dist/qh_RATIOcoplanarapex;
1876 bestmerge= merge;
1877 bestpinched= pinched;
1878 bestvertex= nearest;
1879 }else if (dist < bestdist) {
1880 bestdist= dist;
1881 bestmerge= merge;
1882 bestpinched= pinched;
1883 bestvertex= nearest;
1884 }
1885 }
1886 result= False;
1887 if (bestmerge && bestdist < maxdupdist) {
1888 ratio= bestdist / prevdist;
1889 if (ratio > qh_WIDEpinched) {
1890 if (bestmerge->facet1->mergehorizon || bestmerge->facet2->mergehorizon) { /* e.g., rbox 175 C3,2e-13 t1539182828 | qhull d */
1891 trace1((qh ferr, 1051, "qh_getpinchedmerges: dupridge (MRGdupridge) of coplanar horizon would produce a wide merge (%.0fx) due to pinched vertices v%d and v%d (dist %2.2g) for f%d and f%d. qh_mergecycle_all will merge one or both facets\n",
1892 ratio, bestpinched->id, bestvertex->id, bestdist, bestmerge->facet1->id, bestmerge->facet2->id));
1893 }else {
1894 qh_fprintf(qh ferr, 7081, "qhull precision warning (qh_getpinchedmerges): pinched vertices v%d and v%d (dist %2.2g, %.0fx) would produce a wide merge for f%d and f%d. Will merge dupridge instead\n",
1895 bestpinched->id, bestvertex->id, bestdist, ratio, bestmerge->facet1->id, bestmerge->facet2->id);
1896 }
1897 }else {
1898 if (bestpinched == apex) {
1899 trace2((qh ferr, 2063, "qh_getpinchedmerges: will make the apex a coplanar point. apex p%d(v%d) is the nearest vertex to v%d on dupridge. Dist %2.2g\n",
1900 qh_pointid(apex->point), apex->id, bestvertex->id, bestdist*qh_RATIOcoplanarapex));
1901 qh coplanar_apex= apex->point;
1902 *iscoplanar= True;
1903 result= True;
1904 }else if (qh_setin(bestmerge->facet1->vertices, bestpinched) != qh_setin(bestmerge->facet2->vertices, bestpinched)) { /* pinched in one facet but not the other facet */
1905 trace2((qh ferr, 2064, "qh_getpinchedmerges: will merge new facets to resolve dupridge between f%d and f%d with pinched v%d and v%d\n",
1906 bestmerge->facet1->id, bestmerge->facet2->id, bestpinched->id, bestvertex->id));
1907 qh_appendvertexmerge(bestpinched, bestvertex, MRGsubridge, bestdist, NULL, NULL);
1908 result= True;
1909 }else {
1910 trace2((qh ferr, 2065, "qh_getpinchedmerges: will merge pinched v%d into v%d to resolve dupridge between f%d and f%d\n",
1911 bestpinched->id, bestvertex->id, bestmerge->facet1->id, bestmerge->facet2->id));
1912 qh_appendvertexmerge(bestpinched, bestvertex, MRGsubridge, bestdist, NULL, NULL);
1913 result= True;
1914 }
1915 }
1916 }
1917 /* delete MRGdupridge, qh_mark_dupridges is called a second time in qh_premerge */
1918 while ((merge= (mergeT *)qh_setdellast(qh facet_mergeset)))
1919 qh_memfree(merge, (int)sizeof(mergeT));
1920 return result;
1921 }/* getpinchedmerges */
1922
1923 /*-<a href="qh-merge.htm#TOC"
1924 >-------------------------------</a><a name="hasmerge">-</a>
1925
1926 qh_hasmerge( mergeset, mergetype, facetA, facetB )
1927 True if mergeset has mergetype for facetA and facetB
1928 */
qh_hasmerge(setT * mergeset,mergeType type,facetT * facetA,facetT * facetB)1929 boolT qh_hasmerge(setT *mergeset, mergeType type, facetT *facetA, facetT *facetB) {
1930 mergeT *merge, **mergep;
1931
1932 FOREACHmerge_(mergeset) {
1933 if (merge->mergetype == type) {
1934 if (merge->facet1 == facetA && merge->facet2 == facetB)
1935 return True;
1936 if (merge->facet1 == facetB && merge->facet2 == facetA)
1937 return True;
1938 }
1939 }
1940 return False;
1941 }/* hasmerge */
1942
1943 /*-<a href="qh-merge.htm#TOC"
1944 >-------------------------------</a><a name="hashridge">-</a>
1945
1946 qh_hashridge( hashtable, hashsize, ridge, oldvertex )
1947 add ridge to hashtable without oldvertex
1948
1949 notes:
1950 assumes hashtable is large enough
1951
1952 design:
1953 determine hash value for ridge without oldvertex
1954 find next empty slot for ridge
1955 */
qh_hashridge(setT * hashtable,int hashsize,ridgeT * ridge,vertexT * oldvertex)1956 void qh_hashridge(setT *hashtable, int hashsize, ridgeT *ridge, vertexT *oldvertex) {
1957 int hash;
1958 ridgeT *ridgeA;
1959
1960 hash= qh_gethash(hashsize, ridge->vertices, qh hull_dim-1, 0, oldvertex);
1961 while (True) {
1962 if (!(ridgeA= SETelemt_(hashtable, hash, ridgeT))) {
1963 SETelem_(hashtable, hash)= ridge;
1964 break;
1965 }else if (ridgeA == ridge)
1966 break;
1967 if (++hash == hashsize)
1968 hash= 0;
1969 }
1970 } /* hashridge */
1971
1972
1973 /*-<a href="qh-merge.htm#TOC"
1974 >-------------------------------</a><a name="hashridge_find">-</a>
1975
1976 qh_hashridge_find( hashtable, hashsize, ridge, vertex, oldvertex, hashslot )
1977 returns matching ridge without oldvertex in hashtable
1978 for ridge without vertex
1979 if oldvertex is NULL
1980 matches with any one skip
1981
1982 returns:
1983 matching ridge or NULL
1984 if no match,
1985 if ridge already in table
1986 hashslot= -1
1987 else
1988 hashslot= next NULL index
1989
1990 notes:
1991 assumes hashtable is large enough
1992 can't match ridge to itself
1993
1994 design:
1995 get hash value for ridge without vertex
1996 for each hashslot
1997 return match if ridge matches ridgeA without oldvertex
1998 */
qh_hashridge_find(setT * hashtable,int hashsize,ridgeT * ridge,vertexT * vertex,vertexT * oldvertex,int * hashslot)1999 ridgeT *qh_hashridge_find(setT *hashtable, int hashsize, ridgeT *ridge,
2000 vertexT *vertex, vertexT *oldvertex, int *hashslot) {
2001 int hash;
2002 ridgeT *ridgeA;
2003
2004 *hashslot= 0;
2005 zinc_(Zhashridge);
2006 hash= qh_gethash(hashsize, ridge->vertices, qh hull_dim-1, 0, vertex);
2007 while ((ridgeA= SETelemt_(hashtable, hash, ridgeT))) {
2008 if (ridgeA == ridge)
2009 *hashslot= -1;
2010 else {
2011 zinc_(Zhashridgetest);
2012 if (qh_setequal_except(ridge->vertices, vertex, ridgeA->vertices, oldvertex))
2013 return ridgeA;
2014 }
2015 if (++hash == hashsize)
2016 hash= 0;
2017 }
2018 if (!*hashslot)
2019 *hashslot= hash;
2020 return NULL;
2021 } /* hashridge_find */
2022
2023
2024 /*-<a href="qh-merge.htm#TOC"
2025 >-------------------------------</a><a name="initmergesets">-</a>
2026
2027 qh_initmergesets( )
2028 initialize the merge sets
2029 if 'all', include qh.degen_mergeset
2030
2031 notes:
2032 matches qh_freemergesets
2033 */
qh_initmergesets(void)2034 void qh_initmergesets(void /* qh.facet_mergeset,degen_mergeset,vertex_mergeset */) {
2035
2036 if (qh facet_mergeset || qh degen_mergeset || qh vertex_mergeset) {
2037 qh_fprintf(qh ferr, 6386, "qhull internal error (qh_initmergesets): expecting NULL mergesets. Got qh.facet_mergeset (0x%x), qh.degen_mergeset (0x%x), qh.vertex_mergeset (0x%x)\n",
2038 qh facet_mergeset, qh degen_mergeset, qh vertex_mergeset);
2039 qh_errexit(qh_ERRqhull, NULL, NULL);
2040 }
2041 qh degen_mergeset= qh_settemp(qh TEMPsize);
2042 qh vertex_mergeset= qh_settemp(qh TEMPsize);
2043 qh facet_mergeset= qh_settemp(qh TEMPsize); /* last temporary set for qh_forcedmerges */
2044 } /* initmergesets */
2045
2046 /*-<a href="qh-merge.htm#TOC"
2047 >-------------------------------</a><a name="makeridges">-</a>
2048
2049 qh_makeridges( facet )
2050 creates explicit ridges between simplicial facets
2051
2052 returns:
2053 facet with ridges and without qh_MERGEridge
2054 ->simplicial is False
2055 if facet was tested, new ridges are tested
2056
2057 notes:
2058 allows qh_MERGEridge flag
2059 uses existing ridges
2060 duplicate neighbors ok if ridges already exist (qh_mergecycle_ridges)
2061
2062 see:
2063 qh_mergecycle_ridges()
2064 qh_rename_adjacentvertex for qh_merge_pinchedvertices
2065
2066 design:
2067 look for qh_MERGEridge neighbors
2068 mark neighbors that already have ridges
2069 for each unprocessed neighbor of facet
2070 create a ridge for neighbor and facet
2071 if any qh_MERGEridge neighbors
2072 delete qh_MERGEridge flags (previously processed by qh_mark_dupridges)
2073 */
qh_makeridges(facetT * facet)2074 void qh_makeridges(facetT *facet) {
2075 facetT *neighbor, **neighborp;
2076 ridgeT *ridge, **ridgep;
2077 int neighbor_i, neighbor_n;
2078 boolT toporient, mergeridge= False;
2079
2080 if (!facet->simplicial)
2081 return;
2082 trace4((qh ferr, 4027, "qh_makeridges: make ridges for f%d\n", facet->id));
2083 facet->simplicial= False;
2084 FOREACHneighbor_(facet) {
2085 if (neighbor == qh_MERGEridge)
2086 mergeridge= True;
2087 else
2088 neighbor->seen= False;
2089 }
2090 FOREACHridge_(facet->ridges)
2091 otherfacet_(ridge, facet)->seen= True;
2092 FOREACHneighbor_i_(facet) {
2093 if (neighbor == qh_MERGEridge)
2094 continue; /* fixed by qh_mark_dupridges */
2095 else if (!neighbor->seen) { /* no current ridges */
2096 ridge= qh_newridge();
2097 ridge->vertices= qh_setnew_delnthsorted(facet->vertices, qh hull_dim,
2098 neighbor_i, 0);
2099 toporient= (boolT)(facet->toporient ^ (neighbor_i & 0x1));
2100 if (toporient) {
2101 ridge->top= facet;
2102 ridge->bottom= neighbor;
2103 ridge->simplicialtop= True;
2104 ridge->simplicialbot= neighbor->simplicial;
2105 }else {
2106 ridge->top= neighbor;
2107 ridge->bottom= facet;
2108 ridge->simplicialtop= neighbor->simplicial;
2109 ridge->simplicialbot= True;
2110 }
2111 if (facet->tested && !mergeridge)
2112 ridge->tested= True;
2113 #if 0 /* this also works */
2114 flip= (facet->toporient ^ neighbor->toporient)^(skip1 & 0x1) ^ (skip2 & 0x1);
2115 if (facet->toporient ^ (skip1 & 0x1) ^ flip) {
2116 ridge->top= neighbor;
2117 ridge->bottom= facet;
2118 ridge->simplicialtop= True;
2119 ridge->simplicialbot= neighbor->simplicial;
2120 }else {
2121 ridge->top= facet;
2122 ridge->bottom= neighbor;
2123 ridge->simplicialtop= neighbor->simplicial;
2124 ridge->simplicialbot= True;
2125 }
2126 #endif
2127 qh_setappend(&(facet->ridges), ridge);
2128 trace5((qh ferr, 5005, "makeridges: appended r%d to ridges for f%d. Next is ridges for neighbor f%d\n",
2129 ridge->id, facet->id, neighbor->id));
2130 qh_setappend(&(neighbor->ridges), ridge);
2131 if (qh ridge_id == qh traceridge_id)
2132 qh traceridge= ridge;
2133 }
2134 }
2135 if (mergeridge) {
2136 while (qh_setdel(facet->neighbors, qh_MERGEridge))
2137 ; /* delete each one */
2138 }
2139 } /* makeridges */
2140
2141
2142 /*-<a href="qh-merge.htm#TOC"
2143 >-------------------------------</a><a name="mark_dupridges">-</a>
2144
2145 qh_mark_dupridges( facetlist, allmerges )
2146 add duplicated ridges to qh.facet_mergeset
2147 facet-dupridge is true if it contains a subridge shared by more than one new facet
2148 for each such facet, one has a neighbor marked qh_MERGEridge
2149 allmerges is true if merging dupridges
2150 allmerges is false if merging pinched vertices followed by retry addpoint
2151 qh_mark_dupridges will be called again if pinched vertices not found
2152
2153 returns:
2154 dupridges on qh.facet_mergeset (MRGdupridge)
2155 f.mergeridge and f.mergeridge2 set for facet
2156 f.mergeridge set for neighbor
2157 if allmerges is true
2158 make ridges for facets with dupridges as marked by qh_MERGEridge and both sides facet->dupridge
2159 removes qh_MERGEridge from neighbor sets
2160
2161 notes:
2162 called by qh_premerge and qh_getpinchedmerges
2163 dupridges are due to duplicate subridges
2164 i.e. a subridge occurs in more than two horizon ridges.
2165 i.e., a ridge has more than two neighboring facets
2166 dupridges occur in at least two cases
2167 1) a pinched horizon with nearly adjacent vertices -> merge the vertices (qh_getpinchedmerges)
2168 2) more than one newfacet for a horizon face -> merge coplanar facets (qh_premerge)
2169 qh_matchdupridge previously identified the furthest apart pair of facets to retain
2170 they must have a matching subridge and the same orientation
2171 only way to set facet->mergeridge and mergeridge2
2172 uses qh.visit_id
2173
2174 design:
2175 for all facets on facetlist
2176 if facet contains a dupridge
2177 for each neighbor of facet
2178 if neighbor marked qh_MERGEridge (one side of the merge)
2179 set facet->mergeridge
2180 else
2181 if neighbor contains a dupridge
2182 and the back link is qh_MERGEridge
2183 append dupridge to qh.facet_mergeset
2184 exit if !allmerges for repeating qh_mark_dupridges later
2185 for each dupridge
2186 make ridge sets in preparation for merging
2187 remove qh_MERGEridge from neighbor set
2188 for each dupridge
2189 restore the missing neighbor from the neighbor set that was qh_MERGEridge
2190 add the missing ridge for this neighbor
2191 */
qh_mark_dupridges(facetT * facetlist,boolT allmerges)2192 void qh_mark_dupridges(facetT *facetlist, boolT allmerges) {
2193 facetT *facet, *neighbor, **neighborp;
2194 int nummerge=0;
2195 mergeT *merge, **mergep;
2196
2197 trace4((qh ferr, 4028, "qh_mark_dupridges: identify dupridges in facetlist f%d, allmerges? %d\n",
2198 facetlist->id, allmerges));
2199 FORALLfacet_(facetlist) { /* not necessary for first call */
2200 facet->mergeridge2= False;
2201 facet->mergeridge= False;
2202 }
2203 FORALLfacet_(facetlist) {
2204 if (facet->dupridge) {
2205 FOREACHneighbor_(facet) {
2206 if (neighbor == qh_MERGEridge) {
2207 facet->mergeridge= True;
2208 continue;
2209 }
2210 if (neighbor->dupridge) {
2211 if (!qh_setin(neighbor->neighbors, facet)) { /* i.e., it is qh_MERGEridge, neighbors are distinct */
2212 qh_appendmergeset(facet, neighbor, MRGdupridge, 0.0, 1.0);
2213 facet->mergeridge2= True;
2214 facet->mergeridge= True;
2215 nummerge++;
2216 }else if (qh_setequal(facet->vertices, neighbor->vertices)) { /* neighbors are the same except for horizon and qh_MERGEridge, see QH7085 */
2217 trace3((qh ferr, 3043, "qh_mark_dupridges): dupridge due to duplicate vertices for subridges f%d and f%d\n",
2218 facet->id, neighbor->id));
2219 qh_appendmergeset(facet, neighbor, MRGdupridge, 0.0, 1.0);
2220 facet->mergeridge2= True;
2221 facet->mergeridge= True;
2222 nummerge++;
2223 break; /* same for all neighbors */
2224 }
2225 }
2226 }
2227 }
2228 }
2229 if (!nummerge)
2230 return;
2231 if (!allmerges) {
2232 trace1((qh ferr, 1012, "qh_mark_dupridges: found %d duplicated ridges (MRGdupridge) for qh_getpinchedmerges\n", nummerge));
2233 return;
2234 }
2235 trace1((qh ferr, 1048, "qh_mark_dupridges: found %d duplicated ridges (MRGdupridge) for qh_premerge. Prepare facets for merging\n", nummerge));
2236 /* make ridges in preparation for merging */
2237 FORALLfacet_(facetlist) {
2238 if (facet->mergeridge && !facet->mergeridge2)
2239 qh_makeridges(facet);
2240 }
2241 trace3((qh ferr, 3075, "qh_mark_dupridges: restore missing neighbors and ridges due to qh_MERGEridge\n"));
2242 FOREACHmerge_(qh facet_mergeset) { /* restore the missing neighbors */
2243 if (merge->mergetype == MRGdupridge) { /* only between simplicial facets */
2244 if (merge->facet2->mergeridge2 && qh_setin(merge->facet2->neighbors, merge->facet1)) {
2245 /* Due to duplicate or multiple subridges, e.g., ../eg/qtest.sh t712682 '200 s W1e-13 C1,1e-13 D5' 'd'
2246 merge->facet1: - neighboring facets: f27779 f59186 f59186 f59186 MERGEridge f59186
2247 merge->facet2: - neighboring facets: f27779 f59100 f59100 f59100 f59100 f59100
2248 or, ../eg/qtest.sh 100 '500 s W1e-13 C1,1e-13 D4' 'd'
2249 both facets will be degenerate after merge, consider for special case handling
2250 */
2251 qh_fprintf(qh ferr, 6361, "qhull topological error (qh_mark_dupridges): multiple dupridges for f%d and f%d, including reverse\n",
2252 merge->facet1->id, merge->facet2->id);
2253 qh_errexit2(qh_ERRtopology, merge->facet1, merge->facet2);
2254 }else
2255 qh_setappend(&merge->facet2->neighbors, merge->facet1);
2256 qh_makeridges(merge->facet1); /* and the missing ridges */
2257 }
2258 }
2259 } /* mark_dupridges */
2260
2261 /*-<a href="qh-merge.htm#TOC"
2262 >-------------------------------</a><a name="maybe_duplicateridge">-</a>
2263
2264 qh_maybe_duplicateridge( ridge )
2265 add MRGvertices if neighboring facet has another ridge with the same vertices
2266
2267 returns:
2268 adds rename requests to qh.vertex_mergeset
2269
2270 notes:
2271 called by qh_renamevertex
2272 nop if 2-D
2273 expensive test
2274 Duplicate ridges may lead to new facets with same vertex set (QH7084), will try merging vertices
2275 same as qh_maybe_duplicateridges
2276
2277 design:
2278 for the two neighbors
2279 if non-simplicial
2280 for each ridge with the same first and last vertices (max id and min id)
2281 if the remaining vertices are the same
2282 get the closest pair of vertices
2283 add to vertex_mergeset for merging
2284 */
qh_maybe_duplicateridge(ridgeT * ridgeA)2285 void qh_maybe_duplicateridge(ridgeT *ridgeA) {
2286 ridgeT *ridge, **ridgep;
2287 vertexT *vertex, *pinched;
2288 facetT *neighbor;
2289 coordT dist;
2290 int i, k, last= qh hull_dim-2;
2291
2292 if (qh hull_dim < 3 )
2293 return;
2294
2295 for (neighbor= ridgeA->top, i=0; i<2; neighbor= ridgeA->bottom, i++) {
2296 if (!neighbor->simplicial && neighbor->nummerge > 0) { /* skip degenerate neighbors with both new and old vertices that will be merged */
2297 FOREACHridge_(neighbor->ridges) {
2298 if (ridge != ridgeA && SETfirst_(ridge->vertices) == SETfirst_(ridgeA->vertices)) {
2299 if (SETelem_(ridge->vertices, last) == SETelem_(ridgeA->vertices, last)) {
2300 for (k=1; k<last; k++) {
2301 if (SETelem_(ridge->vertices, k) != SETelem_(ridgeA->vertices, k))
2302 break;
2303 }
2304 if (k == last) {
2305 vertex= qh_findbest_ridgevertex(ridge, &pinched, &dist);
2306 trace2((qh ferr, 2069, "qh_maybe_duplicateridge: will merge v%d into v%d (dist %2.2g) due to duplicate ridges r%d/r%d with the same vertices. mergevertex set\n",
2307 pinched->id, vertex->id, dist, ridgeA->id, ridge->id, ridgeA->top->id, ridgeA->bottom->id, ridge->top->id, ridge->bottom->id));
2308 qh_appendvertexmerge(pinched, vertex, MRGvertices, dist, ridgeA, ridge);
2309 ridge->mergevertex= True; /* disables check for duplicate vertices in qh_checkfacet */
2310 ridgeA->mergevertex= True;
2311 }
2312 }
2313 }
2314 }
2315 }
2316 }
2317 } /* maybe_duplicateridge */
2318
2319 /*-<a href="qh-merge.htm#TOC"
2320 >-------------------------------</a><a name="maybe_duplicateridges">-</a>
2321
2322 qh_maybe_duplicateridges( facet )
2323 if Q15, add MRGvertices if facet has ridges with the same vertices
2324
2325 returns:
2326 adds rename requests to qh.vertex_mergeset
2327
2328 notes:
2329 called at end of qh_mergefacet and qh_mergecycle_all
2330 only enabled if qh.CHECKduplicates ('Q15') and 3-D or more
2331 expensive test, not worth it
2332 same as qh_maybe_duplicateridge
2333
2334 design:
2335 for all ridge pairs in facet
2336 if the same first and last vertices (max id and min id)
2337 if the remaining vertices are the same
2338 get the closest pair of vertices
2339 add to vertex_mergeset for merging
2340 */
qh_maybe_duplicateridges(facetT * facet)2341 void qh_maybe_duplicateridges(facetT *facet) {
2342 facetT *otherfacet;
2343 ridgeT *ridge, *ridge2;
2344 vertexT *vertex, *pinched;
2345 coordT dist;
2346 int ridge_i, ridge_n, i, k, last_v= qh hull_dim-2;
2347
2348 if (qh hull_dim < 3 || !qh CHECKduplicates)
2349 return;
2350
2351 FOREACHridge_i_(facet->ridges) {
2352 otherfacet= otherfacet_(ridge, facet);
2353 if (otherfacet->degenerate || otherfacet->redundant || otherfacet->dupridge || otherfacet->flipped) /* will merge */
2354 continue;
2355 for (i=ridge_i+1; i < ridge_n; i++) {
2356 ridge2= SETelemt_(facet->ridges, i, ridgeT);
2357 otherfacet= otherfacet_(ridge2, facet);
2358 if (otherfacet->degenerate || otherfacet->redundant || otherfacet->dupridge || otherfacet->flipped) /* will merge */
2359 continue;
2360 /* optimize qh_setequal(ridge->vertices, ridge2->vertices) */
2361 if (SETelem_(ridge->vertices, last_v) == SETelem_(ridge2->vertices, last_v)) { /* SETfirst is likely to be the same */
2362 if (SETfirst_(ridge->vertices) == SETfirst_(ridge2->vertices)) {
2363 for (k=1; k<last_v; k++) {
2364 if (SETelem_(ridge->vertices, k) != SETelem_(ridge2->vertices, k))
2365 break;
2366 }
2367 if (k == last_v) {
2368 vertex= qh_findbest_ridgevertex(ridge, &pinched, &dist);
2369 if (ridge->top == ridge2->bottom && ridge->bottom == ridge2->top) {
2370 /* proof that ridges may have opposite orientation */
2371 trace2((qh ferr, 2088, "qh_maybe_duplicateridges: will merge v%d into v%d (dist %2.2g) due to opposite oriented ridges r%d/r%d for f%d and f%d\n",
2372 pinched->id, vertex->id, dist, ridge->id, ridge2->id, ridge->top->id, ridge->bottom->id));
2373 }else {
2374 trace2((qh ferr, 2083, "qh_maybe_duplicateridges: will merge v%d into v%d (dist %2.2g) due to duplicate ridges with the same vertices r%d/r%d in merged facet f%d\n",
2375 pinched->id, vertex->id, dist, ridge->id, ridge2->id, facet->id));
2376 }
2377 qh_appendvertexmerge(pinched, vertex, MRGvertices, dist, ridge, ridge2);
2378 ridge->mergevertex= True; /* disables check for duplicate vertices in qh_checkfacet */
2379 ridge2->mergevertex= True;
2380 }
2381 }
2382 }
2383 }
2384 }
2385 } /* maybe_duplicateridges */
2386
2387 /*-<a href="qh-merge.htm#TOC"
2388 >-------------------------------</a><a name="maydropneighbor">-</a>
2389
2390 qh_maydropneighbor( facet )
2391 drop neighbor relationship if ridge was deleted between a non-simplicial facet and its neighbors
2392
2393 returns:
2394 for deleted ridges
2395 ridges made for simplicial neighbors
2396 neighbor sets updated
2397 appends degenerate facets to qh.facet_mergeset
2398
2399 notes:
2400 called by qh_renamevertex
2401 assumes neighbors do not include qh_MERGEridge (qh_makeridges)
2402 won't cause redundant facets since vertex inclusion is the same
2403 may drop vertex and neighbor if no ridge
2404 uses qh.visit_id
2405
2406 design:
2407 visit all neighbors with ridges
2408 for each unvisited neighbor of facet
2409 delete neighbor and facet from the non-simplicial neighbor sets
2410 if neighbor becomes degenerate
2411 append neighbor to qh.degen_mergeset
2412 if facet is degenerate
2413 append facet to qh.degen_mergeset
2414 */
qh_maydropneighbor(facetT * facet)2415 void qh_maydropneighbor(facetT *facet) {
2416 ridgeT *ridge, **ridgep;
2417 facetT *neighbor, **neighborp;
2418
2419 qh visit_id++;
2420 trace4((qh ferr, 4029, "qh_maydropneighbor: test f%d for no ridges to a neighbor\n",
2421 facet->id));
2422 if (facet->simplicial) {
2423 qh_fprintf(qh ferr, 6278, "qhull internal error (qh_maydropneighbor): not valid for simplicial f%d while adding furthest p%d\n",
2424 facet->id, qh furthest_id);
2425 qh_errexit(qh_ERRqhull, facet, NULL);
2426 }
2427 FOREACHridge_(facet->ridges) {
2428 ridge->top->visitid= qh visit_id;
2429 ridge->bottom->visitid= qh visit_id;
2430 }
2431 FOREACHneighbor_(facet) {
2432 if (neighbor->visible) {
2433 qh_fprintf(qh ferr, 6358, "qhull internal error (qh_maydropneighbor): facet f%d has deleted neighbor f%d (qh.visible_list)\n",
2434 facet->id, neighbor->id);
2435 qh_errexit2(qh_ERRqhull, facet, neighbor);
2436 }
2437 if (neighbor->visitid != qh visit_id) {
2438 trace2((qh ferr, 2104, "qh_maydropneighbor: facets f%d and f%d are no longer neighbors while adding furthest p%d\n",
2439 facet->id, neighbor->id, qh furthest_id));
2440 if (neighbor->simplicial) {
2441 qh_fprintf(qh ferr, 6280, "qhull internal error (qh_maydropneighbor): not valid for simplicial neighbor f%d of f%d while adding furthest p%d\n",
2442 neighbor->id, facet->id, qh furthest_id);
2443 qh_errexit2(qh_ERRqhull, neighbor, facet);
2444 }
2445 zinc_(Zdropneighbor);
2446 qh_setdel(neighbor->neighbors, facet);
2447 if (qh_setsize(neighbor->neighbors) < qh hull_dim) {
2448 zinc_(Zdropdegen);
2449 qh_appendmergeset(neighbor, neighbor, MRGdegen, 0.0, qh_ANGLEnone);
2450 trace2((qh ferr, 2023, "qh_maydropneighbors: f%d is degenerate.\n", neighbor->id));
2451 }
2452 qh_setdel(facet->neighbors, neighbor);
2453 neighborp--; /* repeat, deleted a neighbor */
2454 }
2455 }
2456 if (qh_setsize(facet->neighbors) < qh hull_dim) {
2457 zinc_(Zdropdegen);
2458 qh_appendmergeset(facet, facet, MRGdegen, 0.0, qh_ANGLEnone);
2459 trace2((qh ferr, 2024, "qh_maydropneighbors: f%d is degenerate.\n", facet->id));
2460 }
2461 } /* maydropneighbor */
2462
2463
2464 /*-<a href="qh-merge.htm#TOC"
2465 >-------------------------------</a><a name="merge_degenredundant">-</a>
2466
2467 qh_merge_degenredundant( )
2468 merge all degenerate and redundant facets
2469 qh.degen_mergeset contains merges from qh_test_degen_neighbors, qh_test_redundant_neighbors, and qh_degen_redundant_facet
2470
2471 returns:
2472 number of merges performed
2473 resets facet->degenerate/redundant
2474 if deleted (visible) facet has no neighbors
2475 sets ->f.replace to NULL
2476
2477 notes:
2478 redundant merges happen before degenerate ones
2479 merging and renaming vertices can result in degen/redundant facets
2480 check for coplanar and convex neighbors afterwards
2481
2482 design:
2483 for each merge on qh.degen_mergeset
2484 if redundant merge
2485 if non-redundant facet merged into redundant facet
2486 recheck facet for redundancy
2487 else
2488 merge redundant facet into other facet
2489 */
qh_merge_degenredundant(void)2490 int qh_merge_degenredundant(void) {
2491 int size;
2492 mergeT *merge;
2493 facetT *bestneighbor, *facet1, *facet2, *facet3;
2494 realT dist, mindist, maxdist;
2495 vertexT *vertex, **vertexp;
2496 int nummerges= 0;
2497 mergeType mergetype;
2498 setT *mergedfacets;
2499
2500 trace2((qh ferr, 2095, "qh_merge_degenredundant: merge %d degenerate, redundant, and mirror facets\n",
2501 qh_setsize(qh degen_mergeset)));
2502 mergedfacets= qh_settemp(qh TEMPsize);
2503 while ((merge= (mergeT *)qh_setdellast(qh degen_mergeset))) {
2504 facet1= merge->facet1;
2505 facet2= merge->facet2;
2506 mergetype= merge->mergetype;
2507 qh_memfree(merge, (int)sizeof(mergeT)); /* 'merge' is invalidated */
2508 if (facet1->visible)
2509 continue;
2510 facet1->degenerate= False;
2511 facet1->redundant= False;
2512 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
2513 qhmem.IStracing= qh IStracing= qh TRACElevel;
2514 if (mergetype == MRGredundant) {
2515 zinc_(Zredundant);
2516 facet3= qh_getreplacement(facet2); /* the same facet if !facet2.visible */
2517 if (!facet3) {
2518 qh_fprintf(qh ferr, 6097, "qhull internal error (qh_merge_degenredunant): f%d is redundant but visible f%d has no replacement\n",
2519 facet1->id, getid_(facet2));
2520 qh_errexit2(qh_ERRqhull, facet1, facet2);
2521 }
2522 qh_setunique(&mergedfacets, facet3);
2523 if (facet1 == facet3) {
2524 continue;
2525 }
2526 trace2((qh ferr, 2025, "qh_merge_degenredundant: merge redundant f%d into f%d (arg f%d)\n",
2527 facet1->id, facet3->id, facet2->id));
2528 qh_mergefacet(facet1, facet3, mergetype, NULL, NULL, !qh_MERGEapex);
2529 /* merge distance is already accounted for */
2530 nummerges++;
2531 }else { /* mergetype == MRGdegen or MRGmirror, other merges may have fixed */
2532 if (!(size= qh_setsize(facet1->neighbors))) {
2533 zinc_(Zdelfacetdup);
2534 trace2((qh ferr, 2026, "qh_merge_degenredundant: facet f%d has no neighbors. Deleted\n", facet1->id));
2535 qh_willdelete(facet1, NULL);
2536 FOREACHvertex_(facet1->vertices) {
2537 qh_setdel(vertex->neighbors, facet1);
2538 if (!SETfirst_(vertex->neighbors)) {
2539 zinc_(Zdegenvertex);
2540 trace2((qh ferr, 2027, "qh_merge_degenredundant: deleted v%d because f%d has no neighbors\n",
2541 vertex->id, facet1->id));
2542 vertex->deleted= True;
2543 qh_setappend(&qh del_vertices, vertex);
2544 }
2545 }
2546 nummerges++;
2547 }else if (size < qh hull_dim) {
2548 bestneighbor= qh_findbestneighbor(facet1, &dist, &mindist, &maxdist);
2549 trace2((qh ferr, 2028, "qh_merge_degenredundant: facet f%d has %d neighbors, merge into f%d dist %2.2g\n",
2550 facet1->id, size, bestneighbor->id, dist));
2551 qh_mergefacet(facet1, bestneighbor, mergetype, &mindist, &maxdist, !qh_MERGEapex);
2552 nummerges++;
2553 if (qh PRINTstatistics) {
2554 zinc_(Zdegen);
2555 wadd_(Wdegentot, dist);
2556 wmax_(Wdegenmax, dist);
2557 }
2558 } /* else, another merge fixed the degeneracy and redundancy tested */
2559 }
2560 }
2561 qh_settempfree(&mergedfacets);
2562 return nummerges;
2563 } /* merge_degenredundant */
2564
2565 /*-<a href="qh-merge.htm#TOC"
2566 >-------------------------------</a><a name="merge_nonconvex">-</a>
2567
2568 qh_merge_nonconvex( facet1, facet2, mergetype )
2569 remove non-convex ridge between facet1 into facet2
2570 mergetype gives why the facet's are non-convex
2571
2572 returns:
2573 merges one of the facets into the best neighbor
2574
2575 notes:
2576 mergetype is MRGcoplanar..MRGconvex
2577
2578 design:
2579 if one of the facets is a new facet
2580 prefer merging new facet into old facet
2581 find best neighbors for both facets
2582 merge the nearest facet into its best neighbor
2583 update the statistics
2584 */
qh_merge_nonconvex(facetT * facet1,facetT * facet2,mergeType mergetype)2585 void qh_merge_nonconvex(facetT *facet1, facetT *facet2, mergeType mergetype) {
2586 facetT *bestfacet, *bestneighbor, *neighbor, *merging, *merged;
2587 realT dist, dist2, mindist, mindist2, maxdist, maxdist2;
2588
2589 if (mergetype < MRGcoplanar || mergetype > MRGconcavecoplanar) {
2590 qh_fprintf(qh ferr, 6398, "qhull internal error (qh_merge_nonconvex): expecting mergetype MRGcoplanar..MRGconcavecoplanar. Got merge f%d and f%d type %d\n",
2591 facet1->id, facet2->id, mergetype);
2592 qh_errexit2(qh_ERRqhull, facet1, facet2);
2593 }
2594 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
2595 qhmem.IStracing= qh IStracing= qh TRACElevel;
2596 trace3((qh ferr, 3003, "qh_merge_nonconvex: merge #%d for f%d and f%d type %d\n",
2597 zzval_(Ztotmerge) + 1, facet1->id, facet2->id, mergetype));
2598 /* concave or coplanar */
2599 if (!facet1->newfacet) {
2600 bestfacet= facet2; /* avoid merging old facet if new is ok */
2601 facet2= facet1;
2602 facet1= bestfacet;
2603 }else
2604 bestfacet= facet1;
2605 bestneighbor= qh_findbestneighbor(bestfacet, &dist, &mindist, &maxdist);
2606 neighbor= qh_findbestneighbor(facet2, &dist2, &mindist2, &maxdist2);
2607 if (dist < dist2) {
2608 merging= bestfacet;
2609 merged= bestneighbor;
2610 }else if (qh AVOIDold && !facet2->newfacet
2611 && ((mindist >= -qh MAXcoplanar && maxdist <= qh max_outside)
2612 || dist * 1.5 < dist2)) {
2613 zinc_(Zavoidold);
2614 wadd_(Wavoidoldtot, dist);
2615 wmax_(Wavoidoldmax, dist);
2616 trace2((qh ferr, 2029, "qh_merge_nonconvex: avoid merging old facet f%d dist %2.2g. Use f%d dist %2.2g instead\n",
2617 facet2->id, dist2, facet1->id, dist2));
2618 merging= bestfacet;
2619 merged= bestneighbor;
2620 }else {
2621 merging= facet2;
2622 merged= neighbor;
2623 dist= dist2;
2624 mindist= mindist2;
2625 maxdist= maxdist2;
2626 }
2627 qh_mergefacet(merging, merged, mergetype, &mindist, &maxdist, !qh_MERGEapex);
2628 /* caller merges qh_degenredundant */
2629 if (qh PRINTstatistics) {
2630 if (mergetype == MRGanglecoplanar) {
2631 zinc_(Zacoplanar);
2632 wadd_(Wacoplanartot, dist);
2633 wmax_(Wacoplanarmax, dist);
2634 }else if (mergetype == MRGconcave) {
2635 zinc_(Zconcave);
2636 wadd_(Wconcavetot, dist);
2637 wmax_(Wconcavemax, dist);
2638 }else if (mergetype == MRGconcavecoplanar) {
2639 zinc_(Zconcavecoplanar);
2640 wadd_(Wconcavecoplanartot, dist);
2641 wmax_(Wconcavecoplanarmax, dist);
2642 }else { /* MRGcoplanar */
2643 zinc_(Zcoplanar);
2644 wadd_(Wcoplanartot, dist);
2645 wmax_(Wcoplanarmax, dist);
2646 }
2647 }
2648 } /* merge_nonconvex */
2649
2650 /*-<a href="qh-merge.htm#TOC"
2651 >-------------------------------</a><a name="merge_pinchedvertices">-</a>
2652
2653 qh_merge_pinchedvertices( apex )
2654 merge pinched vertices in qh.vertex_mergeset to avoid qh_forcedmerges of dupridges
2655
2656 notes:
2657 only called by qh_all_vertexmerges
2658 hull_dim >= 3
2659
2660 design:
2661 make vertex neighbors if necessary
2662 for each pinched vertex
2663 determine the ridges for the pinched vertex (make ridges as needed)
2664 merge the pinched vertex into the horizon vertex
2665 merge the degenerate and redundant facets that result
2666 check and resolve new dupridges
2667 */
qh_merge_pinchedvertices(int apexpointid)2668 void qh_merge_pinchedvertices(int apexpointid /* qh.newfacet_list */) {
2669 mergeT *merge, *mergeA, **mergeAp;
2670 vertexT *vertex, *vertex2;
2671 realT dist;
2672 boolT firstmerge= True;
2673
2674 qh_vertexneighbors();
2675 if (qh visible_list || qh newfacet_list || qh newvertex_list) {
2676 qh_fprintf(qh ferr, 6402, "qhull internal error (qh_merge_pinchedvertices): qh.visible_list (f%d), newfacet_list (f%d), or newvertex_list (v%d) not empty\n",
2677 getid_(qh visible_list), getid_(qh newfacet_list), getid_(qh newvertex_list));
2678 qh_errexit(qh_ERRqhull, NULL, NULL);
2679 }
2680 qh visible_list= qh newfacet_list= qh facet_tail;
2681 qh newvertex_list= qh vertex_tail;
2682 qh isRenameVertex= True; /* disable duplicate ridge vertices check in qh_checkfacet */
2683 while ((merge= qh_next_vertexmerge(/* qh.vertex_mergeset */))) { /* only one at a time from qh_getpinchedmerges */
2684 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
2685 qhmem.IStracing= qh IStracing= qh TRACElevel;
2686 if (merge->mergetype == MRGsubridge) {
2687 zzinc_(Zpinchedvertex);
2688 trace1((qh ferr, 1050, "qh_merge_pinchedvertices: merge one of %d pinched vertices before adding apex p%d. Try to resolve duplicate ridges in newfacets\n",
2689 qh_setsize(qh vertex_mergeset)+1, apexpointid));
2690 qh_remove_mergetype(qh vertex_mergeset, MRGsubridge);
2691 }else {
2692 zzinc_(Zpinchduplicate);
2693 if (firstmerge)
2694 trace1((qh ferr, 1056, "qh_merge_pinchedvertices: merge %d pinched vertices from dupridges in merged facets, apex p%d\n",
2695 qh_setsize(qh vertex_mergeset)+1, apexpointid));
2696 firstmerge= False;
2697 }
2698 vertex= merge->vertex1;
2699 vertex2= merge->vertex2;
2700 dist= merge->distance;
2701 qh_memfree(merge, (int)sizeof(mergeT)); /* merge is invalidated */
2702 qh_rename_adjacentvertex(vertex, vertex2, dist);
2703 #ifndef qh_NOtrace
2704 if (qh IStracing >= 2) {
2705 FOREACHmergeA_(qh degen_mergeset) {
2706 if (mergeA->mergetype== MRGdegen) {
2707 qh_fprintf(qh ferr, 2072, "qh_merge_pinchedvertices: merge degenerate f%d into an adjacent facet\n", mergeA->facet1->id);
2708 }else {
2709 qh_fprintf(qh ferr, 2084, "qh_merge_pinchedvertices: merge f%d into f%d mergeType %d\n", mergeA->facet1->id, mergeA->facet2->id, mergeA->mergetype);
2710 }
2711 }
2712 }
2713 #endif
2714 qh_merge_degenredundant(); /* simplicial facets with both old and new vertices */
2715 }
2716 qh isRenameVertex= False;
2717 }/* merge_pinchedvertices */
2718
2719 /*-<a href="qh-merge.htm#TOC"
2720 >-------------------------------</a><a name="merge_twisted">-</a>
2721
2722 qh_merge_twisted( facet1, facet2 )
2723 remove twisted ridge between facet1 into facet2 or report error
2724
2725 returns:
2726 merges one of the facets into the best neighbor
2727
2728 notes:
2729 a twisted ridge has opposite vertices that are convex and concave
2730
2731 design:
2732 find best neighbors for both facets
2733 error if wide merge
2734 merge the nearest facet into its best neighbor
2735 update statistics
2736 */
qh_merge_twisted(facetT * facet1,facetT * facet2)2737 void qh_merge_twisted(facetT *facet1, facetT *facet2) {
2738 facetT *neighbor2, *neighbor, *merging, *merged;
2739 vertexT *bestvertex, *bestpinched;
2740 realT dist, dist2, mindist, mindist2, maxdist, maxdist2, mintwisted, bestdist;
2741
2742 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
2743 qhmem.IStracing= qh IStracing= qh TRACElevel;
2744 trace3((qh ferr, 3050, "qh_merge_twisted: merge #%d for twisted f%d and f%d\n",
2745 zzval_(Ztotmerge) + 1, facet1->id, facet2->id));
2746 /* twisted */
2747 neighbor= qh_findbestneighbor(facet1, &dist, &mindist, &maxdist);
2748 neighbor2= qh_findbestneighbor(facet2, &dist2, &mindist2, &maxdist2);
2749 mintwisted= qh_RATIOtwisted * qh ONEmerge;
2750 maximize_(mintwisted, facet1->maxoutside);
2751 maximize_(mintwisted, facet2->maxoutside);
2752 if (dist > mintwisted && dist2 > mintwisted) {
2753 bestdist= qh_vertex_bestdist2(facet1->vertices, &bestvertex, &bestpinched);
2754 if (bestdist > mintwisted) {
2755 qh_fprintf(qh ferr, 6417, "qhull precision error (qh_merge_twisted): twisted facet f%d does not contain pinched vertices. Too wide to merge into neighbor. mindist %2.2g maxdist %2.2g vertexdist %2.2g maxpinched %2.2g neighbor f%d mindist %2.2g maxdist %2.2g\n",
2756 facet1->id, mindist, maxdist, bestdist, mintwisted, facet2->id, mindist2, maxdist2);
2757 }else {
2758 qh_fprintf(qh ferr, 6418, "qhull precision error (qh_merge_twisted): twisted facet f%d with pinched vertices. Could merge vertices, but too wide to merge into neighbor. mindist %2.2g maxdist %2.2g vertexdist %2.2g neighbor f%d mindist %2.2g maxdist %2.2g\n",
2759 facet1->id, mindist, maxdist, bestdist, facet2->id, mindist2, maxdist2);
2760 }
2761 qh_errexit2(qh_ERRwide, facet1, facet2);
2762 }
2763 if (dist < dist2) {
2764 merging= facet1;
2765 merged= neighbor;
2766 }else {
2767 /* ignores qh.AVOIDold ('Q4') */
2768 merging= facet2;
2769 merged= neighbor2;
2770 dist= dist2;
2771 mindist= mindist2;
2772 maxdist= maxdist2;
2773 }
2774 qh_mergefacet(merging, merged, MRGtwisted, &mindist, &maxdist, !qh_MERGEapex);
2775 /* caller merges qh_degenredundant */
2776 zinc_(Ztwisted);
2777 wadd_(Wtwistedtot, dist);
2778 wmax_(Wtwistedmax, dist);
2779 } /* merge_twisted */
2780
2781 /*-<a href="qh-merge.htm#TOC"
2782 >-------------------------------</a><a name="mergecycle">-</a>
2783
2784 qh_mergecycle( samecycle, newfacet )
2785 merge a cycle of facets starting at samecycle into a newfacet
2786 newfacet is a horizon facet with ->normal
2787 samecycle facets are simplicial from an apex
2788
2789 returns:
2790 initializes vertex neighbors on first merge
2791 samecycle deleted (placed on qh.visible_list)
2792 newfacet at end of qh.facet_list
2793 deleted vertices on qh.del_vertices
2794
2795 notes:
2796 only called by qh_mergecycle_all for multiple, same cycle facets
2797 see qh_mergefacet
2798
2799 design:
2800 make vertex neighbors if necessary
2801 make ridges for newfacet
2802 merge neighbor sets of samecycle into newfacet
2803 merge ridges of samecycle into newfacet
2804 merge vertex neighbors of samecycle into newfacet
2805 make apex of samecycle the apex of newfacet
2806 if newfacet wasn't a new facet
2807 add its vertices to qh.newvertex_list
2808 delete samecycle facets a make newfacet a newfacet
2809 */
qh_mergecycle(facetT * samecycle,facetT * newfacet)2810 void qh_mergecycle(facetT *samecycle, facetT *newfacet) {
2811 int traceonce= False, tracerestore= 0;
2812 vertexT *apex;
2813 #ifndef qh_NOtrace
2814 facetT *same;
2815 #endif
2816
2817 zzinc_(Ztotmerge);
2818 if (qh REPORTfreq2 && qh POSTmerging) {
2819 if (zzval_(Ztotmerge) > qh mergereport + qh REPORTfreq2)
2820 qh_tracemerging();
2821 }
2822 #ifndef qh_NOtrace
2823 if (qh TRACEmerge == zzval_(Ztotmerge))
2824 qhmem.IStracing= qh IStracing= qh TRACElevel;
2825 trace2((qh ferr, 2030, "qh_mergecycle: merge #%d for facets from cycle f%d into coplanar horizon f%d\n",
2826 zzval_(Ztotmerge), samecycle->id, newfacet->id));
2827 if (newfacet == qh tracefacet) {
2828 tracerestore= qh IStracing;
2829 qh IStracing= 4;
2830 qh_fprintf(qh ferr, 8068, "qh_mergecycle: ========= trace merge %d of samecycle %d into trace f%d, furthest is p%d\n",
2831 zzval_(Ztotmerge), samecycle->id, newfacet->id, qh furthest_id);
2832 traceonce= True;
2833 }
2834 if (qh IStracing >=4) {
2835 qh_fprintf(qh ferr, 8069, " same cycle:");
2836 FORALLsame_cycle_(samecycle)
2837 qh_fprintf(qh ferr, 8070, " f%d", same->id);
2838 qh_fprintf(qh ferr, 8071, "\n");
2839 }
2840 if (qh IStracing >=4)
2841 qh_errprint("MERGING CYCLE", samecycle, newfacet, NULL, NULL);
2842 #endif /* !qh_NOtrace */
2843 if (newfacet->tricoplanar) {
2844 if (!qh TRInormals) {
2845 qh_fprintf(qh ferr, 6224, "qhull internal error (qh_mergecycle): does not work for tricoplanar facets. Use option 'Q11'\n");
2846 qh_errexit(qh_ERRqhull, newfacet, NULL);
2847 }
2848 newfacet->tricoplanar= False;
2849 newfacet->keepcentrum= False;
2850 }
2851 if (qh CHECKfrequently)
2852 qh_checkdelridge();
2853 if (!qh VERTEXneighbors)
2854 qh_vertexneighbors();
2855 apex= SETfirstt_(samecycle->vertices, vertexT);
2856 qh_makeridges(newfacet);
2857 qh_mergecycle_neighbors(samecycle, newfacet);
2858 qh_mergecycle_ridges(samecycle, newfacet);
2859 qh_mergecycle_vneighbors(samecycle, newfacet);
2860 if (SETfirstt_(newfacet->vertices, vertexT) != apex)
2861 qh_setaddnth(&newfacet->vertices, 0, apex); /* apex has last id */
2862 if (!newfacet->newfacet)
2863 qh_newvertices(newfacet->vertices);
2864 qh_mergecycle_facets(samecycle, newfacet);
2865 qh_tracemerge(samecycle, newfacet, MRGcoplanarhorizon);
2866 /* check for degen_redundant_neighbors after qh_forcedmerges() */
2867 if (traceonce) {
2868 qh_fprintf(qh ferr, 8072, "qh_mergecycle: end of trace facet\n");
2869 qh IStracing= tracerestore;
2870 }
2871 } /* mergecycle */
2872
2873 /*-<a href="qh-merge.htm#TOC"
2874 >-------------------------------</a><a name="mergecycle_all">-</a>
2875
2876 qh_mergecycle_all( facetlist, wasmerge )
2877 merge all samecycles of coplanar facets into horizon
2878 don't merge facets with ->mergeridge (these already have ->normal)
2879 all facets are simplicial from apex
2880 all facet->cycledone == False
2881
2882 returns:
2883 all newfacets merged into coplanar horizon facets
2884 deleted vertices on qh.del_vertices
2885 sets wasmerge if any merge
2886
2887 notes:
2888 called by qh_premerge
2889 calls qh_mergecycle for multiple, same cycle facets
2890
2891 design:
2892 for each facet on facetlist
2893 skip facets with dupridges and normals
2894 check that facet is in a samecycle (->mergehorizon)
2895 if facet only member of samecycle
2896 sets vertex->delridge for all vertices except apex
2897 merge facet into horizon
2898 else
2899 mark all facets in samecycle
2900 remove facets with dupridges from samecycle
2901 merge samecycle into horizon (deletes facets from facetlist)
2902 */
qh_mergecycle_all(facetT * facetlist,boolT * wasmerge)2903 void qh_mergecycle_all(facetT *facetlist, boolT *wasmerge) {
2904 facetT *facet, *same, *prev, *horizon, *newfacet;
2905 facetT *samecycle= NULL, *nextfacet, *nextsame;
2906 vertexT *apex, *vertex, **vertexp;
2907 int cycles=0, total=0, facets, nummerge, numdegen= 0;
2908
2909 trace2((qh ferr, 2031, "qh_mergecycle_all: merge new facets into coplanar horizon facets. Bulk merge a cycle of facets with the same horizon facet\n"));
2910 for (facet=facetlist; facet && (nextfacet= facet->next); facet= nextfacet) {
2911 if (facet->normal)
2912 continue;
2913 if (!facet->mergehorizon) {
2914 qh_fprintf(qh ferr, 6225, "qhull internal error (qh_mergecycle_all): f%d without normal\n", facet->id);
2915 qh_errexit(qh_ERRqhull, facet, NULL);
2916 }
2917 horizon= SETfirstt_(facet->neighbors, facetT);
2918 if (facet->f.samecycle == facet) {
2919 if (qh TRACEmerge-1 == zzval_(Ztotmerge))
2920 qhmem.IStracing= qh IStracing= qh TRACElevel;
2921 zinc_(Zonehorizon);
2922 /* merge distance done in qh_findhorizon */
2923 apex= SETfirstt_(facet->vertices, vertexT);
2924 FOREACHvertex_(facet->vertices) {
2925 if (vertex != apex)
2926 vertex->delridge= True;
2927 }
2928 horizon->f.newcycle= NULL;
2929 qh_mergefacet(facet, horizon, MRGcoplanarhorizon, NULL, NULL, qh_MERGEapex);
2930 }else {
2931 samecycle= facet;
2932 facets= 0;
2933 prev= facet;
2934 for (same= facet->f.samecycle; same; /* FORALLsame_cycle_(facet) */
2935 same= (same == facet ? NULL :nextsame)) { /* ends at facet */
2936 nextsame= same->f.samecycle;
2937 if (same->cycledone || same->visible)
2938 qh_infiniteloop(same);
2939 same->cycledone= True;
2940 if (same->normal) {
2941 prev->f.samecycle= same->f.samecycle; /* unlink ->mergeridge */
2942 same->f.samecycle= NULL;
2943 }else {
2944 prev= same;
2945 facets++;
2946 }
2947 }
2948 while (nextfacet && nextfacet->cycledone) /* will delete samecycle */
2949 nextfacet= nextfacet->next;
2950 horizon->f.newcycle= NULL;
2951 qh_mergecycle(samecycle, horizon);
2952 nummerge= horizon->nummerge + facets;
2953 if (nummerge > qh_MAXnummerge)
2954 horizon->nummerge= qh_MAXnummerge;
2955 else
2956 horizon->nummerge= (short unsigned int)nummerge; /* limited to 9 bits by qh_MAXnummerge, -Wconversion */
2957 zzinc_(Zcyclehorizon);
2958 total += facets;
2959 zzadd_(Zcyclefacettot, facets);
2960 zmax_(Zcyclefacetmax, facets);
2961 }
2962 cycles++;
2963 }
2964 if (cycles) {
2965 FORALLnew_facets {
2966 /* qh_maybe_duplicateridges postponed since qh_mergecycle_ridges deletes ridges without calling qh_delridge_merge */
2967 if (newfacet->coplanarhorizon) {
2968 qh_test_redundant_neighbors(newfacet);
2969 qh_maybe_duplicateridges(newfacet);
2970 newfacet->coplanarhorizon= False;
2971 }
2972 }
2973 numdegen += qh_merge_degenredundant();
2974 *wasmerge= True;
2975 trace1((qh ferr, 1013, "qh_mergecycle_all: merged %d same cycles or facets into coplanar horizons and %d degenredundant facets\n",
2976 cycles, numdegen));
2977 }
2978 } /* mergecycle_all */
2979
2980 /*-<a href="qh-merge.htm#TOC"
2981 >-------------------------------</a><a name="mergecycle_facets">-</a>
2982
2983 qh_mergecycle_facets( samecycle, newfacet )
2984 finish merge of samecycle into newfacet
2985
2986 returns:
2987 samecycle prepended to visible_list for later deletion and partitioning
2988 each facet->f.replace == newfacet
2989
2990 newfacet moved to end of qh.facet_list
2991 makes newfacet a newfacet (get's facet1->id if it was old)
2992 sets newfacet->newmerge
2993 clears newfacet->center (unless merging into a large facet)
2994 clears newfacet->tested and ridge->tested for facet1
2995
2996 adds neighboring facets to facet_mergeset if redundant or degenerate
2997
2998 design:
2999 make newfacet a new facet and set its flags
3000 move samecycle facets to qh.visible_list for later deletion
3001 unless newfacet is large
3002 remove its centrum
3003 */
qh_mergecycle_facets(facetT * samecycle,facetT * newfacet)3004 void qh_mergecycle_facets(facetT *samecycle, facetT *newfacet) {
3005 facetT *same, *next;
3006
3007 trace4((qh ferr, 4030, "qh_mergecycle_facets: make newfacet new and samecycle deleted\n"));
3008 qh_removefacet(newfacet); /* append as a newfacet to end of qh facet_list */
3009 qh_appendfacet(newfacet);
3010 newfacet->newfacet= True;
3011 newfacet->simplicial= False;
3012 newfacet->newmerge= True;
3013
3014 for (same= samecycle->f.samecycle; same; same= (same == samecycle ? NULL : next)) {
3015 next= same->f.samecycle; /* reused by willdelete */
3016 qh_willdelete(same, newfacet);
3017 }
3018 if (newfacet->center
3019 && qh_setsize(newfacet->vertices) <= qh hull_dim + qh_MAXnewcentrum) {
3020 qh_memfree(newfacet->center, qh normal_size);
3021 newfacet->center= NULL;
3022 }
3023 trace3((qh ferr, 3004, "qh_mergecycle_facets: merged facets from cycle f%d into f%d\n",
3024 samecycle->id, newfacet->id));
3025 } /* mergecycle_facets */
3026
3027 /*-<a href="qh-merge.htm#TOC"
3028 >-------------------------------</a><a name="mergecycle_neighbors">-</a>
3029
3030 qh_mergecycle_neighbors( samecycle, newfacet )
3031 add neighbors for samecycle facets to newfacet
3032
3033 returns:
3034 newfacet with updated neighbors and vice-versa
3035 newfacet has ridges
3036 all neighbors of newfacet marked with qh.visit_id
3037 samecycle facets marked with qh.visit_id-1
3038 ridges updated for simplicial neighbors of samecycle with a ridge
3039
3040 notes:
3041 assumes newfacet not in samecycle
3042 usually, samecycle facets are new, simplicial facets without internal ridges
3043 not so if horizon facet is coplanar to two different samecycles
3044
3045 see:
3046 qh_mergeneighbors()
3047
3048 design:
3049 check samecycle
3050 delete neighbors from newfacet that are also in samecycle
3051 for each neighbor of a facet in samecycle
3052 if neighbor is simplicial
3053 if first visit
3054 move the neighbor relation to newfacet
3055 update facet links for its ridges
3056 else
3057 make ridges for neighbor
3058 remove samecycle reference
3059 else
3060 update neighbor sets
3061 */
qh_mergecycle_neighbors(facetT * samecycle,facetT * newfacet)3062 void qh_mergecycle_neighbors(facetT *samecycle, facetT *newfacet) {
3063 facetT *same, *neighbor, **neighborp;
3064 int delneighbors= 0, newneighbors= 0;
3065 unsigned int samevisitid;
3066 ridgeT *ridge, **ridgep;
3067
3068 samevisitid= ++qh visit_id;
3069 FORALLsame_cycle_(samecycle) {
3070 if (same->visitid == samevisitid || same->visible)
3071 qh_infiniteloop(samecycle);
3072 same->visitid= samevisitid;
3073 }
3074 newfacet->visitid= ++qh visit_id;
3075 trace4((qh ferr, 4031, "qh_mergecycle_neighbors: delete shared neighbors from newfacet\n"));
3076 FOREACHneighbor_(newfacet) {
3077 if (neighbor->visitid == samevisitid) {
3078 SETref_(neighbor)= NULL; /* samecycle neighbors deleted */
3079 delneighbors++;
3080 }else
3081 neighbor->visitid= qh visit_id;
3082 }
3083 qh_setcompact(newfacet->neighbors);
3084
3085 trace4((qh ferr, 4032, "qh_mergecycle_neighbors: update neighbors\n"));
3086 FORALLsame_cycle_(samecycle) {
3087 FOREACHneighbor_(same) {
3088 if (neighbor->visitid == samevisitid)
3089 continue;
3090 if (neighbor->simplicial) {
3091 if (neighbor->visitid != qh visit_id) {
3092 qh_setappend(&newfacet->neighbors, neighbor);
3093 qh_setreplace(neighbor->neighbors, same, newfacet);
3094 newneighbors++;
3095 neighbor->visitid= qh visit_id;
3096 FOREACHridge_(neighbor->ridges) { /* update ridge in case of qh_makeridges */
3097 if (ridge->top == same) {
3098 ridge->top= newfacet;
3099 break;
3100 }else if (ridge->bottom == same) {
3101 ridge->bottom= newfacet;
3102 break;
3103 }
3104 }
3105 }else {
3106 qh_makeridges(neighbor);
3107 qh_setdel(neighbor->neighbors, same);
3108 /* same can't be horizon facet for neighbor */
3109 }
3110 }else { /* non-simplicial neighbor */
3111 qh_setdel(neighbor->neighbors, same);
3112 if (neighbor->visitid != qh visit_id) {
3113 qh_setappend(&neighbor->neighbors, newfacet);
3114 qh_setappend(&newfacet->neighbors, neighbor);
3115 neighbor->visitid= qh visit_id;
3116 newneighbors++;
3117 }
3118 }
3119 }
3120 }
3121 trace2((qh ferr, 2032, "qh_mergecycle_neighbors: deleted %d neighbors and added %d\n",
3122 delneighbors, newneighbors));
3123 } /* mergecycle_neighbors */
3124
3125 /*-<a href="qh-merge.htm#TOC"
3126 >-------------------------------</a><a name="mergecycle_ridges">-</a>
3127
3128 qh_mergecycle_ridges( samecycle, newfacet )
3129 add ridges/neighbors for facets in samecycle to newfacet
3130 all new/old neighbors of newfacet marked with qh.visit_id
3131 facets in samecycle marked with qh.visit_id-1
3132 newfacet marked with qh.visit_id
3133
3134 returns:
3135 newfacet has merged ridges
3136
3137 notes:
3138 ridge already updated for simplicial neighbors of samecycle with a ridge
3139 qh_checkdelridge called by qh_mergecycle
3140
3141 see:
3142 qh_mergeridges()
3143 qh_makeridges()
3144
3145 design:
3146 remove ridges between newfacet and samecycle
3147 for each facet in samecycle
3148 for each ridge in facet
3149 update facet pointers in ridge
3150 skip ridges processed in qh_mergecycle_neighors
3151 free ridges between newfacet and samecycle
3152 free ridges between facets of samecycle (on 2nd visit)
3153 append remaining ridges to newfacet
3154 if simplicial facet
3155 for each neighbor of facet
3156 if simplicial facet
3157 and not samecycle facet or newfacet
3158 make ridge between neighbor and newfacet
3159 */
qh_mergecycle_ridges(facetT * samecycle,facetT * newfacet)3160 void qh_mergecycle_ridges(facetT *samecycle, facetT *newfacet) {
3161 facetT *same, *neighbor= NULL;
3162 int numold=0, numnew=0;
3163 int neighbor_i, neighbor_n;
3164 unsigned int samevisitid;
3165 ridgeT *ridge, **ridgep;
3166 boolT toporient;
3167 void **freelistp; /* used if !qh_NOmem by qh_memfree_() */
3168
3169 trace4((qh ferr, 4033, "qh_mergecycle_ridges: delete shared ridges from newfacet\n"));
3170 samevisitid= qh visit_id -1;
3171 FOREACHridge_(newfacet->ridges) {
3172 neighbor= otherfacet_(ridge, newfacet);
3173 if (neighbor->visitid == samevisitid)
3174 SETref_(ridge)= NULL; /* ridge free'd below */
3175 }
3176 qh_setcompact(newfacet->ridges);
3177
3178 trace4((qh ferr, 4034, "qh_mergecycle_ridges: add ridges to newfacet\n"));
3179 FORALLsame_cycle_(samecycle) {
3180 FOREACHridge_(same->ridges) {
3181 if (ridge->top == same) {
3182 ridge->top= newfacet;
3183 neighbor= ridge->bottom;
3184 }else if (ridge->bottom == same) {
3185 ridge->bottom= newfacet;
3186 neighbor= ridge->top;
3187 }else if (ridge->top == newfacet || ridge->bottom == newfacet) {
3188 qh_setappend(&newfacet->ridges, ridge);
3189 numold++; /* already set by qh_mergecycle_neighbors */
3190 continue;
3191 }else {
3192 qh_fprintf(qh ferr, 6098, "qhull internal error (qh_mergecycle_ridges): bad ridge r%d\n", ridge->id);
3193 qh_errexit(qh_ERRqhull, NULL, ridge);
3194 }
3195 if (neighbor == newfacet) {
3196 if (qh traceridge == ridge)
3197 qh traceridge= NULL;
3198 qh_setfree(&(ridge->vertices));
3199 qh_memfree_(ridge, (int)sizeof(ridgeT), freelistp);
3200 numold++;
3201 }else if (neighbor->visitid == samevisitid) {
3202 qh_setdel(neighbor->ridges, ridge);
3203 if (qh traceridge == ridge)
3204 qh traceridge= NULL;
3205 qh_setfree(&(ridge->vertices));
3206 qh_memfree_(ridge, (int)sizeof(ridgeT), freelistp);
3207 numold++;
3208 }else {
3209 qh_setappend(&newfacet->ridges, ridge);
3210 numold++;
3211 }
3212 }
3213 if (same->ridges)
3214 qh_settruncate(same->ridges, 0);
3215 if (!same->simplicial)
3216 continue;
3217 FOREACHneighbor_i_(same) { /* note: !newfact->simplicial */
3218 if (neighbor->visitid != samevisitid && neighbor->simplicial) {
3219 ridge= qh_newridge();
3220 ridge->vertices= qh_setnew_delnthsorted(same->vertices, qh hull_dim,
3221 neighbor_i, 0);
3222 toporient= (boolT)(same->toporient ^ (neighbor_i & 0x1));
3223 if (toporient) {
3224 ridge->top= newfacet;
3225 ridge->bottom= neighbor;
3226 ridge->simplicialbot= True;
3227 }else {
3228 ridge->top= neighbor;
3229 ridge->bottom= newfacet;
3230 ridge->simplicialtop= True;
3231 }
3232 qh_setappend(&(newfacet->ridges), ridge);
3233 qh_setappend(&(neighbor->ridges), ridge);
3234 if (qh ridge_id == qh traceridge_id)
3235 qh traceridge= ridge;
3236 numnew++;
3237 }
3238 }
3239 }
3240
3241 trace2((qh ferr, 2033, "qh_mergecycle_ridges: found %d old ridges and %d new ones\n",
3242 numold, numnew));
3243 } /* mergecycle_ridges */
3244
3245 /*-<a href="qh-merge.htm#TOC"
3246 >-------------------------------</a><a name="mergecycle_vneighbors">-</a>
3247
3248 qh_mergecycle_vneighbors( samecycle, newfacet )
3249 create vertex neighbors for newfacet from vertices of facets in samecycle
3250 samecycle marked with visitid == qh.visit_id - 1
3251
3252 returns:
3253 newfacet vertices with updated neighbors
3254 marks newfacet with qh.visit_id-1
3255 deletes vertices that are merged away
3256 sets delridge on all vertices (faster here than in mergecycle_ridges)
3257
3258 see:
3259 qh_mergevertex_neighbors()
3260
3261 design:
3262 for each vertex of samecycle facet
3263 set vertex->delridge
3264 delete samecycle facets from vertex neighbors
3265 append newfacet to vertex neighbors
3266 if vertex only in newfacet
3267 delete it from newfacet
3268 add it to qh.del_vertices for later deletion
3269 */
qh_mergecycle_vneighbors(facetT * samecycle,facetT * newfacet)3270 void qh_mergecycle_vneighbors(facetT *samecycle, facetT *newfacet) {
3271 facetT *neighbor, **neighborp;
3272 unsigned int mergeid;
3273 vertexT *vertex, **vertexp, *apex;
3274 setT *vertices;
3275
3276 trace4((qh ferr, 4035, "qh_mergecycle_vneighbors: update vertex neighbors for newfacet\n"));
3277 mergeid= qh visit_id - 1;
3278 newfacet->visitid= mergeid;
3279 vertices= qh_basevertices(samecycle); /* temp */
3280 apex= SETfirstt_(samecycle->vertices, vertexT);
3281 qh_setappend(&vertices, apex);
3282 FOREACHvertex_(vertices) {
3283 vertex->delridge= True;
3284 FOREACHneighbor_(vertex) {
3285 if (neighbor->visitid == mergeid)
3286 SETref_(neighbor)= NULL;
3287 }
3288 qh_setcompact(vertex->neighbors);
3289 qh_setappend(&vertex->neighbors, newfacet);
3290 if (!SETsecond_(vertex->neighbors)) {
3291 zinc_(Zcyclevertex);
3292 trace2((qh ferr, 2034, "qh_mergecycle_vneighbors: deleted v%d when merging cycle f%d into f%d\n",
3293 vertex->id, samecycle->id, newfacet->id));
3294 qh_setdelsorted(newfacet->vertices, vertex);
3295 vertex->deleted= True;
3296 qh_setappend(&qh del_vertices, vertex);
3297 }
3298 }
3299 qh_settempfree(&vertices);
3300 trace3((qh ferr, 3005, "qh_mergecycle_vneighbors: merged vertices from cycle f%d into f%d\n",
3301 samecycle->id, newfacet->id));
3302 } /* mergecycle_vneighbors */
3303
3304 /*-<a href="qh-merge.htm#TOC"
3305 >-------------------------------</a><a name="mergefacet">-</a>
3306
3307 qh_mergefacet( facet1, facet2, mergetype, mindist, maxdist, mergeapex )
3308 merges facet1 into facet2
3309 mergeapex==qh_MERGEapex if merging new facet into coplanar horizon (optimizes qh_mergesimplex)
3310
3311 returns:
3312 qh.max_outside and qh.min_vertex updated
3313 initializes vertex neighbors on first merge
3314
3315 note:
3316 mergetype only used for logging and error reporting
3317
3318 returns:
3319 facet2 contains facet1's vertices, neighbors, and ridges
3320 facet2 moved to end of qh.facet_list
3321 makes facet2 a newfacet
3322 sets facet2->newmerge set
3323 clears facet2->center (unless merging into a large facet)
3324 clears facet2->tested and ridge->tested for facet1
3325
3326 facet1 prepended to visible_list for later deletion and partitioning
3327 facet1->f.replace == facet2
3328
3329 adds neighboring facets to facet_mergeset if redundant or degenerate
3330
3331 notes:
3332 when done, tests facet1 and facet2 for degenerate or redundant neighbors and dupridges
3333 mindist/maxdist may be NULL (only if both NULL)
3334 traces merge if fmax_(maxdist,-mindist) > TRACEdist
3335
3336 see:
3337 qh_mergecycle()
3338
3339 design:
3340 trace merge and check for degenerate simplex
3341 make ridges for both facets
3342 update qh.max_outside, qh.max_vertex, qh.min_vertex
3343 update facet2->maxoutside and keepcentrum
3344 update facet2->nummerge
3345 update tested flags for facet2
3346 if facet1 is simplicial
3347 merge facet1 into facet2
3348 else
3349 merge facet1's neighbors into facet2
3350 merge facet1's ridges into facet2
3351 merge facet1's vertices into facet2
3352 merge facet1's vertex neighbors into facet2
3353 add facet2's vertices to qh.new_vertexlist
3354 move facet2 to end of qh.newfacet_list
3355 unless MRGcoplanarhorizon
3356 test facet2 for redundant neighbors
3357 test facet1 for degenerate neighbors
3358 test for redundant facet2
3359 maybe test for duplicate ridges ('Q15')
3360 move facet1 to qh.visible_list for later deletion
3361 */
qh_mergefacet(facetT * facet1,facetT * facet2,mergeType mergetype,realT * mindist,realT * maxdist,boolT mergeapex)3362 void qh_mergefacet(facetT *facet1, facetT *facet2, mergeType mergetype, realT *mindist, realT *maxdist, boolT mergeapex) {
3363 boolT traceonce= False;
3364 vertexT *vertex, **vertexp;
3365 realT mintwisted, vertexdist;
3366 realT onemerge;
3367 int tracerestore=0, nummerge;
3368 const char *mergename;
3369
3370 if(mergetype > 0 && mergetype <= sizeof(mergetypes))
3371 mergename= mergetypes[mergetype];
3372 else
3373 mergename= mergetypes[MRGnone];
3374 if (facet1->tricoplanar || facet2->tricoplanar) {
3375 if (!qh TRInormals) {
3376 qh_fprintf(qh ferr, 6226, "qhull internal error (qh_mergefacet): merge f%d into f%d for mergetype %d (%s) does not work for tricoplanar facets. Use option 'Q11'\n",
3377 facet1->id, facet2->id, mergetype, mergename);
3378 qh_errexit2(qh_ERRqhull, facet1, facet2);
3379 }
3380 if (facet2->tricoplanar) {
3381 facet2->tricoplanar= False;
3382 facet2->keepcentrum= False;
3383 }
3384 }
3385 zzinc_(Ztotmerge);
3386 if (qh REPORTfreq2 && qh POSTmerging) {
3387 if (zzval_(Ztotmerge) > qh mergereport + qh REPORTfreq2)
3388 qh_tracemerging();
3389 }
3390 #ifndef qh_NOtrace
3391 if (qh build_cnt >= qh RERUN) {
3392 if (mindist && (-*mindist > qh TRACEdist || *maxdist > qh TRACEdist)) {
3393 tracerestore= 0;
3394 qh IStracing= qh TRACElevel;
3395 traceonce= True;
3396 qh_fprintf(qh ferr, 8075, "qh_mergefacet: ========= trace wide merge #%d(%2.2g) for f%d into f%d for mergetype %d (%s), last point was p%d\n",
3397 zzval_(Ztotmerge), fmax_(-*mindist, *maxdist), facet1->id, facet2->id, mergetype, mergename, qh furthest_id);
3398 }else if (facet1 == qh tracefacet || facet2 == qh tracefacet) {
3399 tracerestore= qh IStracing;
3400 qh IStracing= 4;
3401 traceonce= True;
3402 qh_fprintf(qh ferr, 8076, "qh_mergefacet: ========= trace merge #%d for f%d into f%d for mergetype %d (%s), furthest is p%d\n",
3403 zzval_(Ztotmerge), facet1->id, facet2->id, mergetype, mergename, qh furthest_id);
3404 }
3405 }
3406 if (qh IStracing >= 2) {
3407 realT mergemin= -2;
3408 realT mergemax= -2;
3409
3410 if (mindist) {
3411 mergemin= *mindist;
3412 mergemax= *maxdist;
3413 }
3414 qh_fprintf(qh ferr, 2081, "qh_mergefacet: #%d merge f%d into f%d for merge for mergetype %d (%s), mindist= %2.2g, maxdist= %2.2g, max_outside %2.2g\n",
3415 zzval_(Ztotmerge), facet1->id, facet2->id, mergetype, mergename, mergemin, mergemax, qh max_outside);
3416 }
3417 #endif /* !qh_NOtrace */
3418 if(!qh ALLOWwide && mindist) {
3419 mintwisted= qh_WIDEmaxoutside * qh ONEmerge; /* same as qh_merge_twisted and qh_check_maxout (poly2) */
3420 maximize_(mintwisted, facet1->maxoutside);
3421 maximize_(mintwisted, facet2->maxoutside);
3422 if (*maxdist > mintwisted || -*mindist > mintwisted) {
3423 vertexdist= qh_vertex_bestdist(facet1->vertices);
3424 onemerge= qh ONEmerge + qh DISTround;
3425 if (vertexdist > mintwisted) {
3426 qh_fprintf(qh ferr, 6347, "qhull precision error (qh_mergefacet): wide merge for facet f%d into f%d for mergetype %d (%s). maxdist %2.2g (%.1fx) mindist %2.2g (%.1fx) vertexdist %2.2g Allow with 'Q12' (allow-wide)\n",
3427 facet1->id, facet2->id, mergetype, mergename, *maxdist, *maxdist/onemerge, *mindist, -*mindist/onemerge, vertexdist);
3428 }else {
3429 qh_fprintf(qh ferr, 6348, "qhull precision error (qh_mergefacet): wide merge for pinched facet f%d into f%d for mergetype %d (%s). maxdist %2.2g (%.fx) mindist %2.2g (%.1fx) vertexdist %2.2g Allow with 'Q12' (allow-wide)\n",
3430 facet1->id, facet2->id, mergetype, mergename, *maxdist, *maxdist/onemerge, *mindist, -*mindist/onemerge, vertexdist);
3431 }
3432 qh_errexit2(qh_ERRwide, facet1, facet2);
3433 }
3434 }
3435 if (facet1 == facet2 || facet1->visible || facet2->visible) {
3436 qh_fprintf(qh ferr, 6099, "qhull internal error (qh_mergefacet): either f%d and f%d are the same or one is a visible facet, mergetype %d (%s)\n",
3437 facet1->id, facet2->id, mergetype, mergename);
3438 qh_errexit2(qh_ERRqhull, facet1, facet2);
3439 }
3440 if (qh num_facets - qh num_visible <= qh hull_dim + 1) {
3441 qh_fprintf(qh ferr, 6227, "qhull topology error: Only %d facets remain. The input is too degenerate or the convexity constraints are too strong.\n",
3442 qh hull_dim+1);
3443 if (qh hull_dim >= 5 && !qh MERGEexact)
3444 qh_fprintf(qh ferr, 8079, " Option 'Qx' may avoid this problem.\n");
3445 qh_errexit(qh_ERRtopology, NULL, NULL);
3446 }
3447 if (!qh VERTEXneighbors)
3448 qh_vertexneighbors();
3449 qh_makeridges(facet1);
3450 qh_makeridges(facet2);
3451 if (qh IStracing >=4)
3452 qh_errprint("MERGING", facet1, facet2, NULL, NULL);
3453 if (mindist) {
3454 maximize_(qh max_outside, *maxdist);
3455 maximize_(qh max_vertex, *maxdist);
3456 #if qh_MAXoutside
3457 maximize_(facet2->maxoutside, *maxdist);
3458 #endif
3459 minimize_(qh min_vertex, *mindist);
3460 if (!facet2->keepcentrum
3461 && (*maxdist > qh WIDEfacet || *mindist < -qh WIDEfacet)) {
3462 facet2->keepcentrum= True;
3463 zinc_(Zwidefacet);
3464 }
3465 }
3466 nummerge= facet1->nummerge + facet2->nummerge + 1;
3467 if (nummerge >= qh_MAXnummerge)
3468 facet2->nummerge= qh_MAXnummerge;
3469 else
3470 facet2->nummerge= (short unsigned int)nummerge; /* limited to 9 bits by qh_MAXnummerge, -Wconversion */
3471 facet2->newmerge= True;
3472 facet2->dupridge= False;
3473 qh_updatetested(facet1, facet2);
3474 if (qh hull_dim > 2 && qh_setsize(facet1->vertices) == qh hull_dim)
3475 qh_mergesimplex(facet1, facet2, mergeapex);
3476 else {
3477 qh vertex_visit++;
3478 FOREACHvertex_(facet2->vertices)
3479 vertex->visitid= qh vertex_visit;
3480 if (qh hull_dim == 2)
3481 qh_mergefacet2d(facet1, facet2);
3482 else {
3483 qh_mergeneighbors(facet1, facet2);
3484 qh_mergevertices(facet1->vertices, &facet2->vertices);
3485 }
3486 qh_mergeridges(facet1, facet2);
3487 qh_mergevertex_neighbors(facet1, facet2);
3488 if (!facet2->newfacet)
3489 qh_newvertices(facet2->vertices);
3490 }
3491 if (facet2->coplanarhorizon) {
3492 zinc_(Zmergeintocoplanar);
3493 }else if (!facet2->newfacet) {
3494 zinc_(Zmergeintohorizon);
3495 }else if (!facet1->newfacet && facet2->newfacet) {
3496 zinc_(Zmergehorizon);
3497 }else {
3498 zinc_(Zmergenew);
3499 }
3500 qh_removefacet(facet2); /* append as a newfacet to end of qh facet_list */
3501 qh_appendfacet(facet2);
3502 facet2->newfacet= True;
3503 facet2->tested= False;
3504 qh_tracemerge(facet1, facet2, mergetype);
3505 if (traceonce) {
3506 qh_fprintf(qh ferr, 8080, "qh_mergefacet: end of wide tracing\n");
3507 qh IStracing= tracerestore;
3508 }
3509 if (mergetype != MRGcoplanarhorizon) {
3510 trace3((qh ferr, 3076, "qh_mergefacet: check f%d and f%d for redundant and degenerate neighbors\n",
3511 facet1->id, facet2->id));
3512 qh_test_redundant_neighbors(facet2);
3513 qh_test_degen_neighbors(facet1); /* after qh_test_redundant_neighbors since MRGdegen more difficult than MRGredundant
3514 and before qh_willdelete which clears facet1.neighbors */
3515 qh_degen_redundant_facet(facet2); /* may occur in qh_merge_pinchedvertices, e.g., rbox 175 C3,2e-13 D4 t1545228104 | qhull d */
3516 qh_maybe_duplicateridges(facet2);
3517 }
3518 qh_willdelete(facet1, facet2);
3519 } /* mergefacet */
3520
3521
3522 /*-<a href="qh-merge.htm#TOC"
3523 >-------------------------------</a><a name="mergefacet2d">-</a>
3524
3525 qh_mergefacet2d( facet1, facet2 )
3526 in 2d, merges neighbors and vertices of facet1 into facet2
3527
3528 returns:
3529 build ridges for neighbors if necessary
3530 facet2 looks like a simplicial facet except for centrum, ridges
3531 neighbors are opposite the corresponding vertex
3532 maintains orientation of facet2
3533
3534 notes:
3535 qh_mergefacet() retains non-simplicial structures
3536 they are not needed in 2d, but later routines may use them
3537 preserves qh.vertex_visit for qh_mergevertex_neighbors()
3538
3539 design:
3540 get vertices and neighbors
3541 determine new vertices and neighbors
3542 set new vertices and neighbors and adjust orientation
3543 make ridges for new neighbor if needed
3544 */
qh_mergefacet2d(facetT * facet1,facetT * facet2)3545 void qh_mergefacet2d(facetT *facet1, facetT *facet2) {
3546 vertexT *vertex1A, *vertex1B, *vertex2A, *vertex2B, *vertexA, *vertexB;
3547 facetT *neighbor1A, *neighbor1B, *neighbor2A, *neighbor2B, *neighborA, *neighborB;
3548
3549 vertex1A= SETfirstt_(facet1->vertices, vertexT);
3550 vertex1B= SETsecondt_(facet1->vertices, vertexT);
3551 vertex2A= SETfirstt_(facet2->vertices, vertexT);
3552 vertex2B= SETsecondt_(facet2->vertices, vertexT);
3553 neighbor1A= SETfirstt_(facet1->neighbors, facetT);
3554 neighbor1B= SETsecondt_(facet1->neighbors, facetT);
3555 neighbor2A= SETfirstt_(facet2->neighbors, facetT);
3556 neighbor2B= SETsecondt_(facet2->neighbors, facetT);
3557 if (vertex1A == vertex2A) {
3558 vertexA= vertex1B;
3559 vertexB= vertex2B;
3560 neighborA= neighbor2A;
3561 neighborB= neighbor1A;
3562 }else if (vertex1A == vertex2B) {
3563 vertexA= vertex1B;
3564 vertexB= vertex2A;
3565 neighborA= neighbor2B;
3566 neighborB= neighbor1A;
3567 }else if (vertex1B == vertex2A) {
3568 vertexA= vertex1A;
3569 vertexB= vertex2B;
3570 neighborA= neighbor2A;
3571 neighborB= neighbor1B;
3572 }else { /* 1B == 2B */
3573 vertexA= vertex1A;
3574 vertexB= vertex2A;
3575 neighborA= neighbor2B;
3576 neighborB= neighbor1B;
3577 }
3578 /* vertexB always from facet2, neighborB always from facet1 */
3579 if (vertexA->id > vertexB->id) {
3580 SETfirst_(facet2->vertices)= vertexA;
3581 SETsecond_(facet2->vertices)= vertexB;
3582 if (vertexB == vertex2A)
3583 facet2->toporient= !facet2->toporient;
3584 SETfirst_(facet2->neighbors)= neighborA;
3585 SETsecond_(facet2->neighbors)= neighborB;
3586 }else {
3587 SETfirst_(facet2->vertices)= vertexB;
3588 SETsecond_(facet2->vertices)= vertexA;
3589 if (vertexB == vertex2B)
3590 facet2->toporient= !facet2->toporient;
3591 SETfirst_(facet2->neighbors)= neighborB;
3592 SETsecond_(facet2->neighbors)= neighborA;
3593 }
3594 /* qh_makeridges not needed since neighborB is not degenerate */
3595 qh_setreplace(neighborB->neighbors, facet1, facet2);
3596 trace4((qh ferr, 4036, "qh_mergefacet2d: merged v%d and neighbor f%d of f%d into f%d\n",
3597 vertexA->id, neighborB->id, facet1->id, facet2->id));
3598 } /* mergefacet2d */
3599
3600
3601 /*-<a href="qh-merge.htm#TOC"
3602 >-------------------------------</a><a name="mergeneighbors">-</a>
3603
3604 qh_mergeneighbors( facet1, facet2 )
3605 merges the neighbors of facet1 into facet2
3606
3607 notes:
3608 only called by qh_mergefacet
3609 qh.hull_dim >= 3
3610 see qh_mergecycle_neighbors
3611
3612 design:
3613 for each neighbor of facet1
3614 if neighbor is also a neighbor of facet2
3615 if neighbor is simplicial
3616 make ridges for later deletion as a degenerate facet
3617 update its neighbor set
3618 else
3619 move the neighbor relation to facet2
3620 remove the neighbor relation for facet1 and facet2
3621 */
qh_mergeneighbors(facetT * facet1,facetT * facet2)3622 void qh_mergeneighbors(facetT *facet1, facetT *facet2) {
3623 facetT *neighbor, **neighborp;
3624
3625 trace4((qh ferr, 4037, "qh_mergeneighbors: merge neighbors of f%d and f%d\n",
3626 facet1->id, facet2->id));
3627 qh visit_id++;
3628 FOREACHneighbor_(facet2) {
3629 neighbor->visitid= qh visit_id;
3630 }
3631 FOREACHneighbor_(facet1) {
3632 if (neighbor->visitid == qh visit_id) {
3633 if (neighbor->simplicial) /* is degen, needs ridges */
3634 qh_makeridges(neighbor);
3635 if (SETfirstt_(neighbor->neighbors, facetT) != facet1) /*keep newfacet->horizon*/
3636 qh_setdel(neighbor->neighbors, facet1);
3637 else {
3638 qh_setdel(neighbor->neighbors, facet2);
3639 qh_setreplace(neighbor->neighbors, facet1, facet2);
3640 }
3641 }else if (neighbor != facet2) {
3642 qh_setappend(&(facet2->neighbors), neighbor);
3643 qh_setreplace(neighbor->neighbors, facet1, facet2);
3644 }
3645 }
3646 qh_setdel(facet1->neighbors, facet2); /* here for makeridges */
3647 qh_setdel(facet2->neighbors, facet1);
3648 } /* mergeneighbors */
3649
3650
3651 /*-<a href="qh-merge.htm#TOC"
3652 >-------------------------------</a><a name="mergeridges">-</a>
3653
3654 qh_mergeridges( facet1, facet2 )
3655 merges the ridge set of facet1 into facet2
3656
3657 returns:
3658 may delete all ridges for a vertex
3659 sets vertex->delridge on deleted ridges
3660
3661 see:
3662 qh_mergecycle_ridges()
3663
3664 design:
3665 delete ridges between facet1 and facet2
3666 mark (delridge) vertices on these ridges for later testing
3667 for each remaining ridge
3668 rename facet1 to facet2
3669 */
qh_mergeridges(facetT * facet1,facetT * facet2)3670 void qh_mergeridges(facetT *facet1, facetT *facet2) {
3671 ridgeT *ridge, **ridgep;
3672
3673 trace4((qh ferr, 4038, "qh_mergeridges: merge ridges of f%d into f%d\n",
3674 facet1->id, facet2->id));
3675 FOREACHridge_(facet2->ridges) {
3676 if ((ridge->top == facet1) || (ridge->bottom == facet1)) {
3677 /* ridge.nonconvex is irrelevant due to merge */
3678 qh_delridge_merge(ridge); /* expensive in high-d, could rebuild */
3679 ridgep--; /* deleted this ridge, repeat with next ridge*/
3680 }
3681 }
3682 FOREACHridge_(facet1->ridges) {
3683 if (ridge->top == facet1) {
3684 ridge->top= facet2;
3685 ridge->simplicialtop= False;
3686 }else { /* ridge.bottom is facet1 */
3687 ridge->bottom= facet2;
3688 ridge->simplicialbot= False;
3689 }
3690 qh_setappend(&(facet2->ridges), ridge);
3691 }
3692 } /* mergeridges */
3693
3694
3695 /*-<a href="qh-merge.htm#TOC"
3696 >-------------------------------</a><a name="mergesimplex">-</a>
3697
3698 qh_mergesimplex( facet1, facet2, mergeapex )
3699 merge simplicial facet1 into facet2
3700 mergeapex==qh_MERGEapex if merging samecycle into horizon facet
3701 vertex id is latest (most recently created)
3702 facet1 may be contained in facet2
3703 ridges exist for both facets
3704
3705 returns:
3706 facet2 with updated vertices, ridges, neighbors
3707 updated neighbors for facet1's vertices
3708 facet1 not deleted
3709 sets vertex->delridge on deleted ridges
3710
3711 notes:
3712 special case code since this is the most common merge
3713 called from qh_mergefacet()
3714
3715 design:
3716 if qh_MERGEapex
3717 add vertices of facet2 to qh.new_vertexlist if necessary
3718 add apex to facet2
3719 else
3720 for each ridge between facet1 and facet2
3721 set vertex->delridge
3722 determine the apex for facet1 (i.e., vertex to be merged)
3723 unless apex already in facet2
3724 insert apex into vertices for facet2
3725 add vertices of facet2 to qh.new_vertexlist if necessary
3726 add apex to qh.new_vertexlist if necessary
3727 for each vertex of facet1
3728 if apex
3729 rename facet1 to facet2 in its vertex neighbors
3730 else
3731 delete facet1 from vertex neighbors
3732 if only in facet2
3733 add vertex to qh.del_vertices for later deletion
3734 for each ridge of facet1
3735 delete ridges between facet1 and facet2
3736 append other ridges to facet2 after renaming facet to facet2
3737 */
qh_mergesimplex(facetT * facet1,facetT * facet2,boolT mergeapex)3738 void qh_mergesimplex(facetT *facet1, facetT *facet2, boolT mergeapex) {
3739 vertexT *vertex, **vertexp, *opposite;
3740 ridgeT *ridge, **ridgep;
3741 boolT isnew= False;
3742 facetT *neighbor, **neighborp, *otherfacet;
3743
3744 if (mergeapex) {
3745 opposite= SETfirstt_(facet1->vertices, vertexT); /* apex is opposite facet2. It has the last vertex id */
3746 trace4((qh ferr, 4086, "qh_mergesimplex: merge apex v%d of f%d into facet f%d\n",
3747 opposite->id, facet1->id, facet2->id));
3748 if (!facet2->newfacet)
3749 qh_newvertices(facet2->vertices); /* apex, the first vertex, is already new */
3750 if (SETfirstt_(facet2->vertices, vertexT) != opposite) {
3751 qh_setaddnth(&facet2->vertices, 0, opposite);
3752 isnew= True;
3753 }
3754 }else {
3755 zinc_(Zmergesimplex);
3756 FOREACHvertex_(facet1->vertices)
3757 vertex->seen= False;
3758 FOREACHridge_(facet1->ridges) {
3759 if (otherfacet_(ridge, facet1) == facet2) {
3760 FOREACHvertex_(ridge->vertices) {
3761 vertex->seen= True;
3762 vertex->delridge= True;
3763 }
3764 break;
3765 }
3766 }
3767 FOREACHvertex_(facet1->vertices) {
3768 if (!vertex->seen)
3769 break; /* must occur */
3770 }
3771 opposite= vertex;
3772 trace4((qh ferr, 4039, "qh_mergesimplex: merge opposite v%d of f%d into facet f%d\n",
3773 opposite->id, facet1->id, facet2->id));
3774 isnew= qh_addfacetvertex(facet2, opposite);
3775 if (!facet2->newfacet)
3776 qh_newvertices(facet2->vertices);
3777 else if (!opposite->newfacet) {
3778 qh_removevertex(opposite);
3779 qh_appendvertex(opposite);
3780 }
3781 }
3782 trace4((qh ferr, 4040, "qh_mergesimplex: update vertex neighbors of f%d\n",
3783 facet1->id));
3784 FOREACHvertex_(facet1->vertices) {
3785 if (vertex == opposite && isnew)
3786 qh_setreplace(vertex->neighbors, facet1, facet2);
3787 else {
3788 qh_setdel(vertex->neighbors, facet1);
3789 if (!SETsecond_(vertex->neighbors))
3790 qh_mergevertex_del(vertex, facet1, facet2);
3791 }
3792 }
3793 trace4((qh ferr, 4041, "qh_mergesimplex: merge ridges and neighbors of f%d into f%d\n",
3794 facet1->id, facet2->id));
3795 qh visit_id++;
3796 FOREACHneighbor_(facet2)
3797 neighbor->visitid= qh visit_id;
3798 FOREACHridge_(facet1->ridges) {
3799 otherfacet= otherfacet_(ridge, facet1);
3800 if (otherfacet == facet2) {
3801 /* ridge.nonconvex is irrelevant due to merge */
3802 qh_delridge_merge(ridge); /* expensive in high-d, could rebuild */
3803 ridgep--; /* deleted this ridge, repeat with next ridge*/
3804 qh_setdel(facet2->neighbors, facet1); /* a simplicial facet may have duplicate neighbors, need to delete each one */
3805 }else if (otherfacet->dupridge && !qh_setin(otherfacet->neighbors, facet1)) {
3806 qh_fprintf(qh ferr, 6356, "qhull topology error (qh_mergesimplex): f%d is a dupridge of f%d, cannot merge f%d into f%d\n",
3807 facet1->id, otherfacet->id, facet1->id, facet2->id);
3808 qh_errexit2(qh_ERRqhull, facet1, otherfacet);
3809 }else {
3810 trace4((qh ferr, 4059, "qh_mergesimplex: move r%d with f%d to f%d, new neighbor? %d, maybe horizon? %d\n",
3811 ridge->id, otherfacet->id, facet2->id, (otherfacet->visitid != qh visit_id), (SETfirstt_(otherfacet->neighbors, facetT) == facet1)));
3812 qh_setappend(&facet2->ridges, ridge);
3813 if (otherfacet->visitid != qh visit_id) {
3814 qh_setappend(&facet2->neighbors, otherfacet);
3815 qh_setreplace(otherfacet->neighbors, facet1, facet2);
3816 otherfacet->visitid= qh visit_id;
3817 }else {
3818 if (otherfacet->simplicial) /* is degen, needs ridges */
3819 qh_makeridges(otherfacet);
3820 if (SETfirstt_(otherfacet->neighbors, facetT) == facet1) {
3821 /* keep new, otherfacet->neighbors->horizon */
3822 qh_setdel(otherfacet->neighbors, facet2);
3823 qh_setreplace(otherfacet->neighbors, facet1, facet2);
3824 }else {
3825 /* facet2 is already a neighbor of otherfacet, by f.visitid */
3826 qh_setdel(otherfacet->neighbors, facet1);
3827 }
3828 }
3829 if (ridge->top == facet1) { /* wait until after qh_makeridges */
3830 ridge->top= facet2;
3831 ridge->simplicialtop= False;
3832 }else {
3833 ridge->bottom= facet2;
3834 ridge->simplicialbot= False;
3835 }
3836 }
3837 }
3838 trace3((qh ferr, 3006, "qh_mergesimplex: merged simplex f%d v%d into facet f%d\n",
3839 facet1->id, opposite->id, facet2->id));
3840 } /* mergesimplex */
3841
3842 /*-<a href="qh-merge.htm#TOC"
3843 >-------------------------------</a><a name="mergevertex_del">-</a>
3844
3845 qh_mergevertex_del( vertex, facet1, facet2 )
3846 delete a vertex because of merging facet1 into facet2
3847
3848 returns:
3849 deletes vertex from facet2
3850 adds vertex to qh.del_vertices for later deletion
3851 */
qh_mergevertex_del(vertexT * vertex,facetT * facet1,facetT * facet2)3852 void qh_mergevertex_del(vertexT *vertex, facetT *facet1, facetT *facet2) {
3853
3854 zinc_(Zmergevertex);
3855 trace2((qh ferr, 2035, "qh_mergevertex_del: deleted v%d when merging f%d into f%d\n",
3856 vertex->id, facet1->id, facet2->id));
3857 qh_setdelsorted(facet2->vertices, vertex);
3858 vertex->deleted= True;
3859 qh_setappend(&qh del_vertices, vertex);
3860 } /* mergevertex_del */
3861
3862 /*-<a href="qh-merge.htm#TOC"
3863 >-------------------------------</a><a name="mergevertex_neighbors">-</a>
3864
3865 qh_mergevertex_neighbors( facet1, facet2 )
3866 merge the vertex neighbors of facet1 to facet2
3867
3868 returns:
3869 if vertex is current qh.vertex_visit
3870 deletes facet1 from vertex->neighbors
3871 else
3872 renames facet1 to facet2 in vertex->neighbors
3873 deletes vertices if only one neighbor
3874
3875 notes:
3876 assumes vertex neighbor sets are good
3877 */
qh_mergevertex_neighbors(facetT * facet1,facetT * facet2)3878 void qh_mergevertex_neighbors(facetT *facet1, facetT *facet2) {
3879 vertexT *vertex, **vertexp;
3880
3881 trace4((qh ferr, 4042, "qh_mergevertex_neighbors: merge vertex neighborset for f%d into f%d\n",
3882 facet1->id, facet2->id));
3883 if (qh tracevertex) {
3884 qh_fprintf(qh ferr, 8081, "qh_mergevertex_neighbors: of f%d into f%d at furthest p%d f0= %p\n",
3885 facet1->id, facet2->id, qh furthest_id, qh tracevertex->neighbors->e[0].p);
3886 qh_errprint("TRACE", NULL, NULL, NULL, qh tracevertex);
3887 }
3888 FOREACHvertex_(facet1->vertices) {
3889 if (vertex->visitid != qh vertex_visit)
3890 qh_setreplace(vertex->neighbors, facet1, facet2);
3891 else {
3892 qh_setdel(vertex->neighbors, facet1);
3893 if (!SETsecond_(vertex->neighbors))
3894 qh_mergevertex_del(vertex, facet1, facet2);
3895 }
3896 }
3897 if (qh tracevertex)
3898 qh_errprint("TRACE", NULL, NULL, NULL, qh tracevertex);
3899 } /* mergevertex_neighbors */
3900
3901
3902 /*-<a href="qh-merge.htm#TOC"
3903 >-------------------------------</a><a name="mergevertices">-</a>
3904
3905 qh_mergevertices( vertices1, vertices2 )
3906 merges the vertex set of facet1 into facet2
3907
3908 returns:
3909 replaces vertices2 with merged set
3910 preserves vertex_visit for qh_mergevertex_neighbors
3911 updates qh.newvertex_list
3912
3913 design:
3914 create a merged set of both vertices (in inverse id order)
3915 */
qh_mergevertices(setT * vertices1,setT ** vertices2)3916 void qh_mergevertices(setT *vertices1, setT **vertices2) {
3917 int newsize= qh_setsize(vertices1)+qh_setsize(*vertices2) - qh hull_dim + 1;
3918 setT *mergedvertices;
3919 vertexT *vertex, **vertexp, **vertex2= SETaddr_(*vertices2, vertexT);
3920
3921 mergedvertices= qh_settemp(newsize);
3922 FOREACHvertex_(vertices1) {
3923 if (!*vertex2 || vertex->id > (*vertex2)->id)
3924 qh_setappend(&mergedvertices, vertex);
3925 else {
3926 while (*vertex2 && (*vertex2)->id > vertex->id)
3927 qh_setappend(&mergedvertices, *vertex2++);
3928 if (!*vertex2 || (*vertex2)->id < vertex->id)
3929 qh_setappend(&mergedvertices, vertex);
3930 else
3931 qh_setappend(&mergedvertices, *vertex2++);
3932 }
3933 }
3934 while (*vertex2)
3935 qh_setappend(&mergedvertices, *vertex2++);
3936 if (newsize < qh_setsize(mergedvertices)) {
3937 qh_fprintf(qh ferr, 6100, "qhull internal error (qh_mergevertices): facets did not share a ridge\n");
3938 qh_errexit(qh_ERRqhull, NULL, NULL);
3939 }
3940 qh_setfree(vertices2);
3941 *vertices2= mergedvertices;
3942 qh_settemppop();
3943 } /* mergevertices */
3944
3945
3946 /*-<a href="qh-merge.htm#TOC"
3947 >-------------------------------</a><a name="neighbor_intersections">-</a>
3948
3949 qh_neighbor_intersections( vertex )
3950 return intersection of all vertices in vertex->neighbors except for vertex
3951
3952 returns:
3953 returns temporary set of vertices
3954 does not include vertex
3955 NULL if a neighbor is simplicial
3956 NULL if empty set
3957
3958 notes:
3959 only called by qh_redundant_vertex for qh_reducevertices
3960 so f.vertices does not contain extraneous vertices that are not in f.ridges
3961 used for renaming vertices
3962
3963 design:
3964 initialize the intersection set with vertices of the first two neighbors
3965 delete vertex from the intersection
3966 for each remaining neighbor
3967 intersect its vertex set with the intersection set
3968 return NULL if empty
3969 return the intersection set
3970 */
qh_neighbor_intersections(vertexT * vertex)3971 setT *qh_neighbor_intersections(vertexT *vertex) {
3972 facetT *neighbor, **neighborp, *neighborA, *neighborB;
3973 setT *intersect;
3974 int neighbor_i, neighbor_n;
3975
3976 FOREACHneighbor_(vertex) {
3977 if (neighbor->simplicial)
3978 return NULL;
3979 }
3980 neighborA= SETfirstt_(vertex->neighbors, facetT);
3981 neighborB= SETsecondt_(vertex->neighbors, facetT);
3982 zinc_(Zintersectnum);
3983 if (!neighborA)
3984 return NULL;
3985 if (!neighborB)
3986 intersect= qh_setcopy(neighborA->vertices, 0);
3987 else
3988 intersect= qh_vertexintersect_new(neighborA->vertices, neighborB->vertices);
3989 qh_settemppush(intersect);
3990 qh_setdelsorted(intersect, vertex);
3991 FOREACHneighbor_i_(vertex) {
3992 if (neighbor_i >= 2) {
3993 zinc_(Zintersectnum);
3994 qh_vertexintersect(&intersect, neighbor->vertices);
3995 if (!SETfirst_(intersect)) {
3996 zinc_(Zintersectfail);
3997 qh_settempfree(&intersect);
3998 return NULL;
3999 }
4000 }
4001 }
4002 trace3((qh ferr, 3007, "qh_neighbor_intersections: %d vertices in neighbor intersection of v%d\n",
4003 qh_setsize(intersect), vertex->id));
4004 return intersect;
4005 } /* neighbor_intersections */
4006
4007 /*-<a href="qh-merge.htm#TOC"
4008 >-------------------------------</a><a name="neighbor_vertices">-</a>
4009
4010 qh_neighbor_vertices( vertex )
4011 return neighboring vertices for a vertex (not in subridge)
4012 assumes vertices have full vertex->neighbors
4013
4014 returns:
4015 temporary set of vertices
4016
4017 notes:
4018 updates qh.visit_id and qh.vertex_visit
4019 similar to qh_vertexridges
4020
4021 */
qh_neighbor_vertices(vertexT * vertexA,setT * subridge)4022 setT *qh_neighbor_vertices(vertexT *vertexA, setT *subridge) {
4023 facetT *neighbor, **neighborp;
4024 vertexT *vertex, **vertexp;
4025 setT *vertices= qh_settemp(qh TEMPsize);
4026
4027 qh visit_id++;
4028 FOREACHneighbor_(vertexA)
4029 neighbor->visitid= qh visit_id;
4030 qh vertex_visit++;
4031 vertexA->visitid= qh vertex_visit;
4032 FOREACHvertex_(subridge) {
4033 vertex->visitid= qh vertex_visit;
4034 }
4035 FOREACHneighbor_(vertexA) {
4036 if (*neighborp) /* no new ridges in last neighbor */
4037 qh_neighbor_vertices_facet(vertexA, neighbor, &vertices);
4038 }
4039 trace3((qh ferr, 3035, "qh_neighbor_vertices: %d non-subridge, vertex neighbors for v%d\n",
4040 qh_setsize(vertices), vertexA->id));
4041 return vertices;
4042 } /* neighbor_vertices */
4043
4044 /*-<a href="qh-merge.htm#TOC"
4045 >-------------------------------</a><a name="neighbor_vertices_facet">-</a>
4046
4047 qh_neighbor_vertices_facet( vertex, facet, vertices )
4048 add neighboring vertices on ridges for vertex in facet
4049 neighbor->visitid==qh.visit_id if it hasn't been visited
4050 v.visitid==qh.vertex_visit if it is already in vertices
4051
4052 returns:
4053 vertices updated
4054 sets facet->visitid to qh.visit_id-1
4055
4056 notes:
4057 only called by qh_neighbor_vertices
4058 similar to qh_vertexridges_facet
4059
4060 design:
4061 for each ridge of facet
4062 if ridge of visited neighbor (i.e., unprocessed)
4063 if vertex in ridge
4064 append unprocessed vertices of ridge
4065 mark facet processed
4066 */
qh_neighbor_vertices_facet(vertexT * vertexA,facetT * facet,setT ** vertices)4067 void qh_neighbor_vertices_facet(vertexT *vertexA, facetT *facet, setT **vertices) {
4068 ridgeT *ridge, **ridgep;
4069 facetT *neighbor;
4070 vertexT *second, *last, *vertex, **vertexp;
4071 int last_i= qh hull_dim-2, count= 0;
4072 boolT isridge;
4073
4074 if (facet->simplicial) {
4075 FOREACHvertex_(facet->vertices) {
4076 if (vertex->visitid != qh vertex_visit) {
4077 vertex->visitid= qh vertex_visit;
4078 qh_setappend(vertices, vertex);
4079 count++;
4080 }
4081 }
4082 }else {
4083 FOREACHridge_(facet->ridges) {
4084 neighbor= otherfacet_(ridge, facet);
4085 if (neighbor->visitid == qh visit_id) {
4086 isridge= False;
4087 if (SETfirst_(ridge->vertices) == vertexA) {
4088 isridge= True;
4089 }else if (last_i > 2) {
4090 second= SETsecondt_(ridge->vertices, vertexT);
4091 last= SETelemt_(ridge->vertices, last_i, vertexT);
4092 if (second->id >= vertexA->id && last->id <= vertexA->id) { /* vertices inverse sorted by id */
4093 if (second == vertexA || last == vertexA)
4094 isridge= True;
4095 else if (qh_setin(ridge->vertices, vertexA))
4096 isridge= True;
4097 }
4098 }else if (SETelem_(ridge->vertices, last_i) == vertexA) {
4099 isridge= True;
4100 }else if (last_i > 1 && SETsecond_(ridge->vertices) == vertexA) {
4101 isridge= True;
4102 }
4103 if (isridge) {
4104 FOREACHvertex_(ridge->vertices) {
4105 if (vertex->visitid != qh vertex_visit) {
4106 vertex->visitid= qh vertex_visit;
4107 qh_setappend(vertices, vertex);
4108 count++;
4109 }
4110 }
4111 }
4112 }
4113 }
4114 }
4115 facet->visitid= qh visit_id-1;
4116 if (count) {
4117 trace4((qh ferr, 4079, "qh_neighbor_vertices_facet: found %d vertex neighbors for v%d in f%d (simplicial? %d)\n",
4118 count, vertexA->id, facet->id, facet->simplicial));
4119 }
4120 } /* neighbor_vertices_facet */
4121
4122
4123 /*-<a href="qh-merge.htm#TOC"
4124 >-------------------------------</a><a name="newvertices">-</a>
4125
4126 qh_newvertices( vertices )
4127 add vertices to end of qh.vertex_list (marks as new vertices)
4128
4129 returns:
4130 vertices on qh.newvertex_list
4131 vertex->newfacet set
4132 */
qh_newvertices(setT * vertices)4133 void qh_newvertices(setT *vertices) {
4134 vertexT *vertex, **vertexp;
4135
4136 FOREACHvertex_(vertices) {
4137 if (!vertex->newfacet) {
4138 qh_removevertex(vertex);
4139 qh_appendvertex(vertex);
4140 }
4141 }
4142 } /* newvertices */
4143
4144 /*-<a href="qh-merge.htm#TOC"
4145 >-------------------------------</a><a name="next_vertexmerge">-</a>
4146
4147 qh_next_vertexmerge( )
4148 return next vertex merge from qh.vertex_mergeset
4149
4150 returns:
4151 vertex merge either MRGvertices or MRGsubridge
4152 drops merges of deleted vertices
4153
4154 notes:
4155 called from qh_merge_pinchedvertices
4156 */
qh_next_vertexmerge(void)4157 mergeT *qh_next_vertexmerge(void /* qh.vertex_mergeset */) {
4158 mergeT *merge;
4159 int merge_i, merge_n, best_i= -1;
4160 realT bestdist= REALmax;
4161
4162 FOREACHmerge_i_(qh vertex_mergeset) {
4163 if (!merge->vertex1 || !merge->vertex2) {
4164 qh_fprintf(qh ferr, 6299, "qhull internal error (qh_next_vertexmerge): expecting two vertices for vertex merge. Got v%d v%d and optional f%d\n",
4165 getid_(merge->vertex1), getid_(merge->vertex2), getid_(merge->facet1));
4166 qh_errexit(qh_ERRqhull, merge->facet1, NULL);
4167 }
4168 if (merge->vertex1->deleted || merge->vertex2->deleted) {
4169 trace3((qh ferr, 3030, "qh_next_vertexmerge: drop merge of v%d (del? %d) into v%d (del? %d) due to deleted vertex of r%d and r%d\n",
4170 merge->vertex1->id, merge->vertex1->deleted, merge->vertex2->id, merge->vertex2->deleted, getid_(merge->ridge1), getid_(merge->ridge2)));
4171 qh_drop_mergevertex(merge);
4172 qh_setdelnth(qh vertex_mergeset, merge_i);
4173 merge_i--; merge_n--;
4174 qh_memfree(merge, (int)sizeof(mergeT));
4175 }else if (merge->distance < bestdist) {
4176 bestdist= merge->distance;
4177 best_i= merge_i;
4178 }
4179 }
4180 merge= NULL;
4181 if (best_i >= 0) {
4182 merge= SETelemt_(qh vertex_mergeset, best_i, mergeT);
4183 if (bestdist/qh ONEmerge > qh_WIDEpinched) {
4184 if (merge->mergetype==MRGvertices) {
4185 if (merge->ridge1->top == merge->ridge2->bottom && merge->ridge1->bottom == merge->ridge2->top)
4186 qh_fprintf(qh ferr, 6391, "qhull topology error (qh_next_vertexmerge): no nearly adjacent vertices to resolve opposite oriented ridges r%d and r%d in f%d and f%d. Nearest v%d and v%d dist %2.2g (%.1fx)\n",
4187 merge->ridge1->id, merge->ridge2->id, merge->ridge1->top->id, merge->ridge1->bottom->id, merge->vertex1->id, merge->vertex2->id, bestdist, bestdist/qh ONEmerge);
4188 else
4189 qh_fprintf(qh ferr, 6381, "qhull topology error (qh_next_vertexmerge): no nearly adjacent vertices to resolve duplicate ridges r%d and r%d. Nearest v%d and v%d dist %2.2g (%.1fx)\n",
4190 merge->ridge1->id, merge->ridge2->id, merge->vertex1->id, merge->vertex2->id, bestdist, bestdist/qh ONEmerge);
4191 }else {
4192 qh_fprintf(qh ferr, 6208, "qhull topology error (qh_next_vertexmerge): no nearly adjacent vertices to resolve dupridge. Nearest v%d and v%d dist %2.2g (%.1fx)\n",
4193 merge->vertex1->id, merge->vertex2->id, bestdist, bestdist/qh ONEmerge);
4194 }
4195 /* it may be possible to find a different vertex, after other vertex merges have occurred */
4196 qh_errexit(qh_ERRtopology, NULL, merge->ridge1);
4197 }
4198 qh_setdelnth(qh vertex_mergeset, best_i);
4199 }
4200 return merge;
4201 } /* next_vertexmerge */
4202
4203 /*-<a href="qh-merge.htm#TOC"
4204 >-------------------------------</a><a name="opposite_horizonfacet">-</a>
4205
4206 qh_opposite_horizonfacet( merge, opposite )
4207 return horizon facet for one of the merge facets, and its opposite vertex across the ridge
4208 assumes either facet1 or facet2 of merge is 'mergehorizon'
4209 assumes both facets are simplicial facets on qh.new_facetlist
4210
4211 returns:
4212 horizon facet and opposite vertex
4213
4214 notes:
4215 called by qh_getpinchedmerges
4216 */
qh_opposite_horizonfacet(mergeT * merge,vertexT ** opposite)4217 facetT *qh_opposite_horizonfacet(mergeT *merge, vertexT **opposite) {
4218 facetT *facet, *horizon, *otherfacet;
4219 int neighbor_i;
4220
4221 if (!merge->facet1->simplicial || !merge->facet2->simplicial || (!merge->facet1->mergehorizon && !merge->facet2->mergehorizon)) {
4222 qh_fprintf(qh ferr, 6273, "qhull internal error (qh_opposite_horizonfacet): expecting merge of simplicial facets, at least one of which is mergehorizon. Either simplicial or mergehorizon is wrong\n");
4223 qh_errexit2(qh_ERRqhull, merge->facet1, merge->facet2);
4224 }
4225 if (merge->facet1->mergehorizon) {
4226 facet= merge->facet1;
4227 otherfacet= merge->facet2;
4228 }else {
4229 facet= merge->facet2;
4230 otherfacet= merge->facet1;
4231 }
4232 horizon= SETfirstt_(facet->neighbors, facetT);
4233 neighbor_i= qh_setindex(otherfacet->neighbors, facet);
4234 if (neighbor_i==-1)
4235 neighbor_i= qh_setindex(otherfacet->neighbors, qh_MERGEridge);
4236 if (neighbor_i==-1) {
4237 qh_fprintf(qh ferr, 6238, "qhull internal error (qh_opposite_horizonfacet): merge facet f%d not connected to mergehorizon f%d\n",
4238 otherfacet->id, facet->id);
4239 qh_errexit2(qh_ERRqhull, otherfacet, facet);
4240 }
4241 *opposite= SETelemt_(otherfacet->vertices, neighbor_i, vertexT);
4242 return horizon;
4243 } /* opposite_horizonfacet */
4244
4245
4246 /*-<a href="qh-merge.htm#TOC"
4247 >-------------------------------</a><a name="reducevertices">-</a>
4248
4249 qh_reducevertices( )
4250 reduce extra vertices, shared vertices, and redundant vertices
4251 facet->newmerge is set if merged since last call
4252 vertex->delridge is set if vertex was on a deleted ridge
4253 if !qh.MERGEvertices, only removes extra vertices
4254
4255 returns:
4256 True if also merged degen_redundant facets
4257 vertices are renamed if possible
4258 clears facet->newmerge and vertex->delridge
4259
4260 notes:
4261 called by qh_all_merges and qh_postmerge
4262 ignored if 2-d
4263
4264 design:
4265 merge any degenerate or redundant facets
4266 repeat until no more degenerate or redundant facets
4267 for each newly merged facet
4268 remove extra vertices
4269 if qh.MERGEvertices
4270 for each newly merged facet
4271 for each vertex
4272 if vertex was on a deleted ridge
4273 rename vertex if it is shared
4274 for each new, undeleted vertex
4275 remove delridge flag
4276 if vertex is redundant
4277 merge degenerate or redundant facets
4278 */
qh_reducevertices(void)4279 boolT qh_reducevertices(void) {
4280 int numshare=0, numrename= 0;
4281 boolT degenredun= False;
4282 facetT *newfacet;
4283 vertexT *vertex, **vertexp;
4284
4285 if (qh hull_dim == 2)
4286 return False;
4287 trace2((qh ferr, 2101, "qh_reducevertices: reduce extra vertices, shared vertices, and redundant vertices\n"));
4288 if (qh_merge_degenredundant())
4289 degenredun= True;
4290 LABELrestart:
4291 FORALLnew_facets {
4292 if (newfacet->newmerge) {
4293 if (!qh MERGEvertices)
4294 newfacet->newmerge= False;
4295 if (qh_remove_extravertices(newfacet)) {
4296 qh_degen_redundant_facet(newfacet);
4297 if (qh_merge_degenredundant()) {
4298 degenredun= True;
4299 goto LABELrestart;
4300 }
4301 }
4302 }
4303 }
4304 if (!qh MERGEvertices)
4305 return False;
4306 FORALLnew_facets {
4307 if (newfacet->newmerge) {
4308 newfacet->newmerge= False;
4309 FOREACHvertex_(newfacet->vertices) {
4310 if (vertex->delridge) {
4311 if (qh_rename_sharedvertex(vertex, newfacet)) {
4312 numshare++;
4313 if (qh_merge_degenredundant()) {
4314 degenredun= True;
4315 goto LABELrestart;
4316 }
4317 vertexp--; /* repeat since deleted vertex */
4318 }
4319 }
4320 }
4321 }
4322 }
4323 FORALLvertex_(qh newvertex_list) {
4324 if (vertex->delridge && !vertex->deleted) {
4325 vertex->delridge= False;
4326 if (qh hull_dim >= 4 && qh_redundant_vertex(vertex)) {
4327 numrename++;
4328 if (qh_merge_degenredundant()) {
4329 degenredun= True;
4330 goto LABELrestart;
4331 }
4332 }
4333 }
4334 }
4335 trace1((qh ferr, 1014, "qh_reducevertices: renamed %d shared vertices and %d redundant vertices. Degen? %d\n",
4336 numshare, numrename, degenredun));
4337 return degenredun;
4338 } /* reducevertices */
4339
4340 /*-<a href="qh-merge.htm#TOC"
4341 >-------------------------------</a><a name="redundant_vertex">-</a>
4342
4343 qh_redundant_vertex( vertex )
4344 rename a redundant vertex if qh_find_newvertex succeeds
4345 assumes vertices have full vertex->neighbors
4346
4347 returns:
4348 if find a replacement vertex
4349 returns new vertex
4350 qh_renamevertex sets vertex->deleted for redundant vertex
4351
4352 notes:
4353 only called by qh_reducevertices for vertex->delridge and hull_dim >= 4
4354 may add degenerate facets to qh.facet_mergeset
4355 doesn't change vertex->neighbors or create redundant facets
4356
4357 design:
4358 intersect vertices of all facet neighbors of vertex
4359 determine ridges for these vertices
4360 if find a new vertex for vertex among these ridges and vertices
4361 rename vertex to the new vertex
4362 */
qh_redundant_vertex(vertexT * vertex)4363 vertexT *qh_redundant_vertex(vertexT *vertex) {
4364 vertexT *newvertex= NULL;
4365 setT *vertices, *ridges;
4366
4367 trace3((qh ferr, 3008, "qh_redundant_vertex: check if v%d from a deleted ridge can be renamed\n", vertex->id));
4368 if ((vertices= qh_neighbor_intersections(vertex))) {
4369 ridges= qh_vertexridges(vertex, !qh_ALL);
4370 if ((newvertex= qh_find_newvertex(vertex, vertices, ridges))) {
4371 zinc_(Zrenameall);
4372 qh_renamevertex(vertex, newvertex, ridges, NULL, NULL); /* ridges invalidated */
4373 }
4374 qh_settempfree(&ridges);
4375 qh_settempfree(&vertices);
4376 }
4377 return newvertex;
4378 } /* redundant_vertex */
4379
4380 /*-<a href="qh-merge.htm#TOC"
4381 >-------------------------------</a><a name="remove_extravertices">-</a>
4382
4383 qh_remove_extravertices( facet )
4384 remove extra vertices from non-simplicial facets
4385
4386 returns:
4387 returns True if it finds them
4388 deletes facet from vertex neighbors
4389 facet may be redundant (test with qh_degen_redundant)
4390
4391 notes:
4392 called by qh_renamevertex and qh_reducevertices
4393 a merge (qh_reducevertices) or qh_renamevertex may drop all ridges for a vertex in a facet
4394
4395 design:
4396 for each vertex in facet
4397 if vertex not in a ridge (i.e., no longer used)
4398 delete vertex from facet
4399 delete facet from vertex's neighbors
4400 unless vertex in another facet
4401 add vertex to qh.del_vertices for later deletion
4402 */
qh_remove_extravertices(facetT * facet)4403 boolT qh_remove_extravertices(facetT *facet) {
4404 ridgeT *ridge, **ridgep;
4405 vertexT *vertex, **vertexp;
4406 boolT foundrem= False;
4407
4408 if (facet->simplicial) {
4409 return False;
4410 }
4411 trace4((qh ferr, 4043, "qh_remove_extravertices: test non-simplicial f%d for extra vertices\n",
4412 facet->id));
4413 FOREACHvertex_(facet->vertices)
4414 vertex->seen= False;
4415 FOREACHridge_(facet->ridges) {
4416 FOREACHvertex_(ridge->vertices)
4417 vertex->seen= True;
4418 }
4419 FOREACHvertex_(facet->vertices) {
4420 if (!vertex->seen) {
4421 foundrem= True;
4422 zinc_(Zremvertex);
4423 qh_setdelsorted(facet->vertices, vertex);
4424 qh_setdel(vertex->neighbors, facet);
4425 if (!qh_setsize(vertex->neighbors)) {
4426 vertex->deleted= True;
4427 qh_setappend(&qh del_vertices, vertex);
4428 zinc_(Zremvertexdel);
4429 trace2((qh ferr, 2036, "qh_remove_extravertices: v%d deleted because it's lost all ridges\n", vertex->id));
4430 }else
4431 trace3((qh ferr, 3009, "qh_remove_extravertices: v%d removed from f%d because it's lost all ridges\n", vertex->id, facet->id));
4432 vertexp--; /*repeat*/
4433 }
4434 }
4435 return foundrem;
4436 } /* remove_extravertices */
4437
4438 /*-<a href="qh-merge.htm#TOC"
4439 >-------------------------------</a><a name="remove_mergetype">-</a>
4440
4441 qh_remove_mergetype( mergeset, mergetype )
4442 Remove mergetype merges from mergeset
4443
4444 notes:
4445 Does not preserve order
4446 */
qh_remove_mergetype(setT * mergeset,mergeType type)4447 void qh_remove_mergetype(setT *mergeset, mergeType type) {
4448 mergeT *merge;
4449 int merge_i, merge_n;
4450
4451 FOREACHmerge_i_(mergeset) {
4452 if (merge->mergetype == type) {
4453 trace3((qh ferr, 3037, "qh_remove_mergetype: remove merge f%d f%d v%d v%d r%d r%d dist %2.2g type %d",
4454 getid_(merge->facet1), getid_(merge->facet2), getid_(merge->vertex1), getid_(merge->vertex2), getid_(merge->ridge1), getid_(merge->ridge2), merge->distance, type));
4455 qh_setdelnth(mergeset, merge_i);
4456 merge_i--; merge_n--; /* repeat with next merge */
4457 }
4458 }
4459 } /* remove_mergetype */
4460
4461 /*-<a href="qh-merge.htm#TOC"
4462 >-------------------------------</a><a name="rename_adjacentvertex">-</a>
4463
4464 qh_rename_adjacentvertex( oldvertex, newvertex )
4465 renames oldvertex as newvertex. Must be adjacent (i.e., in the same subridge)
4466 no-op if either vertex is deleted
4467
4468 notes:
4469 called from qh_merge_pinchedvertices
4470
4471 design:
4472 for all neighbors of oldvertex
4473 if simplicial, rename oldvertex to newvertex and drop if degenerate
4474 if needed, add oldvertex neighbor to newvertex
4475 determine ridges for oldvertex
4476 rename oldvertex as newvertex in ridges (qh_renamevertex)
4477 */
qh_rename_adjacentvertex(vertexT * oldvertex,vertexT * newvertex,realT dist)4478 void qh_rename_adjacentvertex(vertexT *oldvertex, vertexT *newvertex, realT dist) {
4479 setT *ridges;
4480 facetT *neighbor, **neighborp, *maxfacet= NULL;
4481 ridgeT *ridge, **ridgep;
4482 boolT istrace= False;
4483 int oldsize= qh_setsize(oldvertex->neighbors);
4484 int newsize= qh_setsize(newvertex->neighbors);
4485 coordT maxdist2= -REALmax, dist2;
4486
4487 if (qh IStracing >= 4 || oldvertex->id == qh tracevertex_id || newvertex->id == qh tracevertex_id) {
4488 istrace= True;
4489 }
4490 zzinc_(Ztotmerge);
4491 if (istrace) {
4492 qh_fprintf(qh ferr, 2071, "qh_rename_adjacentvertex: merge #%d rename v%d (%d neighbors) to v%d (%d neighbors) dist %2.2g\n",
4493 zzval_(Ztotmerge), oldvertex->id, oldsize, newvertex->id, newsize, dist);
4494 }
4495 if (oldvertex->deleted || newvertex->deleted) {
4496 if (istrace || qh IStracing >= 2) {
4497 qh_fprintf(qh ferr, 2072, "qh_rename_adjacentvertex: ignore rename. Either v%d (%d) or v%d (%d) is deleted\n",
4498 oldvertex->id, oldvertex->deleted, newvertex->id, newvertex->deleted);
4499 }
4500 return;
4501 }
4502 if (oldsize == 0 || newsize == 0) {
4503 qh_fprintf(qh ferr, 2072, "qhull internal error (qh_rename_adjacentvertex): expecting neighbor facets for v%d and v%d. Got %d and %d neighbors resp.\n",
4504 oldvertex->id, newvertex->id, oldsize, newsize);
4505 qh_errexit(qh_ERRqhull, NULL, NULL);
4506 }
4507 FOREACHneighbor_(oldvertex) {
4508 if (neighbor->simplicial) {
4509 if (qh_setin(neighbor->vertices, newvertex)) {
4510 if (istrace || qh IStracing >= 2) {
4511 qh_fprintf(qh ferr, 2070, "qh_rename_adjacentvertex: simplicial f%d contains old v%d and new v%d. Will be marked degenerate by qh_renamevertex\n",
4512 neighbor->id, oldvertex->id, newvertex->id);
4513 }
4514 qh_makeridges(neighbor); /* no longer simplicial, nummerge==0, skipped by qh_maybe_duplicateridge */
4515 }else {
4516 qh_replacefacetvertex(neighbor, oldvertex, newvertex);
4517 qh_setunique(&newvertex->neighbors, neighbor);
4518 qh_newvertices(neighbor->vertices); /* for qh_update_vertexneighbors of vertex neighbors */
4519 }
4520 }
4521 }
4522 ridges= qh_vertexridges(oldvertex, qh_ALL);
4523 if (istrace) {
4524 FOREACHridge_(ridges) {
4525 qh_printridge(qh ferr, ridge);
4526 }
4527 }
4528 FOREACHneighbor_(oldvertex) {
4529 if (!neighbor->simplicial){
4530 qh_addfacetvertex(neighbor, newvertex);
4531 qh_setunique(&newvertex->neighbors, neighbor);
4532 qh_newvertices(neighbor->vertices); /* for qh_update_vertexneighbors of vertex neighbors */
4533 if (qh newfacet_list == qh facet_tail) {
4534 qh_removefacet(neighbor); /* add a neighbor to newfacet_list so that qh_partitionvisible has a newfacet */
4535 qh_appendfacet(neighbor);
4536 neighbor->newfacet= True;
4537 }
4538 }
4539 }
4540 qh_renamevertex(oldvertex, newvertex, ridges, NULL, NULL); /* ridges invalidated */
4541 if (oldvertex->deleted && !oldvertex->partitioned) {
4542 FOREACHneighbor_(newvertex) {
4543 if (!neighbor->visible) {
4544 qh_distplane(oldvertex->point, neighbor, &dist2);
4545 if (dist2>maxdist2) {
4546 maxdist2= dist2;
4547 maxfacet= neighbor;
4548 }
4549 }
4550 }
4551 trace2((qh ferr, 2096, "qh_rename_adjacentvertex: partition old p%d(v%d) as a coplanar point for furthest f%d dist %2.2g. Maybe repartition later (QH0031)\n",
4552 qh_pointid(oldvertex->point), oldvertex->id, maxfacet->id, maxdist2))
4553 qh_partitioncoplanar(oldvertex->point, maxfacet, NULL, !qh_ALL); /* faster with maxdist2, otherwise duplicates distance tests from maxdist2/dist2 */
4554 oldvertex->partitioned= True;
4555 }
4556 qh_settempfree(&ridges);
4557 } /* rename_adjacentvertex */
4558
4559 /*-<a href="qh-merge.htm#TOC"
4560 >-------------------------------</a><a name="rename_sharedvertex">-</a>
4561
4562 qh_rename_sharedvertex( vertex, facet )
4563 detect and rename if shared vertex in facet
4564 vertices have full ->neighbors
4565
4566 returns:
4567 newvertex or NULL
4568 the vertex may still exist in other facets (i.e., a neighbor was pinched)
4569 does not change facet->neighbors
4570 updates vertex->neighbors
4571
4572 notes:
4573 only called by qh_reducevertices after qh_remove_extravertices
4574 so f.vertices does not contain extraneous vertices
4575 a shared vertex for a facet is only in ridges to one neighbor
4576 this may undo a pinched facet
4577
4578 it does not catch pinches involving multiple facets. These appear
4579 to be difficult to detect, since an exhaustive search is too expensive.
4580
4581 design:
4582 if vertex only has two neighbors
4583 determine the ridges that contain the vertex
4584 determine the vertices shared by both neighbors
4585 if can find a new vertex in this set
4586 rename the vertex to the new vertex
4587 */
qh_rename_sharedvertex(vertexT * vertex,facetT * facet)4588 vertexT *qh_rename_sharedvertex(vertexT *vertex, facetT *facet) {
4589 facetT *neighbor, **neighborp, *neighborA= NULL;
4590 setT *vertices, *ridges;
4591 vertexT *newvertex= NULL;
4592
4593 if (qh_setsize(vertex->neighbors) == 2) {
4594 neighborA= SETfirstt_(vertex->neighbors, facetT);
4595 if (neighborA == facet)
4596 neighborA= SETsecondt_(vertex->neighbors, facetT);
4597 }else if (qh hull_dim == 3)
4598 return NULL;
4599 else {
4600 qh visit_id++;
4601 FOREACHneighbor_(facet)
4602 neighbor->visitid= qh visit_id;
4603 FOREACHneighbor_(vertex) {
4604 if (neighbor->visitid == qh visit_id) {
4605 if (neighborA)
4606 return NULL;
4607 neighborA= neighbor;
4608 }
4609 }
4610 }
4611 if (!neighborA) {
4612 qh_fprintf(qh ferr, 6101, "qhull internal error (qh_rename_sharedvertex): v%d's neighbors not in f%d\n",
4613 vertex->id, facet->id);
4614 qh_errprint("ERRONEOUS", facet, NULL, NULL, vertex);
4615 qh_errexit(qh_ERRqhull, NULL, NULL);
4616 }
4617 if (neighborA) { /* avoid warning */
4618 /* the vertex is shared by facet and neighborA */
4619 ridges= qh_settemp(qh TEMPsize);
4620 neighborA->visitid= ++qh visit_id;
4621 qh_vertexridges_facet(vertex, facet, &ridges);
4622 trace2((qh ferr, 2037, "qh_rename_sharedvertex: p%d(v%d) is shared by f%d(%d ridges) and f%d\n",
4623 qh_pointid(vertex->point), vertex->id, facet->id, qh_setsize(ridges), neighborA->id));
4624 zinc_(Zintersectnum);
4625 vertices= qh_vertexintersect_new(facet->vertices, neighborA->vertices);
4626 qh_setdel(vertices, vertex);
4627 qh_settemppush(vertices);
4628 if ((newvertex= qh_find_newvertex(vertex, vertices, ridges)))
4629 qh_renamevertex(vertex, newvertex, ridges, facet, neighborA); /* ridges invalidated */
4630 qh_settempfree(&vertices);
4631 qh_settempfree(&ridges);
4632 }
4633 return newvertex;
4634 } /* rename_sharedvertex */
4635
4636 /*-<a href="qh-merge.htm#TOC"
4637 >-------------------------------</a><a name="renameridgevertex">-</a>
4638
4639 qh_renameridgevertex( ridge, oldvertex, newvertex )
4640 renames oldvertex as newvertex in ridge
4641
4642 returns:
4643 True if renames oldvertex
4644 False if deleted the ridge
4645
4646 notes:
4647 called by qh_renamevertex
4648 caller sets newvertex->delridge for deleted ridges (qh_reducevertices)
4649
4650 design:
4651 delete oldvertex from ridge
4652 if newvertex already in ridge
4653 copy ridge->noconvex to another ridge if possible
4654 delete the ridge
4655 else
4656 insert newvertex into the ridge
4657 adjust the ridge's orientation
4658 */
qh_renameridgevertex(ridgeT * ridge,vertexT * oldvertex,vertexT * newvertex)4659 boolT qh_renameridgevertex(ridgeT *ridge, vertexT *oldvertex, vertexT *newvertex) {
4660 int nth= 0, oldnth;
4661 facetT *temp;
4662 vertexT *vertex, **vertexp;
4663
4664 oldnth= qh_setindex(ridge->vertices, oldvertex);
4665 if (oldnth < 0) {
4666 qh_fprintf(qh ferr, 6424, "qhull internal error (qh_renameridgevertex): oldvertex v%d not found in r%d. Cannot rename to v%d\n",
4667 oldvertex->id, ridge->id, newvertex->id);
4668 qh_errexit(qh_ERRqhull, NULL, ridge);
4669 }
4670 qh_setdelnthsorted(ridge->vertices, oldnth);
4671 FOREACHvertex_(ridge->vertices) {
4672 if (vertex == newvertex) {
4673 zinc_(Zdelridge);
4674 if (ridge->nonconvex) /* only one ridge has nonconvex set */
4675 qh_copynonconvex(ridge);
4676 trace2((qh ferr, 2038, "qh_renameridgevertex: ridge r%d deleted. It contained both v%d and v%d\n",
4677 ridge->id, oldvertex->id, newvertex->id));
4678 qh_delridge_merge(ridge); /* ridge.vertices deleted */
4679 return False;
4680 }
4681 if (vertex->id < newvertex->id)
4682 break;
4683 nth++;
4684 }
4685 qh_setaddnth(&ridge->vertices, nth, newvertex);
4686 ridge->simplicialtop= False;
4687 ridge->simplicialbot= False;
4688 if (abs(oldnth - nth)%2) {
4689 trace3((qh ferr, 3010, "qh_renameridgevertex: swapped the top and bottom of ridge r%d\n",
4690 ridge->id));
4691 temp= ridge->top;
4692 ridge->top= ridge->bottom;
4693 ridge->bottom= temp;
4694 }
4695 return True;
4696 } /* renameridgevertex */
4697
4698
4699 /*-<a href="qh-merge.htm#TOC"
4700 >-------------------------------</a><a name="renamevertex">-</a>
4701
4702 qh_renamevertex( oldvertex, newvertex, ridges, oldfacet, neighborA )
4703 renames oldvertex as newvertex in ridges of non-simplicial neighbors
4704 set oldfacet/neighborA if oldvertex is shared between two facets (qh_rename_sharedvertex)
4705 otherwise qh_redundant_vertex or qh_rename_adjacentvertex
4706
4707 returns:
4708 if oldfacet and multiple neighbors, oldvertex may still exist afterwards
4709 otherwise sets oldvertex->deleted for later deletion
4710 one or more ridges maybe deleted
4711 ridges is invalidated
4712 merges may be added to degen_mergeset via qh_maydropneighbor or qh_degen_redundant_facet
4713
4714 notes:
4715 qh_rename_sharedvertex can not change neighbors of newvertex (since it's a subset)
4716 qh_redundant_vertex due to vertex->delridge for qh_reducevertices
4717 qh_rename_adjacentvertex for complete renames
4718
4719 design:
4720 for each ridge in ridges
4721 rename oldvertex to newvertex and delete degenerate ridges
4722 if oldfacet not defined
4723 for each non-simplicial neighbor of oldvertex
4724 delete oldvertex from neighbor's vertices
4725 remove extra vertices from neighbor
4726 add oldvertex to qh.del_vertices
4727 else if oldvertex only between oldfacet and neighborA
4728 delete oldvertex from oldfacet and neighborA
4729 add oldvertex to qh.del_vertices
4730 else oldvertex is in oldfacet and neighborA and other facets (i.e., pinched)
4731 delete oldvertex from oldfacet
4732 delete oldfacet from old vertex's neighbors
4733 remove extra vertices (e.g., oldvertex) from neighborA
4734 */
qh_renamevertex(vertexT * oldvertex,vertexT * newvertex,setT * ridges,facetT * oldfacet,facetT * neighborA)4735 void qh_renamevertex(vertexT *oldvertex, vertexT *newvertex, setT *ridges, facetT *oldfacet, facetT *neighborA) {
4736 facetT *neighbor, **neighborp;
4737 ridgeT *ridge, **ridgep;
4738 int topsize, bottomsize;
4739 boolT istrace= False;
4740
4741 #ifndef qh_NOtrace
4742 if (qh IStracing >= 2 || oldvertex->id == qh tracevertex_id ||
4743 newvertex->id == qh tracevertex_id) {
4744 istrace= True;
4745 qh_fprintf(qh ferr, 2086, "qh_renamevertex: rename v%d to v%d in %d ridges with old f%d and neighbor f%d\n",
4746 oldvertex->id, newvertex->id, qh_setsize(ridges), getid_(oldfacet), getid_(neighborA));
4747 }
4748 #endif
4749 FOREACHridge_(ridges) {
4750 if (qh_renameridgevertex(ridge, oldvertex, newvertex)) { /* ridge is deleted if False, invalidating ridges */
4751 topsize= qh_setsize(ridge->top->vertices);
4752 bottomsize= qh_setsize(ridge->bottom->vertices);
4753 if (topsize < qh hull_dim || (topsize == qh hull_dim && !ridge->top->simplicial && qh_setin(ridge->top->vertices, newvertex))) {
4754 trace4((qh ferr, 4070, "qh_renamevertex: ignore duplicate check for r%d. top f%d (size %d) will be degenerate after rename v%d to v%d\n",
4755 ridge->id, ridge->top->id, topsize, oldvertex->id, newvertex->id));
4756 }else if (bottomsize < qh hull_dim || (bottomsize == qh hull_dim && !ridge->bottom->simplicial && qh_setin(ridge->bottom->vertices, newvertex))) {
4757 trace4((qh ferr, 4071, "qh_renamevertex: ignore duplicate check for r%d. bottom f%d (size %d) will be degenerate after rename v%d to v%d\n",
4758 ridge->id, ridge->bottom->id, bottomsize, oldvertex->id, newvertex->id));
4759 }else
4760 qh_maybe_duplicateridge(ridge);
4761 }
4762 }
4763 if (!oldfacet) {
4764 /* stat Zrenameall or Zpinchduplicate */
4765 if (istrace)
4766 qh_fprintf(qh ferr, 2087, "qh_renamevertex: renaming v%d to v%d in several facets for qh_redundant_vertex or MRGsubridge\n",
4767 oldvertex->id, newvertex->id);
4768 FOREACHneighbor_(oldvertex) {
4769 if (neighbor->simplicial) {
4770 qh_degen_redundant_facet(neighbor); /* e.g., rbox 175 C3,2e-13 D4 t1545235541 | qhull d */
4771 }else {
4772 if (istrace)
4773 qh_fprintf(qh ferr, 4080, "qh_renamevertex: rename vertices in non-simplicial neighbor f%d of v%d\n", neighbor->id, oldvertex->id);
4774 qh_maydropneighbor(neighbor);
4775 qh_setdelsorted(neighbor->vertices, oldvertex); /* if degenerate, qh_degen_redundant_facet will add to mergeset */
4776 if (qh_remove_extravertices(neighbor))
4777 neighborp--; /* neighbor deleted from oldvertex neighborset */
4778 qh_degen_redundant_facet(neighbor); /* either direction may be redundant, faster if combine? */
4779 qh_test_redundant_neighbors(neighbor);
4780 qh_test_degen_neighbors(neighbor);
4781 }
4782 }
4783 if (!oldvertex->deleted) {
4784 oldvertex->deleted= True;
4785 qh_setappend(&qh del_vertices, oldvertex);
4786 }
4787 }else if (qh_setsize(oldvertex->neighbors) == 2) {
4788 zinc_(Zrenameshare);
4789 if (istrace)
4790 qh_fprintf(qh ferr, 3039, "qh_renamevertex: renaming v%d to v%d in oldfacet f%d for qh_rename_sharedvertex\n",
4791 oldvertex->id, newvertex->id, oldfacet->id);
4792 FOREACHneighbor_(oldvertex) {
4793 qh_setdelsorted(neighbor->vertices, oldvertex);
4794 qh_degen_redundant_facet(neighbor);
4795 }
4796 oldvertex->deleted= True;
4797 qh_setappend(&qh del_vertices, oldvertex);
4798 }else {
4799 zinc_(Zrenamepinch);
4800 if (istrace || qh IStracing >= 1)
4801 qh_fprintf(qh ferr, 3040, "qh_renamevertex: renaming pinched v%d to v%d between f%d and f%d\n",
4802 oldvertex->id, newvertex->id, oldfacet->id, neighborA->id);
4803 qh_setdelsorted(oldfacet->vertices, oldvertex);
4804 qh_setdel(oldvertex->neighbors, oldfacet);
4805 if (qh_remove_extravertices(neighborA))
4806 qh_degen_redundant_facet(neighborA);
4807 }
4808 if (oldfacet)
4809 qh_degen_redundant_facet(oldfacet);
4810 } /* renamevertex */
4811
4812 /*-<a href="qh-merge.htm#TOC"
4813 >-------------------------------</a><a name="test_appendmerge">-</a>
4814
4815 qh_test_appendmerge( facet, neighbor, simplicial )
4816 test convexity and append to qh.facet_mergeset if non-convex
4817 if pre-merging,
4818 no-op if qh.SKIPconvex, or qh.MERGEexact and coplanar
4819 if simplicial, assumes centrum test is valid (e.g., adjacent, simplicial new facets)
4820
4821 returns:
4822 true if appends facet/neighbor to qh.facet_mergeset
4823 sets facet->center as needed
4824 does not change facet->seen
4825
4826 notes:
4827 called from qh_getmergeset_initial, qh_getmergeset, and qh_test_vneighbors
4828 must be at least as strong as qh_checkconvex (poly2.c)
4829 assumes !f.flipped
4830
4831 design:
4832 exit if qh.SKIPconvex ('Q0') and !qh.POSTmerging
4833 if qh.cos_max ('An') is defined and merging coplanars
4834 if the angle between facet normals is too shallow
4835 append an angle-coplanar merge to qh.mergeset
4836 return True
4837 test convexity of facet and neighbor
4838 */
qh_test_appendmerge(facetT * facet,facetT * neighbor,boolT simplicial)4839 boolT qh_test_appendmerge(facetT *facet, facetT *neighbor, boolT simplicial) {
4840 realT angle= -REALmax;
4841 boolT okangle= False;
4842
4843 if (qh SKIPconvex && !qh POSTmerging)
4844 return False;
4845 if (qh cos_max < REALmax/2 && (!qh MERGEexact || qh POSTmerging)) {
4846 angle= qh_getangle(facet->normal, neighbor->normal);
4847 okangle= True;
4848 zinc_(Zangletests);
4849 if (angle > qh cos_max) {
4850 zinc_(Zcoplanarangle);
4851 qh_appendmergeset(facet, neighbor, MRGanglecoplanar, 0.0, angle);
4852 trace2((qh ferr, 2039, "qh_test_appendmerge: coplanar angle %4.4g between f%d and f%d\n",
4853 angle, facet->id, neighbor->id));
4854 return True;
4855 }
4856 }
4857 if (simplicial || qh hull_dim <= 3)
4858 return qh_test_centrum_merge(facet, neighbor, angle, okangle);
4859 else
4860 return qh_test_nonsimplicial_merge(facet, neighbor, angle, okangle);
4861 } /* test_appendmerge */
4862
4863 /*-<a href="qh-merge.htm#TOC"
4864 >-------------------------------</a><a name="test_centrum_merge">-</a>
4865
4866 qh_test_centrum_merge( facet, neighbor, angle, okangle )
4867 test centrum convexity and append non-convex facets to qh.facet_mergeset
4868 'angle' is angle between facets if okangle is true, otherwise use 0.0
4869
4870 returns:
4871 true if append facet/neighbor to qh.facet_mergeset
4872 sets facet->center as needed
4873 does not change facet->seen
4874
4875 notes:
4876 called from test_appendmerge if adjacent simplicial facets or 2-d/3-d
4877 at least as strict as qh_checkconvex, including qh.DISTround ('En' and 'Rn')
4878
4879 design:
4880 make facet's centrum if needed
4881 if facet's centrum is above the neighbor (qh.centrum_radius)
4882 set isconcave
4883
4884 if facet's centrum is not below the neighbor (-qh.centrum_radius)
4885 set iscoplanar
4886 make neighbor's centrum if needed
4887 if neighbor's centrum is above the facet
4888 set isconcave
4889 else if neighbor's centrum is not below the facet
4890 set iscoplanar
4891 if isconcave or iscoplanar and merging coplanars
4892 get angle if needed (qh.ANGLEmerge 'An')
4893 append concave-coplanar, concave ,or coplanar merge to qh.mergeset
4894 */
qh_test_centrum_merge(facetT * facet,facetT * neighbor,realT angle,boolT okangle)4895 boolT qh_test_centrum_merge(facetT *facet, facetT *neighbor, realT angle, boolT okangle) {
4896 coordT dist, dist2, mergedist;
4897 boolT isconcave= False, iscoplanar= False;
4898
4899 if (!facet->center)
4900 facet->center= qh_getcentrum(facet);
4901 zzinc_(Zcentrumtests);
4902 qh_distplane(facet->center, neighbor, &dist);
4903 if (dist > qh centrum_radius)
4904 isconcave= True;
4905 else if (dist >= -qh centrum_radius)
4906 iscoplanar= True;
4907 if (!neighbor->center)
4908 neighbor->center= qh_getcentrum(neighbor);
4909 zzinc_(Zcentrumtests);
4910 qh_distplane(neighbor->center, facet, &dist2);
4911 if (dist2 > qh centrum_radius)
4912 isconcave= True;
4913 else if (!iscoplanar && dist2 >= -qh centrum_radius)
4914 iscoplanar= True;
4915 if (!isconcave && (!iscoplanar || (qh MERGEexact && !qh POSTmerging)))
4916 return False;
4917 if (!okangle && qh ANGLEmerge) {
4918 angle= qh_getangle(facet->normal, neighbor->normal);
4919 zinc_(Zangletests);
4920 }
4921 if (isconcave && iscoplanar) {
4922 zinc_(Zconcavecoplanarridge);
4923 if (dist > dist2)
4924 qh_appendmergeset(facet, neighbor, MRGconcavecoplanar, dist, angle);
4925 else
4926 qh_appendmergeset(neighbor, facet, MRGconcavecoplanar, dist2, angle);
4927 trace0((qh ferr, 36, "qh_test_centrum_merge: concave f%d to coplanar f%d, dist %4.4g and reverse dist %4.4g, angle %4.4g during p%d\n",
4928 facet->id, neighbor->id, dist, dist2, angle, qh furthest_id));
4929 }else if (isconcave) {
4930 mergedist= fmax_(dist, dist2);
4931 zinc_(Zconcaveridge);
4932 qh_appendmergeset(facet, neighbor, MRGconcave, mergedist, angle);
4933 trace0((qh ferr, 37, "qh_test_centrum_merge: concave f%d to f%d, dist %4.4g and reverse dist %4.4g, angle %4.4g during p%d\n",
4934 facet->id, neighbor->id, dist, dist2, angle, qh furthest_id));
4935 }else /* iscoplanar */ {
4936 mergedist= fmin_(fabs_(dist), fabs_(dist2));
4937 zinc_(Zcoplanarcentrum);
4938 qh_appendmergeset(facet, neighbor, MRGcoplanar, mergedist, angle);
4939 trace2((qh ferr, 2097, "qh_test_centrum_merge: coplanar f%d to f%d dist %4.4g, reverse dist %4.4g angle %4.4g\n",
4940 facet->id, neighbor->id, dist, dist2, angle));
4941 }
4942 return True;
4943 } /* test_centrum_merge */
4944
4945 /*-<a href="qh-merge.htm#TOC"
4946 >-------------------------------</a><a name="test_degen_neighbors">-</a>
4947
4948 qh_test_degen_neighbors( facet )
4949 append degenerate neighbors to qh.degen_mergeset
4950
4951 notes:
4952 called at end of qh_mergefacet() and qh_renamevertex()
4953 call after test_redundant_facet() since MRGredundant is less expensive then MRGdegen
4954 a degenerate facet has fewer than hull_dim neighbors
4955 see: qh_merge_degenredundant()
4956
4957 */
qh_test_degen_neighbors(facetT * facet)4958 void qh_test_degen_neighbors(facetT *facet) {
4959 facetT *neighbor, **neighborp;
4960 int size;
4961
4962 trace4((qh ferr, 4073, "qh_test_degen_neighbors: test for degenerate neighbors of f%d\n", facet->id));
4963 FOREACHneighbor_(facet) {
4964 if (neighbor->visible) {
4965 qh_fprintf(qh ferr, 6359, "qhull internal error (qh_test_degen_neighbors): facet f%d has deleted neighbor f%d (qh.visible_list)\n",
4966 facet->id, neighbor->id);
4967 qh_errexit2(qh_ERRqhull, facet, neighbor);
4968 }
4969 if (neighbor->degenerate || neighbor->redundant || neighbor->dupridge) /* will merge or delete */
4970 continue;
4971 /* merge flipped-degenerate facet before flipped facets */
4972 if ((size= qh_setsize(neighbor->neighbors)) < qh hull_dim) {
4973 qh_appendmergeset(neighbor, neighbor, MRGdegen, 0.0, 1.0);
4974 trace2((qh ferr, 2019, "qh_test_degen_neighbors: f%d is degenerate with %d neighbors. Neighbor of f%d.\n", neighbor->id, size, facet->id));
4975 }
4976 }
4977 } /* test_degen_neighbors */
4978
4979
4980 /*-<a href="qh-merge.htm#TOC"
4981 >-------------------------------</a><a name="test_nonsimplicial_merge">-</a>
4982
4983 qh_test_nonsimplicial_merge( facet, neighbor, angle, okangle )
4984 test centrum and vertex convexity and append non-convex or redundant facets to qh.facet_mergeset
4985 'angle' is angle between facets if okangle is true, otherwise use 0.0
4986 skips coplanar merges if pre-merging with qh.MERGEexact ('Qx')
4987
4988 returns:
4989 true if appends facet/neighbor to qh.facet_mergeset
4990 sets facet->center as needed
4991 does not change facet->seen
4992
4993 notes:
4994 only called from test_appendmerge if a non-simplicial facet and at least 4-d
4995 at least as strict as qh_checkconvex, including qh.DISTround ('En' and 'Rn')
4996 centrums must be < -qh.centrum_radius
4997 tests vertices as well as centrums since a facet may be twisted relative to its neighbor
4998
4999 design:
5000 set precision constants for maxoutside, clearlyconcave, minvertex, and coplanarcentrum
5001 use maxoutside for coplanarcentrum if premerging with 'Qx' and qh_MAXcoplanarcentrum merges
5002 otherwise use qh.centrum_radious for coplanarcentrum
5003 make facet and neighbor centrums if needed
5004 isconcave if a centrum is above neighbor (coplanarcentrum)
5005 iscoplanar if a centrum is not below neighbor (-qh.centrum_radius)
5006 maybeconvex if a centrum is clearly below neighbor (-clearyconvex)
5007 return False if both centrums clearly below neighbor (-clearyconvex)
5008 return MRGconcave if isconcave
5009
5010 facets are neither clearly convex nor clearly concave
5011 test vertices as well as centrums
5012 if maybeconvex
5013 determine mindist and maxdist for vertices of the other facet
5014 maybe MRGredundant
5015 otherwise
5016 determine mindist and maxdist for vertices of either facet
5017 maybe MRGredundant
5018 maybeconvex if a vertex is clearly below neighbor (-clearconvex)
5019
5020 vertices are concave if dist > clearlyconcave
5021 vertices are twisted if dist > maxoutside (isconcave and maybeconvex)
5022 return False if not concave and pre-merge of 'Qx' (qh.MERGEexact)
5023 vertices are coplanar if dist in -minvertex..maxoutside
5024 if !isconcave, vertices are coplanar if dist >= -qh.MAXcoplanar (n*qh.premerge_centrum)
5025
5026 return False if neither concave nor coplanar
5027 return MRGtwisted if isconcave and maybeconvex
5028 return MRGconcavecoplanar if isconcave and isconvex
5029 return MRGconcave if isconcave
5030 return MRGcoplanar if iscoplanar
5031 */
qh_test_nonsimplicial_merge(facetT * facet,facetT * neighbor,realT angle,boolT okangle)5032 boolT qh_test_nonsimplicial_merge(facetT *facet, facetT *neighbor, realT angle, boolT okangle) {
5033 coordT dist, mindist, maxdist, mindist2, maxdist2, dist2, maxoutside, clearlyconcave, minvertex, clearlyconvex, mergedist, coplanarcentrum;
5034 boolT isconcave= False, iscoplanar= False, maybeconvex= False, isredundant= False;
5035 vertexT *maxvertex= NULL, *maxvertex2= NULL;
5036
5037 maxoutside= fmax_(neighbor->maxoutside, qh ONEmerge + qh DISTround);
5038 maxoutside= fmax_(maxoutside, facet->maxoutside);
5039 clearlyconcave= qh_RATIOconcavehorizon * maxoutside;
5040 minvertex= fmax_(-qh min_vertex, qh MAXcoplanar); /* non-negative, not available per facet, not used for iscoplanar */
5041 clearlyconvex= qh_RATIOconvexmerge * minvertex; /* must be convex for MRGtwisted */
5042 if (qh MERGEexact && !qh POSTmerging && (facet->nummerge > qh_MAXcoplanarcentrum || neighbor->nummerge > qh_MAXcoplanarcentrum))
5043 coplanarcentrum= maxoutside;
5044 else
5045 coplanarcentrum= qh centrum_radius;
5046
5047 if (!facet->center)
5048 facet->center= qh_getcentrum(facet);
5049 zzinc_(Zcentrumtests);
5050 qh_distplane(facet->center, neighbor, &dist);
5051 if (dist > coplanarcentrum)
5052 isconcave= True;
5053 else if (dist >= -qh centrum_radius)
5054 iscoplanar= True;
5055 else if (dist < -clearlyconvex)
5056 maybeconvex= True;
5057 if (!neighbor->center)
5058 neighbor->center= qh_getcentrum(neighbor);
5059 zzinc_(Zcentrumtests);
5060 qh_distplane(neighbor->center, facet, &dist2);
5061 if (dist2 > coplanarcentrum)
5062 isconcave= True;
5063 else if (dist2 >= -qh centrum_radius)
5064 iscoplanar= True;
5065 else if (dist2 < -clearlyconvex) {
5066 if (maybeconvex)
5067 return False; /* both centrums clearly convex */
5068 maybeconvex= True;
5069 }
5070 if (isconcave) {
5071 if (!okangle && qh ANGLEmerge) {
5072 angle= qh_getangle(facet->normal, neighbor->normal);
5073 zinc_(Zangletests);
5074 }
5075 mergedist= fmax_(dist, dist2);
5076 zinc_(Zconcaveridge);
5077 qh_appendmergeset(facet, neighbor, MRGconcave, mergedist, angle);
5078 trace0((qh ferr, 18, "qh_test_nonsimplicial_merge: concave centrum for f%d or f%d, dist %4.4g and reverse dist %4.4g, angle %4.4g during p%d\n",
5079 facet->id, neighbor->id, dist, dist2, angle, qh furthest_id));
5080 return True;
5081 }
5082 /* neither clearly convex nor clearly concave, test vertices as well as centrums */
5083 if (maybeconvex) {
5084 if (dist < -clearlyconvex) {
5085 maxdist= dist; /* facet centrum clearly convex, no need to test its vertex distance */
5086 mindist= dist;
5087 maxvertex2= qh_furthestvertex(neighbor, facet, &maxdist2, &mindist2);
5088 if (!maxvertex2) {
5089 qh_appendmergeset(neighbor, facet, MRGredundant, maxdist2, qh_ANGLEnone);
5090 isredundant= True;
5091 }
5092 }else { /* dist2 < -clearlyconvex */
5093 maxdist2= dist2; /* neighbor centrum clearly convex, no need to test its vertex distance */
5094 mindist2= dist2;
5095 maxvertex= qh_furthestvertex(facet, neighbor, &maxdist, &mindist);
5096 if (!maxvertex) {
5097 qh_appendmergeset(facet, neighbor, MRGredundant, maxdist, qh_ANGLEnone);
5098 isredundant= True;
5099 }
5100 }
5101 }else {
5102 maxvertex= qh_furthestvertex(facet, neighbor, &maxdist, &mindist);
5103 if (maxvertex) {
5104 maxvertex2= qh_furthestvertex(neighbor, facet, &maxdist2, &mindist2);
5105 if (!maxvertex2) {
5106 qh_appendmergeset(neighbor, facet, MRGredundant, maxdist2, qh_ANGLEnone);
5107 isredundant= True;
5108 }else if (mindist < -clearlyconvex || mindist2 < -clearlyconvex)
5109 maybeconvex= True;
5110 }else { /* !maxvertex */
5111 qh_appendmergeset(facet, neighbor, MRGredundant, maxdist, qh_ANGLEnone);
5112 isredundant= True;
5113 }
5114 }
5115 if (isredundant) {
5116 zinc_(Zredundantmerge);
5117 return True;
5118 }
5119
5120 if (maxdist > clearlyconcave || maxdist2 > clearlyconcave)
5121 isconcave= True;
5122 else if (maybeconvex) {
5123 if (maxdist > maxoutside || maxdist2 > maxoutside)
5124 isconcave= True; /* MRGtwisted */
5125 }
5126 if (!isconcave && qh MERGEexact && !qh POSTmerging)
5127 return False;
5128 if (isconcave && !iscoplanar) {
5129 if (maxdist < maxoutside && (-qh MAXcoplanar || (maxdist2 < maxoutside && mindist2 >= -qh MAXcoplanar)))
5130 iscoplanar= True; /* MRGconcavecoplanar */
5131 }else if (!iscoplanar) {
5132 if (mindist >= -qh MAXcoplanar || mindist2 >= -qh MAXcoplanar)
5133 iscoplanar= True; /* MRGcoplanar */
5134 }
5135 if (!isconcave && !iscoplanar)
5136 return False;
5137 if (!okangle && qh ANGLEmerge) {
5138 angle= qh_getangle(facet->normal, neighbor->normal);
5139 zinc_(Zangletests);
5140 }
5141 if (isconcave && maybeconvex) {
5142 zinc_(Ztwistedridge);
5143 if (maxdist > maxdist2)
5144 qh_appendmergeset(facet, neighbor, MRGtwisted, maxdist, angle);
5145 else
5146 qh_appendmergeset(neighbor, facet, MRGtwisted, maxdist2, angle);
5147 trace0((qh ferr, 27, "qh_test_nonsimplicial_merge: twisted concave f%d v%d to f%d v%d, dist %4.4g and reverse dist %4.4g, angle %4.4g during p%d\n",
5148 facet->id, getid_(maxvertex), neighbor->id, getid_(maxvertex2), maxdist, maxdist2, angle, qh furthest_id));
5149 }else if (isconcave && iscoplanar) {
5150 zinc_(Zconcavecoplanarridge);
5151 if (maxdist > maxdist2)
5152 qh_appendmergeset(facet, neighbor, MRGconcavecoplanar, maxdist, angle);
5153 else
5154 qh_appendmergeset(neighbor, facet, MRGconcavecoplanar, maxdist2, angle);
5155 trace0((qh ferr, 28, "qh_test_nonsimplicial_merge: concave coplanar f%d v%d to f%d v%d, dist %4.4g and reverse dist %4.4g, angle %4.4g during p%d\n",
5156 facet->id, getid_(maxvertex), neighbor->id, getid_(maxvertex2), maxdist, maxdist2, angle, qh furthest_id));
5157 }else if (isconcave) {
5158 mergedist= fmax_(maxdist, maxdist2);
5159 zinc_(Zconcaveridge);
5160 qh_appendmergeset(facet, neighbor, MRGconcave, mergedist, angle);
5161 trace0((qh ferr, 29, "qh_test_nonsimplicial_merge: concave f%d v%d to f%d v%d, dist %4.4g and reverse dist %4.4g, angle %4.4g during p%d\n",
5162 facet->id, getid_(maxvertex), neighbor->id, getid_(maxvertex2), maxdist, maxdist2, angle, qh furthest_id));
5163 }else /* iscoplanar */ {
5164 mergedist= fmax_(fmax_(maxdist, maxdist2), fmax_(-mindist, -mindist2));
5165 zinc_(Zcoplanarcentrum);
5166 qh_appendmergeset(facet, neighbor, MRGcoplanar, mergedist, angle);
5167 trace2((qh ferr, 2099, "qh_test_nonsimplicial_merge: coplanar f%d v%d to f%d v%d, dist %4.4g and reverse dist %4.4g, angle %4.4g during p%d\n",
5168 facet->id, getid_(maxvertex), neighbor->id, getid_(maxvertex2), maxdist, maxdist2, angle, qh furthest_id));
5169 }
5170 return True;
5171 } /* test_nonsimplicial_merge */
5172
5173 /*-<a href="qh-merge.htm#TOC"
5174 >-------------------------------</a><a name="test_redundant_neighbors">-</a>
5175
5176 qh_test_redundant_neighbors( facet )
5177 append degenerate facet or its redundant neighbors to qh.degen_mergeset
5178
5179 returns:
5180 bumps vertex_visit
5181
5182 notes:
5183 called at end of qh_mergefacet(), qh_mergecycle_all(), and qh_renamevertex
5184 call before qh_test_degen_neighbors (MRGdegen are more likely to cause problems)
5185 a redundant neighbor's vertices is a subset of the facet's vertices
5186 with pinched and flipped facets, a redundant neighbor may have a wildly different normal
5187
5188 see qh_merge_degenredundant() and qh_-_facet()
5189
5190 design:
5191 if facet is degenerate
5192 appends facet to degen_mergeset
5193 else
5194 appends redundant neighbors of facet to degen_mergeset
5195 */
qh_test_redundant_neighbors(facetT * facet)5196 void qh_test_redundant_neighbors(facetT *facet) {
5197 vertexT *vertex, **vertexp;
5198 facetT *neighbor, **neighborp;
5199 int size;
5200
5201 trace4((qh ferr, 4022, "qh_test_redundant_neighbors: test neighbors of f%d vertex_visit %d\n",
5202 facet->id, qh vertex_visit+1));
5203 if ((size= qh_setsize(facet->neighbors)) < qh hull_dim) {
5204 qh_appendmergeset(facet, facet, MRGdegen, 0.0, 1.0);
5205 trace2((qh ferr, 2017, "qh_test_redundant_neighbors: f%d is degenerate with %d neighbors.\n", facet->id, size));
5206 }else {
5207 qh vertex_visit++;
5208 FOREACHvertex_(facet->vertices)
5209 vertex->visitid= qh vertex_visit;
5210 FOREACHneighbor_(facet) {
5211 if (neighbor->visible) {
5212 qh_fprintf(qh ferr, 6360, "qhull internal error (qh_test_redundant_neighbors): facet f%d has deleted neighbor f%d (qh.visible_list)\n",
5213 facet->id, neighbor->id);
5214 qh_errexit2(qh_ERRqhull, facet, neighbor);
5215 }
5216 if (neighbor->degenerate || neighbor->redundant || neighbor->dupridge) /* will merge or delete */
5217 continue;
5218 if (facet->flipped && !neighbor->flipped) /* do not merge non-flipped into flipped */
5219 continue;
5220 /* merge redundant-flipped facet first */
5221 /* uses early out instead of checking vertex count */
5222 FOREACHvertex_(neighbor->vertices) {
5223 if (vertex->visitid != qh vertex_visit)
5224 break;
5225 }
5226 if (!vertex) {
5227 qh_appendmergeset(neighbor, facet, MRGredundant, 0.0, 1.0);
5228 trace2((qh ferr, 2018, "qh_test_redundant_neighbors: f%d is contained in f%d. merge\n", neighbor->id, facet->id));
5229 }
5230 }
5231 }
5232 } /* test_redundant_neighbors */
5233
5234 /*-<a href="qh-merge.htm#TOC"
5235 >-------------------------------</a><a name="test_vneighbors">-</a>
5236
5237 qh_test_vneighbors( )
5238 test vertex neighbors for convexity
5239 tests all facets on qh.newfacet_list
5240
5241 returns:
5242 true if non-convex vneighbors appended to qh.facet_mergeset
5243 initializes vertex neighbors if needed
5244
5245 notes:
5246 called by qh_all_merges from qh_postmerge if qh.TESTvneighbors ('Qv')
5247 assumes all facet neighbors have been tested
5248 this can be expensive
5249 this does not guarantee that a centrum is below all facets
5250 but it is unlikely
5251 uses qh.visit_id
5252
5253 design:
5254 build vertex neighbors if necessary
5255 for all new facets
5256 for all vertices
5257 for each unvisited facet neighbor of the vertex
5258 test new facet and neighbor for convexity
5259 */
qh_test_vneighbors(void)5260 boolT qh_test_vneighbors(void /* qh.newfacet_list */) {
5261 facetT *newfacet, *neighbor, **neighborp;
5262 vertexT *vertex, **vertexp;
5263 int nummerges= 0;
5264
5265 trace1((qh ferr, 1015, "qh_test_vneighbors: testing vertex neighbors for convexity\n"));
5266 if (!qh VERTEXneighbors)
5267 qh_vertexneighbors();
5268 FORALLnew_facets
5269 newfacet->seen= False;
5270 FORALLnew_facets {
5271 newfacet->seen= True;
5272 newfacet->visitid= qh visit_id++;
5273 FOREACHneighbor_(newfacet)
5274 newfacet->visitid= qh visit_id;
5275 FOREACHvertex_(newfacet->vertices) {
5276 FOREACHneighbor_(vertex) {
5277 if (neighbor->seen || neighbor->visitid == qh visit_id)
5278 continue;
5279 if (qh_test_appendmerge(newfacet, neighbor, False)) /* ignores optimization for simplicial ridges */
5280 nummerges++;
5281 }
5282 }
5283 }
5284 zadd_(Ztestvneighbor, nummerges);
5285 trace1((qh ferr, 1016, "qh_test_vneighbors: found %d non-convex, vertex neighbors\n",
5286 nummerges));
5287 return (nummerges > 0);
5288 } /* test_vneighbors */
5289
5290 /*-<a href="qh-merge.htm#TOC"
5291 >-------------------------------</a><a name="tracemerge">-</a>
5292
5293 qh_tracemerge( facet1, facet2 )
5294 print trace message after merge
5295 */
qh_tracemerge(facetT * facet1,facetT * facet2,mergeType mergetype)5296 void qh_tracemerge(facetT *facet1, facetT *facet2, mergeType mergetype) {
5297 boolT waserror= False;
5298 const char *mergename;
5299
5300 #ifndef qh_NOtrace
5301 if(mergetype > 0 && mergetype <= sizeof(mergetypes))
5302 mergename= mergetypes[mergetype];
5303 else
5304 mergename= mergetypes[MRGnone];
5305 if (qh IStracing >= 4)
5306 qh_errprint("MERGED", facet2, NULL, NULL, NULL);
5307 if (facet2 == qh tracefacet || (qh tracevertex && qh tracevertex->newfacet)) {
5308 qh_fprintf(qh ferr, 8085, "qh_tracemerge: trace facet and vertex after merge of f%d into f%d type %d (%s), furthest p%d\n",
5309 facet1->id, facet2->id, mergetype, mergename, qh furthest_id);
5310 if (facet2 != qh tracefacet)
5311 qh_errprint("TRACE", qh tracefacet,
5312 (qh tracevertex && qh tracevertex->neighbors) ?
5313 SETfirstt_(qh tracevertex->neighbors, facetT) : NULL,
5314 NULL, qh tracevertex);
5315 }
5316 if (qh tracevertex) {
5317 if (qh tracevertex->deleted)
5318 qh_fprintf(qh ferr, 8086, "qh_tracemerge: trace vertex deleted at furthest p%d\n",
5319 qh furthest_id);
5320 else
5321 qh_checkvertex(qh tracevertex, qh_ALL, &waserror);
5322 }
5323 if (qh tracefacet && qh tracefacet->normal && !qh tracefacet->visible)
5324 qh_checkfacet(qh tracefacet, True /* newmerge */, &waserror);
5325 #endif /* !qh_NOtrace */
5326 if (qh CHECKfrequently || qh IStracing >= 4) { /* can't check polygon here */
5327 if (qh IStracing >= 4 && qh num_facets < 500) {
5328 qh_printlists();
5329 }
5330 qh_checkfacet(facet2, True /* newmerge */, &waserror);
5331 }
5332 if (waserror)
5333 qh_errexit(qh_ERRqhull, NULL, NULL); /* erroneous facet logged by qh_checkfacet */
5334 } /* tracemerge */
5335
5336 /*-<a href="qh-merge.htm#TOC"
5337 >-------------------------------</a><a name="tracemerging">-</a>
5338
5339 qh_tracemerging( )
5340 print trace message during POSTmerging
5341
5342 returns:
5343 updates qh.mergereport
5344
5345 notes:
5346 called from qh_mergecycle() and qh_mergefacet()
5347
5348 see:
5349 qh_buildtracing()
5350 */
qh_tracemerging(void)5351 void qh_tracemerging(void) {
5352 realT cpu;
5353 int total;
5354 time_t timedata;
5355 struct tm *tp;
5356
5357 qh mergereport= zzval_(Ztotmerge);
5358 time(&timedata);
5359 tp= localtime(&timedata);
5360 cpu= qh_CPUclock;
5361 cpu /= qh_SECticks;
5362 total= zzval_(Ztotmerge) - zzval_(Zcyclehorizon) + zzval_(Zcyclefacettot);
5363 qh_fprintf(qh ferr, 8087, "\n\
5364 At %d:%d:%d & %2.5g CPU secs, qhull has merged %d facets with max_outside %2.2g, min_vertex %2.2g.\n\
5365 The hull contains %d facets and %d vertices.\n",
5366 tp->tm_hour, tp->tm_min, tp->tm_sec, cpu, total, qh max_outside, qh min_vertex,
5367 qh num_facets - qh num_visible,
5368 qh num_vertices-qh_setsize(qh del_vertices));
5369 } /* tracemerging */
5370
5371 /*-<a href="qh-merge.htm#TOC"
5372 >-------------------------------</a><a name="updatetested">-</a>
5373
5374 qh_updatetested( facet1, facet2 )
5375 clear facet2->tested and facet1->ridge->tested for merge
5376
5377 returns:
5378 deletes facet2->center unless it's already large
5379 if so, clears facet2->ridge->tested
5380
5381 notes:
5382 only called by qh_mergefacet
5383
5384 design:
5385 clear facet2->tested
5386 clear ridge->tested for facet1's ridges
5387 if facet2 has a centrum
5388 if facet2 is large
5389 set facet2->keepcentrum
5390 else if facet2 has 3 vertices due to many merges, or not large and post merging
5391 clear facet2->keepcentrum
5392 unless facet2->keepcentrum
5393 clear facet2->center to recompute centrum later
5394 clear ridge->tested for facet2's ridges
5395 */
qh_updatetested(facetT * facet1,facetT * facet2)5396 void qh_updatetested(facetT *facet1, facetT *facet2) {
5397 ridgeT *ridge, **ridgep;
5398 int size;
5399
5400 facet2->tested= False;
5401 FOREACHridge_(facet1->ridges)
5402 ridge->tested= False;
5403 if (!facet2->center)
5404 return;
5405 size= qh_setsize(facet2->vertices);
5406 if (!facet2->keepcentrum) {
5407 if (size > qh hull_dim + qh_MAXnewcentrum) {
5408 facet2->keepcentrum= True;
5409 zinc_(Zwidevertices);
5410 }
5411 }else if (size <= qh hull_dim + qh_MAXnewcentrum) {
5412 /* center and keepcentrum was set */
5413 if (size == qh hull_dim || qh POSTmerging)
5414 facet2->keepcentrum= False; /* if many merges need to recompute centrum */
5415 }
5416 if (!facet2->keepcentrum) {
5417 qh_memfree(facet2->center, qh normal_size);
5418 facet2->center= NULL;
5419 FOREACHridge_(facet2->ridges)
5420 ridge->tested= False;
5421 }
5422 } /* updatetested */
5423
5424 /*-<a href="qh-merge.htm#TOC"
5425 >-------------------------------</a><a name="vertexridges">-</a>
5426
5427 qh_vertexridges( vertex, allneighbors )
5428 return temporary set of ridges adjacent to a vertex
5429 vertex->neighbors defined (qh_vertexneighbors)
5430
5431 notes:
5432 uses qh.visit_id
5433 does not include implicit ridges for simplicial facets
5434 skips last neighbor, unless allneighbors. For new facets, the last neighbor shares ridges with adjacent neighbors
5435 if the last neighbor is not simplicial, it will have ridges for its simplicial neighbors
5436 Use allneighbors when a new cone is attached to an existing convex hull
5437 similar to qh_neighbor_vertices
5438
5439 design:
5440 for each neighbor of vertex
5441 add ridges that include the vertex to ridges
5442 */
qh_vertexridges(vertexT * vertex,boolT allneighbors)5443 setT *qh_vertexridges(vertexT *vertex, boolT allneighbors) {
5444 facetT *neighbor, **neighborp;
5445 setT *ridges= qh_settemp(qh TEMPsize);
5446 int size;
5447
5448 qh visit_id += 2; /* visit_id for vertex neighbors, visit_id-1 for facets of visited ridges */
5449 FOREACHneighbor_(vertex)
5450 neighbor->visitid= qh visit_id;
5451 FOREACHneighbor_(vertex) {
5452 if (*neighborp || allneighbors) /* no new ridges in last neighbor */
5453 qh_vertexridges_facet(vertex, neighbor, &ridges);
5454 }
5455 if (qh PRINTstatistics || qh IStracing) {
5456 size= qh_setsize(ridges);
5457 zinc_(Zvertexridge);
5458 zadd_(Zvertexridgetot, size);
5459 zmax_(Zvertexridgemax, size);
5460 trace3((qh ferr, 3011, "qh_vertexridges: found %d ridges for v%d\n",
5461 size, vertex->id));
5462 }
5463 return ridges;
5464 } /* vertexridges */
5465
5466 /*-<a href="qh-merge.htm#TOC"
5467 >-------------------------------</a><a name="vertexridges_facet">-</a>
5468
5469 qh_vertexridges_facet( vertex, facet, ridges )
5470 add adjacent ridges for vertex in facet
5471 neighbor->visitid==qh.visit_id if it hasn't been visited
5472
5473 returns:
5474 ridges updated
5475 sets facet->visitid to qh.visit_id-1
5476
5477 design:
5478 for each ridge of facet
5479 if ridge of visited neighbor (i.e., unprocessed)
5480 if vertex in ridge
5481 append ridge
5482 mark facet processed
5483 */
qh_vertexridges_facet(vertexT * vertex,facetT * facet,setT ** ridges)5484 void qh_vertexridges_facet(vertexT *vertex, facetT *facet, setT **ridges) {
5485 ridgeT *ridge, **ridgep;
5486 facetT *neighbor;
5487 int last_i= qh hull_dim-2;
5488 vertexT *second, *last;
5489
5490 FOREACHridge_(facet->ridges) {
5491 neighbor= otherfacet_(ridge, facet);
5492 if (neighbor->visitid == qh visit_id) {
5493 if (SETfirst_(ridge->vertices) == vertex) {
5494 qh_setappend(ridges, ridge);
5495 }else if (last_i > 2) {
5496 second= SETsecondt_(ridge->vertices, vertexT);
5497 last= SETelemt_(ridge->vertices, last_i, vertexT);
5498 if (second->id >= vertex->id && last->id <= vertex->id) { /* vertices inverse sorted by id */
5499 if (second == vertex || last == vertex)
5500 qh_setappend(ridges, ridge);
5501 else if (qh_setin(ridge->vertices, vertex))
5502 qh_setappend(ridges, ridge);
5503 }
5504 }else if (SETelem_(ridge->vertices, last_i) == vertex
5505 || (last_i > 1 && SETsecond_(ridge->vertices) == vertex)) {
5506 qh_setappend(ridges, ridge);
5507 }
5508 }
5509 }
5510 facet->visitid= qh visit_id-1;
5511 } /* vertexridges_facet */
5512
5513 /*-<a href="qh-merge.htm#TOC"
5514 >-------------------------------</a><a name="willdelete">-</a>
5515
5516 qh_willdelete( facet, replace )
5517 moves facet to visible list for qh_deletevisible
5518 sets facet->f.replace to replace (may be NULL)
5519 clears f.ridges and f.neighbors -- no longer valid
5520
5521 returns:
5522 bumps qh.num_visible
5523 */
qh_willdelete(facetT * facet,facetT * replace)5524 void qh_willdelete(facetT *facet, facetT *replace) {
5525
5526 trace4((qh ferr, 4081, "qh_willdelete: move f%d to visible list, set its replacement as f%d, and clear f.neighbors and f.ridges\n", facet->id, getid_(replace)));
5527 if (!qh visible_list && qh newfacet_list) {
5528 qh_fprintf(qh ferr, 6378, "qhull internal error (qh_willdelete): expecting qh.visible_list at before qh.newfacet_list f%d. Got NULL\n",
5529 qh newfacet_list->id);
5530 qh_errexit2(qh_ERRqhull, NULL, NULL);
5531 }
5532 qh_removefacet(facet);
5533 qh_prependfacet(facet, &qh visible_list);
5534 qh num_visible++;
5535 facet->visible= True;
5536 facet->f.replace= replace;
5537 if (facet->ridges)
5538 SETfirst_(facet->ridges)= NULL;
5539 if (facet->neighbors)
5540 SETfirst_(facet->neighbors)= NULL;
5541 } /* willdelete */
5542
5543 #else /* qh_NOmerge */
5544
qh_all_vertexmerges(int apexpointid,facetT * facet,facetT ** retryfacet)5545 void qh_all_vertexmerges(int apexpointid, facetT *facet, facetT **retryfacet) {
5546 QHULL_UNUSED(apexpointid)
5547 QHULL_UNUSED(facet)
5548 QHULL_UNUSED(retryfacet)
5549 }
qh_premerge(int apexpointid,realT maxcentrum,realT maxangle)5550 void qh_premerge(int apexpointid, realT maxcentrum, realT maxangle) {
5551 QHULL_UNUSED(apexpointid)
5552 QHULL_UNUSED(maxcentrum)
5553 QHULL_UNUSED(maxangle)
5554 }
qh_postmerge(const char * reason,realT maxcentrum,realT maxangle,boolT vneighbors)5555 void qh_postmerge(const char *reason, realT maxcentrum, realT maxangle,
5556 boolT vneighbors) {
5557 QHULL_UNUSED(reason)
5558 QHULL_UNUSED(maxcentrum)
5559 QHULL_UNUSED(maxangle)
5560 QHULL_UNUSED(vneighbors)
5561 }
qh_checkdelfacet(facetT * facet,setT * mergeset)5562 void qh_checkdelfacet(facetT *facet, setT *mergeset) {
5563 QHULL_UNUSED(facet)
5564 QHULL_UNUSED(mergeset)
5565 }
qh_checkdelridge(void)5566 void qh_checkdelridge(void /* qh.visible_facets, vertex_mergeset */) {
5567 }
qh_checkzero(boolT testall)5568 boolT qh_checkzero(boolT testall) {
5569 QHULL_UNUSED(testall)
5570
5571 return True;
5572 }
qh_freemergesets(void)5573 void qh_freemergesets(void) {
5574 }
qh_initmergesets(void)5575 void qh_initmergesets(void) {
5576 }
qh_merge_pinchedvertices(int apexpointid)5577 void qh_merge_pinchedvertices(int apexpointid /* qh.newfacet_list */) {
5578 QHULL_UNUSED(apexpointid)
5579 }
5580 #endif /* qh_NOmerge */
5581
5582