1 /*<html><pre>  -<a                             href="qh-merge_r.htm#TOC"
2   >-------------------------------</a><a name="TOP">-</a>
3 
4    merge_r.c
5    merges non-convex facets
6 
7    see qh-merge_r.htm and merge_r.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_r.c):
14      qh_partitionvisible(qh, !qh_ALL, &numoutside);  // visible_list, newfacet_list
15      qh_deletevisible();         // qh.visible_list
16      qh_resetlists(qh, 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_r/merge_r.c#12 $$Change: 2712 $
25    $DateTime: 2019/06/28 12:57:00 $$Author: bbarber $
26 */
27 
28 #include "qhull_ra.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_r.htm#TOC"
53   >-------------------------------</a><a name="premerge">-</a>
54 
55   qh_premerge(qh, 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(qhT * qh,int apexpointid,realT maxcentrum,realT maxangle)75 void qh_premerge(qhT *qh, int apexpointid, realT maxcentrum, realT maxangle /* qh.newfacet_list */) {
76   boolT othermerge= False;
77 
78   if (qh->ZEROcentrum && qh_checkzero(qh, !qh_ALL))
79     return;
80   trace2((qh, 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(qh);
84   qh->centrum_radius= maxcentrum;
85   qh->cos_max= maxangle;
86   if (qh->hull_dim >=3) {
87     qh_mark_dupridges(qh, qh->newfacet_list, qh_ALL); /* facet_mergeset */
88     qh_mergecycle_all(qh, qh->newfacet_list, &othermerge);
89     qh_forcedmerges(qh, &othermerge /* qh.facet_mergeset */);
90   }else /* qh.hull_dim == 2 */
91     qh_mergecycle_all(qh, qh->newfacet_list, &othermerge);
92   qh_flippedmerges(qh, qh->newfacet_list, &othermerge);
93   if (!qh->MERGEexact || zzval_(Ztotmerge)) {
94     zinc_(Zpremergetot);
95     qh->POSTmerging= False;
96     qh_getmergeset_initial(qh, qh->newfacet_list);
97     qh_all_merges(qh, othermerge, False);
98   }
99 } /* premerge */
100 
101 /*-<a                             href="qh-merge_r.htm#TOC"
102   >-------------------------------</a><a name="postmerge">-</a>
103 
104   qh_postmerge(qh, 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(qhT * qh,const char * reason,realT maxcentrum,realT maxangle,boolT vneighbors)139 void qh_postmerge(qhT *qh, 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(qh, NULL, NULL);
147     qh_printsummary(qh, qh->ferr);
148     if (qh->PRINTstatistics)
149       qh_printallstatistics(qh, qh->ferr, "reason");
150     qh_fprintf(qh, qh->ferr, 8062, "\n%s with 'C%.2g' and 'A%.2g'\n",
151         reason, maxcentrum, maxangle);
152   }
153   trace2((qh, 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);  /* qh_all_merges did not call qh_reducevertices for v.delridge */
173     }
174     if (!qh->PREmerge && !qh->MERGEexact)
175       qh_flippedmerges(qh, qh->newfacet_list, &othermerges);
176   }
177   qh_getmergeset_initial(qh, qh->newfacet_list);
178   qh_all_merges(qh, 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_r.htm#TOC"
184   >-------------------------------</a><a name="all_merges">-</a>
185 
186   qh_all_merges(qh, 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(qhT * qh,boolT othermerge,boolT vneighbors)223 void qh_all_merges(qhT *qh, 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, 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, qh->facet_mergeset), qh_setsize(qh, qh->degen_mergeset), getid_(qh->newfacet_list), othermerge));
235 
236   while (True) {
237     wasmerge= False;
238     while (qh_setsize(qh, qh->facet_mergeset) > 0 || qh_setsize(qh, qh->degen_mergeset) > 0) {
239       if (qh_setsize(qh, qh->degen_mergeset) > 0) {
240         numdegenredun += qh_merge_degenredundant(qh);
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_(qh, merge, (int)sizeof(mergeT), freelistp);   /* 'merge' is invalid */
251         if (facet1->visible || facet2->visible) {
252           trace3((qh, 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, 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, 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(qh, facet1, facet2);
269         else
270           qh_merge_nonconvex(qh, facet1, facet2, mergetype);
271         numnewmerges++;
272         numdegenredun += qh_merge_degenredundant(qh);
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, 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, 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(qh);  /* otherwise large post merges too slow */
293       }
294       qh_getmergeset(qh, 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(qh)) {
307           qh_getmergeset(qh, qh->newfacet_list); /* facet_mergeset */
308           continue;
309         }
310       }
311     }
312     if (vneighbors && qh_test_vneighbors(qh /* qh.newfacet_list */))
313       continue;
314     break;
315   } /* while (True) */
316   if (wasmerge || othermerge) {
317     trace3((qh, 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, qh->newfacet_list, qh_ALGORITHMfault);
326     /* qh_checkconnect(qh); [this is slow and it changes the facet order] */
327     qh->RANDOMdist= qh->old_randomdist;
328   }
329   trace1((qh, 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(qh);
333 } /* all_merges */
334 
335 /*-<a                             href="qh-merge_r.htm#TOC"
336   >-------------------------------</a><a name="all_vertexmerges">-</a>
337 
338   qh_all_vertexmerges(qh, 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(qhT * qh,int apexpointid,facetT * facet,facetT ** retryfacet)358 void qh_all_vertexmerges(qhT *qh, 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, qh->vertex_mergeset) > 0) {
364     trace1((qh, qh->ferr, 1057, "qh_all_vertexmerges: starting to merge %d vertex merges for apex p%d facet f%d\n",
365             qh_setsize(qh, qh->vertex_mergeset), apexpointid, getid_(facet)));
366     if (qh->IStracing >= 4  && qh->num_facets < 1000)
367       qh_printlists(qh);
368     qh_merge_pinchedvertices(qh, apexpointid /* qh.vertex_mergeset, visible_list, newvertex_list, newfacet_list */);
369     qh_update_vertexneighbors(qh); /* 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, qh->newfacet_list);
372     qh_all_merges(qh, True, False); /* calls qh_reducevertices */
373     if (qh->CHECKfrequently)
374       qh_checkpolygon(qh, qh->facet_list);
375     qh_partitionvisible(qh, !qh_ALL, &numpoints /* qh.visible_list qh.del_vertices*/);
376     if (retryfacet)
377       *retryfacet= qh_getreplacement(qh, *retryfacet);
378     qh_deletevisible(qh /* qh.visible_list  qh.del_vertices*/);
379     qh_resetlists(qh, False, qh_RESETvisible /* qh.visible_list newvertex_list qh.newfacet_list */);
380     if (qh->IStracing >= 4  && qh->num_facets < 1000) {
381       qh_printlists(qh);
382       qh_checkpolygon(qh, qh->facet_list);
383     }
384   }
385 } /* all_vertexmerges */
386 
387 /*-<a                             href="qh-merge_r.htm#TOC"
388   >-------------------------------</a><a name="appendmergeset">-</a>
389 
390   qh_appendmergeset(qh, 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(qhT * qh,facetT * facet,facetT * neighbor,mergeType mergetype,coordT dist,realT angle)414 void qh_appendmergeset(qhT *qh, 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, 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, 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, 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, qh_ERRqhull, NULL, NULL);
434   }
435   if (neighbor->flipped && !facet->flipped) {
436     if (mergetype != MRGdupridge) {
437       qh_fprintf(qh, 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, qh_ERRqhull, NULL, NULL);
440     }else {
441       trace2((qh, 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_(qh, (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, &(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, &(qh->degen_mergeset), merge);
466     else
467       qh_setaddnth(qh, &(qh->degen_mergeset), 0, merge);    /* merged last */
468   }else if (mergetype == MRGredundant) {
469     facet->redundant= True;
470     qh_setappend(qh, &(qh->degen_mergeset), merge);
471   }else /* mergetype == MRGmirror */ {
472     if (facet->redundant || neighbor->redundant) {
473       qh_fprintf(qh, 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, qh_ERRqhull, facet, neighbor);
476     }
477     if (!qh_setequal(facet->vertices, neighbor->vertices)) {
478       qh_fprintf(qh, 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, qh_ERRqhull, facet, neighbor);
481     }
482     facet->redundant= True;
483     neighbor->redundant= True;
484     qh_setappend(qh, &(qh->degen_mergeset), merge);
485   }
486   if (merge->mergetype >= MRGdegen) {
487     trace3((qh, 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, qh->degen_mergeset)));
489   }else {
490     trace3((qh, 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, qh->facet_mergeset)));
492   }
493 } /* appendmergeset */
494 
495 
496 /*-<a                             href="qh-merge_r.htm#TOC"
497   >-------------------------------</a><a name="appendvertexmerge">-</a>
498 
499   qh_appendvertexmerge(qh, 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(qhT * qh,vertexT * vertex,vertexT * destination,mergeType mergetype,realT distance,ridgeT * ridge1,ridgeT * ridge2)510 void qh_appendvertexmerge(qhT *qh, 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, 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, qh_ERRqhull, NULL, NULL);
520   }
521   qh_memalloc_(qh, (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, 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, qh_ERRqhull, NULL, ridge1);
540     }
541   }
542   qh_setappend(qh, &(qh->vertex_mergeset), merge);
543   trace3((qh, 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_r.htm#TOC"
549   >-------------------------------</a><a name="basevertices">-</a>
550 
551   qh_basevertices(qh, 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(qhT * qh,facetT * samecycle)568 setT *qh_basevertices(qhT *qh, facetT *samecycle) {
569   facetT *same;
570   vertexT *apex, *vertex, **vertexp;
571   setT *vertices= qh_settemp(qh, 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(qh, &vertices, vertex);
581         vertex->visitid= qh->vertex_visit;
582         vertex->seen= False;
583       }
584     }
585   }
586   trace4((qh, qh->ferr, 4019, "qh_basevertices: found %d vertices\n",
587          qh_setsize(qh, vertices)));
588   return vertices;
589 } /* basevertices */
590 
591 /*-<a                             href="qh-merge_r.htm#TOC"
592   >-------------------------------</a><a name="check_dupridge">-</a>
593 
594   qh_check_dupridge(qh, 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(qhT * qh,facetT * facet1,realT dist1,facetT * facet2,realT dist2)606 void qh_check_dupridge(qhT *qh, 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(qh, 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, 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, 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, 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, qh->ferr, 8145, "- A bounding box for the input sites may alleviate this error.\n");
636     if (!qh->ALLOWwide)
637       qh_errexit2(qh, qh_ERRwide, facet1, facet2);
638   }
639 } /* check_dupridge */
640 
641 /*-<a                             href="qh-merge_r.htm#TOC"
642   >-------------------------------</a><a name="checkconnect">-</a>
643 
644   qh_checkconnect(qh)
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(qhT * qh)659 void qh_checkconnect(qhT *qh /* qh.newfacet_list */) {
660   facetT *facet, *newfacet, *errfacet= NULL, *neighbor, **neighborp;
661 
662   facet= qh->newfacet_list;
663   qh_removefacet(qh, facet);
664   qh_appendfacet(qh, facet);
665   facet->visitid= ++qh->visit_id;
666   FORALLfacet_(facet) {
667     FOREACHneighbor_(facet) {
668       if (neighbor->visitid != qh->visit_id) {
669         qh_removefacet(qh, neighbor);
670         qh_appendfacet(qh, 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, 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, qh_ERRqhull, errfacet, NULL);
684 } /* checkconnect */
685 
686 /*-<a                             href="qh-merge_r.htm#TOC"
687   >-------------------------------</a><a name="checkdelfacet">-</a>
688 
689   qh_checkdelfacet(qh, facet, mergeset )
690     check that mergeset does not reference facet
691 
692 */
qh_checkdelfacet(qhT * qh,facetT * facet,setT * mergeset)693 void qh_checkdelfacet(qhT *qh, facetT *facet, setT *mergeset) {
694   mergeT *merge, **mergep;
695 
696   FOREACHmerge_(mergeset) {
697     if (merge->facet1 == facet || merge->facet2 == facet) {
698       qh_fprintf(qh, 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, qh_ERRqhull, merge->facet1, merge->facet2);
701     }
702   }
703 } /* checkdelfacet */
704 
705 /*-<a                             href="qh-merge_r.htm#TOC"
706   >-------------------------------</a><a name="checkdelridge">-</a>
707 
708   qh_checkdelridge(qh)
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(qhT * qh)717 void qh_checkdelridge(qhT *qh /* qh.visible_facets, vertex_mergeset */) {
718   facetT *newfacet, *visible;
719   ridgeT *ridge, **ridgep;
720 
721   if (!SETempty_(qh->vertex_mergeset)) {
722     qh_fprintf(qh, 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, qh->vertex_mergeset));
723     qh_errexit(qh, qh_ERRqhull, NULL, NULL);
724   }
725 
726   FORALLnew_facets {
727     FOREACHridge_(newfacet->ridges) {
728       if (ridge->nonconvex) {
729         qh_fprintf(qh, 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, qh_ERRqhull, newfacet, ridge);
732       }
733     }
734   }
735 
736   FORALLvisible_facets {
737     FOREACHridge_(visible->ridges) {
738       if (ridge->nonconvex) {
739         qh_fprintf(qh, 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, qh_ERRqhull, visible, ridge);
742       }
743     }
744   }
745 } /* checkdelridge */
746 
747 
748 /*-<a                             href="qh-merge_r.htm#TOC"
749   >-------------------------------</a><a name="checkzero">-</a>
750 
751   qh_checkzero(qh, 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(qhT * qh,boolT testall)785 boolT qh_checkzero(qhT *qh, 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, 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_(qh, 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(qh, 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(qh, 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, 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, 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, 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, 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_r.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_r.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_r.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_r.htm#TOC"
924   >-------------------------------</a><a name="copynonconvex">-</a>
925 
926   qh_copynonconvex(qh, 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(qhT * qh,ridgeT * atridge)937 void qh_copynonconvex(qhT *qh, 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, 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_r.htm#TOC"
957   >-------------------------------</a><a name="degen_redundant_facet">-</a>
958 
959   qh_degen_redundant_facet(qh, 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(qhT * qh,facetT * facet)973 void qh_degen_redundant_facet(qhT *qh, facetT *facet) {
974   vertexT *vertex, **vertexp;
975   facetT *neighbor, **neighborp;
976 
977   trace3((qh, 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, 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, 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, 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, qh->ferr, 2015, "qh_degen_redundant_facet: f%d is contained in f%d.  merge\n", facet->id, neighbor->id));
1000       qh_appendmergeset(qh, facet, neighbor, MRGredundant, 0.0, 1.0);
1001       return;
1002     }
1003   }
1004   if (qh_setsize(qh, facet->neighbors) < qh->hull_dim) {
1005     qh_appendmergeset(qh, facet, facet, MRGdegen, 0.0, 1.0);
1006     trace2((qh, 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_r.htm#TOC"
1012   >-------------------------------</a><a name="delridge_merge">-</a>
1013 
1014   qh_delridge_merge(qh, ridge )
1015     delete ridge due to a merge
1016 
1017   notes:
1018     only called by merge_r.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(qhT * qh,ridgeT * ridge)1028 void qh_delridge_merge(qhT *qh, ridgeT *ridge) {
1029   vertexT *vertex, **vertexp;
1030   mergeT *merge;
1031   int merge_i, merge_n;
1032 
1033   trace3((qh, 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(qh, ridge);
1037   FOREACHvertex_(ridge->vertices)
1038     vertex->delridge= True;
1039   FOREACHmerge_i_(qh, qh->vertex_mergeset) {
1040     if (merge->ridge1 == ridge || merge->ridge2 == ridge) {
1041       trace3((qh, 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, 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(qh, ridge);
1054 } /* delridge_merge */
1055 
1056 
1057 /*-<a                             href="qh-merge_r.htm#TOC"
1058   >-------------------------------</a><a name="drop_mergevertex">-</a>
1059 
1060   qh_drop_mergevertex(qh, merge )
1061 
1062   clear mergevertex flags for ridges of a vertex merge
1063 */
qh_drop_mergevertex(qhT * qh,mergeT * merge)1064 void qh_drop_mergevertex(qhT *qh, 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, 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_r.htm#TOC"
1077   >-------------------------------</a><a name="find_newvertex">-</a>
1078 
1079   qh_find_newvertex(qh, 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(qhT * qh,vertexT * oldvertex,setT * vertices,setT * ridges)1109 vertexT *qh_find_newvertex(qhT *qh, 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, qh->ferr, 8063, "qh_find_newvertex: find new vertex for v%d from ",
1120              oldvertex->id);
1121     FOREACHvertex_(vertices)
1122       qh_fprintf(qh, qh->ferr, 8064, "v%d ", vertex->id);
1123     FOREACHridge_(ridges)
1124       qh_fprintf(qh, qh->ferr, 8065, "r%d ", ridge->id);
1125     qh_fprintf(qh, 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(qh, vertices, SETindex_(vertices,vertex));
1145       vertexp--; /* repeat since deleted this vertex */
1146     }
1147   }
1148   maxvisit= (unsigned int)qh_setsize(qh, ridges);
1149   maximize_(qh->vertex_visit, maxvisit);
1150   if (!qh_setsize(qh, vertices)) {
1151     trace4((qh, 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(qh, vertices),
1156                 sizeof(vertexT *), qh_comparevisit);
1157   /* can now use qh->vertex_visit */
1158   if (qh->PRINTstatistics) {
1159     size= qh_setsize(qh, vertices);
1160     zinc_(Zintersect);
1161     zadd_(Zintersecttot, size);
1162     zmax_(Zintersectmax, size);
1163   }
1164   hashsize= qh_newhashtable(qh, qh_setsize(qh, ridges));
1165   FOREACHridge_(ridges)
1166     qh_hashridge(qh, qh->hash_table, hashsize, ridge, oldvertex);
1167   FOREACHvertex_(vertices) {
1168     newridges= qh_vertexridges(qh, vertex, !qh_ALL);
1169     FOREACHridge_(newridges) {
1170       if (qh_hashridge_find(qh, qh->hash_table, hashsize, ridge, vertex, oldvertex, &hash)) {
1171         zinc_(Zvertexridge);
1172         break;
1173       }
1174     }
1175     qh_settempfree(qh, &newridges);
1176     if (!ridge)
1177       break;  /* found a rename */
1178   }
1179   if (vertex) {
1180     /* counted in qh_renamevertex */
1181     trace2((qh, 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(qh, vertices), qh_setsize(qh, ridges)));
1183   }else {
1184     zinc_(Zfindfail);
1185     trace0((qh, 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, &qh->hash_table);
1189   return vertex;
1190 } /* find_newvertex */
1191 
1192 /*-<a                             href="qh-geom2_r.htm#TOC"
1193   >-------------------------------</a><a name="findbest_pinchedvertex">-</a>
1194 
1195   qh_findbest_pinchedvertex(qh, 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(qhT * qh,mergeT * merge,vertexT * apex,vertexT ** nearestp,coordT * distp)1219 vertexT *qh_findbest_pinchedvertex(qhT *qh, 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, 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, qh_ERRqhull, merge->facet1, merge->facet2);
1230   }
1231   subridge= qh_vertexintersect_new(qh, merge->facet1->vertices, merge->facet2->vertices); /* new setT.  No error_exit() */
1232   if (qh_setsize(qh, subridge) == qh->hull_dim) { /* duplicate vertices */
1233     bestdist= qh_vertex_bestdist2(qh, 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(qh, subridge) != qh->hull_dim - 2) {
1241       qh_fprintf(qh, 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(qh, subridge));
1243       qh_errexit2(qh, 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(qh, 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(qh, &maybepinched);
1279       }
1280     }
1281   }
1282   *distp= bestdist;
1283   qh_setfree(qh, &subridge); /* qh_err_exit not called since allocated */
1284   if (!bestvertex) {  /* should never happen if qh.hull_dim > 2 */
1285     qh_fprintf(qh, 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, qh_ERRqhull, merge->facet1, merge->facet2);
1287   }
1288   *nearestp= bestvertex;
1289   trace2((qh, 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(qh, bestpinched->point), bestpinched->id, qh_pointid(qh, bestvertex->point), bestvertex->id, bestdist, merge->facet1->id, merge->facet2->id));
1291   return bestpinched;
1292 } /* findbest_pinchedvertex */
1293 
1294 /*-<a                             href="qh-geom2_r.htm#TOC"
1295   >-------------------------------</a><a name="findbest_ridgevertex">-</a>
1296 
1297   qh_findbest_ridgevertex(qh, 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(qhT * qh,ridgeT * ridge,vertexT ** pinchedp,coordT * distp)1308 vertexT *qh_findbest_ridgevertex(qhT *qh, ridgeT *ridge, vertexT **pinchedp, coordT *distp) {
1309   vertexT *bestvertex;
1310 
1311   *distp= qh_vertex_bestdist2(qh, ridge->vertices, &bestvertex, pinchedp);
1312   trace4((qh, 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(qh, (*pinchedp)->point), (*pinchedp)->id, qh_pointid(qh, 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_r.htm#TOC"
1318   >-------------------------------</a><a name="findbest_test">-</a>
1319 
1320   qh_findbest_test(qh, 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(qhT * qh,boolT testcentrum,facetT * facet,facetT * neighbor,facetT ** bestfacet,realT * distp,realT * mindistp,realT * maxdistp)1336 void qh_findbest_test(qhT *qh, 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(qh, 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(qh, 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_r.htm#TOC"
1365   >-------------------------------</a><a name="findbestneighbor">-</a>
1366 
1367   qh_findbestneighbor(qh, 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(qhT * qh,facetT * facet,realT * distp,realT * mindistp,realT * maxdistp)1391 facetT *qh_findbestneighbor(qhT *qh, 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(qh, facet->vertices);
1396 
1397   if(qh->CENTERtype==qh_ASvoronoi){
1398     qh_fprintf(qh, 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, 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(qh, 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(qh, testcentrum, facet, neighbor,
1413                           &bestfacet, distp, mindistp, maxdistp);
1414       }
1415     }
1416   }
1417   if (!bestfacet) {
1418     nonconvex= False;
1419     FOREACHneighbor_(facet)
1420       qh_findbest_test(qh, testcentrum, facet, neighbor,
1421                         &bestfacet, distp, mindistp, maxdistp);
1422   }
1423   if (!bestfacet) {
1424     qh_fprintf(qh, qh->ferr, 6095, "qhull internal error (qh_findbestneighbor): no neighbors for f%d\n", facet->id);
1425     qh_errexit(qh, qh_ERRqhull, facet, NULL);
1426   }
1427   if (testcentrum)
1428     qh_getdistance(qh, facet, bestfacet, mindistp, maxdistp);
1429   trace3((qh, 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_r.htm#TOC"
1436   >-------------------------------</a><a name="flippedmerges">-</a>
1437 
1438   qh_flippedmerges(qh, 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(qhT * qh,facetT * facetlist,boolT * wasmerge)1459 void qh_flippedmerges(qhT *qh, 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, qh->ferr, 4024, "qh_flippedmerges: begin\n"));
1467   FORALLfacet_(facetlist) {
1468     if (facet->flipped && !facet->visible)
1469       qh_appendmergeset(qh, facet, facet, MRGflip, 0.0, 1.0);
1470   }
1471   othermerges= qh_settemppop(qh);
1472   if(othermerges != qh->facet_mergeset) {
1473     qh_fprintf(qh, qh->ferr, 6392, "qhull internal error (qh_flippedmerges): facet_mergeset (%d merges) not at top of tempstack (%d merges)\n",
1474         qh_setsize(qh, qh->facet_mergeset), qh_setsize(qh, othermerges));
1475     qh_errexit(qh, qh_ERRqhull, NULL, NULL);
1476   }
1477   qh->facet_mergeset= qh_settemp(qh, qh->TEMPsize);
1478   qh_settemppush(qh, 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       qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1485     neighbor= qh_findbestneighbor(qh, facet1, &dist, &mindist, &maxdist);
1486     trace0((qh, 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(qh, 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(qh, merge, (int)sizeof(mergeT)); /* invalidates merge and othermerges */
1499     else
1500       qh_setappend(qh, &qh->facet_mergeset, merge);
1501   }
1502   qh_settempfree(qh, &othermerges);
1503   numdegen += qh_merge_degenredundant(qh); /* 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, 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_r.htm#TOC"
1512   >-------------------------------</a><a name="forcedmerges">-</a>
1513 
1514   qh_forcedmerges(qh, 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(qhT * qh,boolT * wasmerge)1542 void qh_forcedmerges(qhT *qh, 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     qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1552   trace3((qh, qh->ferr, 3054, "qh_forcedmerges: merge dupridges\n"));
1553   othermerges= qh_settemppop(qh); /* was facet_mergeset */
1554   if (qh->facet_mergeset != othermerges ) {
1555       qh_fprintf(qh, qh->ferr, 6279, "qhull internal error (qh_forcedmerges): qh_settemppop (size %d) is not qh->facet_mergeset (size %d)\n",
1556           qh_setsize(qh, othermerges), qh_setsize(qh, qh->facet_mergeset));
1557       qh_errexit(qh, qh_ERRqhull, NULL, NULL);
1558   }
1559   qh->facet_mergeset= qh_settemp(qh, qh->TEMPsize);
1560   qh_settemppush(qh, othermerges);
1561   FOREACHmerge_(othermerges) {
1562     if (merge->mergetype != MRGdupridge)
1563         continue;
1564     wasdupridge= True;
1565     if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
1566         qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1567     facet1= qh_getreplacement(qh, merge->facet1);  /* must exist, no qh_merge_degenredunant */
1568     facet2= qh_getreplacement(qh, merge->facet2);  /* previously merged facet, if any */
1569     if (facet1 == facet2)
1570       continue;
1571     if (!qh_setin(facet2->neighbors, facet1)) {
1572       qh_fprintf(qh, 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, qh_ERRqhull, facet1, facet2);
1575     }
1576     dist= qh_getdistance(qh, facet1, facet2, &mindist, &maxdist);
1577     dist2= qh_getdistance(qh, facet2, facet1, &mindist2, &maxdist2);
1578     qh_check_dupridge(qh, 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(qh, merging, merged, merge->mergetype, &mindist, &maxdist, !qh_MERGEapex);
1603     numdegen += qh_merge_degenredundant(qh); /* 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(qh, merge, (int)sizeof(mergeT)); /* invalidates merge and othermerges */
1618     else
1619       qh_setappend(qh, &qh->facet_mergeset, merge);
1620   }
1621   qh_settempfree(qh, &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(qh, newfacet->neighbors) < qh->hull_dim) { /* not tested for MRGdupridge */
1629           qh_appendmergeset(qh, newfacet, newfacet, MRGdegen, 0.0, 1.0);
1630           trace2((qh, 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(qh);
1636   }
1637   if (nummerge || numflip) {
1638     *wasmerge= True;
1639     trace1((qh, 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_r.htm#TOC"
1646   >-------------------------------</a><a name="freemergesets">-</a>
1647 
1648   qh_freemergesets(qh )
1649     free the merge sets
1650 
1651   notes:
1652     matches qh_initmergesets
1653 */
qh_freemergesets(qhT * qh)1654 void qh_freemergesets(qhT *qh) {
1655 
1656   if (!qh->facet_mergeset || !qh->degen_mergeset || !qh->vertex_mergeset) {
1657     qh_fprintf(qh, 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, qh_ERRqhull, NULL, NULL);
1660   }
1661   if (!SETempty_(qh->facet_mergeset) || !SETempty_(qh->degen_mergeset) || !SETempty_(qh->vertex_mergeset)) {
1662     qh_fprintf(qh, 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, qh->facet_mergeset), qh_setsize(qh, qh->degen_mergeset), qh_setsize(qh, qh->vertex_mergeset));
1664     qh_errexit(qh, qh_ERRqhull, NULL, NULL);
1665   }
1666   qh_settempfree(qh, &qh->facet_mergeset);
1667   qh_settempfree(qh, &qh->vertex_mergeset);
1668   qh_settempfree(qh, &qh->degen_mergeset);
1669 } /* freemergesets */
1670 
1671 /*-<a                             href="qh-merge_r.htm#TOC"
1672   >-------------------------------</a><a name="getmergeset">-</a>
1673 
1674   qh_getmergeset(qh, 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(qhT * qh,facetT * facetlist)1698 void qh_getmergeset(qhT *qh, facetT *facetlist) {
1699   facetT *facet, *neighbor, **neighborp;
1700   ridgeT *ridge, **ridgep;
1701   int nummerges;
1702   boolT simplicial;
1703 
1704   nummerges= qh_setsize(qh, qh->facet_mergeset);
1705   trace4((qh, 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(qh, facet, neighbor, simplicial))
1729           ridge->nonconvex= True;
1730         ridge->tested= True;
1731       }
1732     }
1733     facet->tested= True;
1734   }
1735   nummerges= qh_setsize(qh, 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, 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, qh->ferr, 2021, "qh_getmergeset: %d merges found\n", nummerges));
1748 } /* getmergeset */
1749 
1750 
1751 /*-<a                             href="qh-merge_r.htm#TOC"
1752   >-------------------------------</a><a name="getmergeset_initial">-</a>
1753 
1754   qh_getmergeset_initial(qh, 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(qhT * qh,facetT * facetlist)1775 void qh_getmergeset_initial(qhT *qh, 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(qh, 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, 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, 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, qh->ferr, 2022, "qh_getmergeset_initial: %d merges found\n", nummerges));
1816 } /* getmergeset_initial */
1817 
1818 /*-<a                             href="qh-merge_r.htm#TOC"
1819   >-------------------------------</a><a name="getpinchedmerges">-</a>
1820 
1821   qh_getpinchedmerges(qh, 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(qhT * qh,vertexT * apex,coordT maxdupdist,boolT * iscoplanar)1851 boolT qh_getpinchedmerges(qhT *qh, 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, 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(qh, 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, 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, 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, qh_ERRqhull, NULL, NULL);
1871     }
1872     /* dist is distance between vertices */
1873     pinched= qh_findbest_pinchedvertex(qh, 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, 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, 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, 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(qh, 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, 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(qh, bestpinched, bestvertex, MRGsubridge, bestdist, NULL, NULL);
1908         result= True;
1909       }else {
1910         trace2((qh, 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(qh, 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(qh, merge, (int)sizeof(mergeT));
1920   return result;
1921 }/* getpinchedmerges */
1922 
1923 /*-<a                             href="qh-merge_r.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_r.htm#TOC"
1944   >-------------------------------</a><a name="hashridge">-</a>
1945 
1946   qh_hashridge(qh, 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(qhT * qh,setT * hashtable,int hashsize,ridgeT * ridge,vertexT * oldvertex)1956 void qh_hashridge(qhT *qh, setT *hashtable, int hashsize, ridgeT *ridge, vertexT *oldvertex) {
1957   int hash;
1958   ridgeT *ridgeA;
1959 
1960   hash= qh_gethash(qh, 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_r.htm#TOC"
1974   >-------------------------------</a><a name="hashridge_find">-</a>
1975 
1976   qh_hashridge_find(qh, 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(qhT * qh,setT * hashtable,int hashsize,ridgeT * ridge,vertexT * vertex,vertexT * oldvertex,int * hashslot)1999 ridgeT *qh_hashridge_find(qhT *qh, 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(qh, 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_r.htm#TOC"
2025   >-------------------------------</a><a name="initmergesets">-</a>
2026 
2027   qh_initmergesets(qh )
2028     initialize the merge sets
2029     if 'all', include qh.degen_mergeset
2030 
2031   notes:
2032     matches qh_freemergesets
2033 */
qh_initmergesets(qhT * qh)2034 void qh_initmergesets(qhT *qh /* qh.facet_mergeset,degen_mergeset,vertex_mergeset */) {
2035 
2036   if (qh->facet_mergeset || qh->degen_mergeset || qh->vertex_mergeset) {
2037     qh_fprintf(qh, 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, qh_ERRqhull, NULL, NULL);
2040   }
2041   qh->degen_mergeset= qh_settemp(qh, qh->TEMPsize);
2042   qh->vertex_mergeset= qh_settemp(qh, qh->TEMPsize);
2043   qh->facet_mergeset= qh_settemp(qh, qh->TEMPsize); /* last temporary set for qh_forcedmerges */
2044 } /* initmergesets */
2045 
2046 /*-<a                             href="qh-merge_r.htm#TOC"
2047   >-------------------------------</a><a name="makeridges">-</a>
2048 
2049   qh_makeridges(qh, 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(qhT * qh,facetT * facet)2074 void qh_makeridges(qhT *qh, 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, 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_(qh, 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(qh);
2097       ridge->vertices= qh_setnew_delnthsorted(qh, 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(qh, &(facet->ridges), ridge);
2128       trace5((qh, 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(qh, &(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_r.htm#TOC"
2143   >-------------------------------</a><a name="mark_dupridges">-</a>
2144 
2145   qh_mark_dupridges(qh, 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(qhT * qh,facetT * facetlist,boolT allmerges)2192 void qh_mark_dupridges(qhT *qh, facetT *facetlist, boolT allmerges) {
2193   facetT *facet, *neighbor, **neighborp;
2194   int nummerge=0;
2195   mergeT *merge, **mergep;
2196 
2197   trace4((qh, 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(qh, 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, 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(qh, 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, qh->ferr, 1012, "qh_mark_dupridges: found %d duplicated ridges (MRGdupridge) for qh_getpinchedmerges\n", nummerge));
2233     return;
2234   }
2235   trace1((qh, 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(qh, facet);
2240   }
2241   trace3((qh, 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, 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, qh_ERRtopology, merge->facet1, merge->facet2);
2254       }else
2255         qh_setappend(qh, &merge->facet2->neighbors, merge->facet1);
2256       qh_makeridges(qh, merge->facet1);   /* and the missing ridges */
2257     }
2258   }
2259 } /* mark_dupridges */
2260 
2261 /*-<a                             href="qh-merge_r.htm#TOC"
2262   >-------------------------------</a><a name="maybe_duplicateridge">-</a>
2263 
2264   qh_maybe_duplicateridge(qh, 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(qhT * qh,ridgeT * ridgeA)2285 void qh_maybe_duplicateridge(qhT *qh, 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(qh, ridge, &pinched, &dist);
2306               trace2((qh, 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(qh, 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_r.htm#TOC"
2320   >-------------------------------</a><a name="maybe_duplicateridges">-</a>
2321 
2322   qh_maybe_duplicateridges(qh, 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(qhT * qh,facetT * facet)2341 void qh_maybe_duplicateridges(qhT *qh, 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_(qh, 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(qh, ridge, &pinched, &dist);
2369             if (ridge->top == ridge2->bottom && ridge->bottom == ridge2->top) {
2370               /* proof that ridges may have opposite orientation */
2371               trace2((qh, 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, 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(qh, 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_r.htm#TOC"
2388   >-------------------------------</a><a name="maydropneighbor">-</a>
2389 
2390   qh_maydropneighbor(qh, 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(qhT * qh,facetT * facet)2415 void qh_maydropneighbor(qhT *qh, facetT *facet) {
2416   ridgeT *ridge, **ridgep;
2417   facetT *neighbor, **neighborp;
2418 
2419   qh->visit_id++;
2420   trace4((qh, 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, 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, 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, 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, qh_ERRqhull, facet, neighbor);
2436     }
2437     if (neighbor->visitid != qh->visit_id) {
2438       trace2((qh, 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, 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, qh_ERRqhull, neighbor, facet);
2444       }
2445       zinc_(Zdropneighbor);
2446       qh_setdel(neighbor->neighbors, facet);
2447       if (qh_setsize(qh, neighbor->neighbors) < qh->hull_dim) {
2448         zinc_(Zdropdegen);
2449         qh_appendmergeset(qh, neighbor, neighbor, MRGdegen, 0.0, qh_ANGLEnone);
2450         trace2((qh, 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(qh, facet->neighbors) < qh->hull_dim) {
2457     zinc_(Zdropdegen);
2458     qh_appendmergeset(qh, facet, facet, MRGdegen, 0.0, qh_ANGLEnone);
2459     trace2((qh, qh->ferr, 2024, "qh_maydropneighbors: f%d is degenerate.\n", facet->id));
2460   }
2461 } /* maydropneighbor */
2462 
2463 
2464 /*-<a                             href="qh-merge_r.htm#TOC"
2465   >-------------------------------</a><a name="merge_degenredundant">-</a>
2466 
2467   qh_merge_degenredundant(qh)
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(qhT * qh)2490 int qh_merge_degenredundant(qhT *qh) {
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, qh->ferr, 2095, "qh_merge_degenredundant: merge %d degenerate, redundant, and mirror facets\n",
2501     qh_setsize(qh, qh->degen_mergeset)));
2502   mergedfacets= qh_settemp(qh, 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(qh, 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       qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
2514     if (mergetype == MRGredundant) {
2515       zinc_(Zredundant);
2516       facet3= qh_getreplacement(qh, facet2); /* the same facet if !facet2.visible */
2517       if (!facet3) {
2518           qh_fprintf(qh, 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, qh_ERRqhull, facet1, facet2);
2521       }
2522       qh_setunique(qh, &mergedfacets, facet3);
2523       if (facet1 == facet3) {
2524         continue;
2525       }
2526       trace2((qh, 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(qh, 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(qh, facet1->neighbors))) {
2533         zinc_(Zdelfacetdup);
2534         trace2((qh, qh->ferr, 2026, "qh_merge_degenredundant: facet f%d has no neighbors.  Deleted\n", facet1->id));
2535         qh_willdelete(qh, facet1, NULL);
2536         FOREACHvertex_(facet1->vertices) {
2537           qh_setdel(vertex->neighbors, facet1);
2538           if (!SETfirst_(vertex->neighbors)) {
2539             zinc_(Zdegenvertex);
2540             trace2((qh, 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, &qh->del_vertices, vertex);
2544           }
2545         }
2546         nummerges++;
2547       }else if (size < qh->hull_dim) {
2548         bestneighbor= qh_findbestneighbor(qh, facet1, &dist, &mindist, &maxdist);
2549         trace2((qh, 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(qh, 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(qh, &mergedfacets);
2562   return nummerges;
2563 } /* merge_degenredundant */
2564 
2565 /*-<a                             href="qh-merge_r.htm#TOC"
2566   >-------------------------------</a><a name="merge_nonconvex">-</a>
2567 
2568   qh_merge_nonconvex(qh, 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(qhT * qh,facetT * facet1,facetT * facet2,mergeType mergetype)2585 void qh_merge_nonconvex(qhT *qh, 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, 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, qh_ERRqhull, facet1, facet2);
2593   }
2594   if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
2595     qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
2596   trace3((qh, 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(qh, bestfacet, &dist, &mindist, &maxdist);
2606   neighbor= qh_findbestneighbor(qh, 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, 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(qh, 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_r.htm#TOC"
2651   >-------------------------------</a><a name="merge_pinchedvertices">-</a>
2652 
2653   qh_merge_pinchedvertices(qh, 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(qhT * qh,int apexpointid)2668 void qh_merge_pinchedvertices(qhT *qh, 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(qh);
2675   if (qh->visible_list || qh->newfacet_list || qh->newvertex_list) {
2676     qh_fprintf(qh, 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, 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 /* qh.vertex_mergeset */))) { /* only one at a time from qh_getpinchedmerges */
2684     if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
2685       qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
2686     if (merge->mergetype == MRGsubridge) {
2687       zzinc_(Zpinchedvertex);
2688       trace1((qh, 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, qh->vertex_mergeset)+1, apexpointid));
2690       qh_remove_mergetype(qh, qh->vertex_mergeset, MRGsubridge);
2691     }else {
2692       zzinc_(Zpinchduplicate);
2693       if (firstmerge)
2694         trace1((qh, qh->ferr, 1056, "qh_merge_pinchedvertices: merge %d pinched vertices from dupridges in merged facets, apex p%d\n",
2695            qh_setsize(qh, qh->vertex_mergeset)+1, apexpointid));
2696       firstmerge= False;
2697     }
2698     vertex= merge->vertex1;
2699     vertex2= merge->vertex2;
2700     dist= merge->distance;
2701     qh_memfree(qh, merge, (int)sizeof(mergeT)); /* merge is invalidated */
2702     qh_rename_adjacentvertex(qh, 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, qh->ferr, 2072, "qh_merge_pinchedvertices: merge degenerate f%d into an adjacent facet\n", mergeA->facet1->id);
2708         }else {
2709           qh_fprintf(qh, 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(qh); /* simplicial facets with both old and new vertices */
2715   }
2716   qh->isRenameVertex= False;
2717 }/* merge_pinchedvertices */
2718 
2719 /*-<a                             href="qh-merge_r.htm#TOC"
2720   >-------------------------------</a><a name="merge_twisted">-</a>
2721 
2722   qh_merge_twisted(qh, 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(qhT * qh,facetT * facet1,facetT * facet2)2737 void qh_merge_twisted(qhT *qh, 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     qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
2744   trace3((qh, 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(qh, facet1, &dist, &mindist, &maxdist);
2748   neighbor2= qh_findbestneighbor(qh, 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(qh, facet1->vertices, &bestvertex, &bestpinched);
2754     if (bestdist > mintwisted) {
2755       qh_fprintf(qh, 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, 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, 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(qh, 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_r.htm#TOC"
2782   >-------------------------------</a><a name="mergecycle">-</a>
2783 
2784   qh_mergecycle(qh, 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(qhT * qh,facetT * samecycle,facetT * newfacet)2810 void qh_mergecycle(qhT *qh, 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(qh);
2821   }
2822 #ifndef qh_NOtrace
2823   if (qh->TRACEmerge == zzval_(Ztotmerge))
2824     qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
2825   trace2((qh, 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, 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, qh->ferr, 8069, "  same cycle:");
2836     FORALLsame_cycle_(samecycle)
2837       qh_fprintf(qh, qh->ferr, 8070, " f%d", same->id);
2838     qh_fprintf(qh, qh->ferr, 8071, "\n");
2839   }
2840   if (qh->IStracing >=4)
2841     qh_errprint(qh, "MERGING CYCLE", samecycle, newfacet, NULL, NULL);
2842 #endif /* !qh_NOtrace */
2843   if (newfacet->tricoplanar) {
2844     if (!qh->TRInormals) {
2845       qh_fprintf(qh, qh->ferr, 6224, "qhull internal error (qh_mergecycle): does not work for tricoplanar facets.  Use option 'Q11'\n");
2846       qh_errexit(qh, qh_ERRqhull, newfacet, NULL);
2847     }
2848     newfacet->tricoplanar= False;
2849     newfacet->keepcentrum= False;
2850   }
2851   if (qh->CHECKfrequently)
2852     qh_checkdelridge(qh);
2853   if (!qh->VERTEXneighbors)
2854     qh_vertexneighbors(qh);
2855   apex= SETfirstt_(samecycle->vertices, vertexT);
2856   qh_makeridges(qh, newfacet);
2857   qh_mergecycle_neighbors(qh, samecycle, newfacet);
2858   qh_mergecycle_ridges(qh, samecycle, newfacet);
2859   qh_mergecycle_vneighbors(qh, samecycle, newfacet);
2860   if (SETfirstt_(newfacet->vertices, vertexT) != apex)
2861     qh_setaddnth(qh, &newfacet->vertices, 0, apex);  /* apex has last id */
2862   if (!newfacet->newfacet)
2863     qh_newvertices(qh, newfacet->vertices);
2864   qh_mergecycle_facets(qh, samecycle, newfacet);
2865   qh_tracemerge(qh, samecycle, newfacet, MRGcoplanarhorizon);
2866   /* check for degen_redundant_neighbors after qh_forcedmerges() */
2867   if (traceonce) {
2868     qh_fprintf(qh, qh->ferr, 8072, "qh_mergecycle: end of trace facet\n");
2869     qh->IStracing= tracerestore;
2870   }
2871 } /* mergecycle */
2872 
2873 /*-<a                             href="qh-merge_r.htm#TOC"
2874   >-------------------------------</a><a name="mergecycle_all">-</a>
2875 
2876   qh_mergecycle_all(qh, 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(qhT * qh,facetT * facetlist,boolT * wasmerge)2903 void qh_mergecycle_all(qhT *qh, 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, 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, qh->ferr, 6225, "qhull internal error (qh_mergecycle_all): f%d without normal\n", facet->id);
2915       qh_errexit(qh, 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         qh->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(qh, 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(qh, 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(qh, 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(qh, newfacet);
2969         qh_maybe_duplicateridges(qh, newfacet);
2970         newfacet->coplanarhorizon= False;
2971       }
2972     }
2973     numdegen += qh_merge_degenredundant(qh);
2974     *wasmerge= True;
2975     trace1((qh, 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_r.htm#TOC"
2981   >-------------------------------</a><a name="mergecycle_facets">-</a>
2982 
2983   qh_mergecycle_facets(qh, 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(qhT * qh,facetT * samecycle,facetT * newfacet)3004 void qh_mergecycle_facets(qhT *qh, facetT *samecycle, facetT *newfacet) {
3005   facetT *same, *next;
3006 
3007   trace4((qh, qh->ferr, 4030, "qh_mergecycle_facets: make newfacet new and samecycle deleted\n"));
3008   qh_removefacet(qh, newfacet);  /* append as a newfacet to end of qh->facet_list */
3009   qh_appendfacet(qh, 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(qh, same, newfacet);
3017   }
3018   if (newfacet->center
3019       && qh_setsize(qh, newfacet->vertices) <= qh->hull_dim + qh_MAXnewcentrum) {
3020     qh_memfree(qh, newfacet->center, qh->normal_size);
3021     newfacet->center= NULL;
3022   }
3023   trace3((qh, 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_r.htm#TOC"
3028   >-------------------------------</a><a name="mergecycle_neighbors">-</a>
3029 
3030   qh_mergecycle_neighbors(qh, 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(qhT * qh,facetT * samecycle,facetT * newfacet)3062 void qh_mergecycle_neighbors(qhT *qh, 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(qh, samecycle);
3072     same->visitid= samevisitid;
3073   }
3074   newfacet->visitid= ++qh->visit_id;
3075   trace4((qh, 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(qh, newfacet->neighbors);
3084 
3085   trace4((qh, 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(qh, &newfacet->neighbors, neighbor);
3093           qh_setreplace(qh, 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(qh, 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(qh, &neighbor->neighbors, newfacet);
3114           qh_setappend(qh, &newfacet->neighbors, neighbor);
3115           neighbor->visitid= qh->visit_id;
3116           newneighbors++;
3117         }
3118       }
3119     }
3120   }
3121   trace2((qh, 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_r.htm#TOC"
3126   >-------------------------------</a><a name="mergecycle_ridges">-</a>
3127 
3128   qh_mergecycle_ridges(qh, 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(qhT * qh,facetT * samecycle,facetT * newfacet)3160 void qh_mergecycle_ridges(qhT *qh, 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, 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(qh, newfacet->ridges);
3177 
3178   trace4((qh, 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(qh, &newfacet->ridges, ridge);
3189         numold++;  /* already set by qh_mergecycle_neighbors */
3190         continue;
3191       }else {
3192         qh_fprintf(qh, qh->ferr, 6098, "qhull internal error (qh_mergecycle_ridges): bad ridge r%d\n", ridge->id);
3193         qh_errexit(qh, qh_ERRqhull, NULL, ridge);
3194       }
3195       if (neighbor == newfacet) {
3196         if (qh->traceridge == ridge)
3197           qh->traceridge= NULL;
3198         qh_setfree(qh, &(ridge->vertices));
3199         qh_memfree_(qh, 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(qh, &(ridge->vertices));
3206         qh_memfree_(qh, ridge, (int)sizeof(ridgeT), freelistp);
3207         numold++;
3208       }else {
3209         qh_setappend(qh, &newfacet->ridges, ridge);
3210         numold++;
3211       }
3212     }
3213     if (same->ridges)
3214       qh_settruncate(qh, same->ridges, 0);
3215     if (!same->simplicial)
3216       continue;
3217     FOREACHneighbor_i_(qh, same) {       /* note: !newfact->simplicial */
3218       if (neighbor->visitid != samevisitid && neighbor->simplicial) {
3219         ridge= qh_newridge(qh);
3220         ridge->vertices= qh_setnew_delnthsorted(qh, 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(qh, &(newfacet->ridges), ridge);
3233         qh_setappend(qh, &(neighbor->ridges), ridge);
3234         if (qh->ridge_id == qh->traceridge_id)
3235           qh->traceridge= ridge;
3236         numnew++;
3237       }
3238     }
3239   }
3240 
3241   trace2((qh, 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_r.htm#TOC"
3246   >-------------------------------</a><a name="mergecycle_vneighbors">-</a>
3247 
3248   qh_mergecycle_vneighbors(qh, 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(qhT * qh,facetT * samecycle,facetT * newfacet)3270 void qh_mergecycle_vneighbors(qhT *qh, facetT *samecycle, facetT *newfacet) {
3271   facetT *neighbor, **neighborp;
3272   unsigned int mergeid;
3273   vertexT *vertex, **vertexp, *apex;
3274   setT *vertices;
3275 
3276   trace4((qh, 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(qh, samecycle); /* temp */
3280   apex= SETfirstt_(samecycle->vertices, vertexT);
3281   qh_setappend(qh, &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(qh, vertex->neighbors);
3289     qh_setappend(qh, &vertex->neighbors, newfacet);
3290     if (!SETsecond_(vertex->neighbors)) {
3291       zinc_(Zcyclevertex);
3292       trace2((qh, 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, &qh->del_vertices, vertex);
3297     }
3298   }
3299   qh_settempfree(qh, &vertices);
3300   trace3((qh, 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_r.htm#TOC"
3305   >-------------------------------</a><a name="mergefacet">-</a>
3306 
3307   qh_mergefacet(qh, 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(qhT * qh,facetT * facet1,facetT * facet2,mergeType mergetype,realT * mindist,realT * maxdist,boolT mergeapex)3362 void qh_mergefacet(qhT *qh, 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, 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, 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(qh);
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, 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, 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, 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(qh, facet1->vertices);
3424       onemerge= qh->ONEmerge + qh->DISTround;
3425       if (vertexdist > mintwisted) {
3426         qh_fprintf(qh, 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, 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, qh_ERRwide, facet1, facet2);
3433     }
3434   }
3435   if (facet1 == facet2 || facet1->visible || facet2->visible) {
3436     qh_fprintf(qh, 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, qh_ERRqhull, facet1, facet2);
3439   }
3440   if (qh->num_facets - qh->num_visible <= qh->hull_dim + 1) {
3441     qh_fprintf(qh, 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, qh->ferr, 8079, "    Option 'Qx' may avoid this problem.\n");
3445     qh_errexit(qh, qh_ERRtopology, NULL, NULL);
3446   }
3447   if (!qh->VERTEXneighbors)
3448     qh_vertexneighbors(qh);
3449   qh_makeridges(qh, facet1);
3450   qh_makeridges(qh, facet2);
3451   if (qh->IStracing >=4)
3452     qh_errprint(qh, "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(qh, facet1, facet2);
3474   if (qh->hull_dim > 2 && qh_setsize(qh, facet1->vertices) == qh->hull_dim)
3475     qh_mergesimplex(qh, 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(qh, facet1, facet2);
3482     else {
3483       qh_mergeneighbors(qh, facet1, facet2);
3484       qh_mergevertices(qh, facet1->vertices, &facet2->vertices);
3485     }
3486     qh_mergeridges(qh, facet1, facet2);
3487     qh_mergevertex_neighbors(qh, facet1, facet2);
3488     if (!facet2->newfacet)
3489       qh_newvertices(qh, 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(qh, facet2);  /* append as a newfacet to end of qh->facet_list */
3501   qh_appendfacet(qh, facet2);
3502   facet2->newfacet= True;
3503   facet2->tested= False;
3504   qh_tracemerge(qh, facet1, facet2, mergetype);
3505   if (traceonce) {
3506     qh_fprintf(qh, qh->ferr, 8080, "qh_mergefacet: end of wide tracing\n");
3507     qh->IStracing= tracerestore;
3508   }
3509   if (mergetype != MRGcoplanarhorizon) {
3510     trace3((qh, 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(qh, facet2);
3513     qh_test_degen_neighbors(qh, 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(qh, facet2); /* may occur in qh_merge_pinchedvertices, e.g., rbox 175 C3,2e-13 D4 t1545228104 | qhull d */
3516     qh_maybe_duplicateridges(qh, facet2);
3517   }
3518   qh_willdelete(qh, facet1, facet2);
3519 } /* mergefacet */
3520 
3521 
3522 /*-<a                             href="qh-merge_r.htm#TOC"
3523   >-------------------------------</a><a name="mergefacet2d">-</a>
3524 
3525   qh_mergefacet2d(qh, 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(qhT * qh,facetT * facet1,facetT * facet2)3545 void qh_mergefacet2d(qhT *qh, 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(qh, neighborB->neighbors, facet1, facet2);
3596   trace4((qh, 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_r.htm#TOC"
3602   >-------------------------------</a><a name="mergeneighbors">-</a>
3603 
3604   qh_mergeneighbors(qh, 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(qhT * qh,facetT * facet1,facetT * facet2)3622 void qh_mergeneighbors(qhT *qh, facetT *facet1, facetT *facet2) {
3623   facetT *neighbor, **neighborp;
3624 
3625   trace4((qh, 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(qh, 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(qh, neighbor->neighbors, facet1, facet2);
3640       }
3641     }else if (neighbor != facet2) {
3642       qh_setappend(qh, &(facet2->neighbors), neighbor);
3643       qh_setreplace(qh, 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_r.htm#TOC"
3652   >-------------------------------</a><a name="mergeridges">-</a>
3653 
3654   qh_mergeridges(qh, 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(qhT * qh,facetT * facet1,facetT * facet2)3670 void qh_mergeridges(qhT *qh, facetT *facet1, facetT *facet2) {
3671   ridgeT *ridge, **ridgep;
3672 
3673   trace4((qh, 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(qh, 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(qh, &(facet2->ridges), ridge);
3691   }
3692 } /* mergeridges */
3693 
3694 
3695 /*-<a                             href="qh-merge_r.htm#TOC"
3696   >-------------------------------</a><a name="mergesimplex">-</a>
3697 
3698   qh_mergesimplex(qh, 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(qhT * qh,facetT * facet1,facetT * facet2,boolT mergeapex)3738 void qh_mergesimplex(qhT *qh, 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, 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(qh, facet2->vertices);  /* apex, the first vertex, is already new */
3750     if (SETfirstt_(facet2->vertices, vertexT) != opposite) {
3751       qh_setaddnth(qh, &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, 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(qh, facet2, opposite);
3775     if (!facet2->newfacet)
3776       qh_newvertices(qh, facet2->vertices);
3777     else if (!opposite->newfacet) {
3778       qh_removevertex(qh, opposite);
3779       qh_appendvertex(qh, opposite);
3780     }
3781   }
3782   trace4((qh, 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(qh, vertex->neighbors, facet1, facet2);
3787     else {
3788       qh_setdel(vertex->neighbors, facet1);
3789       if (!SETsecond_(vertex->neighbors))
3790         qh_mergevertex_del(qh, vertex, facet1, facet2);
3791     }
3792   }
3793   trace4((qh, 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(qh, 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, 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, qh_ERRqhull, facet1, otherfacet);
3809     }else {
3810       trace4((qh, 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(qh, &facet2->ridges, ridge);
3813       if (otherfacet->visitid != qh->visit_id) {
3814         qh_setappend(qh, &facet2->neighbors, otherfacet);
3815         qh_setreplace(qh, otherfacet->neighbors, facet1, facet2);
3816         otherfacet->visitid= qh->visit_id;
3817       }else {
3818         if (otherfacet->simplicial)    /* is degen, needs ridges */
3819           qh_makeridges(qh, otherfacet);
3820         if (SETfirstt_(otherfacet->neighbors, facetT) == facet1) {
3821           /* keep new, otherfacet->neighbors->horizon */
3822           qh_setdel(otherfacet->neighbors, facet2);
3823           qh_setreplace(qh, 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, 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_r.htm#TOC"
3843   >-------------------------------</a><a name="mergevertex_del">-</a>
3844 
3845   qh_mergevertex_del(qh, 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(qhT * qh,vertexT * vertex,facetT * facet1,facetT * facet2)3852 void qh_mergevertex_del(qhT *qh, vertexT *vertex, facetT *facet1, facetT *facet2) {
3853 
3854   zinc_(Zmergevertex);
3855   trace2((qh, 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, &qh->del_vertices, vertex);
3860 } /* mergevertex_del */
3861 
3862 /*-<a                             href="qh-merge_r.htm#TOC"
3863   >-------------------------------</a><a name="mergevertex_neighbors">-</a>
3864 
3865   qh_mergevertex_neighbors(qh, 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(qhT * qh,facetT * facet1,facetT * facet2)3878 void qh_mergevertex_neighbors(qhT *qh, facetT *facet1, facetT *facet2) {
3879   vertexT *vertex, **vertexp;
3880 
3881   trace4((qh, 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, 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(qh, "TRACE", NULL, NULL, NULL, qh->tracevertex);
3887   }
3888   FOREACHvertex_(facet1->vertices) {
3889     if (vertex->visitid != qh->vertex_visit)
3890       qh_setreplace(qh, vertex->neighbors, facet1, facet2);
3891     else {
3892       qh_setdel(vertex->neighbors, facet1);
3893       if (!SETsecond_(vertex->neighbors))
3894         qh_mergevertex_del(qh, vertex, facet1, facet2);
3895     }
3896   }
3897   if (qh->tracevertex)
3898     qh_errprint(qh, "TRACE", NULL, NULL, NULL, qh->tracevertex);
3899 } /* mergevertex_neighbors */
3900 
3901 
3902 /*-<a                             href="qh-merge_r.htm#TOC"
3903   >-------------------------------</a><a name="mergevertices">-</a>
3904 
3905   qh_mergevertices(qh, 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(qhT * qh,setT * vertices1,setT ** vertices2)3916 void qh_mergevertices(qhT *qh, setT *vertices1, setT **vertices2) {
3917   int newsize= qh_setsize(qh, vertices1)+qh_setsize(qh, *vertices2) - qh->hull_dim + 1;
3918   setT *mergedvertices;
3919   vertexT *vertex, **vertexp, **vertex2= SETaddr_(*vertices2, vertexT);
3920 
3921   mergedvertices= qh_settemp(qh, newsize);
3922   FOREACHvertex_(vertices1) {
3923     if (!*vertex2 || vertex->id > (*vertex2)->id)
3924       qh_setappend(qh, &mergedvertices, vertex);
3925     else {
3926       while (*vertex2 && (*vertex2)->id > vertex->id)
3927         qh_setappend(qh, &mergedvertices, *vertex2++);
3928       if (!*vertex2 || (*vertex2)->id < vertex->id)
3929         qh_setappend(qh, &mergedvertices, vertex);
3930       else
3931         qh_setappend(qh, &mergedvertices, *vertex2++);
3932     }
3933   }
3934   while (*vertex2)
3935     qh_setappend(qh, &mergedvertices, *vertex2++);
3936   if (newsize < qh_setsize(qh, mergedvertices)) {
3937     qh_fprintf(qh, qh->ferr, 6100, "qhull internal error (qh_mergevertices): facets did not share a ridge\n");
3938     qh_errexit(qh, qh_ERRqhull, NULL, NULL);
3939   }
3940   qh_setfree(qh, vertices2);
3941   *vertices2= mergedvertices;
3942   qh_settemppop(qh);
3943 } /* mergevertices */
3944 
3945 
3946 /*-<a                             href="qh-merge_r.htm#TOC"
3947   >-------------------------------</a><a name="neighbor_intersections">-</a>
3948 
3949   qh_neighbor_intersections(qh, 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(qhT * qh,vertexT * vertex)3971 setT *qh_neighbor_intersections(qhT *qh, 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(qh, neighborA->vertices, 0);
3987   else
3988     intersect= qh_vertexintersect_new(qh, neighborA->vertices, neighborB->vertices);
3989   qh_settemppush(qh, intersect);
3990   qh_setdelsorted(intersect, vertex);
3991   FOREACHneighbor_i_(qh, vertex) {
3992     if (neighbor_i >= 2) {
3993       zinc_(Zintersectnum);
3994       qh_vertexintersect(qh, &intersect, neighbor->vertices);
3995       if (!SETfirst_(intersect)) {
3996         zinc_(Zintersectfail);
3997         qh_settempfree(qh, &intersect);
3998         return NULL;
3999       }
4000     }
4001   }
4002   trace3((qh, qh->ferr, 3007, "qh_neighbor_intersections: %d vertices in neighbor intersection of v%d\n",
4003           qh_setsize(qh, intersect), vertex->id));
4004   return intersect;
4005 } /* neighbor_intersections */
4006 
4007 /*-<a                             href="qh-merge_r.htm#TOC"
4008   >-------------------------------</a><a name="neighbor_vertices">-</a>
4009 
4010   qh_neighbor_vertices(qh, 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(qhT * qh,vertexT * vertexA,setT * subridge)4022 setT *qh_neighbor_vertices(qhT *qh, vertexT *vertexA, setT *subridge) {
4023   facetT *neighbor, **neighborp;
4024   vertexT *vertex, **vertexp;
4025   setT *vertices= qh_settemp(qh, 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(qh, vertexA, neighbor, &vertices);
4038   }
4039   trace3((qh, qh->ferr, 3035, "qh_neighbor_vertices: %d non-subridge, vertex neighbors for v%d\n",
4040     qh_setsize(qh, vertices), vertexA->id));
4041   return vertices;
4042 } /* neighbor_vertices */
4043 
4044 /*-<a                             href="qh-merge_r.htm#TOC"
4045   >-------------------------------</a><a name="neighbor_vertices_facet">-</a>
4046 
4047   qh_neighbor_vertices_facet(qh, 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(qhT * qh,vertexT * vertexA,facetT * facet,setT ** vertices)4067 void qh_neighbor_vertices_facet(qhT *qh, 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(qh, 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(qh, vertices, vertex);
4108               count++;
4109             }
4110           }
4111         }
4112       }
4113     }
4114   }
4115   facet->visitid= qh->visit_id-1;
4116   if (count) {
4117     trace4((qh, 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_r.htm#TOC"
4124   >-------------------------------</a><a name="newvertices">-</a>
4125 
4126   qh_newvertices(qh, 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(qhT * qh,setT * vertices)4133 void qh_newvertices(qhT *qh, setT *vertices) {
4134   vertexT *vertex, **vertexp;
4135 
4136   FOREACHvertex_(vertices) {
4137     if (!vertex->newfacet) {
4138       qh_removevertex(qh, vertex);
4139       qh_appendvertex(qh, vertex);
4140     }
4141   }
4142 } /* newvertices */
4143 
4144 /*-<a                             href="qh-merge_r.htm#TOC"
4145   >-------------------------------</a><a name="next_vertexmerge">-</a>
4146 
4147   qh_next_vertexmerge(qh )
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(qhT * qh)4157 mergeT *qh_next_vertexmerge(qhT *qh /* qh.vertex_mergeset */) {
4158   mergeT *merge;
4159   int merge_i, merge_n, best_i= -1;
4160   realT bestdist= REALmax;
4161 
4162   FOREACHmerge_i_(qh, qh->vertex_mergeset) {
4163     if (!merge->vertex1 || !merge->vertex2) {
4164       qh_fprintf(qh, 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, qh_ERRqhull, merge->facet1, NULL);
4167     }
4168     if (merge->vertex1->deleted || merge->vertex2->deleted) {
4169       trace3((qh, 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(qh, merge);
4172       qh_setdelnth(qh, qh->vertex_mergeset, merge_i);
4173       merge_i--; merge_n--;
4174       qh_memfree(qh, 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, 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, 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, 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, qh_ERRtopology, NULL, merge->ridge1);
4197     }
4198     qh_setdelnth(qh, qh->vertex_mergeset, best_i);
4199   }
4200   return merge;
4201 } /* next_vertexmerge */
4202 
4203 /*-<a                             href="qh-merge_r.htm#TOC"
4204   >-------------------------------</a><a name="opposite_horizonfacet">-</a>
4205 
4206   qh_opposite_horizonfacet(qh, 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(qhT * qh,mergeT * merge,vertexT ** opposite)4217 facetT *qh_opposite_horizonfacet(qhT *qh, 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, 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, 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, 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, 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_r.htm#TOC"
4247   >-------------------------------</a><a name="reducevertices">-</a>
4248 
4249   qh_reducevertices(qh)
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(qhT * qh)4279 boolT qh_reducevertices(qhT *qh) {
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, qh->ferr, 2101, "qh_reducevertices: reduce extra vertices, shared vertices, and redundant vertices\n"));
4288   if (qh_merge_degenredundant(qh))
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(qh, newfacet)) {
4296         qh_degen_redundant_facet(qh, newfacet);
4297         if (qh_merge_degenredundant(qh)) {
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(qh, vertex, newfacet)) {
4312             numshare++;
4313             if (qh_merge_degenredundant(qh)) {
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(qh, vertex)) {
4327         numrename++;
4328         if (qh_merge_degenredundant(qh)) {
4329           degenredun= True;
4330           goto LABELrestart;
4331         }
4332       }
4333     }
4334   }
4335   trace1((qh, 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_r.htm#TOC"
4341   >-------------------------------</a><a name="redundant_vertex">-</a>
4342 
4343   qh_redundant_vertex(qh, 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(qhT * qh,vertexT * vertex)4363 vertexT *qh_redundant_vertex(qhT *qh, vertexT *vertex) {
4364   vertexT *newvertex= NULL;
4365   setT *vertices, *ridges;
4366 
4367   trace3((qh, 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(qh, vertex))) {
4369     ridges= qh_vertexridges(qh, vertex, !qh_ALL);
4370     if ((newvertex= qh_find_newvertex(qh, vertex, vertices, ridges))) {
4371       zinc_(Zrenameall);
4372       qh_renamevertex(qh, vertex, newvertex, ridges, NULL, NULL); /* ridges invalidated */
4373     }
4374     qh_settempfree(qh, &ridges);
4375     qh_settempfree(qh, &vertices);
4376   }
4377   return newvertex;
4378 } /* redundant_vertex */
4379 
4380 /*-<a                             href="qh-merge_r.htm#TOC"
4381   >-------------------------------</a><a name="remove_extravertices">-</a>
4382 
4383   qh_remove_extravertices(qh, 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(qhT * qh,facetT * facet)4403 boolT qh_remove_extravertices(qhT *qh, 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, 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(qh, vertex->neighbors)) {
4426         vertex->deleted= True;
4427         qh_setappend(qh, &qh->del_vertices, vertex);
4428         zinc_(Zremvertexdel);
4429         trace2((qh, qh->ferr, 2036, "qh_remove_extravertices: v%d deleted because it's lost all ridges\n", vertex->id));
4430       }else
4431         trace3((qh, 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_r.htm#TOC"
4439   >-------------------------------</a><a name="remove_mergetype">-</a>
4440 
4441   qh_remove_mergetype(qh, mergeset, mergetype )
4442     Remove mergetype merges from mergeset
4443 
4444   notes:
4445     Does not preserve order
4446 */
qh_remove_mergetype(qhT * qh,setT * mergeset,mergeType type)4447 void qh_remove_mergetype(qhT *qh, setT *mergeset, mergeType type) {
4448   mergeT *merge;
4449   int merge_i, merge_n;
4450 
4451   FOREACHmerge_i_(qh, mergeset) {
4452     if (merge->mergetype == type) {
4453         trace3((qh, 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(qh, mergeset, merge_i);
4456         merge_i--; merge_n--;  /* repeat with next merge */
4457     }
4458   }
4459 } /* remove_mergetype */
4460 
4461 /*-<a                             href="qh-merge_r.htm#TOC"
4462   >-------------------------------</a><a name="rename_adjacentvertex">-</a>
4463 
4464   qh_rename_adjacentvertex(qh, 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(qhT * qh,vertexT * oldvertex,vertexT * newvertex,realT dist)4478 void qh_rename_adjacentvertex(qhT *qh, 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(qh, oldvertex->neighbors);
4484   int newsize= qh_setsize(qh, 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, 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, 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, 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, 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, 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(qh, neighbor); /* no longer simplicial, nummerge==0, skipped by qh_maybe_duplicateridge */
4515       }else {
4516         qh_replacefacetvertex(qh, neighbor, oldvertex, newvertex);
4517         qh_setunique(qh, &newvertex->neighbors, neighbor);
4518         qh_newvertices(qh, neighbor->vertices);  /* for qh_update_vertexneighbors of vertex neighbors */
4519       }
4520     }
4521   }
4522   ridges= qh_vertexridges(qh, oldvertex, qh_ALL);
4523   if (istrace) {
4524     FOREACHridge_(ridges) {
4525       qh_printridge(qh, qh->ferr, ridge);
4526     }
4527   }
4528   FOREACHneighbor_(oldvertex) {
4529     if (!neighbor->simplicial){
4530       qh_addfacetvertex(qh, neighbor, newvertex);
4531       qh_setunique(qh, &newvertex->neighbors, neighbor);
4532       qh_newvertices(qh, neighbor->vertices);  /* for qh_update_vertexneighbors of vertex neighbors */
4533       if (qh->newfacet_list == qh->facet_tail) {
4534         qh_removefacet(qh, neighbor);  /* add a neighbor to newfacet_list so that qh_partitionvisible has a newfacet */
4535         qh_appendfacet(qh, neighbor);
4536         neighbor->newfacet= True;
4537       }
4538     }
4539   }
4540   qh_renamevertex(qh, oldvertex, newvertex, ridges, NULL, NULL);  /* ridges invalidated */
4541   if (oldvertex->deleted && !oldvertex->partitioned) {
4542     FOREACHneighbor_(newvertex) {
4543       if (!neighbor->visible) {
4544         qh_distplane(qh, oldvertex->point, neighbor, &dist2);
4545         if (dist2>maxdist2) {
4546           maxdist2= dist2;
4547           maxfacet= neighbor;
4548         }
4549       }
4550     }
4551     trace2((qh, 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(qh, oldvertex->point), oldvertex->id, maxfacet->id, maxdist2))
4553     qh_partitioncoplanar(qh, oldvertex->point, maxfacet, NULL, !qh_ALL);  /* faster with maxdist2, otherwise duplicates distance tests from maxdist2/dist2 */
4554     oldvertex->partitioned= True;
4555   }
4556   qh_settempfree(qh, &ridges);
4557 } /* rename_adjacentvertex */
4558 
4559 /*-<a                             href="qh-merge_r.htm#TOC"
4560   >-------------------------------</a><a name="rename_sharedvertex">-</a>
4561 
4562   qh_rename_sharedvertex(qh, 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(qhT * qh,vertexT * vertex,facetT * facet)4588 vertexT *qh_rename_sharedvertex(qhT *qh, vertexT *vertex, facetT *facet) {
4589   facetT *neighbor, **neighborp, *neighborA= NULL;
4590   setT *vertices, *ridges;
4591   vertexT *newvertex= NULL;
4592 
4593   if (qh_setsize(qh, 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, 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(qh, "ERRONEOUS", facet, NULL, NULL, vertex);
4615     qh_errexit(qh, qh_ERRqhull, NULL, NULL);
4616   }
4617   if (neighborA) { /* avoid warning */
4618     /* the vertex is shared by facet and neighborA */
4619     ridges= qh_settemp(qh, qh->TEMPsize);
4620     neighborA->visitid= ++qh->visit_id;
4621     qh_vertexridges_facet(qh, vertex, facet, &ridges);
4622     trace2((qh, qh->ferr, 2037, "qh_rename_sharedvertex: p%d(v%d) is shared by f%d(%d ridges) and f%d\n",
4623       qh_pointid(qh, vertex->point), vertex->id, facet->id, qh_setsize(qh, ridges), neighborA->id));
4624     zinc_(Zintersectnum);
4625     vertices= qh_vertexintersect_new(qh, facet->vertices, neighborA->vertices);
4626     qh_setdel(vertices, vertex);
4627     qh_settemppush(qh, vertices);
4628     if ((newvertex= qh_find_newvertex(qh, vertex, vertices, ridges)))
4629       qh_renamevertex(qh, vertex, newvertex, ridges, facet, neighborA);  /* ridges invalidated */
4630     qh_settempfree(qh, &vertices);
4631     qh_settempfree(qh, &ridges);
4632   }
4633   return newvertex;
4634 } /* rename_sharedvertex */
4635 
4636 /*-<a                             href="qh-merge_r.htm#TOC"
4637   >-------------------------------</a><a name="renameridgevertex">-</a>
4638 
4639   qh_renameridgevertex(qh, 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(qhT * qh,ridgeT * ridge,vertexT * oldvertex,vertexT * newvertex)4659 boolT qh_renameridgevertex(qhT *qh, 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, 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, qh_ERRqhull, NULL, ridge);
4669   }
4670   qh_setdelnthsorted(qh, 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(qh, ridge);
4676       trace2((qh, 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(qh, ridge); /* ridge.vertices deleted */
4679       return False;
4680     }
4681     if (vertex->id < newvertex->id)
4682       break;
4683     nth++;
4684   }
4685   qh_setaddnth(qh, &ridge->vertices, nth, newvertex);
4686   ridge->simplicialtop= False;
4687   ridge->simplicialbot= False;
4688   if (abs(oldnth - nth)%2) {
4689     trace3((qh, 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_r.htm#TOC"
4700   >-------------------------------</a><a name="renamevertex">-</a>
4701 
4702   qh_renamevertex(qh, 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(qhT * qh,vertexT * oldvertex,vertexT * newvertex,setT * ridges,facetT * oldfacet,facetT * neighborA)4735 void qh_renamevertex(qhT *qh, 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, 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(qh, ridges), getid_(oldfacet), getid_(neighborA));
4747   }
4748 #endif
4749   FOREACHridge_(ridges) {
4750     if (qh_renameridgevertex(qh, ridge, oldvertex, newvertex)) { /* ridge is deleted if False, invalidating ridges */
4751       topsize= qh_setsize(qh, ridge->top->vertices);
4752       bottomsize= qh_setsize(qh, 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, 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, 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(qh, ridge);
4761     }
4762   }
4763   if (!oldfacet) {
4764     /* stat Zrenameall or Zpinchduplicate */
4765     if (istrace)
4766       qh_fprintf(qh, 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(qh, neighbor); /* e.g., rbox 175 C3,2e-13 D4 t1545235541 | qhull d */
4771       }else {
4772         if (istrace)
4773           qh_fprintf(qh, qh->ferr, 4080, "qh_renamevertex: rename vertices in non-simplicial neighbor f%d of v%d\n", neighbor->id, oldvertex->id);
4774         qh_maydropneighbor(qh, neighbor);
4775         qh_setdelsorted(neighbor->vertices, oldvertex); /* if degenerate, qh_degen_redundant_facet will add to mergeset */
4776         if (qh_remove_extravertices(qh, neighbor))
4777           neighborp--; /* neighbor deleted from oldvertex neighborset */
4778         qh_degen_redundant_facet(qh, neighbor); /* either direction may be redundant, faster if combine? */
4779         qh_test_redundant_neighbors(qh, neighbor);
4780         qh_test_degen_neighbors(qh, neighbor);
4781       }
4782     }
4783     if (!oldvertex->deleted) {
4784       oldvertex->deleted= True;
4785       qh_setappend(qh, &qh->del_vertices, oldvertex);
4786     }
4787   }else if (qh_setsize(qh, oldvertex->neighbors) == 2) {
4788     zinc_(Zrenameshare);
4789     if (istrace)
4790       qh_fprintf(qh, 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(qh, neighbor);
4795     }
4796     oldvertex->deleted= True;
4797     qh_setappend(qh, &qh->del_vertices, oldvertex);
4798   }else {
4799     zinc_(Zrenamepinch);
4800     if (istrace || qh->IStracing >= 1)
4801       qh_fprintf(qh, 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(qh, neighborA))
4806       qh_degen_redundant_facet(qh, neighborA);
4807   }
4808   if (oldfacet)
4809     qh_degen_redundant_facet(qh, oldfacet);
4810 } /* renamevertex */
4811 
4812 /*-<a                             href="qh-merge_r.htm#TOC"
4813   >-------------------------------</a><a name="test_appendmerge">-</a>
4814 
4815   qh_test_appendmerge(qh, 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_r.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(qhT * qh,facetT * facet,facetT * neighbor,boolT simplicial)4839 boolT qh_test_appendmerge(qhT *qh, 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(qh, facet->normal, neighbor->normal);
4847     okangle= True;
4848     zinc_(Zangletests);
4849     if (angle > qh->cos_max) {
4850       zinc_(Zcoplanarangle);
4851       qh_appendmergeset(qh, facet, neighbor, MRGanglecoplanar, 0.0, angle);
4852       trace2((qh, 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(qh, facet, neighbor, angle, okangle);
4859   else
4860     return qh_test_nonsimplicial_merge(qh, facet, neighbor, angle, okangle);
4861 } /* test_appendmerge */
4862 
4863 /*-<a                             href="qh-merge_r.htm#TOC"
4864   >-------------------------------</a><a name="test_centrum_merge">-</a>
4865 
4866   qh_test_centrum_merge(qh, 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(qhT * qh,facetT * facet,facetT * neighbor,realT angle,boolT okangle)4895 boolT qh_test_centrum_merge(qhT *qh, 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(qh, facet);
4901   zzinc_(Zcentrumtests);
4902   qh_distplane(qh, 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(qh, neighbor);
4909   zzinc_(Zcentrumtests);
4910   qh_distplane(qh, 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(qh, facet->normal, neighbor->normal);
4919     zinc_(Zangletests);
4920   }
4921   if (isconcave && iscoplanar) {
4922     zinc_(Zconcavecoplanarridge);
4923     if (dist > dist2)
4924       qh_appendmergeset(qh, facet, neighbor, MRGconcavecoplanar, dist, angle);
4925     else
4926       qh_appendmergeset(qh, neighbor, facet, MRGconcavecoplanar, dist2, angle);
4927     trace0((qh, 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(qh, facet, neighbor, MRGconcave, mergedist, angle);
4933     trace0((qh, 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(qh, facet, neighbor, MRGcoplanar, mergedist, angle);
4939     trace2((qh, 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_r.htm#TOC"
4946   >-------------------------------</a><a name="test_degen_neighbors">-</a>
4947 
4948   qh_test_degen_neighbors(qh, 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(qhT * qh,facetT * facet)4958 void qh_test_degen_neighbors(qhT *qh, facetT *facet) {
4959   facetT *neighbor, **neighborp;
4960   int size;
4961 
4962   trace4((qh, 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, 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, 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(qh, neighbor->neighbors)) < qh->hull_dim) {
4973       qh_appendmergeset(qh, neighbor, neighbor, MRGdegen, 0.0, 1.0);
4974       trace2((qh, 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_r.htm#TOC"
4981   >-------------------------------</a><a name="test_nonsimplicial_merge">-</a>
4982 
4983   qh_test_nonsimplicial_merge(qh, 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(qhT * qh,facetT * facet,facetT * neighbor,realT angle,boolT okangle)5032 boolT qh_test_nonsimplicial_merge(qhT *qh, 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(qh, facet);
5049   zzinc_(Zcentrumtests);
5050   qh_distplane(qh, 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(qh, neighbor);
5059   zzinc_(Zcentrumtests);
5060   qh_distplane(qh, 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(qh, facet->normal, neighbor->normal);
5073       zinc_(Zangletests);
5074     }
5075     mergedist= fmax_(dist, dist2);
5076     zinc_(Zconcaveridge);
5077     qh_appendmergeset(qh, facet, neighbor, MRGconcave, mergedist, angle);
5078     trace0((qh, 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(qh, neighbor, facet, &maxdist2, &mindist2);
5088       if (!maxvertex2) {
5089         qh_appendmergeset(qh, 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(qh, facet, neighbor, &maxdist, &mindist);
5096       if (!maxvertex) {
5097         qh_appendmergeset(qh, facet, neighbor, MRGredundant, maxdist, qh_ANGLEnone);
5098         isredundant= True;
5099       }
5100     }
5101   }else {
5102     maxvertex= qh_furthestvertex(qh, facet, neighbor, &maxdist, &mindist);
5103     if (maxvertex) {
5104       maxvertex2= qh_furthestvertex(qh, neighbor, facet, &maxdist2, &mindist2);
5105       if (!maxvertex2) {
5106         qh_appendmergeset(qh, 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(qh, 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(qh, facet->normal, neighbor->normal);
5139     zinc_(Zangletests);
5140   }
5141   if (isconcave && maybeconvex) {
5142     zinc_(Ztwistedridge);
5143     if (maxdist > maxdist2)
5144       qh_appendmergeset(qh, facet, neighbor, MRGtwisted, maxdist, angle);
5145     else
5146       qh_appendmergeset(qh, neighbor, facet, MRGtwisted, maxdist2, angle);
5147     trace0((qh, 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(qh, facet, neighbor, MRGconcavecoplanar, maxdist, angle);
5153     else
5154       qh_appendmergeset(qh, neighbor, facet, MRGconcavecoplanar, maxdist2, angle);
5155     trace0((qh, 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(qh, facet, neighbor, MRGconcave, mergedist, angle);
5161     trace0((qh, 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(qh, facet, neighbor, MRGcoplanar, mergedist, angle);
5167     trace2((qh, 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_r.htm#TOC"
5174   >-------------------------------</a><a name="test_redundant_neighbors">-</a>
5175 
5176   qh_test_redundant_neighbors(qh, 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(qhT * qh,facetT * facet)5196 void qh_test_redundant_neighbors(qhT *qh, facetT *facet) {
5197   vertexT *vertex, **vertexp;
5198   facetT *neighbor, **neighborp;
5199   int size;
5200 
5201   trace4((qh, 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(qh, facet->neighbors)) < qh->hull_dim) {
5204     qh_appendmergeset(qh, facet, facet, MRGdegen, 0.0, 1.0);
5205     trace2((qh, 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, 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, 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(qh, neighbor, facet, MRGredundant, 0.0, 1.0);
5228         trace2((qh, 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_r.htm#TOC"
5235   >-------------------------------</a><a name="test_vneighbors">-</a>
5236 
5237   qh_test_vneighbors(qh)
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(qhT * qh)5260 boolT qh_test_vneighbors(qhT *qh /* qh.newfacet_list */) {
5261   facetT *newfacet, *neighbor, **neighborp;
5262   vertexT *vertex, **vertexp;
5263   int nummerges= 0;
5264 
5265   trace1((qh, qh->ferr, 1015, "qh_test_vneighbors: testing vertex neighbors for convexity\n"));
5266   if (!qh->VERTEXneighbors)
5267     qh_vertexneighbors(qh);
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(qh, newfacet, neighbor, False)) /* ignores optimization for simplicial ridges */
5280           nummerges++;
5281       }
5282     }
5283   }
5284   zadd_(Ztestvneighbor, nummerges);
5285   trace1((qh, 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_r.htm#TOC"
5291   >-------------------------------</a><a name="tracemerge">-</a>
5292 
5293   qh_tracemerge(qh, facet1, facet2 )
5294     print trace message after merge
5295 */
qh_tracemerge(qhT * qh,facetT * facet1,facetT * facet2,mergeType mergetype)5296 void qh_tracemerge(qhT *qh, 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(qh, "MERGED", facet2, NULL, NULL, NULL);
5307   if (facet2 == qh->tracefacet || (qh->tracevertex && qh->tracevertex->newfacet)) {
5308     qh_fprintf(qh, 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(qh, "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, qh->ferr, 8086, "qh_tracemerge: trace vertex deleted at furthest p%d\n",
5319             qh->furthest_id);
5320     else
5321       qh_checkvertex(qh, qh->tracevertex, qh_ALL, &waserror);
5322   }
5323   if (qh->tracefacet && qh->tracefacet->normal && !qh->tracefacet->visible)
5324     qh_checkfacet(qh, 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(qh);
5329     }
5330     qh_checkfacet(qh, facet2, True /* newmerge */, &waserror);
5331   }
5332   if (waserror)
5333     qh_errexit(qh, qh_ERRqhull, NULL, NULL); /* erroneous facet logged by qh_checkfacet */
5334 } /* tracemerge */
5335 
5336 /*-<a                             href="qh-merge_r.htm#TOC"
5337   >-------------------------------</a><a name="tracemerging">-</a>
5338 
5339   qh_tracemerging(qh)
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(qhT * qh)5351 void qh_tracemerging(qhT *qh) {
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, 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, qh->del_vertices));
5369 } /* tracemerging */
5370 
5371 /*-<a                             href="qh-merge_r.htm#TOC"
5372   >-------------------------------</a><a name="updatetested">-</a>
5373 
5374   qh_updatetested(qh, 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(qhT * qh,facetT * facet1,facetT * facet2)5396 void qh_updatetested(qhT *qh, 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(qh, 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(qh, 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_r.htm#TOC"
5425   >-------------------------------</a><a name="vertexridges">-</a>
5426 
5427   qh_vertexridges(qh, 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(qhT * qh,vertexT * vertex,boolT allneighbors)5443 setT *qh_vertexridges(qhT *qh, vertexT *vertex, boolT allneighbors) {
5444   facetT *neighbor, **neighborp;
5445   setT *ridges= qh_settemp(qh, 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(qh, vertex, neighbor, &ridges);
5454   }
5455   if (qh->PRINTstatistics || qh->IStracing) {
5456     size= qh_setsize(qh, ridges);
5457     zinc_(Zvertexridge);
5458     zadd_(Zvertexridgetot, size);
5459     zmax_(Zvertexridgemax, size);
5460     trace3((qh, 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_r.htm#TOC"
5467   >-------------------------------</a><a name="vertexridges_facet">-</a>
5468 
5469   qh_vertexridges_facet(qh, 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(qhT * qh,vertexT * vertex,facetT * facet,setT ** ridges)5484 void qh_vertexridges_facet(qhT *qh, 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(qh, 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(qh, ridges, ridge);
5501           else if (qh_setin(ridge->vertices, vertex))
5502             qh_setappend(qh, ridges, ridge);
5503         }
5504       }else if (SETelem_(ridge->vertices, last_i) == vertex
5505           || (last_i > 1 && SETsecond_(ridge->vertices) == vertex)) {
5506         qh_setappend(qh, ridges, ridge);
5507       }
5508     }
5509   }
5510   facet->visitid= qh->visit_id-1;
5511 } /* vertexridges_facet */
5512 
5513 /*-<a                             href="qh-merge_r.htm#TOC"
5514   >-------------------------------</a><a name="willdelete">-</a>
5515 
5516   qh_willdelete(qh, 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(qhT * qh,facetT * facet,facetT * replace)5524 void qh_willdelete(qhT *qh, facetT *facet, facetT *replace) {
5525 
5526   trace4((qh, 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, 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, qh_ERRqhull, NULL, NULL);
5531   }
5532   qh_removefacet(qh, facet);
5533   qh_prependfacet(qh, 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(qhT * qh,int apexpointid,facetT * facet,facetT ** retryfacet)5545 void qh_all_vertexmerges(qhT *qh, int apexpointid, facetT *facet, facetT **retryfacet) {
5546   QHULL_UNUSED(qh)
5547   QHULL_UNUSED(apexpointid)
5548   QHULL_UNUSED(facet)
5549   QHULL_UNUSED(retryfacet)
5550 }
qh_premerge(qhT * qh,int apexpointid,realT maxcentrum,realT maxangle)5551 void qh_premerge(qhT *qh, int apexpointid, realT maxcentrum, realT maxangle) {
5552   QHULL_UNUSED(qh)
5553   QHULL_UNUSED(apexpointid)
5554   QHULL_UNUSED(maxcentrum)
5555   QHULL_UNUSED(maxangle)
5556 }
qh_postmerge(qhT * qh,const char * reason,realT maxcentrum,realT maxangle,boolT vneighbors)5557 void qh_postmerge(qhT *qh, const char *reason, realT maxcentrum, realT maxangle,
5558                       boolT vneighbors) {
5559   QHULL_UNUSED(qh)
5560   QHULL_UNUSED(reason)
5561   QHULL_UNUSED(maxcentrum)
5562   QHULL_UNUSED(maxangle)
5563   QHULL_UNUSED(vneighbors)
5564 }
qh_checkdelfacet(qhT * qh,facetT * facet,setT * mergeset)5565 void qh_checkdelfacet(qhT *qh, facetT *facet, setT *mergeset) {
5566   QHULL_UNUSED(qh)
5567   QHULL_UNUSED(facet)
5568   QHULL_UNUSED(mergeset)
5569 }
qh_checkdelridge(qhT * qh)5570 void qh_checkdelridge(qhT *qh /* qh.visible_facets, vertex_mergeset */) {
5571   QHULL_UNUSED(qh)
5572 }
qh_checkzero(qhT * qh,boolT testall)5573 boolT qh_checkzero(qhT *qh, boolT testall) {
5574   QHULL_UNUSED(qh)
5575   QHULL_UNUSED(testall)
5576 
5577   return True;
5578 }
qh_freemergesets(qhT * qh)5579 void qh_freemergesets(qhT *qh) {
5580   QHULL_UNUSED(qh)
5581 }
qh_initmergesets(qhT * qh)5582 void qh_initmergesets(qhT *qh) {
5583   QHULL_UNUSED(qh)
5584 }
qh_merge_pinchedvertices(qhT * qh,int apexpointid)5585 void qh_merge_pinchedvertices(qhT *qh, int apexpointid /* qh.newfacet_list */) {
5586   QHULL_UNUSED(qh)
5587   QHULL_UNUSED(apexpointid)
5588 }
5589 #endif /* qh_NOmerge */
5590 
5591