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-2015 C.B. Barber.
24    $Id: //main/2015/qhull/src/libqhull_r/merge_r.c#5 $$Change: 2064 $
25    $DateTime: 2016/01/18 12:36:08 $$Author: bbarber $
26 */
27 
28 #include "qhull_ra.h"
29 
30 #ifndef qh_NOmerge
31 
32 /*===== functions(alphabetical after premerge and postmerge) ======*/
33 
34 /*-<a                             href="qh-merge_r.htm#TOC"
35   >-------------------------------</a><a name="premerge">-</a>
36 
37   qh_premerge(qh, apex, maxcentrum )
38     pre-merge nonconvex facets in qh.newfacet_list for apex
39     maxcentrum defines coplanar and concave (qh_test_appendmerge)
40 
41   returns:
42     deleted facets added to qh.visible_list with facet->visible set
43 
44   notes:
45     uses globals, qh.MERGEexact, qh.PREmerge
46 
47   design:
48     mark duplicate ridges in qh.newfacet_list
49     merge facet cycles in qh.newfacet_list
50     merge duplicate ridges and concave facets in qh.newfacet_list
51     check merged facet cycles for degenerate and redundant facets
52     merge degenerate and redundant facets
53     collect coplanar and concave facets
54     merge concave, coplanar, degenerate, and redundant facets
55 */
qh_premerge(qhT * qh,vertexT * apex,realT maxcentrum,realT maxangle)56 void qh_premerge(qhT *qh, vertexT *apex, realT maxcentrum, realT maxangle) {
57   boolT othermerge= False;
58   facetT *newfacet;
59 
60   if (qh->ZEROcentrum && qh_checkzero(qh, !qh_ALL))
61     return;
62   trace2((qh, qh->ferr, 2008, "qh_premerge: premerge centrum %2.2g angle %2.2g for apex v%d facetlist f%d\n",
63             maxcentrum, maxangle, apex->id, getid_(qh->newfacet_list)));
64   if (qh->IStracing >= 4 && qh->num_facets < 50)
65     qh_printlists(qh);
66   qh->centrum_radius= maxcentrum;
67   qh->cos_max= maxangle;
68   qh->degen_mergeset= qh_settemp(qh, qh->TEMPsize);
69   qh->facet_mergeset= qh_settemp(qh, qh->TEMPsize);
70   if (qh->hull_dim >=3) {
71     qh_mark_dupridges(qh, qh->newfacet_list); /* facet_mergeset */
72     qh_mergecycle_all(qh, qh->newfacet_list, &othermerge);
73     qh_forcedmerges(qh, &othermerge /* qh->facet_mergeset */);
74     FORALLnew_facets {  /* test samecycle merges */
75       if (!newfacet->simplicial && !newfacet->mergeridge)
76         qh_degen_redundant_neighbors(qh, newfacet, NULL);
77     }
78     if (qh_merge_degenredundant(qh))
79       othermerge= True;
80   }else /* qh->hull_dim == 2 */
81     qh_mergecycle_all(qh, qh->newfacet_list, &othermerge);
82   qh_flippedmerges(qh, qh->newfacet_list, &othermerge);
83   if (!qh->MERGEexact || zzval_(Ztotmerge)) {
84     zinc_(Zpremergetot);
85     qh->POSTmerging= False;
86     qh_getmergeset_initial(qh, qh->newfacet_list);
87     qh_all_merges(qh, othermerge, False);
88   }
89   qh_settempfree(qh, &qh->facet_mergeset);
90   qh_settempfree(qh, &qh->degen_mergeset);
91 } /* premerge */
92 
93 /*-<a                             href="qh-merge_r.htm#TOC"
94   >-------------------------------</a><a name="postmerge">-</a>
95 
96   qh_postmerge(qh, reason, maxcentrum, maxangle, vneighbors )
97     post-merge nonconvex facets as defined by maxcentrum and maxangle
98     'reason' is for reporting progress
99     if vneighbors,
100       calls qh_test_vneighbors at end of qh_all_merge
101     if firstmerge,
102       calls qh_reducevertices before qh_getmergeset
103 
104   returns:
105     if first call (qh.visible_list != qh.facet_list),
106       builds qh.facet_newlist, qh.newvertex_list
107     deleted facets added to qh.visible_list with facet->visible
108     qh.visible_list == qh.facet_list
109 
110   notes:
111 
112 
113   design:
114     if first call
115       set qh.visible_list and qh.newfacet_list to qh.facet_list
116       add all facets to qh.newfacet_list
117       mark non-simplicial facets, facet->newmerge
118       set qh.newvertext_list to qh.vertex_list
119       add all vertices to qh.newvertex_list
120       if a pre-merge occured
121         set vertex->delridge {will retest the ridge}
122         if qh.MERGEexact
123           call qh_reducevertices()
124       if no pre-merging
125         merge flipped facets
126     determine non-convex facets
127     merge all non-convex facets
128 */
qh_postmerge(qhT * qh,const char * reason,realT maxcentrum,realT maxangle,boolT vneighbors)129 void qh_postmerge(qhT *qh, const char *reason, realT maxcentrum, realT maxangle,
130                       boolT vneighbors) {
131   facetT *newfacet;
132   boolT othermerges= False;
133   vertexT *vertex;
134 
135   if (qh->REPORTfreq || qh->IStracing) {
136     qh_buildtracing(qh, NULL, NULL);
137     qh_printsummary(qh, qh->ferr);
138     if (qh->PRINTstatistics)
139       qh_printallstatistics(qh, qh->ferr, "reason");
140     qh_fprintf(qh, qh->ferr, 8062, "\n%s with 'C%.2g' and 'A%.2g'\n",
141         reason, maxcentrum, maxangle);
142   }
143   trace2((qh, qh->ferr, 2009, "qh_postmerge: postmerge.  test vneighbors? %d\n",
144             vneighbors));
145   qh->centrum_radius= maxcentrum;
146   qh->cos_max= maxangle;
147   qh->POSTmerging= True;
148   qh->degen_mergeset= qh_settemp(qh, qh->TEMPsize);
149   qh->facet_mergeset= qh_settemp(qh, qh->TEMPsize);
150   if (qh->visible_list != qh->facet_list) {  /* first call */
151     qh->NEWfacets= True;
152     qh->visible_list= qh->newfacet_list= qh->facet_list;
153     FORALLnew_facets {
154       newfacet->newfacet= True;
155        if (!newfacet->simplicial)
156         newfacet->newmerge= True;
157      zinc_(Zpostfacets);
158     }
159     qh->newvertex_list= qh->vertex_list;
160     FORALLvertices
161       vertex->newlist= True;
162     if (qh->VERTEXneighbors) { /* a merge has occurred */
163       FORALLvertices
164         vertex->delridge= True; /* test for redundant, needed? */
165       if (qh->MERGEexact) {
166         if (qh->hull_dim <= qh_DIMreduceBuild)
167           qh_reducevertices(qh); /* was skipped during pre-merging */
168       }
169     }
170     if (!qh->PREmerge && !qh->MERGEexact)
171       qh_flippedmerges(qh, qh->newfacet_list, &othermerges);
172   }
173   qh_getmergeset_initial(qh, qh->newfacet_list);
174   qh_all_merges(qh, False, vneighbors);
175   qh_settempfree(qh, &qh->facet_mergeset);
176   qh_settempfree(qh, &qh->degen_mergeset);
177 } /* post_merge */
178 
179 /*-<a                             href="qh-merge_r.htm#TOC"
180   >-------------------------------</a><a name="all_merges">-</a>
181 
182   qh_all_merges(qh, othermerge, vneighbors )
183     merge all non-convex facets
184 
185     set othermerge if already merged facets (for qh_reducevertices)
186     if vneighbors
187       tests vertex neighbors for convexity at end
188     qh.facet_mergeset lists the non-convex ridges in qh_newfacet_list
189     qh.degen_mergeset is defined
190     if qh.MERGEexact && !qh.POSTmerging,
191       does not merge coplanar facets
192 
193   returns:
194     deleted facets added to qh.visible_list with facet->visible
195     deleted vertices added qh.delvertex_list with vertex->delvertex
196 
197   notes:
198     unless !qh.MERGEindependent,
199       merges facets in independent sets
200     uses qh.newfacet_list as argument since merges call qh_removefacet()
201 
202   design:
203     while merges occur
204       for each merge in qh.facet_mergeset
205         unless one of the facets was already merged in this pass
206           merge the facets
207         test merged facets for additional merges
208         add merges to qh.facet_mergeset
209       if vertices record neighboring facets
210         rename redundant vertices
211           update qh.facet_mergeset
212     if vneighbors ??
213       tests vertex neighbors for convexity at end
214 */
qh_all_merges(qhT * qh,boolT othermerge,boolT vneighbors)215 void qh_all_merges(qhT *qh, boolT othermerge, boolT vneighbors) {
216   facetT *facet1, *facet2;
217   mergeT *merge;
218   boolT wasmerge= True, isreduce;
219   void **freelistp;  /* used if !qh_NOmem by qh_memfree_() */
220   vertexT *vertex;
221   mergeType mergetype;
222   int numcoplanar=0, numconcave=0, numdegenredun= 0, numnewmerges= 0;
223 
224   trace2((qh, qh->ferr, 2010, "qh_all_merges: starting to merge facets beginning from f%d\n",
225             getid_(qh->newfacet_list)));
226   while (True) {
227     wasmerge= False;
228     while (qh_setsize(qh, qh->facet_mergeset)) {
229       while ((merge= (mergeT*)qh_setdellast(qh->facet_mergeset))) {
230         facet1= merge->facet1;
231         facet2= merge->facet2;
232         mergetype= merge->type;
233         qh_memfree_(qh, merge, (int)sizeof(mergeT), freelistp);
234         if (facet1->visible || facet2->visible) /*deleted facet*/
235           continue;
236         if ((facet1->newfacet && !facet1->tested)
237                 || (facet2->newfacet && !facet2->tested)) {
238           if (qh->MERGEindependent && mergetype <= MRGanglecoplanar)
239             continue;      /* perform independent sets of merges */
240         }
241         qh_merge_nonconvex(qh, facet1, facet2, mergetype);
242         numdegenredun += qh_merge_degenredundant(qh);
243         numnewmerges++;
244         wasmerge= True;
245         if (mergetype == MRGconcave)
246           numconcave++;
247         else /* MRGcoplanar or MRGanglecoplanar */
248           numcoplanar++;
249       } /* while setdellast */
250       if (qh->POSTmerging && qh->hull_dim <= qh_DIMreduceBuild
251       && numnewmerges > qh_MAXnewmerges) {
252         numnewmerges= 0;
253         qh_reducevertices(qh);  /* otherwise large post merges too slow */
254       }
255       qh_getmergeset(qh, qh->newfacet_list); /* facet_mergeset */
256     } /* while mergeset */
257     if (qh->VERTEXneighbors) {
258       isreduce= False;
259       if (qh->hull_dim >=4 && qh->POSTmerging) {
260         FORALLvertices
261           vertex->delridge= True;
262         isreduce= True;
263       }
264       if ((wasmerge || othermerge) && (!qh->MERGEexact || qh->POSTmerging)
265           && qh->hull_dim <= qh_DIMreduceBuild) {
266         othermerge= False;
267         isreduce= True;
268       }
269       if (isreduce) {
270         if (qh_reducevertices(qh)) {
271           qh_getmergeset(qh, qh->newfacet_list); /* facet_mergeset */
272           continue;
273         }
274       }
275     }
276     if (vneighbors && qh_test_vneighbors(qh /* qh->newfacet_list */))
277       continue;
278     break;
279   } /* while (True) */
280   if (qh->CHECKfrequently && !qh->MERGEexact) {
281     qh->old_randomdist= qh->RANDOMdist;
282     qh->RANDOMdist= False;
283     qh_checkconvex(qh, qh->newfacet_list, qh_ALGORITHMfault);
284     /* qh_checkconnect(qh); [this is slow and it changes the facet order] */
285     qh->RANDOMdist= qh->old_randomdist;
286   }
287   trace1((qh, qh->ferr, 1009, "qh_all_merges: merged %d coplanar facets %d concave facets and %d degen or redundant facets.\n",
288     numcoplanar, numconcave, numdegenredun));
289   if (qh->IStracing >= 4 && qh->num_facets < 50)
290     qh_printlists(qh);
291 } /* all_merges */
292 
293 
294 /*-<a                             href="qh-merge_r.htm#TOC"
295   >-------------------------------</a><a name="appendmergeset">-</a>
296 
297   qh_appendmergeset(qh, facet, neighbor, mergetype, angle )
298     appends an entry to qh.facet_mergeset or qh.degen_mergeset
299 
300     angle ignored if NULL or !qh.ANGLEmerge
301 
302   returns:
303     merge appended to facet_mergeset or degen_mergeset
304       sets ->degenerate or ->redundant if degen_mergeset
305 
306   see:
307     qh_test_appendmerge()
308 
309   design:
310     allocate merge entry
311     if regular merge
312       append to qh.facet_mergeset
313     else if degenerate merge and qh.facet_mergeset is all degenerate
314       append to qh.degen_mergeset
315     else if degenerate merge
316       prepend to qh.degen_mergeset
317     else if redundant merge
318       append to qh.degen_mergeset
319 */
qh_appendmergeset(qhT * qh,facetT * facet,facetT * neighbor,mergeType mergetype,realT * angle)320 void qh_appendmergeset(qhT *qh, facetT *facet, facetT *neighbor, mergeType mergetype, realT *angle) {
321   mergeT *merge, *lastmerge;
322   void **freelistp; /* used if !qh_NOmem by qh_memalloc_() */
323 
324   if (facet->redundant)
325     return;
326   if (facet->degenerate && mergetype == MRGdegen)
327     return;
328   qh_memalloc_(qh, (int)sizeof(mergeT), freelistp, merge, mergeT);
329   merge->facet1= facet;
330   merge->facet2= neighbor;
331   merge->type= mergetype;
332   if (angle && qh->ANGLEmerge)
333     merge->angle= *angle;
334   if (mergetype < MRGdegen)
335     qh_setappend(qh, &(qh->facet_mergeset), merge);
336   else if (mergetype == MRGdegen) {
337     facet->degenerate= True;
338     if (!(lastmerge= (mergeT*)qh_setlast(qh->degen_mergeset))
339     || lastmerge->type == MRGdegen)
340       qh_setappend(qh, &(qh->degen_mergeset), merge);
341     else
342       qh_setaddnth(qh, &(qh->degen_mergeset), 0, merge);
343   }else if (mergetype == MRGredundant) {
344     facet->redundant= True;
345     qh_setappend(qh, &(qh->degen_mergeset), merge);
346   }else /* mergetype == MRGmirror */ {
347     if (facet->redundant || neighbor->redundant) {
348       qh_fprintf(qh, qh->ferr, 6092, "qhull error (qh_appendmergeset): facet f%d or f%d is already a mirrored facet\n",
349            facet->id, neighbor->id);
350       qh_errexit2(qh, qh_ERRqhull, facet, neighbor);
351     }
352     if (!qh_setequal(facet->vertices, neighbor->vertices)) {
353       qh_fprintf(qh, qh->ferr, 6093, "qhull error (qh_appendmergeset): mirrored facets f%d and f%d do not have the same vertices\n",
354            facet->id, neighbor->id);
355       qh_errexit2(qh, qh_ERRqhull, facet, neighbor);
356     }
357     facet->redundant= True;
358     neighbor->redundant= True;
359     qh_setappend(qh, &(qh->degen_mergeset), merge);
360   }
361 } /* appendmergeset */
362 
363 
364 /*-<a                             href="qh-merge_r.htm#TOC"
365   >-------------------------------</a><a name="basevertices">-</a>
366 
367   qh_basevertices(qh, samecycle )
368     return temporary set of base vertices for samecycle
369     samecycle is first facet in the cycle
370     assumes apex is SETfirst_( samecycle->vertices )
371 
372   returns:
373     vertices(settemp)
374     all ->seen are cleared
375 
376   notes:
377     uses qh_vertex_visit;
378 
379   design:
380     for each facet in samecycle
381       for each unseen vertex in facet->vertices
382         append to result
383 */
qh_basevertices(qhT * qh,facetT * samecycle)384 setT *qh_basevertices(qhT *qh, facetT *samecycle) {
385   facetT *same;
386   vertexT *apex, *vertex, **vertexp;
387   setT *vertices= qh_settemp(qh, qh->TEMPsize);
388 
389   apex= SETfirstt_(samecycle->vertices, vertexT);
390   apex->visitid= ++qh->vertex_visit;
391   FORALLsame_cycle_(samecycle) {
392     if (same->mergeridge)
393       continue;
394     FOREACHvertex_(same->vertices) {
395       if (vertex->visitid != qh->vertex_visit) {
396         qh_setappend(qh, &vertices, vertex);
397         vertex->visitid= qh->vertex_visit;
398         vertex->seen= False;
399       }
400     }
401   }
402   trace4((qh, qh->ferr, 4019, "qh_basevertices: found %d vertices\n",
403          qh_setsize(qh, vertices)));
404   return vertices;
405 } /* basevertices */
406 
407 /*-<a                             href="qh-merge_r.htm#TOC"
408   >-------------------------------</a><a name="checkconnect">-</a>
409 
410   qh_checkconnect(qh)
411     check that new facets are connected
412     new facets are on qh.newfacet_list
413 
414   notes:
415     this is slow and it changes the order of the facets
416     uses qh.visit_id
417 
418   design:
419     move first new facet to end of qh.facet_list
420     for all newly appended facets
421       append unvisited neighbors to end of qh.facet_list
422     for all new facets
423       report error if unvisited
424 */
qh_checkconnect(qhT * qh)425 void qh_checkconnect(qhT *qh /* qh->newfacet_list */) {
426   facetT *facet, *newfacet, *errfacet= NULL, *neighbor, **neighborp;
427 
428   facet= qh->newfacet_list;
429   qh_removefacet(qh, facet);
430   qh_appendfacet(qh, facet);
431   facet->visitid= ++qh->visit_id;
432   FORALLfacet_(facet) {
433     FOREACHneighbor_(facet) {
434       if (neighbor->visitid != qh->visit_id) {
435         qh_removefacet(qh, neighbor);
436         qh_appendfacet(qh, neighbor);
437         neighbor->visitid= qh->visit_id;
438       }
439     }
440   }
441   FORALLnew_facets {
442     if (newfacet->visitid == qh->visit_id)
443       break;
444     qh_fprintf(qh, qh->ferr, 6094, "qhull error: f%d is not attached to the new facets\n",
445          newfacet->id);
446     errfacet= newfacet;
447   }
448   if (errfacet)
449     qh_errexit(qh, qh_ERRqhull, errfacet, NULL);
450 } /* checkconnect */
451 
452 /*-<a                             href="qh-merge_r.htm#TOC"
453   >-------------------------------</a><a name="checkzero">-</a>
454 
455   qh_checkzero(qh, testall )
456     check that facets are clearly convex for qh.DISTround with qh.MERGEexact
457 
458     if testall,
459       test all facets for qh.MERGEexact post-merging
460     else
461       test qh.newfacet_list
462 
463     if qh.MERGEexact,
464       allows coplanar ridges
465       skips convexity test while qh.ZEROall_ok
466 
467   returns:
468     True if all facets !flipped, !dupridge, normal
469          if all horizon facets are simplicial
470          if all vertices are clearly below neighbor
471          if all opposite vertices of horizon are below
472     clears qh.ZEROall_ok if any problems or coplanar facets
473 
474   notes:
475     uses qh.vertex_visit
476     horizon facets may define multiple new facets
477 
478   design:
479     for all facets in qh.newfacet_list or qh.facet_list
480       check for flagged faults (flipped, etc.)
481     for all facets in qh.newfacet_list or qh.facet_list
482       for each neighbor of facet
483         skip horizon facets for qh.newfacet_list
484         test the opposite vertex
485       if qh.newfacet_list
486         test the other vertices in the facet's horizon facet
487 */
qh_checkzero(qhT * qh,boolT testall)488 boolT qh_checkzero(qhT *qh, boolT testall) {
489   facetT *facet, *neighbor, **neighborp;
490   facetT *horizon, *facetlist;
491   int neighbor_i;
492   vertexT *vertex, **vertexp;
493   realT dist;
494 
495   if (testall)
496     facetlist= qh->facet_list;
497   else {
498     facetlist= qh->newfacet_list;
499     FORALLfacet_(facetlist) {
500       horizon= SETfirstt_(facet->neighbors, facetT);
501       if (!horizon->simplicial)
502         goto LABELproblem;
503       if (facet->flipped || facet->dupridge || !facet->normal)
504         goto LABELproblem;
505     }
506     if (qh->MERGEexact && qh->ZEROall_ok) {
507       trace2((qh, qh->ferr, 2011, "qh_checkzero: skip convexity check until first pre-merge\n"));
508       return True;
509     }
510   }
511   FORALLfacet_(facetlist) {
512     qh->vertex_visit++;
513     neighbor_i= 0;
514     horizon= NULL;
515     FOREACHneighbor_(facet) {
516       if (!neighbor_i && !testall) {
517         horizon= neighbor;
518         neighbor_i++;
519         continue; /* horizon facet tested in qh_findhorizon */
520       }
521       vertex= SETelemt_(facet->vertices, neighbor_i++, vertexT);
522       vertex->visitid= qh->vertex_visit;
523       zzinc_(Zdistzero);
524       qh_distplane(qh, vertex->point, neighbor, &dist);
525       if (dist >= -qh->DISTround) {
526         qh->ZEROall_ok= False;
527         if (!qh->MERGEexact || testall || dist > qh->DISTround)
528           goto LABELnonconvex;
529       }
530     }
531     if (!testall && horizon) {
532       FOREACHvertex_(horizon->vertices) {
533         if (vertex->visitid != qh->vertex_visit) {
534           zzinc_(Zdistzero);
535           qh_distplane(qh, vertex->point, facet, &dist);
536           if (dist >= -qh->DISTround) {
537             qh->ZEROall_ok= False;
538             if (!qh->MERGEexact || dist > qh->DISTround)
539               goto LABELnonconvex;
540           }
541           break;
542         }
543       }
544     }
545   }
546   trace2((qh, qh->ferr, 2012, "qh_checkzero: testall %d, facets are %s\n", testall,
547         (qh->MERGEexact && !testall) ?
548            "not concave, flipped, or duplicate ridged" : "clearly convex"));
549   return True;
550 
551  LABELproblem:
552   qh->ZEROall_ok= False;
553   trace2((qh, qh->ferr, 2013, "qh_checkzero: facet f%d needs pre-merging\n",
554        facet->id));
555   return False;
556 
557  LABELnonconvex:
558   trace2((qh, qh->ferr, 2014, "qh_checkzero: facet f%d and f%d are not clearly convex.  v%d dist %.2g\n",
559          facet->id, neighbor->id, vertex->id, dist));
560   return False;
561 } /* checkzero */
562 
563 /*-<a                             href="qh-merge_r.htm#TOC"
564   >-------------------------------</a><a name="compareangle">-</a>
565 
566   qh_compareangle(angle1, angle2 )
567     used by qsort() to order merges by angle
568 */
qh_compareangle(const void * p1,const void * p2)569 int qh_compareangle(const void *p1, const void *p2) {
570   const mergeT *a= *((mergeT *const*)p1), *b= *((mergeT *const*)p2);
571 
572   return((a->angle > b->angle) ? 1 : -1);
573 } /* compareangle */
574 
575 /*-<a                             href="qh-merge_r.htm#TOC"
576   >-------------------------------</a><a name="comparemerge">-</a>
577 
578   qh_comparemerge(merge1, merge2 )
579     used by qsort() to order merges
580 */
qh_comparemerge(const void * p1,const void * p2)581 int qh_comparemerge(const void *p1, const void *p2) {
582   const mergeT *a= *((mergeT *const*)p1), *b= *((mergeT *const*)p2);
583 
584   return(a->type - b->type);
585 } /* comparemerge */
586 
587 /*-<a                             href="qh-merge_r.htm#TOC"
588   >-------------------------------</a><a name="comparevisit">-</a>
589 
590   qh_comparevisit(vertex1, vertex2 )
591     used by qsort() to order vertices by their visitid
592 */
qh_comparevisit(const void * p1,const void * p2)593 int qh_comparevisit(const void *p1, const void *p2) {
594   const vertexT *a= *((vertexT *const*)p1), *b= *((vertexT *const*)p2);
595 
596   return(a->visitid - b->visitid);
597 } /* comparevisit */
598 
599 /*-<a                             href="qh-merge_r.htm#TOC"
600   >-------------------------------</a><a name="copynonconvex">-</a>
601 
602   qh_copynonconvex(qh, atridge )
603     set non-convex flag on other ridges (if any) between same neighbors
604 
605   notes:
606     may be faster if use smaller ridge set
607 
608   design:
609     for each ridge of atridge's top facet
610       if ridge shares the same neighbor
611         set nonconvex flag
612 */
qh_copynonconvex(qhT * qh,ridgeT * atridge)613 void qh_copynonconvex(qhT *qh, ridgeT *atridge) {
614   facetT *facet, *otherfacet;
615   ridgeT *ridge, **ridgep;
616 
617   facet= atridge->top;
618   otherfacet= atridge->bottom;
619   FOREACHridge_(facet->ridges) {
620     if (otherfacet == otherfacet_(ridge, facet) && ridge != atridge) {
621       ridge->nonconvex= True;
622       trace4((qh, qh->ferr, 4020, "qh_copynonconvex: moved nonconvex flag from r%d to r%d\n",
623               atridge->id, ridge->id));
624       break;
625     }
626   }
627 } /* copynonconvex */
628 
629 /*-<a                             href="qh-merge_r.htm#TOC"
630   >-------------------------------</a><a name="degen_redundant_facet">-</a>
631 
632   qh_degen_redundant_facet(qh, facet )
633     check facet for degen. or redundancy
634 
635   notes:
636     bumps vertex_visit
637     called if a facet was redundant but no longer is (qh_merge_degenredundant)
638     qh_appendmergeset() only appends first reference to facet (i.e., redundant)
639 
640   see:
641     qh_degen_redundant_neighbors()
642 
643   design:
644     test for redundant neighbor
645     test for degenerate facet
646 */
qh_degen_redundant_facet(qhT * qh,facetT * facet)647 void qh_degen_redundant_facet(qhT *qh, facetT *facet) {
648   vertexT *vertex, **vertexp;
649   facetT *neighbor, **neighborp;
650 
651   trace4((qh, qh->ferr, 4021, "qh_degen_redundant_facet: test facet f%d for degen/redundant\n",
652           facet->id));
653   FOREACHneighbor_(facet) {
654     qh->vertex_visit++;
655     FOREACHvertex_(neighbor->vertices)
656       vertex->visitid= qh->vertex_visit;
657     FOREACHvertex_(facet->vertices) {
658       if (vertex->visitid != qh->vertex_visit)
659         break;
660     }
661     if (!vertex) {
662       qh_appendmergeset(qh, facet, neighbor, MRGredundant, NULL);
663       trace2((qh, qh->ferr, 2015, "qh_degen_redundant_facet: f%d is contained in f%d.  merge\n", facet->id, neighbor->id));
664       return;
665     }
666   }
667   if (qh_setsize(qh, facet->neighbors) < qh->hull_dim) {
668     qh_appendmergeset(qh, facet, facet, MRGdegen, NULL);
669     trace2((qh, qh->ferr, 2016, "qh_degen_redundant_neighbors: f%d is degenerate.\n", facet->id));
670   }
671 } /* degen_redundant_facet */
672 
673 
674 /*-<a                             href="qh-merge_r.htm#TOC"
675   >-------------------------------</a><a name="degen_redundant_neighbors">-</a>
676 
677   qh_degen_redundant_neighbors(qh, facet, delfacet,  )
678     append degenerate and redundant neighbors to facet_mergeset
679     if delfacet,
680       only checks neighbors of both delfacet and facet
681     also checks current facet for degeneracy
682 
683   notes:
684     bumps vertex_visit
685     called for each qh_mergefacet() and qh_mergecycle()
686     merge and statistics occur in merge_nonconvex
687     qh_appendmergeset() only appends first reference to facet (i.e., redundant)
688       it appends redundant facets after degenerate ones
689 
690     a degenerate facet has fewer than hull_dim neighbors
691     a redundant facet's vertices is a subset of its neighbor's vertices
692     tests for redundant merges first (appendmergeset is nop for others)
693     in a merge, only needs to test neighbors of merged facet
694 
695   see:
696     qh_merge_degenredundant() and qh_degen_redundant_facet()
697 
698   design:
699     test for degenerate facet
700     test for redundant neighbor
701     test for degenerate neighbor
702 */
qh_degen_redundant_neighbors(qhT * qh,facetT * facet,facetT * delfacet)703 void qh_degen_redundant_neighbors(qhT *qh, facetT *facet, facetT *delfacet) {
704   vertexT *vertex, **vertexp;
705   facetT *neighbor, **neighborp;
706   int size;
707 
708   trace4((qh, qh->ferr, 4022, "qh_degen_redundant_neighbors: test neighbors of f%d with delfacet f%d\n",
709           facet->id, getid_(delfacet)));
710   if ((size= qh_setsize(qh, facet->neighbors)) < qh->hull_dim) {
711     qh_appendmergeset(qh, facet, facet, MRGdegen, NULL);
712     trace2((qh, qh->ferr, 2017, "qh_degen_redundant_neighbors: f%d is degenerate with %d neighbors.\n", facet->id, size));
713   }
714   if (!delfacet)
715     delfacet= facet;
716   qh->vertex_visit++;
717   FOREACHvertex_(facet->vertices)
718     vertex->visitid= qh->vertex_visit;
719   FOREACHneighbor_(delfacet) {
720     /* uses early out instead of checking vertex count */
721     if (neighbor == facet)
722       continue;
723     FOREACHvertex_(neighbor->vertices) {
724       if (vertex->visitid != qh->vertex_visit)
725         break;
726     }
727     if (!vertex) {
728       qh_appendmergeset(qh, neighbor, facet, MRGredundant, NULL);
729       trace2((qh, qh->ferr, 2018, "qh_degen_redundant_neighbors: f%d is contained in f%d.  merge\n", neighbor->id, facet->id));
730     }
731   }
732   FOREACHneighbor_(delfacet) {   /* redundant merges occur first */
733     if (neighbor == facet)
734       continue;
735     if ((size= qh_setsize(qh, neighbor->neighbors)) < qh->hull_dim) {
736       qh_appendmergeset(qh, neighbor, neighbor, MRGdegen, NULL);
737       trace2((qh, qh->ferr, 2019, "qh_degen_redundant_neighbors: f%d is degenerate with %d neighbors.  Neighbor of f%d.\n", neighbor->id, size, facet->id));
738     }
739   }
740 } /* degen_redundant_neighbors */
741 
742 
743 /*-<a                             href="qh-merge_r.htm#TOC"
744   >-------------------------------</a><a name="find_newvertex">-</a>
745 
746   qh_find_newvertex(qh, oldvertex, vertices, ridges )
747     locate new vertex for renaming old vertex
748     vertices is a set of possible new vertices
749       vertices sorted by number of deleted ridges
750 
751   returns:
752     newvertex or NULL
753       each ridge includes both vertex and oldvertex
754     vertices sorted by number of deleted ridges
755 
756   notes:
757     modifies vertex->visitid
758     new vertex is in one of the ridges
759     renaming will not cause a duplicate ridge
760     renaming will minimize the number of deleted ridges
761     newvertex may not be adjacent in the dual (though unlikely)
762 
763   design:
764     for each vertex in vertices
765       set vertex->visitid to number of references in ridges
766     remove unvisited vertices
767     set qh.vertex_visit above all possible values
768     sort vertices by number of references in ridges
769     add each ridge to qh.hash_table
770     for each vertex in vertices
771       look for a vertex that would not cause a duplicate ridge after a rename
772 */
qh_find_newvertex(qhT * qh,vertexT * oldvertex,setT * vertices,setT * ridges)773 vertexT *qh_find_newvertex(qhT *qh, vertexT *oldvertex, setT *vertices, setT *ridges) {
774   vertexT *vertex, **vertexp;
775   setT *newridges;
776   ridgeT *ridge, **ridgep;
777   int size, hashsize;
778   int hash;
779 
780 #ifndef qh_NOtrace
781   if (qh->IStracing >= 4) {
782     qh_fprintf(qh, qh->ferr, 8063, "qh_find_newvertex: find new vertex for v%d from ",
783              oldvertex->id);
784     FOREACHvertex_(vertices)
785       qh_fprintf(qh, qh->ferr, 8064, "v%d ", vertex->id);
786     FOREACHridge_(ridges)
787       qh_fprintf(qh, qh->ferr, 8065, "r%d ", ridge->id);
788     qh_fprintf(qh, qh->ferr, 8066, "\n");
789   }
790 #endif
791   FOREACHvertex_(vertices)
792     vertex->visitid= 0;
793   FOREACHridge_(ridges) {
794     FOREACHvertex_(ridge->vertices)
795       vertex->visitid++;
796   }
797   FOREACHvertex_(vertices) {
798     if (!vertex->visitid) {
799       qh_setdelnth(qh, vertices, SETindex_(vertices,vertex));
800       vertexp--; /* repeat since deleted this vertex */
801     }
802   }
803   qh->vertex_visit += (unsigned int)qh_setsize(qh, ridges);
804   if (!qh_setsize(qh, vertices)) {
805     trace4((qh, qh->ferr, 4023, "qh_find_newvertex: vertices not in ridges for v%d\n",
806             oldvertex->id));
807     return NULL;
808   }
809   qsort(SETaddr_(vertices, vertexT), (size_t)qh_setsize(qh, vertices),
810                 sizeof(vertexT *), qh_comparevisit);
811   /* can now use qh->vertex_visit */
812   if (qh->PRINTstatistics) {
813     size= qh_setsize(qh, vertices);
814     zinc_(Zintersect);
815     zadd_(Zintersecttot, size);
816     zmax_(Zintersectmax, size);
817   }
818   hashsize= qh_newhashtable(qh, qh_setsize(qh, ridges));
819   FOREACHridge_(ridges)
820     qh_hashridge(qh, qh->hash_table, hashsize, ridge, oldvertex);
821   FOREACHvertex_(vertices) {
822     newridges= qh_vertexridges(qh, vertex);
823     FOREACHridge_(newridges) {
824       if (qh_hashridge_find(qh, qh->hash_table, hashsize, ridge, vertex, oldvertex, &hash)) {
825         zinc_(Zdupridge);
826         break;
827       }
828     }
829     qh_settempfree(qh, &newridges);
830     if (!ridge)
831       break;  /* found a rename */
832   }
833   if (vertex) {
834     /* counted in qh_renamevertex */
835     trace2((qh, qh->ferr, 2020, "qh_find_newvertex: found v%d for old v%d from %d vertices and %d ridges.\n",
836       vertex->id, oldvertex->id, qh_setsize(qh, vertices), qh_setsize(qh, ridges)));
837   }else {
838     zinc_(Zfindfail);
839     trace0((qh, qh->ferr, 14, "qh_find_newvertex: no vertex for renaming v%d(all duplicated ridges) during p%d\n",
840       oldvertex->id, qh->furthest_id));
841   }
842   qh_setfree(qh, &qh->hash_table);
843   return vertex;
844 } /* find_newvertex */
845 
846 /*-<a                             href="qh-merge_r.htm#TOC"
847   >-------------------------------</a><a name="findbest_test">-</a>
848 
849   qh_findbest_test(qh, testcentrum, facet, neighbor, bestfacet, dist, mindist, maxdist )
850     test neighbor of facet for qh_findbestneighbor()
851     if testcentrum,
852       tests centrum (assumes it is defined)
853     else
854       tests vertices
855 
856   returns:
857     if a better facet (i.e., vertices/centrum of facet closer to neighbor)
858       updates bestfacet, dist, mindist, and maxdist
859 */
qh_findbest_test(qhT * qh,boolT testcentrum,facetT * facet,facetT * neighbor,facetT ** bestfacet,realT * distp,realT * mindistp,realT * maxdistp)860 void qh_findbest_test(qhT *qh, boolT testcentrum, facetT *facet, facetT *neighbor,
861       facetT **bestfacet, realT *distp, realT *mindistp, realT *maxdistp) {
862   realT dist, mindist, maxdist;
863 
864   if (testcentrum) {
865     zzinc_(Zbestdist);
866     qh_distplane(qh, facet->center, neighbor, &dist);
867     dist *= qh->hull_dim; /* estimate furthest vertex */
868     if (dist < 0) {
869       maxdist= 0;
870       mindist= dist;
871       dist= -dist;
872     }else {
873       mindist= 0;
874       maxdist= dist;
875     }
876   }else
877     dist= qh_getdistance(qh, facet, neighbor, &mindist, &maxdist);
878   if (dist < *distp) {
879     *bestfacet= neighbor;
880     *mindistp= mindist;
881     *maxdistp= maxdist;
882     *distp= dist;
883   }
884 } /* findbest_test */
885 
886 /*-<a                             href="qh-merge_r.htm#TOC"
887   >-------------------------------</a><a name="findbestneighbor">-</a>
888 
889   qh_findbestneighbor(qh, facet, dist, mindist, maxdist )
890     finds best neighbor (least dist) of a facet for merging
891 
892   returns:
893     returns min and max distances and their max absolute value
894 
895   notes:
896     error if qh_ASvoronoi
897     avoids merging old into new
898     assumes ridge->nonconvex only set on one ridge between a pair of facets
899     could use an early out predicate but not worth it
900 
901   design:
902     if a large facet
903       will test centrum
904     else
905       will test vertices
906     if a large facet
907       test nonconvex neighbors for best merge
908     else
909       test all neighbors for the best merge
910     if testing centrum
911       get distance information
912 */
qh_findbestneighbor(qhT * qh,facetT * facet,realT * distp,realT * mindistp,realT * maxdistp)913 facetT *qh_findbestneighbor(qhT *qh, facetT *facet, realT *distp, realT *mindistp, realT *maxdistp) {
914   facetT *neighbor, **neighborp, *bestfacet= NULL;
915   ridgeT *ridge, **ridgep;
916   boolT nonconvex= True, testcentrum= False;
917   int size= qh_setsize(qh, facet->vertices);
918 
919   if(qh->CENTERtype==qh_ASvoronoi){
920     qh_fprintf(qh, qh->ferr, 6272, "qhull error: cannot call qh_findbestneighor for f%d while qh.CENTERtype is qh_ASvoronoi\n", facet->id);
921     qh_errexit(qh, qh_ERRqhull, facet, NULL);
922   }
923   *distp= REALmax;
924   if (size > qh_BESTcentrum2 * qh->hull_dim + qh_BESTcentrum) {
925     testcentrum= True;
926     zinc_(Zbestcentrum);
927     if (!facet->center)
928        facet->center= qh_getcentrum(qh, facet);
929   }
930   if (size > qh->hull_dim + qh_BESTnonconvex) {
931     FOREACHridge_(facet->ridges) {
932       if (ridge->nonconvex) {
933         neighbor= otherfacet_(ridge, facet);
934         qh_findbest_test(qh, testcentrum, facet, neighbor,
935                           &bestfacet, distp, mindistp, maxdistp);
936       }
937     }
938   }
939   if (!bestfacet) {
940     nonconvex= False;
941     FOREACHneighbor_(facet)
942       qh_findbest_test(qh, testcentrum, facet, neighbor,
943                         &bestfacet, distp, mindistp, maxdistp);
944   }
945   if (!bestfacet) {
946     qh_fprintf(qh, qh->ferr, 6095, "qhull internal error (qh_findbestneighbor): no neighbors for f%d\n", facet->id);
947     qh_errexit(qh, qh_ERRqhull, facet, NULL);
948   }
949   if (testcentrum)
950     qh_getdistance(qh, facet, bestfacet, mindistp, maxdistp);
951   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",
952      bestfacet->id, facet->id, testcentrum, nonconvex, *distp, *mindistp, *maxdistp));
953   return(bestfacet);
954 } /* findbestneighbor */
955 
956 
957 /*-<a                             href="qh-merge_r.htm#TOC"
958   >-------------------------------</a><a name="flippedmerges">-</a>
959 
960   qh_flippedmerges(qh, facetlist, wasmerge )
961     merge flipped facets into best neighbor
962     assumes qh.facet_mergeset at top of temporary stack
963 
964   returns:
965     no flipped facets on facetlist
966     sets wasmerge if merge occurred
967     degen/redundant merges passed through
968 
969   notes:
970     othermerges not needed since qh.facet_mergeset is empty before & after
971       keep it in case of change
972 
973   design:
974     append flipped facets to qh.facetmergeset
975     for each flipped merge
976       find best neighbor
977       merge facet into neighbor
978       merge degenerate and redundant facets
979     remove flipped merges from qh.facet_mergeset
980 */
qh_flippedmerges(qhT * qh,facetT * facetlist,boolT * wasmerge)981 void qh_flippedmerges(qhT *qh, facetT *facetlist, boolT *wasmerge) {
982   facetT *facet, *neighbor, *facet1;
983   realT dist, mindist, maxdist;
984   mergeT *merge, **mergep;
985   setT *othermerges;
986   int nummerge=0;
987 
988   trace4((qh, qh->ferr, 4024, "qh_flippedmerges: begin\n"));
989   FORALLfacet_(facetlist) {
990     if (facet->flipped && !facet->visible)
991       qh_appendmergeset(qh, facet, facet, MRGflip, NULL);
992   }
993   othermerges= qh_settemppop(qh); /* was facet_mergeset */
994   qh->facet_mergeset= qh_settemp(qh, qh->TEMPsize);
995   qh_settemppush(qh, othermerges);
996   FOREACHmerge_(othermerges) {
997     facet1= merge->facet1;
998     if (merge->type != MRGflip || facet1->visible)
999       continue;
1000     if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
1001       qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1002     neighbor= qh_findbestneighbor(qh, facet1, &dist, &mindist, &maxdist);
1003     trace0((qh, qh->ferr, 15, "qh_flippedmerges: merge flipped f%d into f%d dist %2.2g during p%d\n",
1004       facet1->id, neighbor->id, dist, qh->furthest_id));
1005     qh_mergefacet(qh, facet1, neighbor, &mindist, &maxdist, !qh_MERGEapex);
1006     nummerge++;
1007     if (qh->PRINTstatistics) {
1008       zinc_(Zflipped);
1009       wadd_(Wflippedtot, dist);
1010       wmax_(Wflippedmax, dist);
1011     }
1012     qh_merge_degenredundant(qh);
1013   }
1014   FOREACHmerge_(othermerges) {
1015     if (merge->facet1->visible || merge->facet2->visible)
1016       qh_memfree(qh, merge, (int)sizeof(mergeT));
1017     else
1018       qh_setappend(qh, &qh->facet_mergeset, merge);
1019   }
1020   qh_settempfree(qh, &othermerges);
1021   if (nummerge)
1022     *wasmerge= True;
1023   trace1((qh, qh->ferr, 1010, "qh_flippedmerges: merged %d flipped facets into a good neighbor\n", nummerge));
1024 } /* flippedmerges */
1025 
1026 
1027 /*-<a                             href="qh-merge_r.htm#TOC"
1028   >-------------------------------</a><a name="forcedmerges">-</a>
1029 
1030   qh_forcedmerges(qh, wasmerge )
1031     merge duplicated ridges
1032 
1033   returns:
1034     removes all duplicate ridges on facet_mergeset
1035     wasmerge set if merge
1036     qh.facet_mergeset may include non-forced merges(none for now)
1037     qh.degen_mergeset includes degen/redun merges
1038 
1039   notes:
1040     duplicate ridges occur when the horizon is pinched,
1041         i.e. a subridge occurs in more than two horizon ridges.
1042      could rename vertices that pinch the horizon
1043     assumes qh_merge_degenredundant() has not be called
1044     othermerges isn't needed since facet_mergeset is empty afterwards
1045       keep it in case of change
1046 
1047   design:
1048     for each duplicate ridge
1049       find current facets by chasing f.replace links
1050       check for wide merge due to duplicate ridge
1051       determine best direction for facet
1052       merge one facet into the other
1053       remove duplicate ridges from qh.facet_mergeset
1054 */
qh_forcedmerges(qhT * qh,boolT * wasmerge)1055 void qh_forcedmerges(qhT *qh, boolT *wasmerge) {
1056   facetT *facet1, *facet2;
1057   mergeT *merge, **mergep;
1058   realT dist1, dist2, mindist1, mindist2, maxdist1, maxdist2;
1059   setT *othermerges;
1060   int nummerge=0, numflip=0;
1061 
1062   if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
1063     qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1064   trace4((qh, qh->ferr, 4025, "qh_forcedmerges: begin\n"));
1065   othermerges= qh_settemppop(qh); /* was facet_mergeset */
1066   qh->facet_mergeset= qh_settemp(qh, qh->TEMPsize);
1067   qh_settemppush(qh, othermerges);
1068   FOREACHmerge_(othermerges) {
1069     if (merge->type != MRGridge)
1070         continue;
1071     if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
1072         qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1073     facet1= merge->facet1;
1074     facet2= merge->facet2;
1075     while (facet1->visible)      /* must exist, no qh_merge_degenredunant */
1076       facet1= facet1->f.replace; /* previously merged facet */
1077     while (facet2->visible)
1078       facet2= facet2->f.replace; /* previously merged facet */
1079     if (facet1 == facet2)
1080       continue;
1081     if (!qh_setin(facet2->neighbors, facet1)) {
1082       qh_fprintf(qh, qh->ferr, 6096, "qhull internal error (qh_forcedmerges): f%d and f%d had a duplicate ridge but as f%d and f%d they are no longer neighbors\n",
1083                merge->facet1->id, merge->facet2->id, facet1->id, facet2->id);
1084       qh_errexit2(qh, qh_ERRqhull, facet1, facet2);
1085     }
1086     dist1= qh_getdistance(qh, facet1, facet2, &mindist1, &maxdist1);
1087     dist2= qh_getdistance(qh, facet2, facet1, &mindist2, &maxdist2);
1088     qh_check_dupridge(qh, facet1, dist1, facet2, dist2);
1089     if (dist1 < dist2)
1090       qh_mergefacet(qh, facet1, facet2, &mindist1, &maxdist1, !qh_MERGEapex);
1091     else {
1092       qh_mergefacet(qh, facet2, facet1, &mindist2, &maxdist2, !qh_MERGEapex);
1093       dist1= dist2;
1094       facet1= facet2;
1095     }
1096     if (facet1->flipped) {
1097       zinc_(Zmergeflipdup);
1098       numflip++;
1099     }else
1100       nummerge++;
1101     if (qh->PRINTstatistics) {
1102       zinc_(Zduplicate);
1103       wadd_(Wduplicatetot, dist1);
1104       wmax_(Wduplicatemax, dist1);
1105     }
1106   }
1107   FOREACHmerge_(othermerges) {
1108     if (merge->type == MRGridge)
1109       qh_memfree(qh, merge, (int)sizeof(mergeT));
1110     else
1111       qh_setappend(qh, &qh->facet_mergeset, merge);
1112   }
1113   qh_settempfree(qh, &othermerges);
1114   if (nummerge)
1115     *wasmerge= True;
1116   trace1((qh, qh->ferr, 1011, "qh_forcedmerges: merged %d facets and %d flipped facets across duplicated ridges\n",
1117                 nummerge, numflip));
1118 } /* forcedmerges */
1119 
1120 
1121 /*-<a                             href="qh-merge_r.htm#TOC"
1122   >-------------------------------</a><a name="getmergeset">-</a>
1123 
1124   qh_getmergeset(qh, facetlist )
1125     determines nonconvex facets on facetlist
1126     tests !tested ridges and nonconvex ridges of !tested facets
1127 
1128   returns:
1129     returns sorted qh.facet_mergeset of facet-neighbor pairs to be merged
1130     all ridges tested
1131 
1132   notes:
1133     assumes no nonconvex ridges with both facets tested
1134     uses facet->tested/ridge->tested to prevent duplicate tests
1135     can not limit tests to modified ridges since the centrum changed
1136     uses qh.visit_id
1137 
1138   see:
1139     qh_getmergeset_initial()
1140 
1141   design:
1142     for each facet on facetlist
1143       for each ridge of facet
1144         if untested ridge
1145           test ridge for convexity
1146           if non-convex
1147             append ridge to qh.facet_mergeset
1148     sort qh.facet_mergeset by angle
1149 */
qh_getmergeset(qhT * qh,facetT * facetlist)1150 void qh_getmergeset(qhT *qh, facetT *facetlist) {
1151   facetT *facet, *neighbor, **neighborp;
1152   ridgeT *ridge, **ridgep;
1153   int nummerges;
1154 
1155   nummerges= qh_setsize(qh, qh->facet_mergeset);
1156   trace4((qh, qh->ferr, 4026, "qh_getmergeset: started.\n"));
1157   qh->visit_id++;
1158   FORALLfacet_(facetlist) {
1159     if (facet->tested)
1160       continue;
1161     facet->visitid= qh->visit_id;
1162     facet->tested= True;  /* must be non-simplicial due to merge */
1163     FOREACHneighbor_(facet)
1164       neighbor->seen= False;
1165     FOREACHridge_(facet->ridges) {
1166       if (ridge->tested && !ridge->nonconvex)
1167         continue;
1168       /* if tested & nonconvex, need to append merge */
1169       neighbor= otherfacet_(ridge, facet);
1170       if (neighbor->seen) {
1171         ridge->tested= True;
1172         ridge->nonconvex= False;
1173       }else if (neighbor->visitid != qh->visit_id) {
1174         ridge->tested= True;
1175         ridge->nonconvex= False;
1176         neighbor->seen= True;      /* only one ridge is marked nonconvex */
1177         if (qh_test_appendmerge(qh, facet, neighbor))
1178           ridge->nonconvex= True;
1179       }
1180     }
1181   }
1182   nummerges= qh_setsize(qh, qh->facet_mergeset);
1183   if (qh->ANGLEmerge)
1184     qsort(SETaddr_(qh->facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_compareangle);
1185   else
1186     qsort(SETaddr_(qh->facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_comparemerge);
1187   if (qh->POSTmerging) {
1188     zadd_(Zmergesettot2, nummerges);
1189   }else {
1190     zadd_(Zmergesettot, nummerges);
1191     zmax_(Zmergesetmax, nummerges);
1192   }
1193   trace2((qh, qh->ferr, 2021, "qh_getmergeset: %d merges found\n", nummerges));
1194 } /* getmergeset */
1195 
1196 
1197 /*-<a                             href="qh-merge_r.htm#TOC"
1198   >-------------------------------</a><a name="getmergeset_initial">-</a>
1199 
1200   qh_getmergeset_initial(qh, facetlist )
1201     determine initial qh.facet_mergeset for facets
1202     tests all facet/neighbor pairs on facetlist
1203 
1204   returns:
1205     sorted qh.facet_mergeset with nonconvex ridges
1206     sets facet->tested, ridge->tested, and ridge->nonconvex
1207 
1208   notes:
1209     uses visit_id, assumes ridge->nonconvex is False
1210 
1211   see:
1212     qh_getmergeset()
1213 
1214   design:
1215     for each facet on facetlist
1216       for each untested neighbor of facet
1217         test facet and neighbor for convexity
1218         if non-convex
1219           append merge to qh.facet_mergeset
1220           mark one of the ridges as nonconvex
1221     sort qh.facet_mergeset by angle
1222 */
qh_getmergeset_initial(qhT * qh,facetT * facetlist)1223 void qh_getmergeset_initial(qhT *qh, facetT *facetlist) {
1224   facetT *facet, *neighbor, **neighborp;
1225   ridgeT *ridge, **ridgep;
1226   int nummerges;
1227 
1228   qh->visit_id++;
1229   FORALLfacet_(facetlist) {
1230     facet->visitid= qh->visit_id;
1231     facet->tested= True;
1232     FOREACHneighbor_(facet) {
1233       if (neighbor->visitid != qh->visit_id) {
1234         if (qh_test_appendmerge(qh, facet, neighbor)) {
1235           FOREACHridge_(neighbor->ridges) {
1236             if (facet == otherfacet_(ridge, neighbor)) {
1237               ridge->nonconvex= True;
1238               break;    /* only one ridge is marked nonconvex */
1239             }
1240           }
1241         }
1242       }
1243     }
1244     FOREACHridge_(facet->ridges)
1245       ridge->tested= True;
1246   }
1247   nummerges= qh_setsize(qh, qh->facet_mergeset);
1248   if (qh->ANGLEmerge)
1249     qsort(SETaddr_(qh->facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_compareangle);
1250   else
1251     qsort(SETaddr_(qh->facet_mergeset, mergeT), (size_t)nummerges, sizeof(mergeT *), qh_comparemerge);
1252   if (qh->POSTmerging) {
1253     zadd_(Zmergeinittot2, nummerges);
1254   }else {
1255     zadd_(Zmergeinittot, nummerges);
1256     zmax_(Zmergeinitmax, nummerges);
1257   }
1258   trace2((qh, qh->ferr, 2022, "qh_getmergeset_initial: %d merges found\n", nummerges));
1259 } /* getmergeset_initial */
1260 
1261 
1262 /*-<a                             href="qh-merge_r.htm#TOC"
1263   >-------------------------------</a><a name="hashridge">-</a>
1264 
1265   qh_hashridge(qh, hashtable, hashsize, ridge, oldvertex )
1266     add ridge to hashtable without oldvertex
1267 
1268   notes:
1269     assumes hashtable is large enough
1270 
1271   design:
1272     determine hash value for ridge without oldvertex
1273     find next empty slot for ridge
1274 */
qh_hashridge(qhT * qh,setT * hashtable,int hashsize,ridgeT * ridge,vertexT * oldvertex)1275 void qh_hashridge(qhT *qh, setT *hashtable, int hashsize, ridgeT *ridge, vertexT *oldvertex) {
1276   int hash;
1277   ridgeT *ridgeA;
1278 
1279   hash= qh_gethash(qh, hashsize, ridge->vertices, qh->hull_dim-1, 0, oldvertex);
1280   while (True) {
1281     if (!(ridgeA= SETelemt_(hashtable, hash, ridgeT))) {
1282       SETelem_(hashtable, hash)= ridge;
1283       break;
1284     }else if (ridgeA == ridge)
1285       break;
1286     if (++hash == hashsize)
1287       hash= 0;
1288   }
1289 } /* hashridge */
1290 
1291 
1292 /*-<a                             href="qh-merge_r.htm#TOC"
1293   >-------------------------------</a><a name="hashridge_find">-</a>
1294 
1295   qh_hashridge_find(qh, hashtable, hashsize, ridge, vertex, oldvertex, hashslot )
1296     returns matching ridge without oldvertex in hashtable
1297       for ridge without vertex
1298     if oldvertex is NULL
1299       matches with any one skip
1300 
1301   returns:
1302     matching ridge or NULL
1303     if no match,
1304       if ridge already in   table
1305         hashslot= -1
1306       else
1307         hashslot= next NULL index
1308 
1309   notes:
1310     assumes hashtable is large enough
1311     can't match ridge to itself
1312 
1313   design:
1314     get hash value for ridge without vertex
1315     for each hashslot
1316       return match if ridge matches ridgeA without oldvertex
1317 */
qh_hashridge_find(qhT * qh,setT * hashtable,int hashsize,ridgeT * ridge,vertexT * vertex,vertexT * oldvertex,int * hashslot)1318 ridgeT *qh_hashridge_find(qhT *qh, setT *hashtable, int hashsize, ridgeT *ridge,
1319               vertexT *vertex, vertexT *oldvertex, int *hashslot) {
1320   int hash;
1321   ridgeT *ridgeA;
1322 
1323   *hashslot= 0;
1324   zinc_(Zhashridge);
1325   hash= qh_gethash(qh, hashsize, ridge->vertices, qh->hull_dim-1, 0, vertex);
1326   while ((ridgeA= SETelemt_(hashtable, hash, ridgeT))) {
1327     if (ridgeA == ridge)
1328       *hashslot= -1;
1329     else {
1330       zinc_(Zhashridgetest);
1331       if (qh_setequal_except(ridge->vertices, vertex, ridgeA->vertices, oldvertex))
1332         return ridgeA;
1333     }
1334     if (++hash == hashsize)
1335       hash= 0;
1336   }
1337   if (!*hashslot)
1338     *hashslot= hash;
1339   return NULL;
1340 } /* hashridge_find */
1341 
1342 
1343 /*-<a                             href="qh-merge_r.htm#TOC"
1344   >-------------------------------</a><a name="makeridges">-</a>
1345 
1346   qh_makeridges(qh, facet )
1347     creates explicit ridges between simplicial facets
1348 
1349   returns:
1350     facet with ridges and without qh_MERGEridge
1351     ->simplicial is False
1352 
1353   notes:
1354     allows qh_MERGEridge flag
1355     uses existing ridges
1356     duplicate neighbors ok if ridges already exist (qh_mergecycle_ridges)
1357 
1358   see:
1359     qh_mergecycle_ridges()
1360 
1361   design:
1362     look for qh_MERGEridge neighbors
1363     mark neighbors that already have ridges
1364     for each unprocessed neighbor of facet
1365       create a ridge for neighbor and facet
1366     if any qh_MERGEridge neighbors
1367       delete qh_MERGEridge flags (already handled by qh_mark_dupridges)
1368 */
qh_makeridges(qhT * qh,facetT * facet)1369 void qh_makeridges(qhT *qh, facetT *facet) {
1370   facetT *neighbor, **neighborp;
1371   ridgeT *ridge, **ridgep;
1372   int neighbor_i, neighbor_n;
1373   boolT toporient, mergeridge= False;
1374 
1375   if (!facet->simplicial)
1376     return;
1377   trace4((qh, qh->ferr, 4027, "qh_makeridges: make ridges for f%d\n", facet->id));
1378   facet->simplicial= False;
1379   FOREACHneighbor_(facet) {
1380     if (neighbor == qh_MERGEridge)
1381       mergeridge= True;
1382     else
1383       neighbor->seen= False;
1384   }
1385   FOREACHridge_(facet->ridges)
1386     otherfacet_(ridge, facet)->seen= True;
1387   FOREACHneighbor_i_(qh, facet) {
1388     if (neighbor == qh_MERGEridge)
1389       continue;  /* fixed by qh_mark_dupridges */
1390     else if (!neighbor->seen) {  /* no current ridges */
1391       ridge= qh_newridge(qh);
1392       ridge->vertices= qh_setnew_delnthsorted(qh, facet->vertices, qh->hull_dim,
1393                                                           neighbor_i, 0);
1394       toporient= facet->toporient ^ (neighbor_i & 0x1);
1395       if (toporient) {
1396         ridge->top= facet;
1397         ridge->bottom= neighbor;
1398       }else {
1399         ridge->top= neighbor;
1400         ridge->bottom= facet;
1401       }
1402 #if 0 /* this also works */
1403       flip= (facet->toporient ^ neighbor->toporient)^(skip1 & 0x1) ^ (skip2 & 0x1);
1404       if (facet->toporient ^ (skip1 & 0x1) ^ flip) {
1405         ridge->top= neighbor;
1406         ridge->bottom= facet;
1407       }else {
1408         ridge->top= facet;
1409         ridge->bottom= neighbor;
1410       }
1411 #endif
1412       qh_setappend(qh, &(facet->ridges), ridge);
1413       qh_setappend(qh, &(neighbor->ridges), ridge);
1414     }
1415   }
1416   if (mergeridge) {
1417     while (qh_setdel(facet->neighbors, qh_MERGEridge))
1418       ; /* delete each one */
1419   }
1420 } /* makeridges */
1421 
1422 
1423 /*-<a                             href="qh-merge_r.htm#TOC"
1424   >-------------------------------</a><a name="mark_dupridges">-</a>
1425 
1426   qh_mark_dupridges(qh, facetlist )
1427     add duplicated ridges to qh.facet_mergeset
1428     facet->dupridge is true
1429 
1430   returns:
1431     duplicate ridges on qh.facet_mergeset
1432     ->mergeridge/->mergeridge2 set
1433     duplicate ridges marked by qh_MERGEridge and both sides facet->dupridge
1434     no MERGEridges in neighbor sets
1435 
1436   notes:
1437     duplicate ridges occur when the horizon is pinched,
1438         i.e. a subridge occurs in more than two horizon ridges.
1439     could rename vertices that pinch the horizon (thus removing subridge)
1440     uses qh.visit_id
1441 
1442   design:
1443     for all facets on facetlist
1444       if facet contains a duplicate ridge
1445         for each neighbor of facet
1446           if neighbor marked qh_MERGEridge (one side of the merge)
1447             set facet->mergeridge
1448           else
1449             if neighbor contains a duplicate ridge
1450             and the back link is qh_MERGEridge
1451               append duplicate ridge to qh.facet_mergeset
1452    for each duplicate ridge
1453      make ridge sets in preparation for merging
1454      remove qh_MERGEridge from neighbor set
1455    for each duplicate ridge
1456      restore the missing neighbor from the neighbor set that was qh_MERGEridge
1457      add the missing ridge for this neighbor
1458 */
qh_mark_dupridges(qhT * qh,facetT * facetlist)1459 void qh_mark_dupridges(qhT *qh, facetT *facetlist) {
1460   facetT *facet, *neighbor, **neighborp;
1461   int nummerge=0;
1462   mergeT *merge, **mergep;
1463 
1464 
1465   trace4((qh, qh->ferr, 4028, "qh_mark_dupridges: identify duplicate ridges\n"));
1466   FORALLfacet_(facetlist) {
1467     if (facet->dupridge) {
1468       FOREACHneighbor_(facet) {
1469         if (neighbor == qh_MERGEridge) {
1470           facet->mergeridge= True;
1471           continue;
1472         }
1473         if (neighbor->dupridge
1474         && !qh_setin(neighbor->neighbors, facet)) { /* qh_MERGEridge */
1475           qh_appendmergeset(qh, facet, neighbor, MRGridge, NULL);
1476           facet->mergeridge2= True;
1477           facet->mergeridge= True;
1478           nummerge++;
1479         }
1480       }
1481     }
1482   }
1483   if (!nummerge)
1484     return;
1485   FORALLfacet_(facetlist) {            /* gets rid of qh_MERGEridge */
1486     if (facet->mergeridge && !facet->mergeridge2)
1487       qh_makeridges(qh, facet);
1488   }
1489   FOREACHmerge_(qh->facet_mergeset) {   /* restore the missing neighbors */
1490     if (merge->type == MRGridge) {
1491       qh_setappend(qh, &merge->facet2->neighbors, merge->facet1);
1492       qh_makeridges(qh, merge->facet1);   /* and the missing ridges */
1493     }
1494   }
1495   trace1((qh, qh->ferr, 1012, "qh_mark_dupridges: found %d duplicated ridges\n",
1496                 nummerge));
1497 } /* mark_dupridges */
1498 
1499 /*-<a                             href="qh-merge_r.htm#TOC"
1500   >-------------------------------</a><a name="maydropneighbor">-</a>
1501 
1502   qh_maydropneighbor(qh, facet )
1503     drop neighbor relationship if no ridge between facet and neighbor
1504 
1505   returns:
1506     neighbor sets updated
1507     appends degenerate facets to qh.facet_mergeset
1508 
1509   notes:
1510     won't cause redundant facets since vertex inclusion is the same
1511     may drop vertex and neighbor if no ridge
1512     uses qh.visit_id
1513 
1514   design:
1515     visit all neighbors with ridges
1516     for each unvisited neighbor of facet
1517       delete neighbor and facet from the neighbor sets
1518       if neighbor becomes degenerate
1519         append neighbor to qh.degen_mergeset
1520     if facet is degenerate
1521       append facet to qh.degen_mergeset
1522 */
qh_maydropneighbor(qhT * qh,facetT * facet)1523 void qh_maydropneighbor(qhT *qh, facetT *facet) {
1524   ridgeT *ridge, **ridgep;
1525   realT angledegen= qh_ANGLEdegen;
1526   facetT *neighbor, **neighborp;
1527 
1528   qh->visit_id++;
1529   trace4((qh, qh->ferr, 4029, "qh_maydropneighbor: test f%d for no ridges to a neighbor\n",
1530           facet->id));
1531   FOREACHridge_(facet->ridges) {
1532     ridge->top->visitid= qh->visit_id;
1533     ridge->bottom->visitid= qh->visit_id;
1534   }
1535   FOREACHneighbor_(facet) {
1536     if (neighbor->visitid != qh->visit_id) {
1537       trace0((qh, qh->ferr, 17, "qh_maydropneighbor: facets f%d and f%d are no longer neighbors during p%d\n",
1538             facet->id, neighbor->id, qh->furthest_id));
1539       zinc_(Zdropneighbor);
1540       qh_setdel(facet->neighbors, neighbor);
1541       neighborp--;  /* repeat, deleted a neighbor */
1542       qh_setdel(neighbor->neighbors, facet);
1543       if (qh_setsize(qh, neighbor->neighbors) < qh->hull_dim) {
1544         zinc_(Zdropdegen);
1545         qh_appendmergeset(qh, neighbor, neighbor, MRGdegen, &angledegen);
1546         trace2((qh, qh->ferr, 2023, "qh_maydropneighbors: f%d is degenerate.\n", neighbor->id));
1547       }
1548     }
1549   }
1550   if (qh_setsize(qh, facet->neighbors) < qh->hull_dim) {
1551     zinc_(Zdropdegen);
1552     qh_appendmergeset(qh, facet, facet, MRGdegen, &angledegen);
1553     trace2((qh, qh->ferr, 2024, "qh_maydropneighbors: f%d is degenerate.\n", facet->id));
1554   }
1555 } /* maydropneighbor */
1556 
1557 
1558 /*-<a                             href="qh-merge_r.htm#TOC"
1559   >-------------------------------</a><a name="merge_degenredundant">-</a>
1560 
1561   qh_merge_degenredundant(qh)
1562     merge all degenerate and redundant facets
1563     qh.degen_mergeset contains merges from qh_degen_redundant_neighbors()
1564 
1565   returns:
1566     number of merges performed
1567     resets facet->degenerate/redundant
1568     if deleted (visible) facet has no neighbors
1569       sets ->f.replace to NULL
1570 
1571   notes:
1572     redundant merges happen before degenerate ones
1573     merging and renaming vertices can result in degen/redundant facets
1574 
1575   design:
1576     for each merge on qh.degen_mergeset
1577       if redundant merge
1578         if non-redundant facet merged into redundant facet
1579           recheck facet for redundancy
1580         else
1581           merge redundant facet into other facet
1582 */
qh_merge_degenredundant(qhT * qh)1583 int qh_merge_degenredundant(qhT *qh) {
1584   int size;
1585   mergeT *merge;
1586   facetT *bestneighbor, *facet1, *facet2;
1587   realT dist, mindist, maxdist;
1588   vertexT *vertex, **vertexp;
1589   int nummerges= 0;
1590   mergeType mergetype;
1591 
1592   while ((merge= (mergeT*)qh_setdellast(qh->degen_mergeset))) {
1593     facet1= merge->facet1;
1594     facet2= merge->facet2;
1595     mergetype= merge->type;
1596     qh_memfree(qh, merge, (int)sizeof(mergeT));
1597     if (facet1->visible)
1598       continue;
1599     facet1->degenerate= False;
1600     facet1->redundant= False;
1601     if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
1602       qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1603     if (mergetype == MRGredundant) {
1604       zinc_(Zneighbor);
1605       while (facet2->visible) {
1606         if (!facet2->f.replace) {
1607           qh_fprintf(qh, qh->ferr, 6097, "qhull internal error (qh_merge_degenredunant): f%d redundant but f%d has no replacement\n",
1608                facet1->id, facet2->id);
1609           qh_errexit2(qh, qh_ERRqhull, facet1, facet2);
1610         }
1611         facet2= facet2->f.replace;
1612       }
1613       if (facet1 == facet2) {
1614         qh_degen_redundant_facet(qh, facet1); /* in case of others */
1615         continue;
1616       }
1617       trace2((qh, qh->ferr, 2025, "qh_merge_degenredundant: facet f%d is contained in f%d, will merge\n",
1618             facet1->id, facet2->id));
1619       qh_mergefacet(qh, facet1, facet2, NULL, NULL, !qh_MERGEapex);
1620       /* merge distance is already accounted for */
1621       nummerges++;
1622     }else {  /* mergetype == MRGdegen, other merges may have fixed */
1623       if (!(size= qh_setsize(qh, facet1->neighbors))) {
1624         zinc_(Zdelfacetdup);
1625         trace2((qh, qh->ferr, 2026, "qh_merge_degenredundant: facet f%d has no neighbors.  Deleted\n", facet1->id));
1626         qh_willdelete(qh, facet1, NULL);
1627         FOREACHvertex_(facet1->vertices) {
1628           qh_setdel(vertex->neighbors, facet1);
1629           if (!SETfirst_(vertex->neighbors)) {
1630             zinc_(Zdegenvertex);
1631             trace2((qh, qh->ferr, 2027, "qh_merge_degenredundant: deleted v%d because f%d has no neighbors\n",
1632                  vertex->id, facet1->id));
1633             vertex->deleted= True;
1634             qh_setappend(qh, &qh->del_vertices, vertex);
1635           }
1636         }
1637         nummerges++;
1638       }else if (size < qh->hull_dim) {
1639         bestneighbor= qh_findbestneighbor(qh, facet1, &dist, &mindist, &maxdist);
1640         trace2((qh, qh->ferr, 2028, "qh_merge_degenredundant: facet f%d has %d neighbors, merge into f%d dist %2.2g\n",
1641               facet1->id, size, bestneighbor->id, dist));
1642         qh_mergefacet(qh, facet1, bestneighbor, &mindist, &maxdist, !qh_MERGEapex);
1643         nummerges++;
1644         if (qh->PRINTstatistics) {
1645           zinc_(Zdegen);
1646           wadd_(Wdegentot, dist);
1647           wmax_(Wdegenmax, dist);
1648         }
1649       } /* else, another merge fixed the degeneracy and redundancy tested */
1650     }
1651   }
1652   return nummerges;
1653 } /* merge_degenredundant */
1654 
1655 /*-<a                             href="qh-merge_r.htm#TOC"
1656   >-------------------------------</a><a name="merge_nonconvex">-</a>
1657 
1658   qh_merge_nonconvex(qh, facet1, facet2, mergetype )
1659     remove non-convex ridge between facet1 into facet2
1660     mergetype gives why the facet's are non-convex
1661 
1662   returns:
1663     merges one of the facets into the best neighbor
1664 
1665   design:
1666     if one of the facets is a new facet
1667       prefer merging new facet into old facet
1668     find best neighbors for both facets
1669     merge the nearest facet into its best neighbor
1670     update the statistics
1671 */
qh_merge_nonconvex(qhT * qh,facetT * facet1,facetT * facet2,mergeType mergetype)1672 void qh_merge_nonconvex(qhT *qh, facetT *facet1, facetT *facet2, mergeType mergetype) {
1673   facetT *bestfacet, *bestneighbor, *neighbor;
1674   realT dist, dist2, mindist, mindist2, maxdist, maxdist2;
1675 
1676   if (qh->TRACEmerge-1 == zzval_(Ztotmerge))
1677     qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1678   trace3((qh, qh->ferr, 3003, "qh_merge_nonconvex: merge #%d for f%d and f%d type %d\n",
1679       zzval_(Ztotmerge) + 1, facet1->id, facet2->id, mergetype));
1680   /* concave or coplanar */
1681   if (!facet1->newfacet) {
1682     bestfacet= facet2;   /* avoid merging old facet if new is ok */
1683     facet2= facet1;
1684     facet1= bestfacet;
1685   }else
1686     bestfacet= facet1;
1687   bestneighbor= qh_findbestneighbor(qh, bestfacet, &dist, &mindist, &maxdist);
1688   neighbor= qh_findbestneighbor(qh, facet2, &dist2, &mindist2, &maxdist2);
1689   if (dist < dist2) {
1690     qh_mergefacet(qh, bestfacet, bestneighbor, &mindist, &maxdist, !qh_MERGEapex);
1691   }else if (qh->AVOIDold && !facet2->newfacet
1692   && ((mindist >= -qh->MAXcoplanar && maxdist <= qh->max_outside)
1693        || dist * 1.5 < dist2)) {
1694     zinc_(Zavoidold);
1695     wadd_(Wavoidoldtot, dist);
1696     wmax_(Wavoidoldmax, dist);
1697     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",
1698            facet2->id, dist2, facet1->id, dist2));
1699     qh_mergefacet(qh, bestfacet, bestneighbor, &mindist, &maxdist, !qh_MERGEapex);
1700   }else {
1701     qh_mergefacet(qh, facet2, neighbor, &mindist2, &maxdist2, !qh_MERGEapex);
1702     dist= dist2;
1703   }
1704   if (qh->PRINTstatistics) {
1705     if (mergetype == MRGanglecoplanar) {
1706       zinc_(Zacoplanar);
1707       wadd_(Wacoplanartot, dist);
1708       wmax_(Wacoplanarmax, dist);
1709     }else if (mergetype == MRGconcave) {
1710       zinc_(Zconcave);
1711       wadd_(Wconcavetot, dist);
1712       wmax_(Wconcavemax, dist);
1713     }else { /* MRGcoplanar */
1714       zinc_(Zcoplanar);
1715       wadd_(Wcoplanartot, dist);
1716       wmax_(Wcoplanarmax, dist);
1717     }
1718   }
1719 } /* merge_nonconvex */
1720 
1721 /*-<a                             href="qh-merge_r.htm#TOC"
1722   >-------------------------------</a><a name="mergecycle">-</a>
1723 
1724   qh_mergecycle(qh, samecycle, newfacet )
1725     merge a cycle of facets starting at samecycle into a newfacet
1726     newfacet is a horizon facet with ->normal
1727     samecycle facets are simplicial from an apex
1728 
1729   returns:
1730     initializes vertex neighbors on first merge
1731     samecycle deleted (placed on qh.visible_list)
1732     newfacet at end of qh.facet_list
1733     deleted vertices on qh.del_vertices
1734 
1735   see:
1736     qh_mergefacet()
1737     called by qh_mergecycle_all() for multiple, same cycle facets
1738 
1739   design:
1740     make vertex neighbors if necessary
1741     make ridges for newfacet
1742     merge neighbor sets of samecycle into newfacet
1743     merge ridges of samecycle into newfacet
1744     merge vertex neighbors of samecycle into newfacet
1745     make apex of samecycle the apex of newfacet
1746     if newfacet wasn't a new facet
1747       add its vertices to qh.newvertex_list
1748     delete samecycle facets a make newfacet a newfacet
1749 */
qh_mergecycle(qhT * qh,facetT * samecycle,facetT * newfacet)1750 void qh_mergecycle(qhT *qh, facetT *samecycle, facetT *newfacet) {
1751   int traceonce= False, tracerestore= 0;
1752   vertexT *apex;
1753 #ifndef qh_NOtrace
1754   facetT *same;
1755 #endif
1756 
1757   if (newfacet->tricoplanar) {
1758     if (!qh->TRInormals) {
1759       qh_fprintf(qh, qh->ferr, 6224, "Qhull internal error (qh_mergecycle): does not work for tricoplanar facets.  Use option 'Q11'\n");
1760       qh_errexit(qh, qh_ERRqhull, newfacet, NULL);
1761     }
1762     newfacet->tricoplanar= False;
1763     newfacet->keepcentrum= False;
1764   }
1765   if (!qh->VERTEXneighbors)
1766     qh_vertexneighbors(qh);
1767   zzinc_(Ztotmerge);
1768   if (qh->REPORTfreq2 && qh->POSTmerging) {
1769     if (zzval_(Ztotmerge) > qh->mergereport + qh->REPORTfreq2)
1770       qh_tracemerging(qh);
1771   }
1772 #ifndef qh_NOtrace
1773   if (qh->TRACEmerge == zzval_(Ztotmerge))
1774     qh->qhmem.IStracing= qh->IStracing= qh->TRACElevel;
1775   trace2((qh, qh->ferr, 2030, "qh_mergecycle: merge #%d for facets from cycle f%d into coplanar horizon f%d\n",
1776         zzval_(Ztotmerge), samecycle->id, newfacet->id));
1777   if (newfacet == qh->tracefacet) {
1778     tracerestore= qh->IStracing;
1779     qh->IStracing= 4;
1780     qh_fprintf(qh, qh->ferr, 8068, "qh_mergecycle: ========= trace merge %d of samecycle %d into trace f%d, furthest is p%d\n",
1781                zzval_(Ztotmerge), samecycle->id, newfacet->id,  qh->furthest_id);
1782     traceonce= True;
1783   }
1784   if (qh->IStracing >=4) {
1785     qh_fprintf(qh, qh->ferr, 8069, "  same cycle:");
1786     FORALLsame_cycle_(samecycle)
1787       qh_fprintf(qh, qh->ferr, 8070, " f%d", same->id);
1788     qh_fprintf(qh, qh->ferr, 8071, "\n");
1789   }
1790   if (qh->IStracing >=4)
1791     qh_errprint(qh, "MERGING CYCLE", samecycle, newfacet, NULL, NULL);
1792 #endif /* !qh_NOtrace */
1793   apex= SETfirstt_(samecycle->vertices, vertexT);
1794   qh_makeridges(qh, newfacet);
1795   qh_mergecycle_neighbors(qh, samecycle, newfacet);
1796   qh_mergecycle_ridges(qh, samecycle, newfacet);
1797   qh_mergecycle_vneighbors(qh, samecycle, newfacet);
1798   if (SETfirstt_(newfacet->vertices, vertexT) != apex)
1799     qh_setaddnth(qh, &newfacet->vertices, 0, apex);  /* apex has last id */
1800   if (!newfacet->newfacet)
1801     qh_newvertices(qh, newfacet->vertices);
1802   qh_mergecycle_facets(qh, samecycle, newfacet);
1803   qh_tracemerge(qh, samecycle, newfacet);
1804   /* check for degen_redundant_neighbors after qh_forcedmerges() */
1805   if (traceonce) {
1806     qh_fprintf(qh, qh->ferr, 8072, "qh_mergecycle: end of trace facet\n");
1807     qh->IStracing= tracerestore;
1808   }
1809 } /* mergecycle */
1810 
1811 /*-<a                             href="qh-merge_r.htm#TOC"
1812   >-------------------------------</a><a name="mergecycle_all">-</a>
1813 
1814   qh_mergecycle_all(qh, facetlist, wasmerge )
1815     merge all samecycles of coplanar facets into horizon
1816     don't merge facets with ->mergeridge (these already have ->normal)
1817     all facets are simplicial from apex
1818     all facet->cycledone == False
1819 
1820   returns:
1821     all newfacets merged into coplanar horizon facets
1822     deleted vertices on  qh.del_vertices
1823     sets wasmerge if any merge
1824 
1825   see:
1826     calls qh_mergecycle for multiple, same cycle facets
1827 
1828   design:
1829     for each facet on facetlist
1830       skip facets with duplicate ridges and normals
1831       check that facet is in a samecycle (->mergehorizon)
1832       if facet only member of samecycle
1833         sets vertex->delridge for all vertices except apex
1834         merge facet into horizon
1835       else
1836         mark all facets in samecycle
1837         remove facets with duplicate ridges from samecycle
1838         merge samecycle into horizon (deletes facets from facetlist)
1839 */
qh_mergecycle_all(qhT * qh,facetT * facetlist,boolT * wasmerge)1840 void qh_mergecycle_all(qhT *qh, facetT *facetlist, boolT *wasmerge) {
1841   facetT *facet, *same, *prev, *horizon;
1842   facetT *samecycle= NULL, *nextfacet, *nextsame;
1843   vertexT *apex, *vertex, **vertexp;
1844   int cycles=0, total=0, facets, nummerge;
1845 
1846   trace2((qh, qh->ferr, 2031, "qh_mergecycle_all: begin\n"));
1847   for (facet= facetlist; facet && (nextfacet= facet->next); facet= nextfacet) {
1848     if (facet->normal)
1849       continue;
1850     if (!facet->mergehorizon) {
1851       qh_fprintf(qh, qh->ferr, 6225, "Qhull internal error (qh_mergecycle_all): f%d without normal\n", facet->id);
1852       qh_errexit(qh, qh_ERRqhull, facet, NULL);
1853     }
1854     horizon= SETfirstt_(facet->neighbors, facetT);
1855     if (facet->f.samecycle == facet) {
1856       zinc_(Zonehorizon);
1857       /* merge distance done in qh_findhorizon */
1858       apex= SETfirstt_(facet->vertices, vertexT);
1859       FOREACHvertex_(facet->vertices) {
1860         if (vertex != apex)
1861           vertex->delridge= True;
1862       }
1863       horizon->f.newcycle= NULL;
1864       qh_mergefacet(qh, facet, horizon, NULL, NULL, qh_MERGEapex);
1865     }else {
1866       samecycle= facet;
1867       facets= 0;
1868       prev= facet;
1869       for (same= facet->f.samecycle; same;  /* FORALLsame_cycle_(facet) */
1870            same= (same == facet ? NULL :nextsame)) { /* ends at facet */
1871         nextsame= same->f.samecycle;
1872         if (same->cycledone || same->visible)
1873           qh_infiniteloop(qh, same);
1874         same->cycledone= True;
1875         if (same->normal) {
1876           prev->f.samecycle= same->f.samecycle; /* unlink ->mergeridge */
1877           same->f.samecycle= NULL;
1878         }else {
1879           prev= same;
1880           facets++;
1881         }
1882       }
1883       while (nextfacet && nextfacet->cycledone)  /* will delete samecycle */
1884         nextfacet= nextfacet->next;
1885       horizon->f.newcycle= NULL;
1886       qh_mergecycle(qh, samecycle, horizon);
1887       nummerge= horizon->nummerge + facets;
1888       if (nummerge > qh_MAXnummerge)
1889         horizon->nummerge= qh_MAXnummerge;
1890       else
1891         horizon->nummerge= (short unsigned int)nummerge;
1892       zzinc_(Zcyclehorizon);
1893       total += facets;
1894       zzadd_(Zcyclefacettot, facets);
1895       zmax_(Zcyclefacetmax, facets);
1896     }
1897     cycles++;
1898   }
1899   if (cycles)
1900     *wasmerge= True;
1901   trace1((qh, qh->ferr, 1013, "qh_mergecycle_all: merged %d same cycles or facets into coplanar horizons\n", cycles));
1902 } /* mergecycle_all */
1903 
1904 /*-<a                             href="qh-merge_r.htm#TOC"
1905   >-------------------------------</a><a name="mergecycle_facets">-</a>
1906 
1907   qh_mergecycle_facets(qh, samecycle, newfacet )
1908     finish merge of samecycle into newfacet
1909 
1910   returns:
1911     samecycle prepended to visible_list for later deletion and partitioning
1912       each facet->f.replace == newfacet
1913 
1914     newfacet moved to end of qh.facet_list
1915       makes newfacet a newfacet (get's facet1->id if it was old)
1916       sets newfacet->newmerge
1917       clears newfacet->center (unless merging into a large facet)
1918       clears newfacet->tested and ridge->tested for facet1
1919 
1920     adds neighboring facets to facet_mergeset if redundant or degenerate
1921 
1922   design:
1923     make newfacet a new facet and set its flags
1924     move samecycle facets to qh.visible_list for later deletion
1925     unless newfacet is large
1926       remove its centrum
1927 */
qh_mergecycle_facets(qhT * qh,facetT * samecycle,facetT * newfacet)1928 void qh_mergecycle_facets(qhT *qh, facetT *samecycle, facetT *newfacet) {
1929   facetT *same, *next;
1930 
1931   trace4((qh, qh->ferr, 4030, "qh_mergecycle_facets: make newfacet new and samecycle deleted\n"));
1932   qh_removefacet(qh, newfacet);  /* append as a newfacet to end of qh->facet_list */
1933   qh_appendfacet(qh, newfacet);
1934   newfacet->newfacet= True;
1935   newfacet->simplicial= False;
1936   newfacet->newmerge= True;
1937 
1938   for (same= samecycle->f.samecycle; same; same= (same == samecycle ?  NULL : next)) {
1939     next= same->f.samecycle;  /* reused by willdelete */
1940     qh_willdelete(qh, same, newfacet);
1941   }
1942   if (newfacet->center
1943       && qh_setsize(qh, newfacet->vertices) <= qh->hull_dim + qh_MAXnewcentrum) {
1944     qh_memfree(qh, newfacet->center, qh->normal_size);
1945     newfacet->center= NULL;
1946   }
1947   trace3((qh, qh->ferr, 3004, "qh_mergecycle_facets: merged facets from cycle f%d into f%d\n",
1948              samecycle->id, newfacet->id));
1949 } /* mergecycle_facets */
1950 
1951 /*-<a                             href="qh-merge_r.htm#TOC"
1952   >-------------------------------</a><a name="mergecycle_neighbors">-</a>
1953 
1954   qh_mergecycle_neighbors(qh, samecycle, newfacet )
1955     add neighbors for samecycle facets to newfacet
1956 
1957   returns:
1958     newfacet with updated neighbors and vice-versa
1959     newfacet has ridges
1960     all neighbors of newfacet marked with qh.visit_id
1961     samecycle facets marked with qh.visit_id-1
1962     ridges updated for simplicial neighbors of samecycle with a ridge
1963 
1964   notes:
1965     assumes newfacet not in samecycle
1966     usually, samecycle facets are new, simplicial facets without internal ridges
1967       not so if horizon facet is coplanar to two different samecycles
1968 
1969   see:
1970     qh_mergeneighbors()
1971 
1972   design:
1973     check samecycle
1974     delete neighbors from newfacet that are also in samecycle
1975     for each neighbor of a facet in samecycle
1976       if neighbor is simplicial
1977         if first visit
1978           move the neighbor relation to newfacet
1979           update facet links for its ridges
1980         else
1981           make ridges for neighbor
1982           remove samecycle reference
1983       else
1984         update neighbor sets
1985 */
qh_mergecycle_neighbors(qhT * qh,facetT * samecycle,facetT * newfacet)1986 void qh_mergecycle_neighbors(qhT *qh, facetT *samecycle, facetT *newfacet) {
1987   facetT *same, *neighbor, **neighborp;
1988   int delneighbors= 0, newneighbors= 0;
1989   unsigned int samevisitid;
1990   ridgeT *ridge, **ridgep;
1991 
1992   samevisitid= ++qh->visit_id;
1993   FORALLsame_cycle_(samecycle) {
1994     if (same->visitid == samevisitid || same->visible)
1995       qh_infiniteloop(qh, samecycle);
1996     same->visitid= samevisitid;
1997   }
1998   newfacet->visitid= ++qh->visit_id;
1999   trace4((qh, qh->ferr, 4031, "qh_mergecycle_neighbors: delete shared neighbors from newfacet\n"));
2000   FOREACHneighbor_(newfacet) {
2001     if (neighbor->visitid == samevisitid) {
2002       SETref_(neighbor)= NULL;  /* samecycle neighbors deleted */
2003       delneighbors++;
2004     }else
2005       neighbor->visitid= qh->visit_id;
2006   }
2007   qh_setcompact(qh, newfacet->neighbors);
2008 
2009   trace4((qh, qh->ferr, 4032, "qh_mergecycle_neighbors: update neighbors\n"));
2010   FORALLsame_cycle_(samecycle) {
2011     FOREACHneighbor_(same) {
2012       if (neighbor->visitid == samevisitid)
2013         continue;
2014       if (neighbor->simplicial) {
2015         if (neighbor->visitid != qh->visit_id) {
2016           qh_setappend(qh, &newfacet->neighbors, neighbor);
2017           qh_setreplace(qh, neighbor->neighbors, same, newfacet);
2018           newneighbors++;
2019           neighbor->visitid= qh->visit_id;
2020           FOREACHridge_(neighbor->ridges) { /* update ridge in case of qh_makeridges */
2021             if (ridge->top == same) {
2022               ridge->top= newfacet;
2023               break;
2024             }else if (ridge->bottom == same) {
2025               ridge->bottom= newfacet;
2026               break;
2027             }
2028           }
2029         }else {
2030           qh_makeridges(qh, neighbor);
2031           qh_setdel(neighbor->neighbors, same);
2032           /* same can't be horizon facet for neighbor */
2033         }
2034       }else { /* non-simplicial neighbor */
2035         qh_setdel(neighbor->neighbors, same);
2036         if (neighbor->visitid != qh->visit_id) {
2037           qh_setappend(qh, &neighbor->neighbors, newfacet);
2038           qh_setappend(qh, &newfacet->neighbors, neighbor);
2039           neighbor->visitid= qh->visit_id;
2040           newneighbors++;
2041         }
2042       }
2043     }
2044   }
2045   trace2((qh, qh->ferr, 2032, "qh_mergecycle_neighbors: deleted %d neighbors and added %d\n",
2046              delneighbors, newneighbors));
2047 } /* mergecycle_neighbors */
2048 
2049 /*-<a                             href="qh-merge_r.htm#TOC"
2050   >-------------------------------</a><a name="mergecycle_ridges">-</a>
2051 
2052   qh_mergecycle_ridges(qh, samecycle, newfacet )
2053     add ridges/neighbors for facets in samecycle to newfacet
2054     all new/old neighbors of newfacet marked with qh.visit_id
2055     facets in samecycle marked with qh.visit_id-1
2056     newfacet marked with qh.visit_id
2057 
2058   returns:
2059     newfacet has merged ridges
2060 
2061   notes:
2062     ridge already updated for simplicial neighbors of samecycle with a ridge
2063 
2064   see:
2065     qh_mergeridges()
2066     qh_makeridges()
2067 
2068   design:
2069     remove ridges between newfacet and samecycle
2070     for each facet in samecycle
2071       for each ridge in facet
2072         update facet pointers in ridge
2073         skip ridges processed in qh_mergecycle_neighors
2074         free ridges between newfacet and samecycle
2075         free ridges between facets of samecycle (on 2nd visit)
2076         append remaining ridges to newfacet
2077       if simpilicial facet
2078         for each neighbor of facet
2079           if simplicial facet
2080           and not samecycle facet or newfacet
2081             make ridge between neighbor and newfacet
2082 */
qh_mergecycle_ridges(qhT * qh,facetT * samecycle,facetT * newfacet)2083 void qh_mergecycle_ridges(qhT *qh, facetT *samecycle, facetT *newfacet) {
2084   facetT *same, *neighbor= NULL;
2085   int numold=0, numnew=0;
2086   int neighbor_i, neighbor_n;
2087   unsigned int samevisitid;
2088   ridgeT *ridge, **ridgep;
2089   boolT toporient;
2090   void **freelistp; /* used if !qh_NOmem by qh_memfree_() */
2091 
2092   trace4((qh, qh->ferr, 4033, "qh_mergecycle_ridges: delete shared ridges from newfacet\n"));
2093   samevisitid= qh->visit_id -1;
2094   FOREACHridge_(newfacet->ridges) {
2095     neighbor= otherfacet_(ridge, newfacet);
2096     if (neighbor->visitid == samevisitid)
2097       SETref_(ridge)= NULL; /* ridge free'd below */
2098   }
2099   qh_setcompact(qh, newfacet->ridges);
2100 
2101   trace4((qh, qh->ferr, 4034, "qh_mergecycle_ridges: add ridges to newfacet\n"));
2102   FORALLsame_cycle_(samecycle) {
2103     FOREACHridge_(same->ridges) {
2104       if (ridge->top == same) {
2105         ridge->top= newfacet;
2106         neighbor= ridge->bottom;
2107       }else if (ridge->bottom == same) {
2108         ridge->bottom= newfacet;
2109         neighbor= ridge->top;
2110       }else if (ridge->top == newfacet || ridge->bottom == newfacet) {
2111         qh_setappend(qh, &newfacet->ridges, ridge);
2112         numold++;  /* already set by qh_mergecycle_neighbors */
2113         continue;
2114       }else {
2115         qh_fprintf(qh, qh->ferr, 6098, "qhull internal error (qh_mergecycle_ridges): bad ridge r%d\n", ridge->id);
2116         qh_errexit(qh, qh_ERRqhull, NULL, ridge);
2117       }
2118       if (neighbor == newfacet) {
2119         qh_setfree(qh, &(ridge->vertices));
2120         qh_memfree_(qh, ridge, (int)sizeof(ridgeT), freelistp);
2121         numold++;
2122       }else if (neighbor->visitid == samevisitid) {
2123         qh_setdel(neighbor->ridges, ridge);
2124         qh_setfree(qh, &(ridge->vertices));
2125         qh_memfree_(qh, ridge, (int)sizeof(ridgeT), freelistp);
2126         numold++;
2127       }else {
2128         qh_setappend(qh, &newfacet->ridges, ridge);
2129         numold++;
2130       }
2131     }
2132     if (same->ridges)
2133       qh_settruncate(qh, same->ridges, 0);
2134     if (!same->simplicial)
2135       continue;
2136     FOREACHneighbor_i_(qh, same) {       /* note: !newfact->simplicial */
2137       if (neighbor->visitid != samevisitid && neighbor->simplicial) {
2138         ridge= qh_newridge(qh);
2139         ridge->vertices= qh_setnew_delnthsorted(qh, same->vertices, qh->hull_dim,
2140                                                           neighbor_i, 0);
2141         toporient= same->toporient ^ (neighbor_i & 0x1);
2142         if (toporient) {
2143           ridge->top= newfacet;
2144           ridge->bottom= neighbor;
2145         }else {
2146           ridge->top= neighbor;
2147           ridge->bottom= newfacet;
2148         }
2149         qh_setappend(qh, &(newfacet->ridges), ridge);
2150         qh_setappend(qh, &(neighbor->ridges), ridge);
2151         numnew++;
2152       }
2153     }
2154   }
2155 
2156   trace2((qh, qh->ferr, 2033, "qh_mergecycle_ridges: found %d old ridges and %d new ones\n",
2157              numold, numnew));
2158 } /* mergecycle_ridges */
2159 
2160 /*-<a                             href="qh-merge_r.htm#TOC"
2161   >-------------------------------</a><a name="mergecycle_vneighbors">-</a>
2162 
2163   qh_mergecycle_vneighbors(qh, samecycle, newfacet )
2164     create vertex neighbors for newfacet from vertices of facets in samecycle
2165     samecycle marked with visitid == qh.visit_id - 1
2166 
2167   returns:
2168     newfacet vertices with updated neighbors
2169     marks newfacet with qh.visit_id-1
2170     deletes vertices that are merged away
2171     sets delridge on all vertices (faster here than in mergecycle_ridges)
2172 
2173   see:
2174     qh_mergevertex_neighbors()
2175 
2176   design:
2177     for each vertex of samecycle facet
2178       set vertex->delridge
2179       delete samecycle facets from vertex neighbors
2180       append newfacet to vertex neighbors
2181       if vertex only in newfacet
2182         delete it from newfacet
2183         add it to qh.del_vertices for later deletion
2184 */
qh_mergecycle_vneighbors(qhT * qh,facetT * samecycle,facetT * newfacet)2185 void qh_mergecycle_vneighbors(qhT *qh, facetT *samecycle, facetT *newfacet) {
2186   facetT *neighbor, **neighborp;
2187   unsigned int mergeid;
2188   vertexT *vertex, **vertexp, *apex;
2189   setT *vertices;
2190 
2191   trace4((qh, qh->ferr, 4035, "qh_mergecycle_vneighbors: update vertex neighbors for newfacet\n"));
2192   mergeid= qh->visit_id - 1;
2193   newfacet->visitid= mergeid;
2194   vertices= qh_basevertices(qh, samecycle); /* temp */
2195   apex= SETfirstt_(samecycle->vertices, vertexT);
2196   qh_setappend(qh, &vertices, apex);
2197   FOREACHvertex_(vertices) {
2198     vertex->delridge= True;
2199     FOREACHneighbor_(vertex) {
2200       if (neighbor->visitid == mergeid)
2201         SETref_(neighbor)= NULL;
2202     }
2203     qh_setcompact(qh, vertex->neighbors);
2204     qh_setappend(qh, &vertex->neighbors, newfacet);
2205     if (!SETsecond_(vertex->neighbors)) {
2206       zinc_(Zcyclevertex);
2207       trace2((qh, qh->ferr, 2034, "qh_mergecycle_vneighbors: deleted v%d when merging cycle f%d into f%d\n",
2208         vertex->id, samecycle->id, newfacet->id));
2209       qh_setdelsorted(newfacet->vertices, vertex);
2210       vertex->deleted= True;
2211       qh_setappend(qh, &qh->del_vertices, vertex);
2212     }
2213   }
2214   qh_settempfree(qh, &vertices);
2215   trace3((qh, qh->ferr, 3005, "qh_mergecycle_vneighbors: merged vertices from cycle f%d into f%d\n",
2216              samecycle->id, newfacet->id));
2217 } /* mergecycle_vneighbors */
2218 
2219 /*-<a                             href="qh-merge_r.htm#TOC"
2220   >-------------------------------</a><a name="mergefacet">-</a>
2221 
2222   qh_mergefacet(qh, facet1, facet2, mindist, maxdist, mergeapex )
2223     merges facet1 into facet2
2224     mergeapex==qh_MERGEapex if merging new facet into coplanar horizon
2225 
2226   returns:
2227     qh.max_outside and qh.min_vertex updated
2228     initializes vertex neighbors on first merge
2229 
2230   returns:
2231     facet2 contains facet1's vertices, neighbors, and ridges
2232       facet2 moved to end of qh.facet_list
2233       makes facet2 a newfacet
2234       sets facet2->newmerge set
2235       clears facet2->center (unless merging into a large facet)
2236       clears facet2->tested and ridge->tested for facet1
2237 
2238     facet1 prepended to visible_list for later deletion and partitioning
2239       facet1->f.replace == facet2
2240 
2241     adds neighboring facets to facet_mergeset if redundant or degenerate
2242 
2243   notes:
2244     mindist/maxdist may be NULL (only if both NULL)
2245     traces merge if fmax_(maxdist,-mindist) > TRACEdist
2246 
2247   see:
2248     qh_mergecycle()
2249 
2250   design:
2251     trace merge and check for degenerate simplex
2252     make ridges for both facets
2253     update qh.max_outside, qh.max_vertex, qh.min_vertex
2254     update facet2->maxoutside and keepcentrum
2255     update facet2->nummerge
2256     update tested flags for facet2
2257     if facet1 is simplicial
2258       merge facet1 into facet2
2259     else
2260       merge facet1's neighbors into facet2
2261       merge facet1's ridges into facet2
2262       merge facet1's vertices into facet2
2263       merge facet1's vertex neighbors into facet2
2264       add facet2's vertices to qh.new_vertexlist
2265       unless qh_MERGEapex
2266         test facet2 for degenerate or redundant neighbors
2267       move facet1 to qh.visible_list for later deletion
2268       move facet2 to end of qh.newfacet_list
2269 */
qh_mergefacet(qhT * qh,facetT * facet1,facetT * facet2,realT * mindist,realT * maxdist,boolT mergeapex)2270 void qh_mergefacet(qhT *qh, facetT *facet1, facetT *facet2, realT *mindist, realT *maxdist, boolT mergeapex) {
2271   boolT traceonce= False;
2272   vertexT *vertex, **vertexp;
2273   int tracerestore=0, nummerge;
2274 
2275   if (facet1->tricoplanar || facet2->tricoplanar) {
2276     if (!qh->TRInormals) {
2277       qh_fprintf(qh, qh->ferr, 6226, "Qhull internal error (qh_mergefacet): does not work for tricoplanar facets.  Use option 'Q11'\n");
2278       qh_errexit2(qh, qh_ERRqhull, facet1, facet2);
2279     }
2280     if (facet2->tricoplanar) {
2281       facet2->tricoplanar= False;
2282       facet2->keepcentrum= False;
2283     }
2284   }
2285   zzinc_(Ztotmerge);
2286   if (qh->REPORTfreq2 && qh->POSTmerging) {
2287     if (zzval_(Ztotmerge) > qh->mergereport + qh->REPORTfreq2)
2288       qh_tracemerging(qh);
2289   }
2290 #ifndef qh_NOtrace
2291   if (qh->build_cnt >= qh->RERUN) {
2292     if (mindist && (-*mindist > qh->TRACEdist || *maxdist > qh->TRACEdist)) {
2293       tracerestore= 0;
2294       qh->IStracing= qh->TRACElevel;
2295       traceonce= True;
2296       qh_fprintf(qh, qh->ferr, 8075, "qh_mergefacet: ========= trace wide merge #%d(%2.2g) for f%d into f%d, last point was p%d\n", zzval_(Ztotmerge),
2297              fmax_(-*mindist, *maxdist), facet1->id, facet2->id, qh->furthest_id);
2298     }else if (facet1 == qh->tracefacet || facet2 == qh->tracefacet) {
2299       tracerestore= qh->IStracing;
2300       qh->IStracing= 4;
2301       traceonce= True;
2302       qh_fprintf(qh, qh->ferr, 8076, "qh_mergefacet: ========= trace merge #%d involving f%d, furthest is p%d\n",
2303                  zzval_(Ztotmerge), qh->tracefacet_id,  qh->furthest_id);
2304     }
2305   }
2306   if (qh->IStracing >= 2) {
2307     realT mergemin= -2;
2308     realT mergemax= -2;
2309 
2310     if (mindist) {
2311       mergemin= *mindist;
2312       mergemax= *maxdist;
2313     }
2314     qh_fprintf(qh, qh->ferr, 8077, "qh_mergefacet: #%d merge f%d into f%d, mindist= %2.2g, maxdist= %2.2g\n",
2315     zzval_(Ztotmerge), facet1->id, facet2->id, mergemin, mergemax);
2316   }
2317 #endif /* !qh_NOtrace */
2318   if (facet1 == facet2 || facet1->visible || facet2->visible) {
2319     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\n",
2320              facet1->id, facet2->id);
2321     qh_errexit2(qh, qh_ERRqhull, facet1, facet2);
2322   }
2323   if (qh->num_facets - qh->num_visible <= qh->hull_dim + 1) {
2324     qh_fprintf(qh, qh->ferr, 6227, "\n\
2325 qhull precision error: Only %d facets remain.  Can not merge another\n\
2326 pair.  The input is too degenerate or the convexity constraints are\n\
2327 too strong.\n", qh->hull_dim+1);
2328     if (qh->hull_dim >= 5 && !qh->MERGEexact)
2329       qh_fprintf(qh, qh->ferr, 8079, "Option 'Qx' may avoid this problem.\n");
2330     qh_errexit(qh, qh_ERRprec, NULL, NULL);
2331   }
2332   if (!qh->VERTEXneighbors)
2333     qh_vertexneighbors(qh);
2334   qh_makeridges(qh, facet1);
2335   qh_makeridges(qh, facet2);
2336   if (qh->IStracing >=4)
2337     qh_errprint(qh, "MERGING", facet1, facet2, NULL, NULL);
2338   if (mindist) {
2339     maximize_(qh->max_outside, *maxdist);
2340     maximize_(qh->max_vertex, *maxdist);
2341 #if qh_MAXoutside
2342     maximize_(facet2->maxoutside, *maxdist);
2343 #endif
2344     minimize_(qh->min_vertex, *mindist);
2345     if (!facet2->keepcentrum
2346     && (*maxdist > qh->WIDEfacet || *mindist < -qh->WIDEfacet)) {
2347       facet2->keepcentrum= True;
2348       zinc_(Zwidefacet);
2349     }
2350   }
2351   nummerge= facet1->nummerge + facet2->nummerge + 1;
2352   if (nummerge >= qh_MAXnummerge)
2353     facet2->nummerge= qh_MAXnummerge;
2354   else
2355     facet2->nummerge= (short unsigned int)nummerge;
2356   facet2->newmerge= True;
2357   facet2->dupridge= False;
2358   qh_updatetested(qh, facet1, facet2);
2359   if (qh->hull_dim > 2 && qh_setsize(qh, facet1->vertices) == qh->hull_dim)
2360     qh_mergesimplex(qh, facet1, facet2, mergeapex);
2361   else {
2362     qh->vertex_visit++;
2363     FOREACHvertex_(facet2->vertices)
2364       vertex->visitid= qh->vertex_visit;
2365     if (qh->hull_dim == 2)
2366       qh_mergefacet2d(qh, facet1, facet2);
2367     else {
2368       qh_mergeneighbors(qh, facet1, facet2);
2369       qh_mergevertices(qh, facet1->vertices, &facet2->vertices);
2370     }
2371     qh_mergeridges(qh, facet1, facet2);
2372     qh_mergevertex_neighbors(qh, facet1, facet2);
2373     if (!facet2->newfacet)
2374       qh_newvertices(qh, facet2->vertices);
2375   }
2376   if (!mergeapex)
2377     qh_degen_redundant_neighbors(qh, facet2, facet1);
2378   if (facet2->coplanar || !facet2->newfacet) {
2379     zinc_(Zmergeintohorizon);
2380   }else if (!facet1->newfacet && facet2->newfacet) {
2381     zinc_(Zmergehorizon);
2382   }else {
2383     zinc_(Zmergenew);
2384   }
2385   qh_willdelete(qh, facet1, facet2);
2386   qh_removefacet(qh, facet2);  /* append as a newfacet to end of qh->facet_list */
2387   qh_appendfacet(qh, facet2);
2388   facet2->newfacet= True;
2389   facet2->tested= False;
2390   qh_tracemerge(qh, facet1, facet2);
2391   if (traceonce) {
2392     qh_fprintf(qh, qh->ferr, 8080, "qh_mergefacet: end of wide tracing\n");
2393     qh->IStracing= tracerestore;
2394   }
2395 } /* mergefacet */
2396 
2397 
2398 /*-<a                             href="qh-merge_r.htm#TOC"
2399   >-------------------------------</a><a name="mergefacet2d">-</a>
2400 
2401   qh_mergefacet2d(qh, facet1, facet2 )
2402     in 2d, merges neighbors and vertices of facet1 into facet2
2403 
2404   returns:
2405     build ridges for neighbors if necessary
2406     facet2 looks like a simplicial facet except for centrum, ridges
2407       neighbors are opposite the corresponding vertex
2408       maintains orientation of facet2
2409 
2410   notes:
2411     qh_mergefacet() retains non-simplicial structures
2412       they are not needed in 2d, but later routines may use them
2413     preserves qh.vertex_visit for qh_mergevertex_neighbors()
2414 
2415   design:
2416     get vertices and neighbors
2417     determine new vertices and neighbors
2418     set new vertices and neighbors and adjust orientation
2419     make ridges for new neighbor if needed
2420 */
qh_mergefacet2d(qhT * qh,facetT * facet1,facetT * facet2)2421 void qh_mergefacet2d(qhT *qh, facetT *facet1, facetT *facet2) {
2422   vertexT *vertex1A, *vertex1B, *vertex2A, *vertex2B, *vertexA, *vertexB;
2423   facetT *neighbor1A, *neighbor1B, *neighbor2A, *neighbor2B, *neighborA, *neighborB;
2424 
2425   vertex1A= SETfirstt_(facet1->vertices, vertexT);
2426   vertex1B= SETsecondt_(facet1->vertices, vertexT);
2427   vertex2A= SETfirstt_(facet2->vertices, vertexT);
2428   vertex2B= SETsecondt_(facet2->vertices, vertexT);
2429   neighbor1A= SETfirstt_(facet1->neighbors, facetT);
2430   neighbor1B= SETsecondt_(facet1->neighbors, facetT);
2431   neighbor2A= SETfirstt_(facet2->neighbors, facetT);
2432   neighbor2B= SETsecondt_(facet2->neighbors, facetT);
2433   if (vertex1A == vertex2A) {
2434     vertexA= vertex1B;
2435     vertexB= vertex2B;
2436     neighborA= neighbor2A;
2437     neighborB= neighbor1A;
2438   }else if (vertex1A == vertex2B) {
2439     vertexA= vertex1B;
2440     vertexB= vertex2A;
2441     neighborA= neighbor2B;
2442     neighborB= neighbor1A;
2443   }else if (vertex1B == vertex2A) {
2444     vertexA= vertex1A;
2445     vertexB= vertex2B;
2446     neighborA= neighbor2A;
2447     neighborB= neighbor1B;
2448   }else { /* 1B == 2B */
2449     vertexA= vertex1A;
2450     vertexB= vertex2A;
2451     neighborA= neighbor2B;
2452     neighborB= neighbor1B;
2453   }
2454   /* vertexB always from facet2, neighborB always from facet1 */
2455   if (vertexA->id > vertexB->id) {
2456     SETfirst_(facet2->vertices)= vertexA;
2457     SETsecond_(facet2->vertices)= vertexB;
2458     if (vertexB == vertex2A)
2459       facet2->toporient= !facet2->toporient;
2460     SETfirst_(facet2->neighbors)= neighborA;
2461     SETsecond_(facet2->neighbors)= neighborB;
2462   }else {
2463     SETfirst_(facet2->vertices)= vertexB;
2464     SETsecond_(facet2->vertices)= vertexA;
2465     if (vertexB == vertex2B)
2466       facet2->toporient= !facet2->toporient;
2467     SETfirst_(facet2->neighbors)= neighborB;
2468     SETsecond_(facet2->neighbors)= neighborA;
2469   }
2470   qh_makeridges(qh, neighborB);
2471   qh_setreplace(qh, neighborB->neighbors, facet1, facet2);
2472   trace4((qh, qh->ferr, 4036, "qh_mergefacet2d: merged v%d and neighbor f%d of f%d into f%d\n",
2473        vertexA->id, neighborB->id, facet1->id, facet2->id));
2474 } /* mergefacet2d */
2475 
2476 
2477 /*-<a                             href="qh-merge_r.htm#TOC"
2478   >-------------------------------</a><a name="mergeneighbors">-</a>
2479 
2480   qh_mergeneighbors(qh, facet1, facet2 )
2481     merges the neighbors of facet1 into facet2
2482 
2483   see:
2484     qh_mergecycle_neighbors()
2485 
2486   design:
2487     for each neighbor of facet1
2488       if neighbor is also a neighbor of facet2
2489         if neighbor is simpilicial
2490           make ridges for later deletion as a degenerate facet
2491         update its neighbor set
2492       else
2493         move the neighbor relation to facet2
2494     remove the neighbor relation for facet1 and facet2
2495 */
qh_mergeneighbors(qhT * qh,facetT * facet1,facetT * facet2)2496 void qh_mergeneighbors(qhT *qh, facetT *facet1, facetT *facet2) {
2497   facetT *neighbor, **neighborp;
2498 
2499   trace4((qh, qh->ferr, 4037, "qh_mergeneighbors: merge neighbors of f%d and f%d\n",
2500           facet1->id, facet2->id));
2501   qh->visit_id++;
2502   FOREACHneighbor_(facet2) {
2503     neighbor->visitid= qh->visit_id;
2504   }
2505   FOREACHneighbor_(facet1) {
2506     if (neighbor->visitid == qh->visit_id) {
2507       if (neighbor->simplicial)    /* is degen, needs ridges */
2508         qh_makeridges(qh, neighbor);
2509       if (SETfirstt_(neighbor->neighbors, facetT) != facet1) /*keep newfacet->horizon*/
2510         qh_setdel(neighbor->neighbors, facet1);
2511       else {
2512         qh_setdel(neighbor->neighbors, facet2);
2513         qh_setreplace(qh, neighbor->neighbors, facet1, facet2);
2514       }
2515     }else if (neighbor != facet2) {
2516       qh_setappend(qh, &(facet2->neighbors), neighbor);
2517       qh_setreplace(qh, neighbor->neighbors, facet1, facet2);
2518     }
2519   }
2520   qh_setdel(facet1->neighbors, facet2);  /* here for makeridges */
2521   qh_setdel(facet2->neighbors, facet1);
2522 } /* mergeneighbors */
2523 
2524 
2525 /*-<a                             href="qh-merge_r.htm#TOC"
2526   >-------------------------------</a><a name="mergeridges">-</a>
2527 
2528   qh_mergeridges(qh, facet1, facet2 )
2529     merges the ridge set of facet1 into facet2
2530 
2531   returns:
2532     may delete all ridges for a vertex
2533     sets vertex->delridge on deleted ridges
2534 
2535   see:
2536     qh_mergecycle_ridges()
2537 
2538   design:
2539     delete ridges between facet1 and facet2
2540       mark (delridge) vertices on these ridges for later testing
2541     for each remaining ridge
2542       rename facet1 to facet2
2543 */
qh_mergeridges(qhT * qh,facetT * facet1,facetT * facet2)2544 void qh_mergeridges(qhT *qh, facetT *facet1, facetT *facet2) {
2545   ridgeT *ridge, **ridgep;
2546   vertexT *vertex, **vertexp;
2547 
2548   trace4((qh, qh->ferr, 4038, "qh_mergeridges: merge ridges of f%d and f%d\n",
2549           facet1->id, facet2->id));
2550   FOREACHridge_(facet2->ridges) {
2551     if ((ridge->top == facet1) || (ridge->bottom == facet1)) {
2552       FOREACHvertex_(ridge->vertices)
2553         vertex->delridge= True;
2554       qh_delridge(qh, ridge);  /* expensive in high-d, could rebuild */
2555       ridgep--; /*repeat*/
2556     }
2557   }
2558   FOREACHridge_(facet1->ridges) {
2559     if (ridge->top == facet1)
2560       ridge->top= facet2;
2561     else
2562       ridge->bottom= facet2;
2563     qh_setappend(qh, &(facet2->ridges), ridge);
2564   }
2565 } /* mergeridges */
2566 
2567 
2568 /*-<a                             href="qh-merge_r.htm#TOC"
2569   >-------------------------------</a><a name="mergesimplex">-</a>
2570 
2571   qh_mergesimplex(qh, facet1, facet2, mergeapex )
2572     merge simplicial facet1 into facet2
2573     mergeapex==qh_MERGEapex if merging samecycle into horizon facet
2574       vertex id is latest (most recently created)
2575     facet1 may be contained in facet2
2576     ridges exist for both facets
2577 
2578   returns:
2579     facet2 with updated vertices, ridges, neighbors
2580     updated neighbors for facet1's vertices
2581     facet1 not deleted
2582     sets vertex->delridge on deleted ridges
2583 
2584   notes:
2585     special case code since this is the most common merge
2586     called from qh_mergefacet()
2587 
2588   design:
2589     if qh_MERGEapex
2590       add vertices of facet2 to qh.new_vertexlist if necessary
2591       add apex to facet2
2592     else
2593       for each ridge between facet1 and facet2
2594         set vertex->delridge
2595       determine the apex for facet1 (i.e., vertex to be merged)
2596       unless apex already in facet2
2597         insert apex into vertices for facet2
2598       add vertices of facet2 to qh.new_vertexlist if necessary
2599       add apex to qh.new_vertexlist if necessary
2600       for each vertex of facet1
2601         if apex
2602           rename facet1 to facet2 in its vertex neighbors
2603         else
2604           delete facet1 from vertex neighors
2605           if only in facet2
2606             add vertex to qh.del_vertices for later deletion
2607       for each ridge of facet1
2608         delete ridges between facet1 and facet2
2609         append other ridges to facet2 after renaming facet to facet2
2610 */
qh_mergesimplex(qhT * qh,facetT * facet1,facetT * facet2,boolT mergeapex)2611 void qh_mergesimplex(qhT *qh, facetT *facet1, facetT *facet2, boolT mergeapex) {
2612   vertexT *vertex, **vertexp, *apex;
2613   ridgeT *ridge, **ridgep;
2614   boolT issubset= False;
2615   int vertex_i= -1, vertex_n;
2616   facetT *neighbor, **neighborp, *otherfacet;
2617 
2618   if (mergeapex) {
2619     if (!facet2->newfacet)
2620       qh_newvertices(qh, facet2->vertices);  /* apex is new */
2621     apex= SETfirstt_(facet1->vertices, vertexT);
2622     if (SETfirstt_(facet2->vertices, vertexT) != apex)
2623       qh_setaddnth(qh, &facet2->vertices, 0, apex);  /* apex has last id */
2624     else
2625       issubset= True;
2626   }else {
2627     zinc_(Zmergesimplex);
2628     FOREACHvertex_(facet1->vertices)
2629       vertex->seen= False;
2630     FOREACHridge_(facet1->ridges) {
2631       if (otherfacet_(ridge, facet1) == facet2) {
2632         FOREACHvertex_(ridge->vertices) {
2633           vertex->seen= True;
2634           vertex->delridge= True;
2635         }
2636         break;
2637       }
2638     }
2639     FOREACHvertex_(facet1->vertices) {
2640       if (!vertex->seen)
2641         break;  /* must occur */
2642     }
2643     apex= vertex;
2644     trace4((qh, qh->ferr, 4039, "qh_mergesimplex: merge apex v%d of f%d into facet f%d\n",
2645           apex->id, facet1->id, facet2->id));
2646     FOREACHvertex_i_(qh, facet2->vertices) {
2647       if (vertex->id < apex->id) {
2648         break;
2649       }else if (vertex->id == apex->id) {
2650         issubset= True;
2651         break;
2652       }
2653     }
2654     if (!issubset)
2655       qh_setaddnth(qh, &facet2->vertices, vertex_i, apex);
2656     if (!facet2->newfacet)
2657       qh_newvertices(qh, facet2->vertices);
2658     else if (!apex->newlist) {
2659       qh_removevertex(qh, apex);
2660       qh_appendvertex(qh, apex);
2661     }
2662   }
2663   trace4((qh, qh->ferr, 4040, "qh_mergesimplex: update vertex neighbors of f%d\n",
2664           facet1->id));
2665   FOREACHvertex_(facet1->vertices) {
2666     if (vertex == apex && !issubset)
2667       qh_setreplace(qh, vertex->neighbors, facet1, facet2);
2668     else {
2669       qh_setdel(vertex->neighbors, facet1);
2670       if (!SETsecond_(vertex->neighbors))
2671         qh_mergevertex_del(qh, vertex, facet1, facet2);
2672     }
2673   }
2674   trace4((qh, qh->ferr, 4041, "qh_mergesimplex: merge ridges and neighbors of f%d into f%d\n",
2675           facet1->id, facet2->id));
2676   qh->visit_id++;
2677   FOREACHneighbor_(facet2)
2678     neighbor->visitid= qh->visit_id;
2679   FOREACHridge_(facet1->ridges) {
2680     otherfacet= otherfacet_(ridge, facet1);
2681     if (otherfacet == facet2) {
2682       qh_setdel(facet2->ridges, ridge);
2683       qh_setfree(qh, &(ridge->vertices));
2684       qh_memfree(qh, ridge, (int)sizeof(ridgeT));
2685       qh_setdel(facet2->neighbors, facet1);
2686     }else {
2687       qh_setappend(qh, &facet2->ridges, ridge);
2688       if (otherfacet->visitid != qh->visit_id) {
2689         qh_setappend(qh, &facet2->neighbors, otherfacet);
2690         qh_setreplace(qh, otherfacet->neighbors, facet1, facet2);
2691         otherfacet->visitid= qh->visit_id;
2692       }else {
2693         if (otherfacet->simplicial)    /* is degen, needs ridges */
2694           qh_makeridges(qh, otherfacet);
2695         if (SETfirstt_(otherfacet->neighbors, facetT) != facet1)
2696           qh_setdel(otherfacet->neighbors, facet1);
2697         else {   /*keep newfacet->neighbors->horizon*/
2698           qh_setdel(otherfacet->neighbors, facet2);
2699           qh_setreplace(qh, otherfacet->neighbors, facet1, facet2);
2700         }
2701       }
2702       if (ridge->top == facet1) /* wait until after qh_makeridges */
2703         ridge->top= facet2;
2704       else
2705         ridge->bottom= facet2;
2706     }
2707   }
2708   SETfirst_(facet1->ridges)= NULL; /* it will be deleted */
2709   trace3((qh, qh->ferr, 3006, "qh_mergesimplex: merged simplex f%d apex v%d into facet f%d\n",
2710           facet1->id, getid_(apex), facet2->id));
2711 } /* mergesimplex */
2712 
2713 /*-<a                             href="qh-merge_r.htm#TOC"
2714   >-------------------------------</a><a name="mergevertex_del">-</a>
2715 
2716   qh_mergevertex_del(qh, vertex, facet1, facet2 )
2717     delete a vertex because of merging facet1 into facet2
2718 
2719   returns:
2720     deletes vertex from facet2
2721     adds vertex to qh.del_vertices for later deletion
2722 */
qh_mergevertex_del(qhT * qh,vertexT * vertex,facetT * facet1,facetT * facet2)2723 void qh_mergevertex_del(qhT *qh, vertexT *vertex, facetT *facet1, facetT *facet2) {
2724 
2725   zinc_(Zmergevertex);
2726   trace2((qh, qh->ferr, 2035, "qh_mergevertex_del: deleted v%d when merging f%d into f%d\n",
2727           vertex->id, facet1->id, facet2->id));
2728   qh_setdelsorted(facet2->vertices, vertex);
2729   vertex->deleted= True;
2730   qh_setappend(qh, &qh->del_vertices, vertex);
2731 } /* mergevertex_del */
2732 
2733 /*-<a                             href="qh-merge_r.htm#TOC"
2734   >-------------------------------</a><a name="mergevertex_neighbors">-</a>
2735 
2736   qh_mergevertex_neighbors(qh, facet1, facet2 )
2737     merge the vertex neighbors of facet1 to facet2
2738 
2739   returns:
2740     if vertex is current qh.vertex_visit
2741       deletes facet1 from vertex->neighbors
2742     else
2743       renames facet1 to facet2 in vertex->neighbors
2744     deletes vertices if only one neighbor
2745 
2746   notes:
2747     assumes vertex neighbor sets are good
2748 */
qh_mergevertex_neighbors(qhT * qh,facetT * facet1,facetT * facet2)2749 void qh_mergevertex_neighbors(qhT *qh, facetT *facet1, facetT *facet2) {
2750   vertexT *vertex, **vertexp;
2751 
2752   trace4((qh, qh->ferr, 4042, "qh_mergevertex_neighbors: merge vertex neighbors of f%d and f%d\n",
2753           facet1->id, facet2->id));
2754   if (qh->tracevertex) {
2755     qh_fprintf(qh, qh->ferr, 8081, "qh_mergevertex_neighbors: of f%d and f%d at furthest p%d f0= %p\n",
2756              facet1->id, facet2->id, qh->furthest_id, qh->tracevertex->neighbors->e[0].p);
2757     qh_errprint(qh, "TRACE", NULL, NULL, NULL, qh->tracevertex);
2758   }
2759   FOREACHvertex_(facet1->vertices) {
2760     if (vertex->visitid != qh->vertex_visit)
2761       qh_setreplace(qh, vertex->neighbors, facet1, facet2);
2762     else {
2763       qh_setdel(vertex->neighbors, facet1);
2764       if (!SETsecond_(vertex->neighbors))
2765         qh_mergevertex_del(qh, vertex, facet1, facet2);
2766     }
2767   }
2768   if (qh->tracevertex)
2769     qh_errprint(qh, "TRACE", NULL, NULL, NULL, qh->tracevertex);
2770 } /* mergevertex_neighbors */
2771 
2772 
2773 /*-<a                             href="qh-merge_r.htm#TOC"
2774   >-------------------------------</a><a name="mergevertices">-</a>
2775 
2776   qh_mergevertices(qh, vertices1, vertices2 )
2777     merges the vertex set of facet1 into facet2
2778 
2779   returns:
2780     replaces vertices2 with merged set
2781     preserves vertex_visit for qh_mergevertex_neighbors
2782     updates qh.newvertex_list
2783 
2784   design:
2785     create a merged set of both vertices (in inverse id order)
2786 */
qh_mergevertices(qhT * qh,setT * vertices1,setT ** vertices2)2787 void qh_mergevertices(qhT *qh, setT *vertices1, setT **vertices2) {
2788   int newsize= qh_setsize(qh, vertices1)+qh_setsize(qh, *vertices2) - qh->hull_dim + 1;
2789   setT *mergedvertices;
2790   vertexT *vertex, **vertexp, **vertex2= SETaddr_(*vertices2, vertexT);
2791 
2792   mergedvertices= qh_settemp(qh, newsize);
2793   FOREACHvertex_(vertices1) {
2794     if (!*vertex2 || vertex->id > (*vertex2)->id)
2795       qh_setappend(qh, &mergedvertices, vertex);
2796     else {
2797       while (*vertex2 && (*vertex2)->id > vertex->id)
2798         qh_setappend(qh, &mergedvertices, *vertex2++);
2799       if (!*vertex2 || (*vertex2)->id < vertex->id)
2800         qh_setappend(qh, &mergedvertices, vertex);
2801       else
2802         qh_setappend(qh, &mergedvertices, *vertex2++);
2803     }
2804   }
2805   while (*vertex2)
2806     qh_setappend(qh, &mergedvertices, *vertex2++);
2807   if (newsize < qh_setsize(qh, mergedvertices)) {
2808     qh_fprintf(qh, qh->ferr, 6100, "qhull internal error (qh_mergevertices): facets did not share a ridge\n");
2809     qh_errexit(qh, qh_ERRqhull, NULL, NULL);
2810   }
2811   qh_setfree(qh, vertices2);
2812   *vertices2= mergedvertices;
2813   qh_settemppop(qh);
2814 } /* mergevertices */
2815 
2816 
2817 /*-<a                             href="qh-merge_r.htm#TOC"
2818   >-------------------------------</a><a name="neighbor_intersections">-</a>
2819 
2820   qh_neighbor_intersections(qh, vertex )
2821     return intersection of all vertices in vertex->neighbors except for vertex
2822 
2823   returns:
2824     returns temporary set of vertices
2825     does not include vertex
2826     NULL if a neighbor is simplicial
2827     NULL if empty set
2828 
2829   notes:
2830     used for renaming vertices
2831 
2832   design:
2833     initialize the intersection set with vertices of the first two neighbors
2834     delete vertex from the intersection
2835     for each remaining neighbor
2836       intersect its vertex set with the intersection set
2837       return NULL if empty
2838     return the intersection set
2839 */
qh_neighbor_intersections(qhT * qh,vertexT * vertex)2840 setT *qh_neighbor_intersections(qhT *qh, vertexT *vertex) {
2841   facetT *neighbor, **neighborp, *neighborA, *neighborB;
2842   setT *intersect;
2843   int neighbor_i, neighbor_n;
2844 
2845   FOREACHneighbor_(vertex) {
2846     if (neighbor->simplicial)
2847       return NULL;
2848   }
2849   neighborA= SETfirstt_(vertex->neighbors, facetT);
2850   neighborB= SETsecondt_(vertex->neighbors, facetT);
2851   zinc_(Zintersectnum);
2852   if (!neighborA)
2853     return NULL;
2854   if (!neighborB)
2855     intersect= qh_setcopy(qh, neighborA->vertices, 0);
2856   else
2857     intersect= qh_vertexintersect_new(qh, neighborA->vertices, neighborB->vertices);
2858   qh_settemppush(qh, intersect);
2859   qh_setdelsorted(intersect, vertex);
2860   FOREACHneighbor_i_(qh, vertex) {
2861     if (neighbor_i >= 2) {
2862       zinc_(Zintersectnum);
2863       qh_vertexintersect(qh, &intersect, neighbor->vertices);
2864       if (!SETfirst_(intersect)) {
2865         zinc_(Zintersectfail);
2866         qh_settempfree(qh, &intersect);
2867         return NULL;
2868       }
2869     }
2870   }
2871   trace3((qh, qh->ferr, 3007, "qh_neighbor_intersections: %d vertices in neighbor intersection of v%d\n",
2872           qh_setsize(qh, intersect), vertex->id));
2873   return intersect;
2874 } /* neighbor_intersections */
2875 
2876 /*-<a                             href="qh-merge_r.htm#TOC"
2877   >-------------------------------</a><a name="newvertices">-</a>
2878 
2879   qh_newvertices(qh, vertices )
2880     add vertices to end of qh.vertex_list (marks as new vertices)
2881 
2882   returns:
2883     vertices on qh.newvertex_list
2884     vertex->newlist set
2885 */
qh_newvertices(qhT * qh,setT * vertices)2886 void qh_newvertices(qhT *qh, setT *vertices) {
2887   vertexT *vertex, **vertexp;
2888 
2889   FOREACHvertex_(vertices) {
2890     if (!vertex->newlist) {
2891       qh_removevertex(qh, vertex);
2892       qh_appendvertex(qh, vertex);
2893     }
2894   }
2895 } /* newvertices */
2896 
2897 /*-<a                             href="qh-merge_r.htm#TOC"
2898   >-------------------------------</a><a name="reducevertices">-</a>
2899 
2900   qh_reducevertices(qh)
2901     reduce extra vertices, shared vertices, and redundant vertices
2902     facet->newmerge is set if merged since last call
2903     if !qh.MERGEvertices, only removes extra vertices
2904 
2905   returns:
2906     True if also merged degen_redundant facets
2907     vertices are renamed if possible
2908     clears facet->newmerge and vertex->delridge
2909 
2910   notes:
2911     ignored if 2-d
2912 
2913   design:
2914     merge any degenerate or redundant facets
2915     for each newly merged facet
2916       remove extra vertices
2917     if qh.MERGEvertices
2918       for each newly merged facet
2919         for each vertex
2920           if vertex was on a deleted ridge
2921             rename vertex if it is shared
2922       remove delridge flag from new vertices
2923 */
qh_reducevertices(qhT * qh)2924 boolT qh_reducevertices(qhT *qh) {
2925   int numshare=0, numrename= 0;
2926   boolT degenredun= False;
2927   facetT *newfacet;
2928   vertexT *vertex, **vertexp;
2929 
2930   if (qh->hull_dim == 2)
2931     return False;
2932   if (qh_merge_degenredundant(qh))
2933     degenredun= True;
2934  LABELrestart:
2935   FORALLnew_facets {
2936     if (newfacet->newmerge) {
2937       if (!qh->MERGEvertices)
2938         newfacet->newmerge= False;
2939       qh_remove_extravertices(qh, newfacet);
2940     }
2941   }
2942   if (!qh->MERGEvertices)
2943     return False;
2944   FORALLnew_facets {
2945     if (newfacet->newmerge) {
2946       newfacet->newmerge= False;
2947       FOREACHvertex_(newfacet->vertices) {
2948         if (vertex->delridge) {
2949           if (qh_rename_sharedvertex(qh, vertex, newfacet)) {
2950             numshare++;
2951             vertexp--; /* repeat since deleted vertex */
2952           }
2953         }
2954       }
2955     }
2956   }
2957   FORALLvertex_(qh->newvertex_list) {
2958     if (vertex->delridge && !vertex->deleted) {
2959       vertex->delridge= False;
2960       if (qh->hull_dim >= 4 && qh_redundant_vertex(qh, vertex)) {
2961         numrename++;
2962         if (qh_merge_degenredundant(qh)) {
2963           degenredun= True;
2964           goto LABELrestart;
2965         }
2966       }
2967     }
2968   }
2969   trace1((qh, qh->ferr, 1014, "qh_reducevertices: renamed %d shared vertices and %d redundant vertices. Degen? %d\n",
2970           numshare, numrename, degenredun));
2971   return degenredun;
2972 } /* reducevertices */
2973 
2974 /*-<a                             href="qh-merge_r.htm#TOC"
2975   >-------------------------------</a><a name="redundant_vertex">-</a>
2976 
2977   qh_redundant_vertex(qh, vertex )
2978     detect and rename a redundant vertex
2979     vertices have full vertex->neighbors
2980 
2981   returns:
2982     returns true if find a redundant vertex
2983       deletes vertex(vertex->deleted)
2984 
2985   notes:
2986     only needed if vertex->delridge and hull_dim >= 4
2987     may add degenerate facets to qh.facet_mergeset
2988     doesn't change vertex->neighbors or create redundant facets
2989 
2990   design:
2991     intersect vertices of all facet neighbors of vertex
2992     determine ridges for these vertices
2993     if find a new vertex for vertex amoung these ridges and vertices
2994       rename vertex to the new vertex
2995 */
qh_redundant_vertex(qhT * qh,vertexT * vertex)2996 vertexT *qh_redundant_vertex(qhT *qh, vertexT *vertex) {
2997   vertexT *newvertex= NULL;
2998   setT *vertices, *ridges;
2999 
3000   trace3((qh, qh->ferr, 3008, "qh_redundant_vertex: check if v%d can be renamed\n", vertex->id));
3001   if ((vertices= qh_neighbor_intersections(qh, vertex))) {
3002     ridges= qh_vertexridges(qh, vertex);
3003     if ((newvertex= qh_find_newvertex(qh, vertex, vertices, ridges)))
3004       qh_renamevertex(qh, vertex, newvertex, ridges, NULL, NULL);
3005     qh_settempfree(qh, &ridges);
3006     qh_settempfree(qh, &vertices);
3007   }
3008   return newvertex;
3009 } /* redundant_vertex */
3010 
3011 /*-<a                             href="qh-merge_r.htm#TOC"
3012   >-------------------------------</a><a name="remove_extravertices">-</a>
3013 
3014   qh_remove_extravertices(qh, facet )
3015     remove extra vertices from non-simplicial facets
3016 
3017   returns:
3018     returns True if it finds them
3019 
3020   design:
3021     for each vertex in facet
3022       if vertex not in a ridge (i.e., no longer used)
3023         delete vertex from facet
3024         delete facet from vertice's neighbors
3025         unless vertex in another facet
3026           add vertex to qh.del_vertices for later deletion
3027 */
qh_remove_extravertices(qhT * qh,facetT * facet)3028 boolT qh_remove_extravertices(qhT *qh, facetT *facet) {
3029   ridgeT *ridge, **ridgep;
3030   vertexT *vertex, **vertexp;
3031   boolT foundrem= False;
3032 
3033   trace4((qh, qh->ferr, 4043, "qh_remove_extravertices: test f%d for extra vertices\n",
3034           facet->id));
3035   FOREACHvertex_(facet->vertices)
3036     vertex->seen= False;
3037   FOREACHridge_(facet->ridges) {
3038     FOREACHvertex_(ridge->vertices)
3039       vertex->seen= True;
3040   }
3041   FOREACHvertex_(facet->vertices) {
3042     if (!vertex->seen) {
3043       foundrem= True;
3044       zinc_(Zremvertex);
3045       qh_setdelsorted(facet->vertices, vertex);
3046       qh_setdel(vertex->neighbors, facet);
3047       if (!qh_setsize(qh, vertex->neighbors)) {
3048         vertex->deleted= True;
3049         qh_setappend(qh, &qh->del_vertices, vertex);
3050         zinc_(Zremvertexdel);
3051         trace2((qh, qh->ferr, 2036, "qh_remove_extravertices: v%d deleted because it's lost all ridges\n", vertex->id));
3052       }else
3053         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));
3054       vertexp--; /*repeat*/
3055     }
3056   }
3057   return foundrem;
3058 } /* remove_extravertices */
3059 
3060 /*-<a                             href="qh-merge_r.htm#TOC"
3061   >-------------------------------</a><a name="rename_sharedvertex">-</a>
3062 
3063   qh_rename_sharedvertex(qh, vertex, facet )
3064     detect and rename if shared vertex in facet
3065     vertices have full ->neighbors
3066 
3067   returns:
3068     newvertex or NULL
3069     the vertex may still exist in other facets (i.e., a neighbor was pinched)
3070     does not change facet->neighbors
3071     updates vertex->neighbors
3072 
3073   notes:
3074     a shared vertex for a facet is only in ridges to one neighbor
3075     this may undo a pinched facet
3076 
3077     it does not catch pinches involving multiple facets.  These appear
3078       to be difficult to detect, since an exhaustive search is too expensive.
3079 
3080   design:
3081     if vertex only has two neighbors
3082       determine the ridges that contain the vertex
3083       determine the vertices shared by both neighbors
3084       if can find a new vertex in this set
3085         rename the vertex to the new vertex
3086 */
qh_rename_sharedvertex(qhT * qh,vertexT * vertex,facetT * facet)3087 vertexT *qh_rename_sharedvertex(qhT *qh, vertexT *vertex, facetT *facet) {
3088   facetT *neighbor, **neighborp, *neighborA= NULL;
3089   setT *vertices, *ridges;
3090   vertexT *newvertex;
3091 
3092   if (qh_setsize(qh, vertex->neighbors) == 2) {
3093     neighborA= SETfirstt_(vertex->neighbors, facetT);
3094     if (neighborA == facet)
3095       neighborA= SETsecondt_(vertex->neighbors, facetT);
3096   }else if (qh->hull_dim == 3)
3097     return NULL;
3098   else {
3099     qh->visit_id++;
3100     FOREACHneighbor_(facet)
3101       neighbor->visitid= qh->visit_id;
3102     FOREACHneighbor_(vertex) {
3103       if (neighbor->visitid == qh->visit_id) {
3104         if (neighborA)
3105           return NULL;
3106         neighborA= neighbor;
3107       }
3108     }
3109     if (!neighborA) {
3110       qh_fprintf(qh, qh->ferr, 6101, "qhull internal error (qh_rename_sharedvertex): v%d's neighbors not in f%d\n",
3111         vertex->id, facet->id);
3112       qh_errprint(qh, "ERRONEOUS", facet, NULL, NULL, vertex);
3113       qh_errexit(qh, qh_ERRqhull, NULL, NULL);
3114     }
3115   }
3116   /* the vertex is shared by facet and neighborA */
3117   ridges= qh_settemp(qh, qh->TEMPsize);
3118   neighborA->visitid= ++qh->visit_id;
3119   qh_vertexridges_facet(qh, vertex, facet, &ridges);
3120   trace2((qh, qh->ferr, 2037, "qh_rename_sharedvertex: p%d(v%d) is shared by f%d(%d ridges) and f%d\n",
3121     qh_pointid(qh, vertex->point), vertex->id, facet->id, qh_setsize(qh, ridges), neighborA->id));
3122   zinc_(Zintersectnum);
3123   vertices= qh_vertexintersect_new(qh, facet->vertices, neighborA->vertices);
3124   qh_setdel(vertices, vertex);
3125   qh_settemppush(qh, vertices);
3126   if ((newvertex= qh_find_newvertex(qh, vertex, vertices, ridges)))
3127     qh_renamevertex(qh, vertex, newvertex, ridges, facet, neighborA);
3128   qh_settempfree(qh, &vertices);
3129   qh_settempfree(qh, &ridges);
3130   return newvertex;
3131 } /* rename_sharedvertex */
3132 
3133 /*-<a                             href="qh-merge_r.htm#TOC"
3134   >-------------------------------</a><a name="renameridgevertex">-</a>
3135 
3136   qh_renameridgevertex(qh, ridge, oldvertex, newvertex )
3137     renames oldvertex as newvertex in ridge
3138 
3139   returns:
3140 
3141   design:
3142     delete oldvertex from ridge
3143     if newvertex already in ridge
3144       copy ridge->noconvex to another ridge if possible
3145       delete the ridge
3146     else
3147       insert newvertex into the ridge
3148       adjust the ridge's orientation
3149 */
qh_renameridgevertex(qhT * qh,ridgeT * ridge,vertexT * oldvertex,vertexT * newvertex)3150 void qh_renameridgevertex(qhT *qh, ridgeT *ridge, vertexT *oldvertex, vertexT *newvertex) {
3151   int nth= 0, oldnth;
3152   facetT *temp;
3153   vertexT *vertex, **vertexp;
3154 
3155   oldnth= qh_setindex(ridge->vertices, oldvertex);
3156   qh_setdelnthsorted(qh, ridge->vertices, oldnth);
3157   FOREACHvertex_(ridge->vertices) {
3158     if (vertex == newvertex) {
3159       zinc_(Zdelridge);
3160       if (ridge->nonconvex) /* only one ridge has nonconvex set */
3161         qh_copynonconvex(qh, ridge);
3162       trace2((qh, qh->ferr, 2038, "qh_renameridgevertex: ridge r%d deleted.  It contained both v%d and v%d\n",
3163         ridge->id, oldvertex->id, newvertex->id));
3164       qh_delridge(qh, ridge);
3165       return;
3166     }
3167     if (vertex->id < newvertex->id)
3168       break;
3169     nth++;
3170   }
3171   qh_setaddnth(qh, &ridge->vertices, nth, newvertex);
3172   if (abs(oldnth - nth)%2) {
3173     trace3((qh, qh->ferr, 3010, "qh_renameridgevertex: swapped the top and bottom of ridge r%d\n",
3174             ridge->id));
3175     temp= ridge->top;
3176     ridge->top= ridge->bottom;
3177     ridge->bottom= temp;
3178   }
3179 } /* renameridgevertex */
3180 
3181 
3182 /*-<a                             href="qh-merge_r.htm#TOC"
3183   >-------------------------------</a><a name="renamevertex">-</a>
3184 
3185   qh_renamevertex(qh, oldvertex, newvertex, ridges, oldfacet, neighborA )
3186     renames oldvertex as newvertex in ridges
3187     gives oldfacet/neighborA if oldvertex is shared between two facets
3188 
3189   returns:
3190     oldvertex may still exist afterwards
3191 
3192 
3193   notes:
3194     can not change neighbors of newvertex (since it's a subset)
3195 
3196   design:
3197     for each ridge in ridges
3198       rename oldvertex to newvertex and delete degenerate ridges
3199     if oldfacet not defined
3200       for each neighbor of oldvertex
3201         delete oldvertex from neighbor's vertices
3202         remove extra vertices from neighbor
3203       add oldvertex to qh.del_vertices
3204     else if oldvertex only between oldfacet and neighborA
3205       delete oldvertex from oldfacet and neighborA
3206       add oldvertex to qh.del_vertices
3207     else oldvertex is in oldfacet and neighborA and other facets (i.e., pinched)
3208       delete oldvertex from oldfacet
3209       delete oldfacet from oldvertice's neighbors
3210       remove extra vertices (e.g., oldvertex) from neighborA
3211 */
qh_renamevertex(qhT * qh,vertexT * oldvertex,vertexT * newvertex,setT * ridges,facetT * oldfacet,facetT * neighborA)3212 void qh_renamevertex(qhT *qh, vertexT *oldvertex, vertexT *newvertex, setT *ridges, facetT *oldfacet, facetT *neighborA) {
3213   facetT *neighbor, **neighborp;
3214   ridgeT *ridge, **ridgep;
3215   boolT istrace= False;
3216 
3217   if (qh->IStracing >= 2 || oldvertex->id == qh->tracevertex_id ||
3218         newvertex->id == qh->tracevertex_id)
3219     istrace= True;
3220   FOREACHridge_(ridges)
3221     qh_renameridgevertex(qh, ridge, oldvertex, newvertex);
3222   if (!oldfacet) {
3223     zinc_(Zrenameall);
3224     if (istrace)
3225       qh_fprintf(qh, qh->ferr, 8082, "qh_renamevertex: renamed v%d to v%d in several facets\n",
3226                oldvertex->id, newvertex->id);
3227     FOREACHneighbor_(oldvertex) {
3228       qh_maydropneighbor(qh, neighbor);
3229       qh_setdelsorted(neighbor->vertices, oldvertex);
3230       if (qh_remove_extravertices(qh, neighbor))
3231         neighborp--; /* neighbor may be deleted */
3232     }
3233     if (!oldvertex->deleted) {
3234       oldvertex->deleted= True;
3235       qh_setappend(qh, &qh->del_vertices, oldvertex);
3236     }
3237   }else if (qh_setsize(qh, oldvertex->neighbors) == 2) {
3238     zinc_(Zrenameshare);
3239     if (istrace)
3240       qh_fprintf(qh, qh->ferr, 8083, "qh_renamevertex: renamed v%d to v%d in oldfacet f%d\n",
3241                oldvertex->id, newvertex->id, oldfacet->id);
3242     FOREACHneighbor_(oldvertex)
3243       qh_setdelsorted(neighbor->vertices, oldvertex);
3244     oldvertex->deleted= True;
3245     qh_setappend(qh, &qh->del_vertices, oldvertex);
3246   }else {
3247     zinc_(Zrenamepinch);
3248     if (istrace || qh->IStracing)
3249       qh_fprintf(qh, qh->ferr, 8084, "qh_renamevertex: renamed pinched v%d to v%d between f%d and f%d\n",
3250                oldvertex->id, newvertex->id, oldfacet->id, neighborA->id);
3251     qh_setdelsorted(oldfacet->vertices, oldvertex);
3252     qh_setdel(oldvertex->neighbors, oldfacet);
3253     qh_remove_extravertices(qh, neighborA);
3254   }
3255 } /* renamevertex */
3256 
3257 
3258 /*-<a                             href="qh-merge_r.htm#TOC"
3259   >-------------------------------</a><a name="test_appendmerge">-</a>
3260 
3261   qh_test_appendmerge(qh, facet, neighbor )
3262     tests facet/neighbor for convexity
3263     appends to mergeset if non-convex
3264     if pre-merging,
3265       nop if qh.SKIPconvex, or qh.MERGEexact and coplanar
3266 
3267   returns:
3268     true if appends facet/neighbor to mergeset
3269     sets facet->center as needed
3270     does not change facet->seen
3271 
3272   design:
3273     if qh.cos_max is defined
3274       if the angle between facet normals is too shallow
3275         append an angle-coplanar merge to qh.mergeset
3276         return True
3277     make facet's centrum if needed
3278     if facet's centrum is above the neighbor
3279       set isconcave
3280     else
3281       if facet's centrum is not below the neighbor
3282         set iscoplanar
3283       make neighbor's centrum if needed
3284       if neighbor's centrum is above the facet
3285         set isconcave
3286       else if neighbor's centrum is not below the facet
3287         set iscoplanar
3288    if isconcave or iscoplanar
3289      get angle if needed
3290      append concave or coplanar merge to qh.mergeset
3291 */
qh_test_appendmerge(qhT * qh,facetT * facet,facetT * neighbor)3292 boolT qh_test_appendmerge(qhT *qh, facetT *facet, facetT *neighbor) {
3293   realT dist, dist2= -REALmax, angle= -REALmax;
3294   boolT isconcave= False, iscoplanar= False, okangle= False;
3295 
3296   if (qh->SKIPconvex && !qh->POSTmerging)
3297     return False;
3298   if ((!qh->MERGEexact || qh->POSTmerging) && qh->cos_max < REALmax/2) {
3299     angle= qh_getangle(qh, facet->normal, neighbor->normal);
3300     zinc_(Zangletests);
3301     if (angle > qh->cos_max) {
3302       zinc_(Zcoplanarangle);
3303       qh_appendmergeset(qh, facet, neighbor, MRGanglecoplanar, &angle);
3304       trace2((qh, qh->ferr, 2039, "qh_test_appendmerge: coplanar angle %4.4g between f%d and f%d\n",
3305          angle, facet->id, neighbor->id));
3306       return True;
3307     }else
3308       okangle= True;
3309   }
3310   if (!facet->center)
3311     facet->center= qh_getcentrum(qh, facet);
3312   zzinc_(Zcentrumtests);
3313   qh_distplane(qh, facet->center, neighbor, &dist);
3314   if (dist > qh->centrum_radius)
3315     isconcave= True;
3316   else {
3317     if (dist > -qh->centrum_radius)
3318       iscoplanar= True;
3319     if (!neighbor->center)
3320       neighbor->center= qh_getcentrum(qh, neighbor);
3321     zzinc_(Zcentrumtests);
3322     qh_distplane(qh, neighbor->center, facet, &dist2);
3323     if (dist2 > qh->centrum_radius)
3324       isconcave= True;
3325     else if (!iscoplanar && dist2 > -qh->centrum_radius)
3326       iscoplanar= True;
3327   }
3328   if (!isconcave && (!iscoplanar || (qh->MERGEexact && !qh->POSTmerging)))
3329     return False;
3330   if (!okangle && qh->ANGLEmerge) {
3331     angle= qh_getangle(qh, facet->normal, neighbor->normal);
3332     zinc_(Zangletests);
3333   }
3334   if (isconcave) {
3335     zinc_(Zconcaveridge);
3336     if (qh->ANGLEmerge)
3337       angle += qh_ANGLEconcave + 0.5;
3338     qh_appendmergeset(qh, facet, neighbor, MRGconcave, &angle);
3339     trace0((qh, qh->ferr, 18, "qh_test_appendmerge: concave f%d to f%d dist %4.4g and reverse dist %4.4g angle %4.4g during p%d\n",
3340            facet->id, neighbor->id, dist, dist2, angle, qh->furthest_id));
3341   }else /* iscoplanar */ {
3342     zinc_(Zcoplanarcentrum);
3343     qh_appendmergeset(qh, facet, neighbor, MRGcoplanar, &angle);
3344     trace2((qh, qh->ferr, 2040, "qh_test_appendmerge: coplanar f%d to f%d dist %4.4g, reverse dist %4.4g angle %4.4g\n",
3345               facet->id, neighbor->id, dist, dist2, angle));
3346   }
3347   return True;
3348 } /* test_appendmerge */
3349 
3350 /*-<a                             href="qh-merge_r.htm#TOC"
3351   >-------------------------------</a><a name="test_vneighbors">-</a>
3352 
3353   qh_test_vneighbors(qh)
3354     test vertex neighbors for convexity
3355     tests all facets on qh.newfacet_list
3356 
3357   returns:
3358     true if non-convex vneighbors appended to qh.facet_mergeset
3359     initializes vertex neighbors if needed
3360 
3361   notes:
3362     assumes all facet neighbors have been tested
3363     this can be expensive
3364     this does not guarantee that a centrum is below all facets
3365       but it is unlikely
3366     uses qh.visit_id
3367 
3368   design:
3369     build vertex neighbors if necessary
3370     for all new facets
3371       for all vertices
3372         for each unvisited facet neighbor of the vertex
3373           test new facet and neighbor for convexity
3374 */
qh_test_vneighbors(qhT * qh)3375 boolT qh_test_vneighbors(qhT *qh /* qh->newfacet_list */) {
3376   facetT *newfacet, *neighbor, **neighborp;
3377   vertexT *vertex, **vertexp;
3378   int nummerges= 0;
3379 
3380   trace1((qh, qh->ferr, 1015, "qh_test_vneighbors: testing vertex neighbors for convexity\n"));
3381   if (!qh->VERTEXneighbors)
3382     qh_vertexneighbors(qh);
3383   FORALLnew_facets
3384     newfacet->seen= False;
3385   FORALLnew_facets {
3386     newfacet->seen= True;
3387     newfacet->visitid= qh->visit_id++;
3388     FOREACHneighbor_(newfacet)
3389       newfacet->visitid= qh->visit_id;
3390     FOREACHvertex_(newfacet->vertices) {
3391       FOREACHneighbor_(vertex) {
3392         if (neighbor->seen || neighbor->visitid == qh->visit_id)
3393           continue;
3394         if (qh_test_appendmerge(qh, newfacet, neighbor))
3395           nummerges++;
3396       }
3397     }
3398   }
3399   zadd_(Ztestvneighbor, nummerges);
3400   trace1((qh, qh->ferr, 1016, "qh_test_vneighbors: found %d non-convex, vertex neighbors\n",
3401            nummerges));
3402   return (nummerges > 0);
3403 } /* test_vneighbors */
3404 
3405 /*-<a                             href="qh-merge_r.htm#TOC"
3406   >-------------------------------</a><a name="tracemerge">-</a>
3407 
3408   qh_tracemerge(qh, facet1, facet2 )
3409     print trace message after merge
3410 */
qh_tracemerge(qhT * qh,facetT * facet1,facetT * facet2)3411 void qh_tracemerge(qhT *qh, facetT *facet1, facetT *facet2) {
3412   boolT waserror= False;
3413 
3414 #ifndef qh_NOtrace
3415   if (qh->IStracing >= 4)
3416     qh_errprint(qh, "MERGED", facet2, NULL, NULL, NULL);
3417   if (facet2 == qh->tracefacet || (qh->tracevertex && qh->tracevertex->newlist)) {
3418     qh_fprintf(qh, qh->ferr, 8085, "qh_tracemerge: trace facet and vertex after merge of f%d and f%d, furthest p%d\n", facet1->id, facet2->id, qh->furthest_id);
3419     if (facet2 != qh->tracefacet)
3420       qh_errprint(qh, "TRACE", qh->tracefacet,
3421         (qh->tracevertex && qh->tracevertex->neighbors) ?
3422            SETfirstt_(qh->tracevertex->neighbors, facetT) : NULL,
3423         NULL, qh->tracevertex);
3424   }
3425   if (qh->tracevertex) {
3426     if (qh->tracevertex->deleted)
3427       qh_fprintf(qh, qh->ferr, 8086, "qh_tracemerge: trace vertex deleted at furthest p%d\n",
3428             qh->furthest_id);
3429     else
3430       qh_checkvertex(qh, qh->tracevertex);
3431   }
3432   if (qh->tracefacet) {
3433     qh_checkfacet(qh, qh->tracefacet, True, &waserror);
3434     if (waserror)
3435       qh_errexit(qh, qh_ERRqhull, qh->tracefacet, NULL);
3436   }
3437 #endif /* !qh_NOtrace */
3438   if (qh->CHECKfrequently || qh->IStracing >= 4) { /* can't check polygon here */
3439     qh_checkfacet(qh, facet2, True, &waserror);
3440     if (waserror)
3441       qh_errexit(qh, qh_ERRqhull, NULL, NULL);
3442   }
3443 } /* tracemerge */
3444 
3445 /*-<a                             href="qh-merge_r.htm#TOC"
3446   >-------------------------------</a><a name="tracemerging">-</a>
3447 
3448   qh_tracemerging(qh)
3449     print trace message during POSTmerging
3450 
3451   returns:
3452     updates qh.mergereport
3453 
3454   notes:
3455     called from qh_mergecycle() and qh_mergefacet()
3456 
3457   see:
3458     qh_buildtracing()
3459 */
qh_tracemerging(qhT * qh)3460 void qh_tracemerging(qhT *qh) {
3461   realT cpu;
3462   int total;
3463   time_t timedata;
3464   struct tm *tp;
3465 
3466   qh->mergereport= zzval_(Ztotmerge);
3467   time(&timedata);
3468   tp= localtime(&timedata);
3469   cpu= qh_CPUclock;
3470   cpu /= qh_SECticks;
3471   total= zzval_(Ztotmerge) - zzval_(Zcyclehorizon) + zzval_(Zcyclefacettot);
3472   qh_fprintf(qh, qh->ferr, 8087, "\n\
3473 At %d:%d:%d & %2.5g CPU secs, qhull has merged %d facets.  The hull\n\
3474   contains %d facets and %d vertices.\n",
3475       tp->tm_hour, tp->tm_min, tp->tm_sec, cpu,
3476       total, qh->num_facets - qh->num_visible,
3477       qh->num_vertices-qh_setsize(qh, qh->del_vertices));
3478 } /* tracemerging */
3479 
3480 /*-<a                             href="qh-merge_r.htm#TOC"
3481   >-------------------------------</a><a name="updatetested">-</a>
3482 
3483   qh_updatetested(qh, facet1, facet2 )
3484     clear facet2->tested and facet1->ridge->tested for merge
3485 
3486   returns:
3487     deletes facet2->center unless it's already large
3488       if so, clears facet2->ridge->tested
3489 
3490   design:
3491     clear facet2->tested
3492     clear ridge->tested for facet1's ridges
3493     if facet2 has a centrum
3494       if facet2 is large
3495         set facet2->keepcentrum
3496       else if facet2 has 3 vertices due to many merges, or not large and post merging
3497         clear facet2->keepcentrum
3498       unless facet2->keepcentrum
3499         clear facet2->center to recompute centrum later
3500         clear ridge->tested for facet2's ridges
3501 */
qh_updatetested(qhT * qh,facetT * facet1,facetT * facet2)3502 void qh_updatetested(qhT *qh, facetT *facet1, facetT *facet2) {
3503   ridgeT *ridge, **ridgep;
3504   int size;
3505 
3506   facet2->tested= False;
3507   FOREACHridge_(facet1->ridges)
3508     ridge->tested= False;
3509   if (!facet2->center)
3510     return;
3511   size= qh_setsize(qh, facet2->vertices);
3512   if (!facet2->keepcentrum) {
3513     if (size > qh->hull_dim + qh_MAXnewcentrum) {
3514       facet2->keepcentrum= True;
3515       zinc_(Zwidevertices);
3516     }
3517   }else if (size <= qh->hull_dim + qh_MAXnewcentrum) {
3518     /* center and keepcentrum was set */
3519     if (size == qh->hull_dim || qh->POSTmerging)
3520       facet2->keepcentrum= False; /* if many merges need to recompute centrum */
3521   }
3522   if (!facet2->keepcentrum) {
3523     qh_memfree(qh, facet2->center, qh->normal_size);
3524     facet2->center= NULL;
3525     FOREACHridge_(facet2->ridges)
3526       ridge->tested= False;
3527   }
3528 } /* updatetested */
3529 
3530 /*-<a                             href="qh-merge_r.htm#TOC"
3531   >-------------------------------</a><a name="vertexridges">-</a>
3532 
3533   qh_vertexridges(qh, vertex )
3534     return temporary set of ridges adjacent to a vertex
3535     vertex->neighbors defined
3536 
3537   ntoes:
3538     uses qh.visit_id
3539     does not include implicit ridges for simplicial facets
3540 
3541   design:
3542     for each neighbor of vertex
3543       add ridges that include the vertex to ridges
3544 */
qh_vertexridges(qhT * qh,vertexT * vertex)3545 setT *qh_vertexridges(qhT *qh, vertexT *vertex) {
3546   facetT *neighbor, **neighborp;
3547   setT *ridges= qh_settemp(qh, qh->TEMPsize);
3548   int size;
3549 
3550   qh->visit_id++;
3551   FOREACHneighbor_(vertex)
3552     neighbor->visitid= qh->visit_id;
3553   FOREACHneighbor_(vertex) {
3554     if (*neighborp)   /* no new ridges in last neighbor */
3555       qh_vertexridges_facet(qh, vertex, neighbor, &ridges);
3556   }
3557   if (qh->PRINTstatistics || qh->IStracing) {
3558     size= qh_setsize(qh, ridges);
3559     zinc_(Zvertexridge);
3560     zadd_(Zvertexridgetot, size);
3561     zmax_(Zvertexridgemax, size);
3562     trace3((qh, qh->ferr, 3011, "qh_vertexridges: found %d ridges for v%d\n",
3563              size, vertex->id));
3564   }
3565   return ridges;
3566 } /* vertexridges */
3567 
3568 /*-<a                             href="qh-merge_r.htm#TOC"
3569   >-------------------------------</a><a name="vertexridges_facet">-</a>
3570 
3571   qh_vertexridges_facet(qh, vertex, facet, ridges )
3572     add adjacent ridges for vertex in facet
3573     neighbor->visitid==qh.visit_id if it hasn't been visited
3574 
3575   returns:
3576     ridges updated
3577     sets facet->visitid to qh.visit_id-1
3578 
3579   design:
3580     for each ridge of facet
3581       if ridge of visited neighbor (i.e., unprocessed)
3582         if vertex in ridge
3583           append ridge to vertex
3584     mark facet processed
3585 */
qh_vertexridges_facet(qhT * qh,vertexT * vertex,facetT * facet,setT ** ridges)3586 void qh_vertexridges_facet(qhT *qh, vertexT *vertex, facetT *facet, setT **ridges) {
3587   ridgeT *ridge, **ridgep;
3588   facetT *neighbor;
3589 
3590   FOREACHridge_(facet->ridges) {
3591     neighbor= otherfacet_(ridge, facet);
3592     if (neighbor->visitid == qh->visit_id
3593     && qh_setin(ridge->vertices, vertex))
3594       qh_setappend(qh, ridges, ridge);
3595   }
3596   facet->visitid= qh->visit_id-1;
3597 } /* vertexridges_facet */
3598 
3599 /*-<a                             href="qh-merge_r.htm#TOC"
3600   >-------------------------------</a><a name="willdelete">-</a>
3601 
3602   qh_willdelete(qh, facet, replace )
3603     moves facet to visible list
3604     sets facet->f.replace to replace (may be NULL)
3605 
3606   returns:
3607     bumps qh.num_visible
3608 */
qh_willdelete(qhT * qh,facetT * facet,facetT * replace)3609 void qh_willdelete(qhT *qh, facetT *facet, facetT *replace) {
3610 
3611   qh_removefacet(qh, facet);
3612   qh_prependfacet(qh, facet, &qh->visible_list);
3613   qh->num_visible++;
3614   facet->visible= True;
3615   facet->f.replace= replace;
3616 } /* willdelete */
3617 
3618 #else /* qh_NOmerge */
qh_premerge(qhT * qh,vertexT * apex,realT maxcentrum,realT maxangle)3619 void qh_premerge(qhT *qh, vertexT *apex, realT maxcentrum, realT maxangle) {
3620 }
qh_postmerge(qhT * qh,const char * reason,realT maxcentrum,realT maxangle,boolT vneighbors)3621 void qh_postmerge(qhT *qh, const char *reason, realT maxcentrum, realT maxangle,
3622                       boolT vneighbors) {
3623 }
qh_checkzero(qhT * qh,boolT testall)3624 boolT qh_checkzero(qhT *qh, boolT testall) {
3625    }
3626 #endif /* qh_NOmerge */
3627 
3628