Scippy

SCIP

Solving Constraint Integer Programs

heur_local.c
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /* */
3 /* This file is part of the program and library */
4 /* SCIP --- Solving Constraint Integer Programs */
5 /* */
6 /* Copyright (C) 2002-2021 Konrad-Zuse-Zentrum */
7 /* fuer Informationstechnik Berlin */
8 /* */
9 /* SCIP is distributed under the terms of the ZIB Academic License. */
10 /* */
11 /* You should have received a copy of the ZIB Academic License */
12 /* along with SCIP; see the file COPYING. If not visit scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file heur_local.c
17  * @brief Improvement heuristic for Steiner problems
18  * @author Daniel Rehfeldt
19  *
20  * This file implements several local heuristics, including vertex insertion, key-path exchange and key-vertex elimination,
21  * ("Fast Local Search for Steiner Trees in Graphs" by Uchoa and Werneck). Other heuristics are for PCSTP and MWCSP.
22  *
23  * A list of all interface methods can be found in heur_local.h.
24  *
25  */
26 
27 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
28 //#define SCIP_DEBUG
29 #include <assert.h>
30 #include <string.h>
31 #include "heur_local.h"
32 #include "heur_tm.h"
33 #include "probdata_stp.h"
34 #include "cons_stp.h"
35 #include "solstp.h"
36 #include "stpvector.h"
37 
38 
39 /* @note if heuristic is running in root node timing is changed there to (SCIP_HEURTIMING_DURINGLPLOOP |
40  * SCIP_HEURTIMING_BEFORENODE), see SCIP_DECL_HEURINITSOL callback
41  */
42 
43 #define HEUR_NAME "local"
44 #define HEUR_DESC "improvement heuristic for STP"
45 #define HEUR_DISPCHAR '-'
46 #define HEUR_PRIORITY 100
47 #define HEUR_FREQ 1
48 #define HEUR_FREQOFS 0
49 #define HEUR_MAXDEPTH -1
50 #define HEUR_TIMING (SCIP_HEURTIMING_BEFORENODE | SCIP_HEURTIMING_DURINGLPLOOP | SCIP_HEURTIMING_AFTERLPLOOP | SCIP_HEURTIMING_AFTERNODE)
51 
52 #define HEUR_USESSUBSCIP FALSE /**< does the heuristic use a secondary SCIP instance? */
53 
54 #define DEFAULT_DURINGROOT TRUE
55 #define DEFAULT_MAXFREQLOC FALSE
56 #define DEFAULT_MAXNBESTSOLS 50
57 #define DEFAULT_NBESTSOLS 15
58 #define DEFAULT_NBESTSOLS_HARD 50
59 #define DEFAULT_NELITESOLS 3
60 #define DEFAULT_MINNBESTSOLS 10
61 #define DEFAULT_MINNBESTSOLS_HARD 30
62 #define DEFAULT_RANDSEED 1492 /**< random seed */
63 #define LOCAL_MAXRESTARTS 10
64 #define LOCAL_MAXRESTARTS_FAST 3
65 
66 
67 #define GREEDY_MAXRESTARTS 3 /**< Max number of restarts for greedy PC/MW heuristic if improving solution has been found. */
68 #define GREEDY_EXTENSIONS_MW 6 /**< Number of extensions for greedy MW heuristic. MUST BE HIGHER THAN GREEDY_EXTENSIONS */
69 #define GREEDY_EXTENSIONS 5 /**< Number of extensions for greedy PC heuristic. */
70 
71 
72 /*
73  * Data structures
74  */
75 
76 /** primal heuristic data */
77 struct SCIP_HeurData
78 {
79  int nfails; /**< number of fails */
80  int maxnsols; /**< maximal number of best solutions to improve */
81  int nbestsols; /**< number of best solutions to improve */
82  int nbestsols_min; /**< minimum number of best solutions to improve */
83  int* lastsolindices; /**< indices of a number of best solutions already tried */
84  SCIP_Bool maxfreq; /**< should the heuristic be called with maximum frequency? */
85  SCIP_Bool duringroot; /**< should the heuristic be called during the root node? */
86  SCIP_RANDNUMGEN* randnumgen; /**< random number generator */
87 };
88 
89 /** Voronoi data for local heuristic */
91 {
92  PATH* vnoi_path; /**< path */
93  int* vnoi_base; /**< base*/
94  SCIP_Real* memvdist; /**< distance */
95  int* memvbase; /**< base*/
96  int* meminedges; /**< in-edge */
97  int* vnoi_nodestate; /**< node state */
98  int nmems; /**< number of memorized elements */
99  int nkpnodes; /**< number of key path nodes */
100 } VNOILOC;
101 
102 
103 /** Connectivity data */
104 typedef struct connectivity_data
105 {
106  STP_Vectype(int)* blists_start; /**< boundary lists starts (on nodes) */
107  STP_Vectype(int)* lvledges_start; /**< horizontal edges starts (on nodes) */
108  PHNODE** pheap_boundpaths; /**< boundary paths (on nodes) */
109  int* pheap_sizes; /**< size (on nodes) */
110  PHNODE* pheap_elements; /**< elements, one per edge! */
111  int* boundaryedges; /**< current boundary edge */
112  UF* uf; /**< union find */
113  int nboundaryedges; /**< number of boundary edges */
114 } CONN;
115 
116 
117 /** Key-paths data */
119 {
120  int* const kpnodes; /**< key path nodes */
121  int* const kpedges; /**< key path edges */
122  SCIP_Real kpcost; /**< cost of key paths */
123  int nkpnodes; /**< number of key path nodes */
124  int nkpedges; /**< number of key path edges */
125  int rootpathstart; /**< start of key path towards root component */
126  int kptailnode; /**< needed for single path */
127 } KPATHS;
128 
129 
130 /** Solution tree data */
131 typedef struct solution_tree_data
132 {
133  STP_Bool* const solNodes; /**< Steiner tree nodes */
134  LCNODE* const linkcutNodes; /**< Steiner tree nodes */
135  int* const solEdges; /**< array indicating whether an arc is part of the solution (CONNECTED/UNKNOWN) */
136  STP_Bool* const nodeIsPinned; /**< of size nodes */
137  STP_Bool* const nodeIsScanned; /**< of size nodes */
138  int* const newedges; /**< marks new edges of the tree */
139 } SOLTREE;
140 
141 
142 /** Super graph data */
143 typedef struct supergraph_data
144 {
145  int* const supernodes; /**< super nodes */
146  STP_Bool* const nodeIsLowSupernode; /**< marks the current super-vertices
147  * (except for the one representing the root-component) */
148  PATH* mst; /**< MST */
149  SCIP_Real mstcost; /**< cost of MST */
150  int nsupernodes; /**< number of super nodes */
151 } SGRAPH;
152 
153 
154 /** Prize-collecting/maximum-weight connected subgraph data */
155 typedef struct pcmw_data
156 {
157  SCIP_Real* const prize_biased; /**< prize */
158  SCIP_Real* const edgecost_biased; /**< fully biased edge costs (PC only) */
159  SCIP_Real* const edgecost_org; /**< original, fully unbiased, edge costs (PC only) */
160  STP_Bool* const prizemark; /**< marked? */
161  int* const prizemarklist; /**< list of all marked */
162  int nprizemarks; /**< size of list */
163  SCIP_Bool isActive; /**< actively used (i.e. Pc/Mw graph)? */
164  int solroot; /**< root of Steiner tree solution */
165 } PCMW;
166 
167 
168 /** local insertion heuristic data */
169 typedef struct insertion_data
170 {
171  int* chainStarts; /**< pointers to starts of current chains (nInsertions many) */
172  int* chainEnds; /**< pointers to ends of current chains (nInsertions many) */
173  const SCIP_Real* edgecosts; /**< the edge costs (original for PC) */
174  int* const solDegreeNonTerm; /**< degree of node [v] in current solution; (pseudo) terminals are marked as UNKNOWN */
175  int* const addedEdges; /**< added edges */
176  int* const cutedgesStart; /**< cut edges for the chains */
177  int* const cutedgesEnd; /**< cut edges for the chains */
178  STP_Bool* const solNodes; /**< solution nodes */
179  SCIP_Bool* const nodeIsBlocked; /**< is node [v] blocked? */
180  int* const blockedList; /**< list of currently blocked nodes */
181  int blockedListSize; /**< size of list */
182  int nInsertions; /**< number of insertions */
183  int insertionVertex; /**< vertex to be inserted */
184 } INSERT;
185 
186 
187 /*
188  * Local methods
189  */
190 
191 
192 /** prune the solution? */
193 static inline
195  SCIP_SOL* initialsol /**< SCIP data structure */
196 )
197 {
198  assert(initialsol);
199 
200  if( SCIPsolGetHeur(initialsol) == NULL )
201  {
202  return TRUE;
203  }
204 
205  if( strcmp(SCIPheurGetName(SCIPsolGetHeur(initialsol)), "rec") == 0 )
206  {
207  return FALSE;
208  }
209 
210  if( strcmp(SCIPheurGetName(SCIPsolGetHeur(initialsol)), "TM") == 0 )
211  {
212  return FALSE;
213  }
214 
215  return TRUE;
216 }
217 
218 /** prune the solution? */
219 static inline
221  SCIP* scip, /**< SCIP data structure */
222  const GRAPH* graph, /**< graph data structure */
223  int* results /**< Steiner tree edges (in/out) */
224 )
225 {
226  STP_Bool* steinertree;
227  const int nnodes = graph->knots;
228 
229 #ifndef NDEBUG
230  const int nedges = graph->edges;
231  const SCIP_Real initialobj = solstp_getObjBounded(graph, results, 0.0, nedges);
232 #endif
233 
234  assert(solstp_isValid(scip, graph, results));
235 
236  SCIP_CALL( SCIPallocBufferArray(scip, &steinertree, nnodes) );
237 
238  solstp_setVertexFromEdge(graph, results, steinertree);
239 
240  SCIP_CALL( solstp_prune(scip, graph, results, steinertree) );
241 
242  SCIPfreeBufferArray(scip, &steinertree);
243 
244 #ifndef NDEBUG
245  {
246  const SCIP_Real initialobj_pruned = solstp_getObjBounded(graph, results, 0.0, nedges);
247 
248  assert(SCIPisLE(scip, initialobj_pruned, initialobj));
249  }
250 #endif
251 
252  return SCIP_OKAY;
253 }
254 
255 
256 /** fills solution array 'results' */
257 static inline
259  SCIP* scip, /**< SCIP data structure */
260  const GRAPH* graph, /**< graph data structure */
261  SCIP_SOL* initialsol, /**< SCIP data structure */
262  int* results /**< Steiner tree edges (out) */
263 )
264 {
265  SCIP_Real* xvals = SCIPprobdataGetXval(scip, initialsol);
266  const int nedges = graph_get_nEdges(graph);
267 
268  assert(results);
269  assert(xvals);
270 
271  /* set solution array */
272  for( int e = 0; e < nedges; e++ )
273  {
274  if( SCIPisEQ(scip, xvals[e], 1.0) )
275  {
276  results[e] = CONNECT;
277  }
278  else
279  {
280  assert(SCIPisEQ(scip, xvals[e], 0.0));
281 
282  results[e] = UNKNOWN;
283  }
284  }
285 }
286 
287 
288 /** initializes */
289 static
291  SCIP* scip, /**< SCIP data structure */
292  SCIP_HEURDATA* heurdata /**< heuristic's data */
293 )
294 {
295  assert(heurdata->nbestsols == -1 && heurdata->nbestsols_min == -1);
296 
298  {
299  heurdata->nbestsols = DEFAULT_NBESTSOLS_HARD;
300  heurdata->nbestsols_min = DEFAULT_MINNBESTSOLS_HARD;
301  }
302  else
303  {
304  heurdata->nbestsols = DEFAULT_NBESTSOLS;
305  heurdata->nbestsols_min = DEFAULT_MINNBESTSOLS;
306  }
307 }
308 
309 
310 /** tries to add solution stored in 'results' to SCIP solution pool */
311 static inline
313  SCIP* scip, /**< SCIP data structure */
314  SCIP_SOL** sols, /**< SCIP solutions */
315  SCIP_HEUR* heur, /**< heuristic */
316  const GRAPH* graph, /**< graph data structure */
317  SCIP_Real initialsol_obj, /**< objective */
318  SCIP_SOL* initialsol, /**< SCIP data structure */
319  int* results, /**< Steiner tree edges (out) */
320  SCIP_RESULT* result /**< pointer */
321 )
322 {
323  SCIP_Real* newsol_vals;
324  const int nedges = graph_get_nEdges(graph);
325  SCIP_Bool feasible = FALSE;
326  SCIP_HEURDATA* heurdata = SCIPheurGetData(heur);
327 
328  assert(heurdata);
329  assert(*result == SCIP_DIDNOTFIND);
330 
331  SCIP_CALL( SCIPallocBufferArray(scip, &newsol_vals, nedges) );
332 
333  for( int v = 0; v < nedges; v++ )
334  newsol_vals[v] = (results[v] == CONNECT) ? 1.0 : 0.0;
335 
336  SCIP_CALL( SCIPStpValidateSol(scip, graph, newsol_vals, FALSE, &feasible) );
337 
338  assert(feasible);
339 
340  /* is new solution feasible? */
341  if( feasible )
342  {
343  const SCIP_Real newsol_obj = solstp_getObjBounded(graph, results, 0.0, nedges);
344  const SCIP_Bool solIsBetter = SCIPisLT(scip, newsol_obj, initialsol_obj);
345  const SCIP_Bool solIsBetter_org
346  = SCIPisLT(scip, newsol_obj, SCIPgetSolOrigObj(scip, initialsol) - SCIPprobdataGetOffset(scip));
347 
348  assert(SCIPisLE(scip, newsol_obj, initialsol_obj));
349 
350 #ifndef NDEBUG
351  {
352  SCIP_Real pobj = 0.0;
353 
354  for( int v = 0; v < nedges; v++ )
355  pobj += graph->cost[v] * newsol_vals[v];
356 
357  assert(SCIPisEQ(scip, pobj, newsol_obj));
358  }
359 #endif
360 
361  // todo if solIsBetter || solIsBetter_org
362 
363  /* has solution been improved? */
364  if( solIsBetter || solIsBetter_org )
365  {
366  SCIP_SOL* const bestsol = sols[0];
367  SCIP_Bool success = FALSE;
368 
369  SCIP_CALL( SCIPprobdataAddNewSol(scip, newsol_vals, heur, &success) );
370 
371  if( success )
372  {
373  *result = SCIP_FOUNDSOL;
374 
375  if( heurdata->nbestsols < heurdata->maxnsols && SCIPisGT(scip, SCIPgetSolOrigObj(scip, bestsol) - SCIPprobdataGetOffset(scip), newsol_obj) )
376  {
377  heurdata->nfails = 0;
378  heurdata->nbestsols++;
379  }
380  SCIPdebugMessage("success in local: old: %f new: %f \n", (SCIPgetSolOrigObj(scip, bestsol) - SCIPprobdataGetOffset(scip)), newsol_obj);
381  }
382  }
383  }
384 
385  SCIPfreeBufferArray(scip, &newsol_vals);
386 
387  if( *result != SCIP_FOUNDSOL )
388  {
389  heurdata->nfails++;
390  if( heurdata->nbestsols > heurdata->nbestsols_min && heurdata->nfails > 1 )
391  heurdata->nbestsols--;
392 
393  SCIPdebugMessage("fail! %d \n", heurdata->nbestsols);
394  }
395 
396  return SCIP_OKAY;
397 }
398 
399 /** can the problem class be used? */
400 static inline
402  const GRAPH* graph /**< graph data structure */
403 )
404 {
405  const int stp_type = graph->stp_type;
406  if( stp_type != STP_SPG && stp_type != STP_RSMT && stp_type != STP_OARSMT && stp_type != STP_PCSPG
407  && stp_type != STP_RPCSPG && stp_type != STP_GSTP && stp_type != STP_MWCSP && stp_type != STP_RMWCSP )
408  return FALSE;
409 
410  return TRUE;
411 }
412 
413 /** submethod for local extend */
414 static
416  SCIP* scip, /**< SCIP data structure */
417  const GRAPH* graph, /**< graph data structure */
418  const PATH* path, /**< shortest data structure array */
419  int i, /**< node */
420  int greedyextensions, /**< greedyextensions */
421  int* nextensions, /**< nextensions */
422  GNODE* candidates, /**< candidates */
423  SCIP_PQUEUE* pqueue /**< pqueue */
424  )
425 {
426  assert(Is_pseudoTerm(graph->term[i]));
427  assert(!graph_pc_knotIsFixedTerm(graph, i));
428 
429  if( *nextensions < greedyextensions )
430  {
431  candidates[*nextensions].dist = graph->prize[i] - path[i].dist;
432  candidates[*nextensions].number = i;
433 
434  SCIP_CALL( SCIPpqueueInsert(pqueue, &(candidates[(*nextensions)++])) );
435  }
436  else
437  {
438  /* get candidate vertex of minimum value */
439  GNODE* min = (GNODE*) SCIPpqueueFirst(pqueue);
440  if( SCIPisLT(scip, min->dist, graph->prize[i] - path[i].dist) )
441  {
442  min = (GNODE*) SCIPpqueueRemove(pqueue);
443  min->dist = graph->prize[i] - path[i].dist;
444  min->number = i;
445  SCIP_CALL( SCIPpqueueInsert(pqueue, min) );
446  }
447  }
448 
449  return SCIP_OKAY;
450 }
451 
452 
453 /** perform local vertex insertion heuristic on given Steiner tree */
454 static
456  SCIP* scip, /**< SCIP data structure */
457  const GRAPH* graph, /**< graph data structure */
458  const int* solEdges, /**< Steiner tree edges */
459  LCNODE* linkcutNodes, /**< Steiner tree nodes */
460  STP_Bool* solNodes /**< Steiner tree nodes */
461  )
462 {
463  const int nnodes = graph->knots;
464  const int nedges = graph->edges;
465 
466  assert(solstp_isValid(scip, graph, solEdges));
467 
468  for( int i = 0; i < nnodes; i++ )
469  {
470  solNodes[i] = FALSE;
471  SCIPlinkcuttreeInitNode(&linkcutNodes[i]);
472  }
473 
474  /* create a link-cut tree representing the current Steiner tree */
475  for( int e = 0; e < nedges; e++ )
476  {
477  if( solEdges[e] == CONNECT )
478  SCIPlinkcuttreeLink(linkcutNodes, graph->head[e], graph->tail[e], flipedge(e));
479  }
480 
481  /* mark current Steiner tree nodes */
482  for( int e = 0; e < nedges; e++ )
483  {
484  if( solEdges[e] == CONNECT )
485  {
486  solNodes[graph->tail[e]] = TRUE;
487  solNodes[graph->head[e]] = TRUE;
488  }
489  }
490 
491  solNodes[graph->source] = TRUE;
492 
493  assert(linkcutNodes[graph->source].edge == -1);
494 }
495 
496 
497 /** is given Steiner tree a trivial solution (i.e. contains only one vertex?) */
498 static
500  const GRAPH* graph, /**< graph data structure */
501  const int* solEdges /**< Steiner tree edges */
502 )
503 {
504  const int root = graph->source;
505  SCIP_Bool isTrivial = TRUE;
506 
507  assert(graph_pc_isPcMw(graph));
508  assert(graph->extended);
509 
510  if( graph_pc_isRootedPcMw(graph) )
511  {
512  for( int e = graph->outbeg[root]; e != EAT_LAST; e = graph->oeat[e] )
513  {
514  if( solEdges[e] )
515  {
516  const int head = graph->head[e];
517  if( graph_pc_knotIsFixedTerm(graph, head) || !Is_term(graph->term[head]) )
518  {
519  isTrivial = FALSE;
520  break;
521  }
522  }
523  }
524  }
525  else
526  {
527  isTrivial = FALSE;
528  }
529 
530 
531  if( isTrivial )
532  {
533  SCIPdebugMessage("trivial solution given \n");
534  }
535 
536  return isTrivial;
537 }
538 
539 
540 /** computes lowest common ancestors for all pairs {vbase(v), vbase(u)} such that {u,w} is a boundary edge,
541  * first call should be with u := root */
542 static
544  SCIP* scip,
545  const GRAPH* graph,
546  int u,
547  UF* uf,
548  STP_Bool* nodesmark,
549  const int* steineredges,
550  STP_Vectype(int)* lcalists,
551  PHNODE** boundpaths,
552  const int* heapsize,
553  const int* vbase
554  )
555 {
556  int* uboundpaths; /* boundary-paths (each one represented by its boundary edge) having node 'u' as an endpoint */
557  uf->parent[u] = u;
558 
559  for( int oedge = graph->outbeg[u]; oedge != EAT_LAST; oedge = graph->oeat[oedge] )
560  {
561  const int v = graph->head[oedge];
562  if( steineredges[oedge] == CONNECT )
563  {
564  SCIP_CALL( lca(scip, graph, v, uf, nodesmark, steineredges, lcalists, boundpaths, heapsize, vbase) );
565  SCIPStpunionfindUnion(uf, u, v, FALSE);
566  uf->parent[SCIPStpunionfindFind(uf, u)] = u;
567  }
568  }
569 
570  nodesmark[u] = TRUE;
571 
572  /* iterate through all boundary-paths having one endpoint in the voronoi region of node 'u' */
573  SCIP_CALL( SCIPpairheapBuffarr(scip, boundpaths[u], heapsize[u], &uboundpaths) );
574 
575  for( int i = 0; i < heapsize[u]; i++ )
576  {
577  const int oedge = uboundpaths[i];
578  const int v = vbase[graph->head[oedge]];
579  if( nodesmark[v] )
580  {
581  const int ancestor = uf->parent[SCIPStpunionfindFind(uf, v)];
582 
583  /* if the ancestor of 'u' and 'v' is one of the two, the boundary-edge is already in boundpaths[u] */
584  if( ancestor != u && ancestor != v)
585  {
586  StpVecPushBack(scip, lcalists[ancestor], oedge);
587  }
588  }
589  }
590 
591  /* free the boundary-paths array */
592  SCIPfreeBufferArrayNull(scip, &uboundpaths);
593 
594  return SCIP_OKAY;
595 }
596 
597 
598 /** computes lowest common ancestors for all pairs {vbase(v), vbase(u)} such that {u,w} is a boundary edge */
599 static
601  SCIP* scip, /**< SCIP data structure */
602  const GRAPH* graph, /**< graph data structure */
603  const VNOILOC* vnoiData, /**< Voronoi data */
604  const SOLTREE* soltreeData, /**< solution tree data */
605  CONN* connectData /**< data */
606 )
607 {
608  PHNODE** const boundpaths = connectData->pheap_boundpaths;
609  STP_Vectype(int)* lvledges_start = connectData->lvledges_start;
610  const int* const solEdges = soltreeData->solEdges;
611  int* const pheapsize = connectData->pheap_sizes;
612  const int* const vnoibase = vnoiData->vnoi_base;
613  STP_Bool* nodesmark;
614  const int nnodes = graph->knots;
615 
616  assert(SCIPStpunionfindIsClear(scip, connectData->uf));
617 
618  SCIP_CALL( SCIPallocBufferArray(scip, &nodesmark, nnodes) );
619 
620  for( int i = 0; i < nnodes; ++i )
621  nodesmark[i] = FALSE;
622 
623  SCIP_CALL( lca(scip, graph, graph->source, connectData->uf, nodesmark, solEdges, lvledges_start, boundpaths, pheapsize, vnoibase) );
624 
625  SCIPfreeBufferArray(scip, &nodesmark);
626 
627  return SCIP_OKAY;
628 }
629 
630 
631 /** recursive method for a DFS ordering of graph 'g' */
632 static
633 void dfsorder(
634  const GRAPH* graph,
635  const int* edges,
636  int node,
637  int* counter,
638  int* dfst
639  )
640 {
641  int oedge = graph->outbeg[node];
642 
643  while( oedge >= 0 )
644  {
645  if( edges[oedge] >= 0 )
646  {
647  dfsorder(graph, edges, graph->head[oedge], counter, dfst);
648  }
649  oedge = graph->oeat[oedge];
650  }
651 
652  dfst[*counter] = node;
653  (*counter)++;
654 }
655 
656 
657 /** helper method for pcmwGetNewEdgePrize */
658 static inline
660  const GRAPH* graph,
661  const STP_Bool* steinertree,
662  const SCIP_Real* prize_biased,
663  const int* graphmark,
664  int node,
665  STP_Bool* prizemark,
666  int* prizemarklist,
667  int* prizemarkcount
668  )
669 {
670  SCIP_Real prizesum = 0.0;
671  assert(graph_pc_isPcMw(graph));
672 
673  if( graphmark[node] && !steinertree[node] && Is_pseudoTerm(graph->term[node]) && !prizemark[node] )
674  {
675  prizesum += prize_biased[node];
676  prizemark[node] = TRUE;
677  prizemarklist[(*prizemarkcount)++] = node;
678  }
679 
680  return prizesum;
681 }
682 
683 
684 /** gets root of solution for Pc/Mw or -1 otherwise */
685 static
687  const GRAPH* graph,
688  const int* solEdges
689 )
690 {
691  int solroot = -1;
692 
693  if( graph_pc_isPcMw(graph) && !graph_pc_isRootedPcMw(graph) )
694  {
695  const int root = graph->source;
696 
697  for( int e = graph->outbeg[root]; e != EAT_LAST; e = graph->oeat[e] )
698  {
699  const int head = graph->head[e];
700  if( solEdges[e] == CONNECT && !Is_term(graph->term[head]) )
701  {
702  assert(graph->cost[e] == 0.0);
703  solroot = head;
704  break;
705  }
706  }
707  assert(solroot >= 0);
708  }
709 
710  return solroot;
711 }
712 
713 
714 /** get prize not already counted that is associated to edge, not including solution nodes or forbidden ones */
715 static
717  const GRAPH* graph, /**< graph */
718  const STP_Bool* steinertree, /**< Steiner tree */
719  int edge, /**< the edge */
720  PCMW* pcmwData /**< data */
721  )
722 {
723  SCIP_Real prizesum = 0.0;
724 
725  assert(graph && steinertree && pcmwData);
726 
727  if( pcmwData->isActive )
728  {
729  const int mhead = graph->head[edge];
730  const int mtail = graph->tail[edge];
731  STP_Bool* const prizemark = pcmwData->prizemark;
732  int* const prizemarklist = pcmwData->prizemarklist;
733 
734  assert(pcmwData->prize_biased);
735  assert(graph_pc_isPcMw(graph));
736 
737  prizesum += getNewPrizeNode(graph, steinertree, pcmwData->prize_biased, graph->mark, mhead,
738  prizemark, prizemarklist, &(pcmwData->nprizemarks));
739  prizesum += getNewPrizeNode(graph, steinertree, pcmwData->prize_biased, graph->mark, mtail,
740  prizemark, prizemarklist, &(pcmwData->nprizemarks));
741  }
742 
743  return prizesum;
744 }
745 
746 
747 #ifndef NDEBUG
748 /** is the data clean? */
749 static
751  const GRAPH* graph, /**< graph data structure */
752  const PCMW* pcmwData /**< data */
753 )
754 {
755  if( !pcmwData->isActive )
756  return TRUE;
757 
758  if( pcmwData->nprizemarks != 0 )
759  {
760  SCIPdebugMessage("nprizemarks != 0 \n");
761  return FALSE;
762  }
763 
764  for( int k = 0; k < graph->knots; k++ )
765  {
766  if( pcmwData->prizemark[k] )
767  {
768  SCIPdebugMessage("entry %d is dirty \n", k);
769  return FALSE;
770  }
771  }
772 
773  return TRUE;
774 }
775 #endif
776 
777 
778 /** clean the data */
779 static
781  const GRAPH* graph, /**< graph data structure */
782  PCMW* pcmwData /**< data */
783 )
784 {
785  if( pcmwData->isActive )
786  {
787  const int* const prizemarklist = pcmwData->prizemarklist;
788  STP_Bool* const prizemark = pcmwData->prizemark;
789  const int prizemarkcount = pcmwData->nprizemarks;
790 
791 
792  for( int pi = 0; pi < prizemarkcount; pi++ )
793  prizemark[prizemarklist[pi]] = FALSE;
794 
795  pcmwData->nprizemarks = 0;
796 
797  assert(pcmwDataIsClean(graph, pcmwData));
798  }
799 }
800 
801 
802 /** updates for PC/MW: graph->mark and pinned */
803 static
805  SCIP* scip, /**< SCIP data structure */
806  GRAPH* graph, /**< graph data structure */
807  SOLTREE* soltreeData, /**< solution tree data */
808  PCMW* pcmwData /**< data */
809  )
810 {
811  SCIP_Bool* reachable_nodes;
812  STP_Bool* const pinned = soltreeData->nodeIsPinned;
813  const int root = graph->source;
814  const int nnodes = graph->knots;
815  int* const graphmark = graph->mark;
816  const int solstartnode = graph_pc_isRootedPcMw(graph)? graph->source : pcmwData->solroot;
817 
818  assert(graph->extended);
819  assert(graph_pc_isPcMw(graph));
820  assert(solstartnode >= 0 && solstartnode < graph->knots);
821  assert(soltreeData->solNodes[solstartnode]);
822 
823  /* remove unconnected vertices */
824  SCIP_CALL( SCIPallocBufferArray(scip, &reachable_nodes, nnodes) );
825 
826  SCIP_CALL( graph_trail_costAware(scip, graph, solstartnode, reachable_nodes) );
827 
828  for( int k = 0; k < nnodes; k++ )
829  graphmark[k] = reachable_nodes[k];
830 
831  SCIPfreeBufferArray(scip, &reachable_nodes);
832 
833  /* unmark and pin pseudo terminals */
834  for( int e = graph->outbeg[root]; e != EAT_LAST; e = graph->oeat[e] )
835  {
836  const int k = graph->head[e];
837  if( Is_term(graph->term[k]) )
838  {
839  if( !graph_pc_knotIsFixedTerm(graph, k) )
840  {
841  const int pterm = graph->head[graph->term2edge[k]];
842 
843  assert(Is_pseudoTerm(graph->term[pterm]));
844 
845  graphmark[k] = FALSE;
846  pinned[pterm] = TRUE;
847  }
848  }
849  }
850 
851  if( !graph_pc_isRootedPcMw(graph) )
852  {
853  graphmark[root] = FALSE;
854  soltreeData->solNodes[root] = FALSE;
855  }
856 
857  return SCIP_OKAY;
858 }
859 
860 
861 /** checks whether node is crucial, i.e. a terminal or a vertex with degree at least 3 (w.r.t. the Steiner tree) */
862 static
864  const GRAPH* graph, /**< graph data structure */
865  const int* steineredges, /**< Steiner edges */
866  int node /**< node to check */
867 )
868 {
869  assert(graph != NULL);
870  assert(steineredges != NULL);
871  assert(node >= 0 && node < graph->knots);
872 
873  if( !Is_term(graph->term[node]) )
874  {
875  int counter = 0;
876 
877  for( int e = graph->outbeg[node]; e != EAT_LAST; e = graph->oeat[e] )
878  {
879  /* check if the adjacent node is in the Steiner tree */
880  if( steineredges[e] == CONNECT || steineredges[flipedge(e)] == CONNECT )
881  counter++;
882  }
883 
884  /* is node of degree smaller 3 in the Steiner tree? */
885  if( counter < 3 )
886  {
887  return FALSE;
888  }
889  }
890 
891  return TRUE;
892 }
893 
894 
895 /** gets unbiased edge cost */
896 static inline
898  const GRAPH* graph, /**< graph data structure */
899  const PCMW* pcmwData, /**< data */
900  int edge /**< edge */
901  )
902 {
903  SCIP_Real cost;
904 
905  assert(graph && pcmwData);
906  assert(edge >= 0 && edge < graph->edges);
907 
908  if( pcmwData->isActive )
909  {
910  assert(graph->extended);
911 
912  if( !graph_pc_isPc(graph) )
913  {
914  assert(graph->stp_type == STP_MWCSP || graph->stp_type == STP_RMWCSP);
915 
916  cost = 0.0;
917  }
918  else
919  {
920  assert(pcmwData->edgecost_org);
921 
922  cost = pcmwData->edgecost_org[edge];
923  }
924  }
925  else
926  {
927  cost = graph->cost[edge];
928  }
929 
930  return cost;
931 }
932 
933 
934 /** Get cost of shortest path along boundary edge.
935  * For Pc/Mw: Will consider biased edge costs, but not the actual prizes of proper terminals! */
936 static
938  const GRAPH* graph, /**< graph data structure */
939  const VNOILOC* vnoiData, /**< data */
940  const PCMW* pcmwData, /**< data */
941  int boundaryedge /**< boundary edge*/
942  )
943 {
944  const PATH* const vnoipath = vnoiData->vnoi_path;
945  SCIP_Real pathcost;
946  const int tail = graph->tail[boundaryedge];
947  const int head = graph->head[boundaryedge];
948 
949  assert(boundaryedge >= 0);
950  assert(vnoiData->vnoi_base[tail] != vnoiData->vnoi_base[head]);
951  assert(graph_pc_isPcMw(graph) == pcmwData->isActive);
952 
953  pathcost = vnoipath[tail].dist + vnoipath[head].dist;
954 
955  if( pcmwData->isActive )
956  {
957  /* either vnoipath[head] goes into head, or head is the end of the entire boundary path;
958  * either way we don't want to subtract the prize! */
959  pathcost += getEdgeCostUnbiased(graph, pcmwData, boundaryedge);
960  }
961  else
962  {
963  pathcost += graph->cost[boundaryedge];
964  }
965 
966  assert(pathcost >= 0.0);
967 
968  return pathcost;
969 }
970 
971 
972 /** For Pc/Mw: Consider biased edge costs, and the actual prizes of proper terminals! */
973 static
975  const GRAPH* graph, /**< graph data structure */
976  const VNOILOC* vnoiData, /**< data */
977  const SOLTREE* soltreeData, /**< solution tree data */
978  int boundaryedge, /**< boundary edge */
979  PCMW* pcmwData /**< data */
980  )
981 {
982  SCIP_Real pathcost;
983 
984  assert(boundaryedge >= 0);
985  assert(vnoiData->vnoi_base[graph->tail[boundaryedge]] != vnoiData->vnoi_base[graph->head[boundaryedge]]);
986  assert(graph_pc_isPcMw(graph) == pcmwData->isActive);
987 
988  pathcost = getBoundaryPathCost(graph, vnoiData, pcmwData, boundaryedge);
989 
990  if( pcmwData->isActive )
991  {
992  const STP_Bool* const solNodes = soltreeData->solNodes;
993  const PATH* const vnoipath = vnoiData->vnoi_path;
994  const int* const vnoibase = vnoiData->vnoi_base;
995 
996 #ifdef SCIP_DEBUG
997  printf("boundary path edge ");
998  graph_edge_printInfo(graph, boundaryedge);
999 #endif
1000  pathcost -= pcmwGetNewEdgePrize(graph, solNodes, boundaryedge, pcmwData);
1001 
1002  for( int node = graph->tail[boundaryedge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1003  {
1004 #ifdef SCIP_DEBUG
1005  printf("boundary path edge ");
1006  graph_edge_printInfo(graph, vnoipath[node].edge);
1007 #endif
1008  pathcost -= pcmwGetNewEdgePrize(graph, solNodes, vnoipath[node].edge, pcmwData);
1009  }
1010 
1011  for( int node = graph->head[boundaryedge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1012  {
1013 #ifdef SCIP_DEBUG
1014  printf("boundary path edge ");
1015  graph_edge_printInfo(graph, vnoipath[node].edge);
1016 #endif
1017 
1018  pathcost -= pcmwGetNewEdgePrize(graph, solNodes, vnoipath[node].edge, pcmwData);
1019  }
1020 
1021  if( pathcost < 0.0 )
1022  pathcost = 0.0;
1023  }
1024 
1025  assert(pathcost >= 0.0);
1026 
1027  return pathcost;
1028 }
1029 
1030 
1031 /** returns edge costs (biased for Pc/Mw) */
1032 static inline
1034  const GRAPH* graph, /**< graph data structure */
1035  const PCMW* pcmwData /**< data */
1036 )
1037 {
1038  assert( graph_pc_isPcMw(graph) == pcmwData->isActive);
1039 
1040  return (pcmwData->isActive? pcmwData->edgecost_biased : graph->cost);
1041 }
1042 
1043 
1044 /** update for key-vertex elimination: Save current boundary edges */
1045 static
1047  SCIP* scip, /**< SCIP data structure */
1048  const GRAPH* graph, /**< graph data structure */
1049  const VNOILOC* vnoiData, /**< Voronoi data */
1050  const SGRAPH* supergraphData, /**< super-graph*/
1051  int crucnode, /**< node to eliminate */
1052  CONN* connectData /**< data */
1053 )
1054 {
1055  PHNODE** const boundpaths = connectData->pheap_boundpaths;
1056  STP_Vectype(int)* lvledges_start = connectData->lvledges_start;
1057  STP_Vectype(int) lvledges_curr = NULL;
1058  UF* const uf = connectData->uf;
1059  int* const pheapsize = connectData->pheap_sizes;
1060  const int* const supernodes = supergraphData->supernodes;
1061  const int* const vnoibase = vnoiData->vnoi_base;
1062  const STP_Bool* const isLowSupernode = supergraphData->nodeIsLowSupernode;
1063  int* const boundaryedges = connectData->boundaryedges;
1064  const int* const graphmark = graph->mark;
1065  int nboundaryedges = 0;
1066 
1067  connectData->nboundaryedges = -1;
1068 
1069  /* add vertical boundary-paths between the child components and the root-component (w.r.t. node 'crucnode') */
1070  for( int k = 0; k < supergraphData->nsupernodes - 1; k++ )
1071  {
1072  const int supernode = supernodes[k];
1073  int edge = UNKNOWN;
1074 
1075  while( boundpaths[supernode] != NULL )
1076  {
1077  int node;
1078  SCIP_Real edgecost;
1079 
1080  SCIP_CALL( SCIPpairheapDeletemin(scip, &edge, &edgecost, &boundpaths[supernode], &pheapsize[supernode]) );
1081 
1082  node = (vnoibase[graph->head[edge]] == UNKNOWN)? UNKNOWN : SCIPStpunionfindFind(uf, vnoibase[graph->head[edge]]);
1083 
1084  /* check whether edge 'edge' represents a boundary-path having an endpoint in the kth-component and in the root-component respectively */
1085  if( node != UNKNOWN && !isLowSupernode[node] && graphmark[node] )
1086  {
1087  boundaryedges[nboundaryedges++] = edge;
1088  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[supernode],
1089  connectData->pheap_elements, edge, edgecost, &pheapsize[supernode]) );
1090  break;
1091  }
1092  }
1093  }
1094 
1095  lvledges_curr = lvledges_start[crucnode];
1096 
1097  /* add horizontal boundary-paths (between the child-components) */
1098  for( int i = 0; i < StpVecGetSize(lvledges_curr); i++ )
1099  {
1100  const int edge = lvledges_curr[i];
1101  const int basetail = vnoibase[graph->tail[edge]];
1102  const int basehead = vnoibase[graph->head[edge]];
1103  const int node = (basehead == UNKNOWN)? UNKNOWN : SCIPStpunionfindFind(uf, basehead);
1104  const int adjnode = (basetail == UNKNOWN)? UNKNOWN : SCIPStpunionfindFind(uf, basetail);
1105 
1106  /* check whether the current boundary-path connects two child components */
1107  if( node != UNKNOWN && isLowSupernode[node] && adjnode != UNKNOWN && isLowSupernode[adjnode] )
1108  {
1109  assert(graphmark[node] && graphmark[adjnode]);
1110  boundaryedges[nboundaryedges++] = edge;
1111  }
1112  }
1113 
1114  connectData->nboundaryedges = nboundaryedges;
1115 
1116  return SCIP_OKAY;
1117 }
1118 
1119 
1120 /** initialize */
1121 static
1123  SCIP* scip, /**< SCIP data structure */
1124  const GRAPH* graph, /**< graph data structure */
1125  const VNOILOC* vnoiData, /**< Voronoi data */
1126  const SOLTREE* soltreeData, /**< solution tree data */
1127  const PCMW* pcmwData, /**< data */
1128  CONN* connectData /**< data */
1129 )
1130 {
1131  PHNODE** const boundpaths = connectData->pheap_boundpaths;
1132  STP_Vectype(int)* blists_start = connectData->blists_start;
1133  STP_Vectype(int)* lvledges_start = connectData->lvledges_start;
1134  int* const pheapsize = connectData->pheap_sizes;
1135  const int* const vnoibase = vnoiData->vnoi_base;
1136  const int* const graphmark = graph->mark;
1137  const int nnodes = graph->knots;
1138  const int nedges = graph->edges;
1139 
1140  assert(connectData->nboundaryedges == 0);
1141  assert(connectData->boundaryedges);
1142 
1143  for( int k = 0; k < nnodes; k++ )
1144  assert(graph->path_state[k] == CONNECT || !graph->mark[k]);
1145 
1146  for( int k = 0; k < nnodes; ++k )
1147  {
1148  if( blists_start[k] )
1149  StpVecClear(blists_start[k]);
1150  if( lvledges_start[k] )
1151  StpVecClear(lvledges_start[k]);
1152  }
1153 
1154  BMSclearMemoryArray(pheapsize, nnodes);
1155  BMSclearMemoryArray(boundpaths, nnodes);
1156 
1157  /* link all nodes to their (respective) Voronoi base */
1158  for( int k = 0; k < nnodes; ++k )
1159  {
1160  assert(pheapsize[k] == 0 && boundpaths[k] == NULL);
1161 
1162  if( !graphmark[k] )
1163  continue;
1164 
1165  StpVecPushBack(scip, blists_start[vnoibase[k]], k);
1166  }
1167 
1168  /* for each node, store all of its outgoing boundary-edges in a (respective) heap*/
1169  for( int e = 0; e < nedges; e += 2 )
1170  {
1171  if( graph->oeat[e] != EAT_FREE )
1172  {
1173  const int tail = graph->tail[e];
1174  const int head = graph->head[e];
1175 
1176  /* is edge 'e' a boundary-edge? */
1177  if( vnoibase[tail] != vnoibase[head] && graphmark[tail] && graphmark[head] )
1178  {
1179  const SCIP_Real edgecost = getBoundaryPathCost(graph, vnoiData, pcmwData, e);
1180 
1181  assert(vnoibase[tail] != UNKNOWN && vnoibase[head] != UNKNOWN);
1182  assert(SCIPisGE(scip, edgecost, 0.0));
1183 
1184  /* add the boundary-edge 'e' and its reversed to the corresponding heaps */
1185  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[vnoibase[tail]],
1186  connectData->pheap_elements, e, edgecost, &(pheapsize[vnoibase[tail]])) );
1187  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[vnoibase[head]],
1188  connectData->pheap_elements, flipedge(e), edgecost, &(pheapsize[vnoibase[head]])) );
1189  }
1190  }
1191  }
1192 
1193  SCIP_CALL( getLowestCommonAncestors(scip, graph, vnoiData, soltreeData, connectData) );
1194 
1195  return SCIP_OKAY;
1196 }
1197 
1198 
1199 /** get key path above given crucial node */
1200 static
1202  SCIP* scip, /**< SCIP data structure */
1203  int crucnode, /**< crucial node to start from */
1204  const GRAPH* graph, /**< graph data structure */
1205  const SOLTREE* soltreeData, /**< solution tree data */
1206  const PCMW* pcmwData, /**< data */
1207  CONN* connectData, /**< data */
1208  KPATHS* keypathsData /**< key paths */
1209 )
1210 {
1211  int* const kpnodes = keypathsData->kpnodes;
1212  const int* const solEdges = soltreeData->solEdges;
1213  const LCNODE* const linkcutNodes = soltreeData->linkcutNodes;
1214  const STP_Bool* const solNodes = soltreeData->solNodes;
1215  const STP_Bool* const pinned = soltreeData->nodeIsPinned;
1216  PHNODE** const boundpaths = connectData->pheap_boundpaths;
1217  const SCIP_Real* const edgecosts = getEdgeCosts(graph, pcmwData);
1218  int* const pheapsize = connectData->pheap_sizes;
1219  const int* const graphmark = graph->mark;
1220  int nkpnodes = 0;
1221  int firstedge2root = -1;
1222  int kptailnode = -1;
1223  SCIP_Real kpcost = -FARAWAY;
1224 
1225  if( Is_term(graph->term[crucnode]) || pinned[crucnode] )
1226  {
1227  /* update the data structures; basically merge the whole key path */
1228  for( int edge = graph->outbeg[crucnode]; edge != EAT_LAST; edge = graph->oeat[edge] )
1229  {
1230  int adjnode = graph->head[edge];
1231 
1232  /* check whether edge 'edge' leads to an ancestor of terminal 'crucnode' */
1233  if( solEdges[edge] == CONNECT && solNodes[adjnode] && graphmark[adjnode] )
1234  {
1235  assert( SCIPStpunionfindFind(connectData->uf, adjnode) != crucnode);
1236  assert(soltreeData->nodeIsScanned[adjnode]);
1237 
1238  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[adjnode], &pheapsize[crucnode], &pheapsize[adjnode]);
1239 
1240  /* update the union-find data structure */
1241  SCIPStpunionfindUnion(connectData->uf, crucnode, adjnode, FALSE);
1242 
1243  /* move along the key-path until its end (i.e. until a crucial node is reached) */
1244  while( !nodeIsCrucial(graph, solEdges, adjnode) && !pinned[adjnode] )
1245  {
1246  int e;
1247  for( e = graph->outbeg[adjnode]; e != EAT_LAST; e = graph->oeat[e] )
1248  {
1249  if( solEdges[e] == CONNECT )
1250  break;
1251 
1252  assert(UNKNOWN == solEdges[e]);
1253  }
1254 
1255  /* assert that each leaf of the Steiner tree is a terminal */
1256  assert(e != EAT_LAST);
1257  adjnode = graph->head[e];
1258 
1259  if( !solNodes[adjnode] || !graphmark[adjnode] )
1260  break;
1261 
1262  assert(soltreeData->nodeIsScanned[adjnode]);
1263  assert(SCIPStpunionfindFind(connectData->uf, adjnode) != crucnode);
1264 
1265  /* update the union-find data structure */
1266  SCIPStpunionfindUnion(connectData->uf, crucnode, adjnode, FALSE);
1267 
1268  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[adjnode], &pheapsize[crucnode], &pheapsize[adjnode]);
1269  }
1270  }
1271  }
1272  }
1273 
1274 #ifndef NDEBUG
1275  if( SCIPisGE(scip, graph->cost[linkcutNodes[crucnode].edge], FARAWAY)
1276  || SCIPisGE(scip, graph->cost[flipedge(linkcutNodes[crucnode].edge)], FARAWAY) )
1277  {
1278  assert(graph_pc_isPcMw(graph));
1279  assert(graph->head[linkcutNodes[crucnode].edge] == graph->source);
1280  }
1281 #endif
1282 
1283  /* traverse the (unique) key-path containing the parent of the current crucial node 'crucnode' */
1284 
1285  firstedge2root = linkcutNodes[crucnode].edge;
1286  kptailnode = graph->head[firstedge2root];
1287 
1288  if( pcmwData->isActive )
1289  {
1290  /* we don't want to count the prize of 'crucnode' */
1291  kpcost = getEdgeCostUnbiased(graph, pcmwData, flipedge(firstedge2root));
1292  }
1293  else
1294  {
1295  kpcost = edgecosts[firstedge2root];
1296  }
1297 
1298  assert(graph->tail[firstedge2root] == crucnode);
1299  assert(soltreeData->solEdges[flipedge(firstedge2root)] == CONNECT);
1300 
1301 #ifdef SCIP_DEBUG
1302  printf("key path edge ");
1303  graph_edge_printInfo(graph, firstedge2root);
1304 #endif
1305 
1306  while( !nodeIsCrucial(graph, solEdges, kptailnode) && !pinned[kptailnode] )
1307  {
1308  const int kpedge2root = linkcutNodes[kptailnode].edge;
1309  kpcost += edgecosts[flipedge(kpedge2root)];
1310 
1311 #ifdef SCIP_DEBUG
1312  printf("key path edge ");
1313  graph_edge_printInfo(graph, kpedge2root);
1314 #endif
1315 
1316  kpnodes[nkpnodes++] = kptailnode;
1317  kptailnode = graph->head[kpedge2root];
1318  }
1319 
1320  keypathsData->kpcost = kpcost;
1321  keypathsData->kptailnode = kptailnode;
1322  keypathsData->nkpnodes = nkpnodes;
1323 }
1324 
1325 
1326 /** unmarks key path nodes */
1327 static
1329  const GRAPH* graph, /**< graph data structure */
1330  const KPATHS* keypathsData, /**< key paths */
1331  SOLTREE* soltreeData /**< solution tree data */
1332 )
1333 {
1334  const int* const kpnodes = keypathsData->kpnodes;
1335  STP_Bool* const solNodes = soltreeData->solNodes;
1336  const int nkpnodes = keypathsData->nkpnodes;
1337 
1338  for( int k = 0; k < nkpnodes; k++ )
1339  {
1340  const int node = kpnodes[k];
1341 
1342  assert(node >= 0 && node < graph->knots);
1343  assert(solNodes[node]);
1344 
1345  solNodes[node] = FALSE;
1346  }
1347 }
1348 
1349 
1350 /** (re-)marks key path nodes */
1351 static
1353  const GRAPH* graph, /**< graph data structure */
1354  const KPATHS* keypathsData, /**< key paths */
1355  SOLTREE* soltreeData /**< solution tree data */
1356 )
1357 {
1358  const int* const kpnodes = keypathsData->kpnodes;
1359  STP_Bool* const solNodes = soltreeData->solNodes;
1360  const int nkpnodes = keypathsData->nkpnodes;
1361 
1362  for( int k = 0; k < nkpnodes; k++ )
1363  {
1364  const int node = kpnodes[k];
1365 
1366  assert(node >= 0 && node < graph->knots);
1367  assert(!solNodes[node]);
1368 
1369  solNodes[node] = TRUE;
1370  }
1371 }
1372 
1373 
1374 /** exchanges key path */
1375 static
1377  SCIP* scip, /**< SCIP data structure */
1378  GRAPH* graph, /**< graph data structure */
1379  const CONN* connectData, /**< data */
1380  const VNOILOC* vnoiData, /**< data */
1381  const KPATHS* keypathsData, /**< key paths */
1382  const int* dfstree, /**< DFS tree */
1383  const STP_Bool* scanned, /**< array to mark which nodes have been scanned */
1384  int dfstree_pos, /**< current position in dfs tree */
1385  int boundedge_new, /**< new Voronoi boundary edge */
1386  SOLTREE* soltreeData /**< solution tree data */
1387 )
1388 {
1389  UF* const uf = connectData->uf;
1390  PHNODE** const boundpaths = connectData->pheap_boundpaths;
1391  int* const pheapsize = connectData->pheap_sizes;
1392  const PATH* const vnoipath = vnoiData->vnoi_path;
1393  const int* const vnoibase = vnoiData->vnoi_base;
1394  const int* const kpnodes = keypathsData->kpnodes;
1395  STP_Bool* const pinned = soltreeData->nodeIsPinned;
1396  const LCNODE* const linkcutNodes = soltreeData->linkcutNodes;
1397  int* const solEdges = soltreeData->solEdges;
1398  STP_Bool* const solNodes = soltreeData->solNodes;
1399  const int nkpnodes = keypathsData->nkpnodes;
1400  const int crucnode = dfstree[dfstree_pos];
1401  int* const graphmark = graph->mark;
1402  int newpathend = -1;
1403  int newedge = boundedge_new;
1404  int node = SCIPStpunionfindFind(uf, vnoibase[graph->head[newedge]]);
1405 
1406  /* remove old keypath */
1407  assert(solEdges[flipedge(linkcutNodes[crucnode].edge)] != UNKNOWN);
1408 
1409  solEdges[flipedge(linkcutNodes[crucnode].edge)] = UNKNOWN;
1410  solNodes[crucnode] = FALSE;
1411  graphmark[crucnode] = FALSE;
1412 
1413  for( int k = 0; k < nkpnodes; k++ )
1414  {
1415  const int keypathnode = kpnodes[k];
1416  assert(solEdges[flipedge(linkcutNodes[keypathnode].edge)] != UNKNOWN);
1417 
1418  solEdges[flipedge(linkcutNodes[keypathnode].edge)] = UNKNOWN;
1419  solNodes[keypathnode] = FALSE;
1420  graphmark[keypathnode] = FALSE;
1421  }
1422 
1423  assert(graphmark[keypathsData->kptailnode]);
1424 
1425  if( node == crucnode )
1426  newedge = flipedge(newedge);
1427 
1428  for( node = graph->tail[newedge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1429  {
1430  graphmark[node] = FALSE;
1431 
1432  solEdges[flipedge(vnoipath[node].edge)] = CONNECT;
1433  solEdges[vnoipath[node].edge] = UNKNOWN;
1434  }
1435 
1436  for( node = graph->head[newedge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1437  {
1438  graphmark[node] = FALSE;
1439 
1440  solEdges[vnoipath[node].edge] = CONNECT;
1441  }
1442 
1443  solEdges[flipedge(newedge)] = CONNECT;
1444 
1445  newpathend = vnoibase[graph->tail[newedge]];
1446  assert(node == vnoibase[graph->head[newedge]] );
1447 
1448  /* flip all edges on the ST path between the endnode of the new key-path and the current crucial node */
1449  assert(SCIPStpunionfindFind(uf, newpathend) == crucnode);
1450 
1451  for( int k = newpathend; k != crucnode; k = graph->head[linkcutNodes[k].edge] )
1452  {
1453  assert(graphmark[k]);
1454  assert( solEdges[flipedge(linkcutNodes[k].edge)] != -1);
1455 
1456  solEdges[flipedge(linkcutNodes[k].edge)] = UNKNOWN;
1457  solEdges[linkcutNodes[k].edge] = CONNECT;
1458  }
1459 
1460  for( int k = 0; k < dfstree_pos; k++ )
1461  {
1462  if( crucnode == SCIPStpunionfindFind(uf, dfstree[k]) )
1463  {
1464  graphmark[dfstree[k]] = FALSE;
1465  solNodes[dfstree[k]] = FALSE;
1466  }
1467  }
1468 
1469  /* update union find */
1470  if( !Is_term(graph->term[node]) && scanned[node] && !pinned[node] && SCIPStpunionfindFind(uf, node) == node )
1471  {
1472  for( int edge = graph->outbeg[node]; edge != EAT_LAST; edge = graph->oeat[edge] )
1473  {
1474  int adjnode = graph->head[edge];
1475 
1476  /* check whether edge 'edge' leads to an ancestor of terminal 'node' */
1477  if( solEdges[edge] == CONNECT && solNodes[adjnode] && graphmark[adjnode] && SCIPStpunionfindFind(uf, adjnode) != node )
1478  {
1479  assert(scanned[adjnode]);
1480  /* meld the heaps */
1481  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[adjnode], &pheapsize[node], &pheapsize[adjnode]);
1482 
1483  /* update the union-find data structure */
1484  SCIPStpunionfindUnion(uf, node, adjnode, FALSE);
1485 
1486  /* move along the key-path until its end (i.e. until a crucial node is reached) */
1487  while( !nodeIsCrucial(graph, solEdges, adjnode) && !pinned[adjnode] )
1488  {
1489  int e;
1490  for( e = graph->outbeg[adjnode]; e != EAT_LAST; e = graph->oeat[e] )
1491  {
1492  if( solEdges[e] != -1 )
1493  break;
1494  }
1495 
1496  /* assert that each leaf of the ST is a terminal */
1497  assert( e != EAT_LAST );
1498  adjnode = graph->head[e];
1499 
1500  if( !solNodes[adjnode] )
1501  break;
1502 
1503  assert(scanned[adjnode]);
1504  assert(SCIPStpunionfindFind(uf, adjnode) != node);
1505 
1506  /* update the union-find data structure */
1507  SCIPStpunionfindUnion(uf, node, adjnode, FALSE);
1508 
1509  /* meld the heaps */
1510  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[adjnode], &pheapsize[node], &pheapsize[adjnode]);
1511  }
1512  }
1513  }
1514 
1515  }
1516  pinned[node] = TRUE;
1517 
1518  return SCIP_OKAY;
1519 }
1520 
1521 
1522 
1523 /** exchanges key-paths star */
1524 static
1526  SCIP* scip, /**< SCIP data structure */
1527  const GRAPH* graph, /**< graph data structure */
1528  const CONN* connectData, /**< data */
1529  const VNOILOC* vnoiData, /**< data */
1530  const KPATHS* keypathsData, /**< key paths */
1531  const SGRAPH* supergraphData, /**< super-graph */
1532  const int* dfstree, /**< DFS tree */
1533  const STP_Bool* scanned, /**< array to mark which nodes have been scanned */
1534  int dfstree_pos, /**< current position in dfs tree */
1535  SOLTREE* soltreeData /**< solution tree data */
1536 )
1537 {
1538  const PATH* mst = supergraphData->mst;
1539  UF* const uf = connectData->uf;
1540  const STP_Bool* const isSupernode = supergraphData->nodeIsLowSupernode;
1541  PHNODE** const boundpaths = connectData->pheap_boundpaths;
1542  int* const pheapsize = connectData->pheap_sizes;
1543  const PATH* const vnoipath = vnoiData->vnoi_path;
1544  const int* const vnoibase = vnoiData->vnoi_base;
1545  const int* const boundedges = connectData->boundaryedges;
1546  const int* const kpnodes = keypathsData->kpnodes;
1547  const int* const kpedges = keypathsData->kpedges;
1548  STP_Bool* const pinned = soltreeData->nodeIsPinned;
1549  const LCNODE* const linkcutNodes = soltreeData->linkcutNodes;
1550  int* const solEdges = soltreeData->solEdges;
1551  int* const graphmark = graph->mark;
1552  STP_Bool* const solNodes = soltreeData->solNodes;
1553  const int nkpnodes = keypathsData->nkpnodes;
1554  const int nkpedges = keypathsData->nkpedges;
1555  const int nsupernodes = supergraphData->nsupernodes;
1556 
1557  /* unmark the original edges spanning the supergraph */
1558  for( int e = 0; e < nkpedges; e++ )
1559  {
1560  assert(solEdges[kpedges[e]] == CONNECT);
1561  solEdges[kpedges[e]] = UNKNOWN;
1562  }
1563 
1564  /* mark all Steiner tree nodes except for those belonging to the root-component as forbidden */
1565  for( int k = keypathsData->rootpathstart; k < nkpnodes; k++ )
1566  {
1567  graphmark[kpnodes[k]] = FALSE;
1568  solNodes[kpnodes[k]] = FALSE;
1569  }
1570 
1571  for( int k = 0; k < dfstree_pos; k++ )
1572  {
1573  const int node = SCIPStpunionfindFind(uf, dfstree[k]);
1574  if( isSupernode[node] || node == dfstree[dfstree_pos] )
1575  {
1576  graphmark[dfstree[k]] = FALSE;
1577  solNodes[dfstree[k]] = FALSE;
1578  }
1579  }
1580 
1581  /* add the new edges reconnecting the (super-) components */
1582  for( int l = 0; l < nsupernodes - 1; l++ )
1583  {
1584  const int edge = (mst[l].edge % 2 == 0)? boundedges[mst[l].edge / 2] : flipedge(boundedges[mst[l].edge / 2]);
1585 
1586  /* change the orientation within the target-component if necessary */
1587  if( !isSupernode[vnoibase[graph->head[edge]]] )
1588  {
1589  int node = vnoibase[graph->head[edge]];
1590  const int nodebase = SCIPStpunionfindFind(uf, node);
1591  assert(isSupernode[nodebase]);
1592 
1593  while( node != nodebase )
1594  {
1595  /* the Steiner tree edge pointing towards the root */
1596  const int e = linkcutNodes[node].edge;
1597 
1598  assert(solEdges[e] == UNKNOWN && solEdges[flipedge(e)] == CONNECT );
1599 
1600  solEdges[e] = CONNECT;
1601  solEdges[flipedge(e)] = UNKNOWN;
1602  node = graph->head[e];
1603  }
1604  }
1605 
1606  /* is the vbase of the current boundary-edge tail in the root-component? */
1607  if( !isSupernode[SCIPStpunionfindFind(uf, vnoibase[graph->tail[edge]])] )
1608  {
1609  int node;
1610  solEdges[edge] = CONNECT;
1611 
1612  for( node = graph->tail[edge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1613  {
1614  graphmark[node] = FALSE;
1615 
1616  if( solEdges[flipedge(vnoipath[node].edge)] == CONNECT )
1617  solEdges[flipedge(vnoipath[node].edge)] = UNKNOWN;
1618 
1619  solEdges[vnoipath[node].edge] = CONNECT;
1620  }
1621 
1622  assert(!isSupernode[node] && vnoibase[node] == node);
1623  assert(graphmark[node]);
1624 
1625  /* is the pinned node its own component identifier? */
1626  if( !Is_term(graph->term[node]) && scanned[node] && !pinned[node] && SCIPStpunionfindFind(uf, node) == node )
1627  {
1628  graphmark[graph->head[edge]] = FALSE;
1629 
1630  for( int oedge = graph->outbeg[node]; oedge != EAT_LAST; oedge = graph->oeat[oedge] )
1631  {
1632  int head = graph->head[oedge];
1633 
1634  /* check whether edge 'edge' leads to an ancestor of 'node' */
1635  if( solEdges[oedge] == CONNECT && graphmark[head] && solNodes[head] && SCIPStpunionfindFind(uf, head) != node )
1636  {
1637  assert(scanned[head]);
1638 
1639  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[head], &pheapsize[node], &pheapsize[head]);
1640 
1641  /* update the union-find data structure */
1642  SCIPStpunionfindUnion(uf, node, head, FALSE);
1643 
1644  /* move along the key-path until its end (i.e. until a crucial node is reached) */
1645  while( !nodeIsCrucial(graph, solEdges, head) && !pinned[head] )
1646  {
1647  int e;
1648  for( e = graph->outbeg[head]; e != EAT_LAST; e = graph->oeat[e] )
1649  {
1650  if( solEdges[e] != -1 )
1651  break;
1652  }
1653 
1654  /* assert that each leaf of the ST is a terminal */
1655  assert( e != EAT_LAST );
1656  head = graph->head[e];
1657 
1658  if( !solNodes[head] )
1659  break;
1660 
1661  assert(scanned[head]);
1662  assert(SCIPStpunionfindFind(uf, head) != node);
1663 
1664  /* update the union-find data structure */
1665  SCIPStpunionfindUnion(uf, node, head, FALSE);
1666 
1667  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[head], &pheapsize[node], &pheapsize[head]);
1668  }
1669  }
1670  }
1671  }
1672 
1673  /* mark the start node (lying in the root-component of the Steiner tree) of the current boundary-path as pinned,
1674  * so that it may not be removed later on */
1675  pinned[node] = TRUE;
1676 
1677  for( node = graph->head[edge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1678  {
1679  graphmark[node] = FALSE;
1680  if( solEdges[vnoipath[node].edge] == CONNECT )
1681  solEdges[vnoipath[node].edge] = UNKNOWN;
1682 
1683  solEdges[flipedge(vnoipath[node].edge)] = CONNECT;
1684  }
1685  }
1686  else /* the vbase of the current boundary-edge tail is NOT in the root-component */
1687  {
1688  solEdges[edge] = CONNECT;
1689 
1690  for( int node = graph->tail[edge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1691  {
1692  graphmark[node] = FALSE;
1693  if( solEdges[vnoipath[node].edge] != CONNECT && solEdges[flipedge(vnoipath[node].edge)] != CONNECT )
1694  {
1695  solEdges[vnoipath[node].edge] = CONNECT;
1696  }
1697  }
1698 
1699  for( int node = graph->head[edge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1700  {
1701  graphmark[node] = FALSE;
1702 
1703  solEdges[flipedge(vnoipath[node].edge)] = CONNECT;
1704  solEdges[vnoipath[node].edge] = UNKNOWN;
1705  }
1706  }
1707  }
1708 
1709  for( int k = 0; k < nkpnodes; k++ )
1710  {
1711  assert(graphmark[kpnodes[k]] == FALSE);
1712  assert(solNodes[kpnodes[k]] == FALSE);
1713  }
1714 
1715  assert(!graphmark[dfstree[dfstree_pos]]);
1716 
1717  return SCIP_OKAY;
1718 }
1719 
1720 
1721 /** compute cost of alternative key path */
1722 static
1724  SCIP* scip, /**< SCIP data structure */
1725  const GRAPH* graph, /**< graph data structure */
1726  const VNOILOC* vnoiData, /**< data */
1727  const SOLTREE* soltreeData, /**< solution tree data */
1728  SCIP_Real edgecost_old_in, /**< edge cost of old edge */
1729  int boundedge_old, /**< Voronoi boundary edge */
1730  PCMW* pcmwData, /**< data */
1731  int* boundedge_new /**< new Voronoi boundary edge (in/out) */
1732 )
1733 {
1734  SCIP_Real edgecost = FARAWAY;
1735  SCIP_Real edgecost_old = FARAWAY;
1736  SCIP_Real edgecost_new = FARAWAY;
1737  int newedge = *boundedge_new;
1738 
1739  assert(pcmwDataIsClean(graph, pcmwData));
1740 
1741  if( boundedge_old != UNKNOWN )
1742  {
1743  SCIPdebugMessage("get replace path for old edge: \n");
1744 
1745  edgecost_old = getBoundaryPathCostPrized(graph, vnoiData, soltreeData, boundedge_old, pcmwData);
1746 
1747  assert(SCIPisLE(scip, edgecost_old, edgecost_old_in));
1748  assert(SCIPisLT(scip, edgecost_old_in, FARAWAY));
1749 
1750  pcmwDataClean(graph, pcmwData);
1751  }
1752 
1753  if( newedge != UNKNOWN )
1754  {
1755  SCIPdebugMessage("get replace path for new edge: \n");
1756 
1757  edgecost_new = getBoundaryPathCostPrized(graph, vnoiData, soltreeData, newedge, pcmwData);
1758 
1759  pcmwDataClean(graph, pcmwData);
1760  }
1761 
1762  if( SCIPisLT(scip, edgecost_old, edgecost_new) )
1763  {
1764  edgecost = edgecost_old;
1765  *boundedge_new = boundedge_old;
1766  }
1767  else
1768  {
1769  edgecost = edgecost_new;
1770  *boundedge_new = newedge;
1771  }
1772 
1773  assert(SCIPisLT(scip, edgecost, FARAWAY) || (graph_pc_isPcMw(graph) && graph->source == graph->head[*boundedge_new] ) );
1774  assert(*boundedge_new != UNKNOWN);
1775 
1776  return edgecost;
1777 }
1778 
1779 
1780 /** restore supergraph data */
1781 static
1783  const GRAPH* graph, /**< graph data structure */
1784  SGRAPH* supergraphData /**< super-graph*/
1785 )
1786 {
1787  const int* const supernodes = supergraphData->supernodes;
1788  STP_Bool* const nodeIsLowSupernode = supergraphData->nodeIsLowSupernode;
1789 
1790  SCIPfreeMemoryArray(scip, &(supergraphData->mst));
1791 
1792  /* unmark the descendant supervertices */
1793  for( int k = 0; k < supergraphData->nsupernodes - 1; k++ )
1794  nodeIsLowSupernode[supernodes[k]] = FALSE;
1795 
1796 #ifndef NDEBUG
1797  for( int k = 0; k < graph->knots; k++ )
1798  assert(!nodeIsLowSupernode[k]);
1799 #endif
1800 }
1801 
1802 /** compute minimum-spanning tree */
1803 static
1805  SCIP* scip, /**< SCIP data structure */
1806  const GRAPH* graph, /**< graph data structure */
1807  const CONN* connectData, /**< data */
1808  const VNOILOC* vnoiData, /**< data */
1809  const KPATHS* keypathsData, /**< key paths */
1810  int crucnode, /**< node to eliminate */
1811  SOLTREE* soltreeData, /**< solution tree data */
1812  PCMW* pcmwData, /**< data */
1813  SGRAPH* supergraphData /**< super-graph*/
1814 )
1815 {
1816  PATH* mst = NULL;
1817  GRAPH* supergraph = NULL;
1818  UF* const uf = connectData->uf;
1819  const int* const supernodes = supergraphData->supernodes;
1820  int* supernodesid = NULL;
1821  const SCIP_Real* const edgecosts = getEdgeCosts(graph, pcmwData);
1822  const STP_Bool* const isLowerSupernode = supergraphData->nodeIsLowSupernode;
1823  STP_Bool* const solNodes = soltreeData->solNodes;
1824  int* const newedges = soltreeData->newedges;
1825  const PATH* const vnoipath = vnoiData->vnoi_path;
1826  const int* const vnoibase = vnoiData->vnoi_base;
1827  const int* const boundaryedges = connectData->boundaryedges;
1828  SCIP_Real mstcost = 0.0;
1829  const int nboundaryedges = connectData->nboundaryedges;
1830  const int nnodes = graph->knots;
1831  const int nsupernodes = supergraphData->nsupernodes;
1832  /* the (super-) vertex representing the current root-component of the Steiner tree */
1833  const int superroot = supernodes[nsupernodes - 1];
1834 #ifdef SCIP_DEBUG
1835  SCIP_Real mstcost_supergraph = 0.0;
1836  const SCIP_Bool pcmw = graph_pc_isPcMw(graph);
1837 #endif
1838  assert(nboundaryedges > 0);
1839  assert(superroot >= 0);
1840  assert(!supergraphData->mst);
1841 
1842  SCIP_CALL( SCIPallocBufferArray(scip, &supernodesid, nnodes) );
1843 #ifndef NDEBUG
1844  for( int k = 0; k < nnodes; k++ )
1845  supernodesid[k] = -1;
1846 #endif
1847 
1848  /* create a supergraph, having the endpoints of the key-paths incident to the current crucial node as (super-) vertices */
1849  SCIP_CALL( graph_init(scip, &supergraph, nsupernodes, nboundaryedges * 2, 1) );
1850  supergraph->stp_type = STP_SPG;
1851 
1852  for( int k = 0; k < nsupernodes; k++ )
1853  {
1854  supernodesid[supernodes[k]] = k;
1855  graph_knot_add(supergraph, graph->term[supernodes[k]]);
1856  }
1857 
1858  if( pcmwData->isActive )
1859  {
1860  soltreeUnmarkKpNodes(graph, keypathsData, soltreeData);
1861  }
1862 
1863  assert(pcmwDataIsClean(graph, pcmwData));
1864 
1865  SCIPdebugMessage("build super-graph: \n");
1866 
1867  /* add edges to the supergraph */
1868  for( int k = 0; k < nboundaryedges; k++ )
1869  {
1870  SCIP_Real edgecost;
1871  const int edge = boundaryedges[k];
1872  int node = SCIPStpunionfindFind(uf, vnoibase[graph->tail[edge]]);
1873  int adjnode = SCIPStpunionfindFind(uf, vnoibase[graph->head[edge]]);
1874 
1875  /* if node 'node' or 'adjnode' belongs to the root-component, take the (temporary) root-component identifier instead */
1876  node = (isLowerSupernode[node]? node : superroot);
1877  adjnode = (isLowerSupernode[adjnode]? adjnode : superroot);
1878 
1879  /* compute the cost of the boundary-path pertaining to the boundary-edge 'edge' */
1880  edgecost = getBoundaryPathCostPrized(graph, vnoiData, soltreeData, edge, pcmwData);
1881  pcmwDataClean(graph, pcmwData);
1882 
1883  graph_edge_add(scip, supergraph, supernodesid[node], supernodesid[adjnode], edgecost, edgecost);
1884  }
1885 
1886  /* compute a MST on the supergraph */
1887  SCIP_CALL( SCIPallocMemoryArray(scip, &(supergraphData->mst), nsupernodes) );
1888  mst = supergraphData->mst;
1889  SCIP_CALL( graph_path_init(scip, supergraph) );
1890  graph_path_exec(scip, supergraph, MST_MODE, nsupernodes - 1, supergraph->cost, mst);
1891 
1892  assert(pcmwDataIsClean(graph, pcmwData));
1893 
1894  /* compute the cost of the MST */
1895 
1896  mstcost = 0.0;
1897 
1898  for( int l = 0; l < nsupernodes - 1; l++ )
1899  {
1900  /* compute the edge in the original graph corresponding to the current MST edge */
1901  const int edge = (mst[l].edge % 2 == 0)? boundaryedges[mst[l].edge / 2] : flipedge(boundaryedges[mst[l].edge / 2]);
1902 
1903  mstcost += getEdgeCostUnbiased(graph, pcmwData, edge);
1904  mstcost -= pcmwGetNewEdgePrize(graph, solNodes, edge, pcmwData);
1905 
1906 #ifdef SCIP_DEBUG
1907  printf("key vertex mst edge ");
1908  graph_edge_printInfo(graph, edge);
1909  if( pcmw ) printf("nodeweights: %f -> %f \n", graph->prize[graph->tail[edge]], graph->prize[graph->head[edge]]);
1910 
1911  mstcost_supergraph += supergraph->cost[mst[l].edge];
1912 #endif
1913 
1914  assert(newedges[edge] != crucnode && newedges[flipedge(edge)] != crucnode);
1915 
1916  /* mark the edge (in the original graph) as visited */
1917  newedges[edge] = crucnode;
1918 
1919  /* traverse along the boundary-path belonging to the boundary-edge 'edge' */
1920  for( int node = graph->tail[edge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1921  {
1922  const int e = vnoipath[node].edge;
1923 
1924  /* if edge 'e' and its reversed have not been visited yet */
1925  if( newedges[e] != crucnode && newedges[flipedge(e)] != crucnode )
1926  {
1927  newedges[e] = crucnode;
1928  mstcost += edgecosts[e];
1929  mstcost -= pcmwGetNewEdgePrize(graph, solNodes, e, pcmwData);
1930 
1931 #ifdef SCIP_DEBUG
1932  printf("key vertex mst edge ");
1933  graph_edge_printInfo(graph, e);
1934  if( pcmw ) printf("nodeweights: %f -> %f \n", graph->prize[graph->tail[edge]], graph->prize[graph->head[edge]]);
1935 #endif
1936  }
1937  }
1938 
1939  for( int node = graph->head[edge]; node != vnoibase[node]; node = graph->tail[vnoipath[node].edge] )
1940  {
1941  const int e = vnoipath[node].edge;
1942  const int erev = flipedge(e);
1943 
1944  /* if edge 'e' and its reversed have not been visited yet */
1945  if( newedges[e] != crucnode && newedges[erev] != crucnode )
1946  {
1947  newedges[erev] = crucnode;
1948 
1949  mstcost += edgecosts[e];
1950  mstcost -= pcmwGetNewEdgePrize(graph, solNodes, e, pcmwData);
1951 
1952 #ifdef SCIP_DEBUG
1953  printf("key vertex mst edge ");
1954  graph_edge_printInfo(graph, erev);
1955  if( pcmw ) printf("nodeweights: %f -> %f \n", graph->prize[graph->tail[edge]], graph->prize[graph->head[edge]]);
1956 #endif
1957  }
1958  }
1959  }
1960 #ifdef SCIP_DEBUG
1961  SCIPdebugMessage("mstcost=%f supergraph mstcost=%f \n", mstcost, mstcost_supergraph);
1962 #endif
1963 
1964  pcmwDataClean(graph, pcmwData);
1965 
1966  supergraphData->mstcost = mstcost;
1967 
1968  if( pcmwData->isActive )
1969  {
1970  soltreeMarkKpNodes(graph, keypathsData, soltreeData);
1971  }
1972 
1973  SCIPfreeBufferArray(scip, &supernodesid);
1974  graph_path_exit(scip, supergraph);
1975  graph_free(scip, &supergraph, TRUE);
1976 
1977  return SCIP_OKAY;
1978 }
1979 
1980 
1981 /** get the key-paths star starting from 'keyvertex' */
1982 static
1984  int keyvertex, /**< key vertex to start from */
1985  const GRAPH* graph, /**< graph data structure */
1986  const CONN* connectData, /**< data */
1987  const SOLTREE* soltreeData, /**< solution tree data */
1988  const PCMW* pcmwData, /**< data */
1989  KPATHS* keypathsData, /**< key paths */
1990  SGRAPH* supergraphData, /**< super-graph*/
1991  SCIP_Bool* success /**< success? */
1992 )
1993 {
1994  int* const kpnodes = keypathsData->kpnodes;
1995  int* const kpedges = keypathsData->kpedges;
1996  const int* const solEdges = soltreeData->solEdges;
1997  int* const supernodes = supergraphData->supernodes;
1998  STP_Bool* isLowSupernode = supergraphData->nodeIsLowSupernode;
1999  const STP_Bool* const solNodes = soltreeData->solNodes;
2000  const STP_Bool* const pinned = soltreeData->nodeIsPinned;
2001  int edgeFromRoot = UNKNOWN;
2002  int nkpnodes = 0;
2003  int nkpedges = 0;
2004  int inedgescount = 0;
2005  int nsupernodes = 0;
2006  const SCIP_Real* const edgecosts = getEdgeCosts(graph, pcmwData);
2007 
2008  assert(!pinned[keyvertex] && nodeIsCrucial(graph, solEdges, keyvertex) && !Is_term(graph->term[keyvertex]));
2009 
2010  keypathsData->kpcost = 0.0;
2011  keypathsData->rootpathstart = -1;
2012  keypathsData->nkpedges = -1;
2013  keypathsData->nkpnodes = -1;
2014  supergraphData->nsupernodes = -1;
2015  *success = TRUE;
2016 
2017  /* find all key-paths starting in node 'keyvertex' */
2018  for( int edge = graph->outbeg[keyvertex]; edge != EAT_LAST; edge = graph->oeat[edge] )
2019  {
2020  assert(solNodes[graph->tail[edge]]);
2021 
2022  /* check whether the outgoing edge is in the Steiner tree */
2023  if( (solEdges[edge] == CONNECT && solNodes[graph->head[edge]])
2024  || solEdges[flipedge(edge)] == CONNECT )
2025  {
2026  const int edge_rev = flipedge(edge);
2027 
2028  /* we want to have the edges going into 'keyvertex' to get the prizes right */
2029  keypathsData->kpcost += edgecosts[edge_rev];
2030  inedgescount++;
2031 
2032 #ifdef SCIP_DEBUG
2033  printf("key vertex start edge ");
2034  graph_edge_printInfo(graph, edge_rev);
2035 #endif
2036 
2037  /* does current edge lead to the Steiner tree root? */
2038  if( solEdges[edge_rev] == CONNECT )
2039  {
2040  assert(edgeFromRoot == UNKNOWN);
2041 
2042  edgeFromRoot = edge_rev;
2043  kpedges[nkpedges++] = edgeFromRoot;
2044 
2045  assert(edge == soltreeData->linkcutNodes[keyvertex].edge);
2046  }
2047  else
2048  {
2049  int adjnode = graph->head[edge];
2050  int e = edge;
2051 #ifndef NDEBUG
2052  const int nkpedges_org = nkpedges;
2053  assert(solEdges[flipedge(edge)] == UNKNOWN);
2054 #endif
2055  kpedges[nkpedges++] = e;
2056 
2057  /* move along the key-path until its end (i.e. a crucial or pinned node) is reached */
2058  while( !pinned[adjnode] && !nodeIsCrucial(graph, solEdges, adjnode) && solNodes[adjnode] )
2059  {
2060  /* update the union-find data structure */
2061  SCIPStpunionfindUnion(connectData->uf, keyvertex, adjnode, FALSE);
2062 
2063  kpnodes[nkpnodes++] = adjnode;
2064 
2065  for( e = graph->outbeg[adjnode]; e != EAT_LAST; e = graph->oeat[e] )
2066  {
2067  if( solEdges[e] == CONNECT )
2068  {
2069  /* we sum the cost of the edges pointing toward the key vertex */
2070  keypathsData->kpcost += edgecosts[flipedge(e)];
2071  kpedges[nkpedges++] = e;
2072 #ifdef SCIP_DEBUG
2073  printf("key vertex edge ");
2074  graph_edge_printInfo(graph, flipedge(e));
2075 #endif
2076  break;
2077  }
2078  }
2079 
2080  /* assert that each leaf of the ST is a terminal
2081  * todo check why this fails, the following should never actually happen */
2082  if( e == EAT_LAST )
2083  {
2084  *success = FALSE;
2085  goto TERMINATE;
2086  }
2087 
2088  assert(e != EAT_LAST);
2089  adjnode = graph->head[e];
2090  }
2091 
2092  /* does the last node on the path belong to a removed component? */
2093  if( !solNodes[adjnode] )
2094  {
2095  /* assert that e is not the edge incident to the key vertex */
2096  assert(e != edge && e != edge_rev);
2097  assert(nkpedges - nkpedges_org >= 2);
2098 
2099  keypathsData->kpcost -= edgecosts[flipedge(e)];
2100 #ifdef SCIP_DEBUG
2101  printf("key vertex remove edge ");
2102  graph_edge_printInfo(graph, flipedge(e));
2103 #endif
2104  nkpedges--;
2105  adjnode = graph->tail[e];
2106  if( adjnode != keyvertex )
2107  {
2108  supernodes[nsupernodes++] = adjnode;
2109  isLowSupernode[adjnode] = TRUE;
2110  }
2111  }
2112  else
2113  {
2114  supernodes[nsupernodes++] = adjnode;
2115  isLowSupernode[adjnode] = TRUE;
2116  }
2117 
2118  assert(nkpedges > nkpedges_org);
2119  } /* edge does not go to root */
2120  } /* edge is in solution */
2121  } /* find all (unique) key-paths starting in node 'crucnode' */
2122 
2123  assert(inedgescount >= 0);
2124 
2125  if( inedgescount > 1 && pcmwData->isActive )
2126  {
2127  assert(graph_pc_isPcMw(graph));
2128 
2129  if( !graph_pc_isPc(graph) || graph_pc_knotIsNonLeafTerm(graph, keyvertex) )
2130  {
2131  /* we have implicitly subtracted the prize of the key vertex several times, need to revert */
2132  keypathsData->kpcost += (SCIP_Real)(inedgescount - 1) * graph->prize[keyvertex];
2133  }
2134  else
2135  {
2136  assert(0.0 == graph->prize[keyvertex]);
2137  }
2138  }
2139 
2140  /* traverse the key-path leading to the root-component */
2141  keypathsData->rootpathstart = nkpnodes;
2142  if( edgeFromRoot != UNKNOWN )
2143  {
2144  /* begin with the edge starting in the root-component of key vertex */
2145  int tail = graph->tail[edgeFromRoot];
2146 
2147  while( !pinned[tail] && !nodeIsCrucial(graph, solEdges, tail) && solNodes[tail] )
2148  {
2149  int e;
2150 
2151  kpnodes[nkpnodes++] = tail;
2152 
2153  for( e = graph->inpbeg[tail]; e != EAT_LAST; e = graph->ieat[e] )
2154  {
2155  if( solEdges[e] == CONNECT )
2156  {
2157  assert(solNodes[graph->tail[e]]);
2158  keypathsData->kpcost += edgecosts[e];
2159 #ifdef SCIP_DEBUG
2160  printf("key vertex (root) edge ");
2161  graph_edge_printInfo(graph, e);
2162 #endif
2163  kpedges[nkpedges++] = e;
2164  break;
2165  }
2166  }
2167 
2168  assert(e != EAT_LAST);
2169  tail = graph->tail[e];
2170  }
2171 
2172  supernodes[nsupernodes++] = tail;
2173  }
2174 
2175  /* the last of the key-path nodes to be stored is the current key-node */
2176  kpnodes[nkpnodes++] = keyvertex;
2177 
2178  TERMINATE:
2179 
2180  keypathsData->nkpedges = nkpedges;
2181  keypathsData->nkpnodes = nkpnodes;
2182  supergraphData->nsupernodes = nsupernodes;
2183 }
2184 
2185 
2186 /** preprocessing for Voronoi repair method */
2187 static
2189  SCIP* scip, /**< SCIP data structure */
2190  const GRAPH* graph, /**< graph data structure */
2191  const KPATHS* keypathsData, /**< key paths */
2192  const CONN* connectData, /**< base lists */
2193  const PCMW* pcmwData, /**< data */
2194  VNOILOC* vnoiData, /**< data */
2195  int* nheapelems /**< to store */
2196 )
2197 {
2198  STP_Vectype(int)* blists_start = connectData->blists_start;
2199  PATH* const vnoipath = vnoiData->vnoi_path;
2200  const int* const kpnodes = keypathsData->kpnodes;
2201  int* const vnoibase = vnoiData->vnoi_base;
2202  int* const state = vnoiData->vnoi_nodestate;
2203  const int* const graphmark = graph->mark;
2204  const int nkpnodes = keypathsData->nkpnodes;
2205  const SCIP_Real* const edgecosts = getEdgeCosts(graph, pcmwData);
2206  int count = 0;
2207 
2208  assert(nheapelems);
2209 
2210  for( int k = 0; k < nkpnodes; k++ )
2211  {
2212  STP_Vectype(int) blists_curr = blists_start[kpnodes[k]];
2213  const int size = StpVecGetSize(blists_curr);
2214 
2215  assert(blists_curr != NULL);
2216 
2217  for( int i = 0; i < size; i++ )
2218  {
2219  const int node = blists_curr[i];
2220  assert(graph_knot_isInRange(graph, node));
2221 
2222  for( int edge_in = graph->inpbeg[node]; edge_in != EAT_LAST; edge_in = graph->ieat[edge_in] )
2223  {
2224  const int tail = graph->tail[edge_in];
2225 
2226  /* is tail node not in C and not forbidden? */
2227  if( state[tail] == CONNECT && graphmark[tail] && graphmark[vnoibase[tail]] )
2228  {
2229  const SCIP_Real newdist = vnoipath[tail].dist + edgecosts[edge_in];
2230 
2231  /* is the new path via tail shorter? */
2232  if( SCIPisGT(scip, vnoipath[node].dist, newdist) )
2233  {
2234  vnoipath[node].dist = newdist;
2235  vnoipath[node].edge = edge_in;
2236  vnoibase[node] = vnoibase[tail];
2237  }
2238  }
2239  }
2240 
2241  if( vnoibase[node] != UNKNOWN )
2242  {
2243  graph_pathHeapAdd(vnoipath, node, graph->path_heap, state, &count);
2244  }
2245  }
2246  }
2247 
2248  assert(nkpnodes == 0 || count > 0);
2249 
2250  *nheapelems = count;
2251 }
2252 
2253 /** restore data */
2254 static
2256  const CONN* connectData, /**< base lists */
2257  const KPATHS* keypathsData, /**< key paths */
2258  VNOILOC* vnoiData /**< data */
2259 )
2260 {
2261  STP_Vectype(int)* blists_start = connectData->blists_start;
2262  PATH* vnoipath = vnoiData->vnoi_path;
2263  int* memvbase = vnoiData->memvbase;
2264  int* meminedges = vnoiData->meminedges;
2265  int* vnoibase = vnoiData->vnoi_base;
2266  const int* kpnodes = keypathsData->kpnodes;
2267  SCIP_Real* memvdist = vnoiData->memvdist;
2268  const int nkpnodes = keypathsData->nkpnodes;
2269  int memcount = 0;
2270 
2271  for( int k = 0; k < nkpnodes; k++ )
2272  {
2273  /* restore data of all nodes having the current (internal) key-path node as their voronoi base */
2274  STP_Vectype(int) blists_curr = blists_start[kpnodes[k]];
2275  const int size = StpVecGetSize(blists_curr);
2276 
2277  for( int i = 0; i < size; i++ )
2278  {
2279  const int node = blists_curr[i];
2280  vnoibase[node] = memvbase[memcount];
2281  vnoipath[node].dist = memvdist[memcount];
2282  vnoipath[node].edge = meminedges[memcount];
2283  memcount++;
2284  }
2285  }
2286 
2287  assert(memcount == vnoiData->nmems);
2288  assert(vnoiData->nkpnodes == nkpnodes);
2289 }
2290 
2291 /** reset data */
2292 static
2294  const CONN* connectData, /**< base lists */
2295  const KPATHS* keypathsData, /**< key paths */
2296  const int* graphmark, /**< graph mark */
2297  VNOILOC* vnoiData /**< data */
2298 )
2299 {
2300  STP_Vectype(int)* blists_start = connectData->blists_start;
2301  PATH* vnoipath = vnoiData->vnoi_path;
2302  int* memvbase = vnoiData->memvbase;
2303  int* meminedges = vnoiData->meminedges;
2304  int* state = vnoiData->vnoi_nodestate;
2305  int* vnoibase = vnoiData->vnoi_base;
2306  const int* kpnodes = keypathsData->kpnodes;
2307  SCIP_Real* memvdist = vnoiData->memvdist;
2308  const int nkpnodes = keypathsData->nkpnodes;
2309  int nresnodes = 0;
2310 
2311  /* reset all nodes (referred to as 'C') whose bases are internal nodes of the current key-paths */
2312  for( int k = 0; k < nkpnodes; k++ )
2313  {
2314  /* reset all nodes having the current (internal) key-path node as their Voronoi base */
2315  STP_Vectype(int) blists_curr = blists_start[kpnodes[k]];
2316  const int size = StpVecGetSize(blists_curr);
2317 
2318  for( int i = 0; i < size; i++ )
2319  {
2320  const int node = blists_curr[i];
2321 
2322  assert(graphmark[node]);
2323 
2324  /* store data */
2325  memvbase[nresnodes] = vnoibase[node];
2326  memvdist[nresnodes] = vnoipath[node].dist;
2327  meminedges[nresnodes] = vnoipath[node].edge;
2328  nresnodes++;
2329 
2330  /* reset data */
2331  vnoibase[node] = UNKNOWN;
2332  vnoipath[node].dist = FARAWAY;
2333  vnoipath[node].edge = UNKNOWN;
2334  state[node] = UNKNOWN;
2335  }
2336  }
2337 
2338  vnoiData->nmems = nresnodes;
2339  vnoiData->nkpnodes = nkpnodes;
2340 }
2341 
2342 
2343 #ifndef NDEBUG
2344 static
2346  SCIP* scip, /**< SCIP data structure */
2347  const GRAPH* graph, /**< graph data structure */
2348  const int* solDegree, /**< solution degree */
2349  const LCNODE* linkcutNodes /**< Steiner tree nodes */
2350  )
2351 {
2352  int* solDegreeNew;
2353  SCIP_Bool isValid;
2354  const int nnodes = graph->knots;
2355 
2356  SCIP_CALL_ABORT( SCIPallocBufferArray(scip, &solDegreeNew, nnodes) );
2357 
2358  BMSclearMemoryArray(solDegreeNew, nnodes);
2359 
2360  for( int i = 0; i < nnodes; ++i )
2361  {
2362  const int edge = linkcutNodes[i].edge;
2363 
2364  if( edge == -1 )
2365  continue;
2366 
2367  assert(edge >= 0 && edge < graph->edges);
2368  assert(graph->tail[edge] == i);
2369 
2370  solDegreeNew[graph->tail[edge]]++;
2371  solDegreeNew[graph->head[edge]]++;
2372  }
2373 
2374  isValid = TRUE;
2375 
2376  for( int i = 0; i < nnodes; ++i )
2377  {
2378  if( Is_term(graph->term[i]) || Is_pseudoTerm(graph->term[i]) )
2379  {
2380  assert(UNKNOWN == solDegree[i]);
2381  continue;
2382  }
2383 
2384  if( solDegree[i] != solDegreeNew[i] )
2385  {
2386 #ifdef SCIP_DEBUG
2387  graph_knot_printInfo(graph, i);
2388  SCIPdebugMessage("%d != %d (old/new degree) \n", solDegree[i], solDegreeNew[i]);
2389 #endif
2390  isValid = FALSE;
2391  break;
2392  }
2393  }
2394 
2395  SCIPfreeBufferArray(scip, &solDegreeNew);
2396 
2397  return isValid;
2398 }
2399 
2400 
2401 static
2403  SCIP* scip, /**< SCIP data structure */
2404  const GRAPH* graph, /**< graph data structure */
2405  const STP_Bool* solNodes, /**< Steiner tree nodes */
2406  const LCNODE* linkcutNodes /**< Steiner tree nodes */
2407  )
2408 {
2409  SCIP_Bool isValid;
2410  STP_Bool* solNodes_new;
2411  const int nnodes = graph->knots;
2412 
2413  SCIP_CALL_ABORT( SCIPallocBufferArray(scip, &solNodes_new, nnodes) );
2414 
2415  for( int i = 0; i < nnodes; ++i )
2416  solNodes_new[i] = FALSE;
2417 
2418  for( int i = 0; i < nnodes; ++i )
2419  {
2420  const int edge = linkcutNodes[i].edge;
2421 
2422  if( edge >= 0 )
2423  {
2424  assert(graph->tail[edge] == i);
2425 
2426  solNodes_new[graph->head[edge]] = TRUE;
2427  solNodes_new[graph->tail[edge]] = TRUE;
2428  }
2429  }
2430 
2431  isValid = TRUE;
2432 
2433  for( int i = 0; i < nnodes; ++i )
2434  {
2435  if( solNodes_new[i] != solNodes[i] )
2436  {
2437  isValid = FALSE;
2438  break;
2439  }
2440  }
2441 
2442  SCIPfreeBufferArray(scip, &solNodes_new);
2443 
2444  return isValid;
2445 }
2446 #endif
2447 
2448 
2449 /** reduces solution degree */
2450 static inline
2452  const GRAPH* graph, /**< graph data structure */
2453  int node, /**< replacement edge */
2454  INSERT* insertData /**< insertion data */
2455 )
2456 {
2457  int* const solDegreeNonTerm = insertData->solDegreeNonTerm;
2458 
2459  assert(node >= 0);
2460  assert(insertData->solNodes[node]);
2461 
2462  if( Is_term(graph->term[node]) || Is_pseudoTerm(graph->term[node]) )
2463  {
2464  assert(solDegreeNonTerm[node] == UNKNOWN);
2465  }
2466  else
2467  {
2468  solDegreeNonTerm[node] -= 1;
2469 
2470  assert(solDegreeNonTerm[node] >= 0);
2471  }
2472 }
2473 
2474 
2475 /** increases solution degree */
2476 static inline
2478  const GRAPH* graph, /**< graph data structure */
2479  int node, /**< replacement edge */
2480  INSERT* insertData /**< insertion data */
2481 )
2482 {
2483  int* const solDegreeNonTerm = insertData->solDegreeNonTerm;
2484 
2485  assert(node >= 0);
2486  assert(insertData->solNodes[node]);
2487 
2488  if( Is_term(graph->term[node]) || Is_pseudoTerm(graph->term[node]) )
2489  {
2490  assert(solDegreeNonTerm[node] == UNKNOWN);
2491  }
2492  else
2493  {
2494  solDegreeNonTerm[node] += 1;
2495  }
2496 }
2497 
2498 /** subroutine for insertion heuristic */
2499 static
2501  SCIP_HEURDATA* heurdata, /**< heuristic data */
2502  const GRAPH* graph, /**< graph data structure */
2503  const int* solEdges, /**< Steiner tree edges */
2504  int* solDegreeNonTerm, /**< degree */
2505  SCIP_Bool* nodeIsBlocked, /**< blocked */
2506  int* vertices /**< vertices permuted */
2507  )
2508 {
2509  const int nnodes = graph->knots;
2510  const int nedges = graph->edges;
2511 
2512  assert(solDegreeNonTerm);
2513 
2514  BMSclearMemoryArray(solDegreeNonTerm, nnodes);
2515 
2516  for( int e = 0; e < nedges; e++ )
2517  {
2518  if( solEdges[e] == CONNECT )
2519  {
2520  solDegreeNonTerm[graph->tail[e]]++;
2521  solDegreeNonTerm[graph->head[e]]++;
2522  }
2523  }
2524 
2525  for( int i = 0; i < nnodes; ++i )
2526  {
2527  vertices[i] = i;
2528  nodeIsBlocked[i] = FALSE;
2529 
2530  if( Is_term(graph->term[i]) || Is_pseudoTerm(graph->term[i]) )
2531  solDegreeNonTerm[i] = UNKNOWN;
2532  }
2533 
2534  SCIPrandomPermuteIntArray(heurdata->randnumgen, vertices, 0, nnodes);
2535 }
2536 
2537 /** subroutine for insertion heuristic: prepare insertion of new vertex */
2538 static
2540  SCIP* scip, /**< SCIP data structure */
2541  const GRAPH* graph, /**< graph data structure */
2542  int v_insert, /**< the vertex to add */
2543  int initialEdge, /**< first edge from node to tree */
2544  LCNODE* linkcutNodes, /**< Steiner tree nodes */
2545  INSERT* insertData, /**< insertion data (in/out) */
2546  SCIP_Real* diff /**< the initial diff (out) */
2547  )
2548 {
2549  const int candHeadInitial = graph->head[initialEdge];
2550  const int stp_type = graph->stp_type;
2551  const SCIP_Bool mw = (stp_type == STP_MWCSP || stp_type == STP_RMWCSP);
2552 
2553  assert(graph->tail[initialEdge] == v_insert);
2554  assert(!Is_term(graph->term[v_insert]));
2555 
2556  SCIPlinkcuttreeLink(linkcutNodes, v_insert, candHeadInitial, initialEdge);
2557 
2558  SCIPdebugMessage("try to insert vertex %d \n", v_insert);
2559 
2560  if( mw )
2561  {
2562  assert(!SCIPisPositive(scip, graph->prize[v_insert]));
2563  *diff = -graph->prize[v_insert];
2564  }
2565  else
2566  {
2567  const SCIP_Bool pc = graph_pc_isPc(graph);
2568 
2569  *diff = insertData->edgecosts[initialEdge];
2570 
2571  if( pc )
2572  {
2573  *diff -= graph->prize[v_insert];
2574  }
2575  }
2576 
2577  assert(insertData->solDegreeNonTerm[v_insert] == 0
2578  || (insertData->solDegreeNonTerm[v_insert] == UNKNOWN && graph_pc_isPcMw(graph)));
2579  assert(!insertData->nodeIsBlocked[v_insert]);
2580 
2581  insertData->nInsertions = 0;
2582  insertData->insertionVertex = v_insert;
2583  insertData->solNodes[v_insert] = TRUE;
2584  insertData->nodeIsBlocked[v_insert] = TRUE;
2585  insertionIncrementSolDegree(graph, v_insert, insertData);
2586  insertionIncrementSolDegree(graph, candHeadInitial, insertData);
2587 }
2588 
2589 
2590 /** remove (blocked) chains from tree */
2591 static
2593  const GRAPH* graph, /**< graph data structure */
2594  int v_insert, /**< the vertex to insert */
2595  LCNODE* linkcutNodes, /**< Steiner tree nodes */
2596  INSERT* insertData /**< insertion data */
2597 )
2598 {
2599  SCIP_Bool* const nodeIsBlocked = insertData->nodeIsBlocked;
2600  const int* const blockedList = insertData->blockedList;
2601  STP_Bool* const solNodes = insertData->solNodes;
2602  const int size = insertData->blockedListSize;
2603  int* const solDegreeNonTerm = insertData->solDegreeNonTerm;
2604 
2605  for( int i = 0; i < size; i++ )
2606  {
2607  const int node = blockedList[i];
2608 
2609  assert(node >= 0 && node < graph->knots);
2610  assert(nodeIsBlocked[node]);
2611  assert(solNodes[node]);
2612 
2613  nodeIsBlocked[node] = FALSE;
2614  solNodes[node] = FALSE;
2615 
2616  if( !(Is_term(graph->term[node]) || Is_pseudoTerm(graph->term[node])) )
2617  {
2618  assert(2 == solDegreeNonTerm[node]);
2619 
2620  solDegreeNonTerm[node] = 0;
2621  }
2622  else
2623  {
2624  assert(UNKNOWN == solDegreeNonTerm[node]);
2625  }
2626 
2627  SCIPlinkcuttreeInitNode(&linkcutNodes[node]);
2628  }
2629 
2630 
2631  assert(nodeIsBlocked[v_insert]);
2632  nodeIsBlocked[v_insert] = FALSE;
2633 
2634  for( int k = 0; k < graph->knots; k++ )
2635  assert(!nodeIsBlocked[k]);
2636 
2637  insertData->blockedListSize = 0;
2638 }
2639 
2640 
2641 /** remove (blocked) chains from tree */
2642 static inline
2644  INSERT* insertData /**< insertion data */
2645 )
2646 {
2647  SCIP_Bool* const nodeIsBlocked = insertData->nodeIsBlocked;
2648  const int* const blockedList = insertData->blockedList;
2649  const int size = insertData->blockedListSize;
2650 
2651  for( int i = 0; i < size; i++ )
2652  {
2653  const int node = blockedList[i];
2654 
2655  assert(node >= 0);
2656  assert(nodeIsBlocked[node]);
2657  assert(2 == insertData->solDegreeNonTerm[node]);
2658 
2659  nodeIsBlocked[node] = FALSE;
2660  }
2661 
2662  assert(nodeIsBlocked[insertData->insertionVertex]);
2663  nodeIsBlocked[insertData->insertionVertex] = FALSE;
2664 }
2665 
2666 
2667 /** mark inner chain nodes as blocked */
2668 static inline
2670  const GRAPH* graph, /**< graph data structure */
2671  const LCNODE* lctree, /**< tree */
2672  int chainfirst_index, /**< first chain entry */
2673  int chainlast_index, /**< last chain entry (inside) */
2674  INSERT* insertData /**< insertion data */
2675 )
2676 {
2677  if( chainfirst_index != chainlast_index )
2678  {
2679  SCIP_Bool* const nodeIsBlocked = insertData->nodeIsBlocked;
2680  int* const blockedList = insertData->blockedList;
2681 
2682  for( int node = chainfirst_index; node != chainlast_index; node = lctree[node].parent )
2683  {
2684  int head;
2685 
2686  assert(node != -1);
2687  assert(lctree[node].edge >= 0);
2688 
2689  head = graph->head[lctree[node].edge];
2690 
2691  assert(!nodeIsBlocked[head]);
2692 
2693  nodeIsBlocked[head] = TRUE;
2694  blockedList[insertData->blockedListSize++] = head;
2695  }
2696  }
2697 }
2698 
2699 
2700 /** reinsert (blocked) chains to tree */
2701 static
2703  const GRAPH* graph, /**< graph data structure */
2704  const int* insertCands, /**< insertion candidates */
2705  int lcvertex_insert, /**< the vertex tested for insertion */
2706  LCNODE* linkcutNodes, /**< Steiner tree nodes */
2707  INSERT* insertData /**< insertion data */
2708 )
2709 {
2710  int* chainEnds = insertData->chainEnds;
2711  const int* const addedEdges = insertData->addedEdges;
2712  const int* const cutedgesStart = insertData->cutedgesStart;
2713  const int* const cutedgesEnd = insertData->cutedgesEnd;
2714  const int* const graphTail = graph->tail;
2715  const int* const graphHead = graph->head;
2716  const int v = insertData->insertionVertex;
2717 
2718  SCIPlinkcuttreeEvert(linkcutNodes, lcvertex_insert);
2719 
2720  insertionResetBlockedNodes(insertData);
2721 
2722  /* reinsert each chain and remove the */
2723  for( int k = insertData->nInsertions - 1; k >= 0; k-- )
2724  {
2725  const int cutedgeFirst = cutedgesStart[k];
2726  const int cutedgeLast = cutedgesEnd[k];
2727  const int firstNode = graphTail[cutedgeFirst];
2728  const int secondNode = graphHead[cutedgeFirst];
2729  const int lastNode = chainEnds[k];
2730 
2731  assert(cutedgeLast >= 0);
2732  assert(firstNode == insertData->chainStarts[k]);
2733 
2734  /* remove the newly added edge */
2735  SCIPlinkcuttreeCutNode(&linkcutNodes[graphHead[addedEdges[k]]]);
2736  insertionDecrementSolDegree(graph, graphHead[addedEdges[k]], insertData);
2737 
2738  /* re-link the tail of the chain */
2739  SCIPlinkcuttreeEvert(linkcutNodes, firstNode);
2740  SCIPlinkcuttreeLink(linkcutNodes, firstNode, secondNode, cutedgeFirst);
2741 
2742  /* re-link the head of the chain (if not already done) */
2743  if( lastNode != firstNode )
2744  SCIPlinkcuttreeLink(linkcutNodes, lastNode, graphHead[cutedgeLast], cutedgeLast);
2745 
2746  /* reset solution degrees of chain border vertices */
2747  insertionIncrementSolDegree(graph, graphTail[cutedgeFirst], insertData);
2748  insertionIncrementSolDegree(graph, graphHead[cutedgeLast], insertData);
2749  }
2750 
2751  /* finally, cut the edge added first (if it had been cut during the insertion process, it would have been restored above) */
2752  SCIPlinkcuttreeEvert(linkcutNodes, lcvertex_insert);
2753  SCIPlinkcuttreeCutNode(&linkcutNodes[graphHead[insertCands[0]]]);
2754  insertionDecrementSolDegree(graph, graphHead[insertCands[0]], insertData);
2755 
2756  if( Is_term(graph->term[v]) || Is_pseudoTerm(graph->term[v]) )
2757  {
2758  assert(graph_pc_isPcMw(graph));
2759  insertData->solDegreeNonTerm[v] = UNKNOWN;
2760  }
2761  else
2762  {
2763  insertData->solDegreeNonTerm[v] = 0;
2764  }
2765 
2766  assert(insertData->solNodes[v]);
2767  insertData->solNodes[v] = FALSE;
2768 
2769  for( int k = 0; k < graph->knots; k++ )
2770  assert(!insertData->nodeIsBlocked[k]);
2771 
2772  insertData->blockedListSize = 0;
2773 }
2774 
2775 
2776 /** replaces chain by a single edge */
2777 static
2779  const GRAPH* graph, /**< graph data structure */
2780  int newedge, /**< replacement edge */
2781  LCNODE* lctree, /**< tree */
2782  INSERT* insertData, /**< insertion data */
2783  int v_lc, /**< current vertex */
2784  int headCurr_lc, /**< head of newedge */
2785  int chainfirst_index, /**< first chain entry */
2786  int chainlast_index /**< last chain entry (inside) */
2787 )
2788 {
2789  LCNODE* chainfirst = &lctree[chainfirst_index];
2790  LCNODE* chainlast = &lctree[chainlast_index];
2791  int* const chainStarts = insertData->chainStarts;
2792  int* const chainEnds = insertData->chainEnds;
2793  int* const cutedgesStart = insertData->cutedgesStart;
2794  int* const cutedgesEnd = insertData->cutedgesEnd;
2795  int* const addedEdges = insertData->addedEdges;
2796  const int nInsertions = insertData->nInsertions;
2797  const int newhead = graph->head[newedge];
2798  const int v_insert = insertData->insertionVertex;
2799 
2800  assert(chainlast->edge >= 0);
2801  assert(graph->tail[newedge] == v_insert);
2802 
2803  cutedgesStart[nInsertions] = chainfirst->edge;
2804  cutedgesEnd[nInsertions] = chainlast->edge;
2805  chainStarts[nInsertions] = chainfirst_index;
2806  chainEnds[nInsertions] = chainlast_index;
2807  addedEdges[nInsertions] = newedge;
2808 
2809  insertionIncrementSolDegree(graph, v_insert, insertData);
2810  insertionIncrementSolDegree(graph, newhead, insertData);
2811 
2812  /* decrease the degree of the chain border vertices */
2813  insertionDecrementSolDegree(graph, graph->tail[chainfirst->edge], insertData);
2814  insertionDecrementSolDegree(graph, graph->head[chainlast->edge], insertData);
2815 
2816  SCIPdebugMessage("remove chain %d->%d \n", graph->tail[chainfirst->edge], graph->head[chainlast->edge]);
2817 
2818  insertionBlockChain(graph, lctree, chainfirst_index, chainlast_index, insertData);
2819 
2820  SCIPlinkcuttreeCutNode(chainfirst);
2821  SCIPlinkcuttreeCutNode(chainlast);
2822 
2823  SCIPlinkcuttreeLink(lctree, v_lc, headCurr_lc, newedge);
2824  assert(lctree[v_lc].edge == newedge);
2825 
2826  insertData->nInsertions++;
2827 }
2828 
2829 
2830 /** perform local vertex insertion heuristic on given Steiner tree */
2831 static
2833  SCIP* scip, /**< SCIP data structure */
2834  SCIP_HEURDATA* heurdata, /**< heuristic data */
2835  const GRAPH* graph, /**< graph data structure */
2836  const STP_Bool* solNodes, /**< Steiner tree nodes */
2837  int vertex, /**< vertex to be inserted */
2838  int* insertcands, /**< candidates */
2839  int* ncands /**< number of candidates */
2840  )
2841 {
2842  int insertcount = 0;
2843  const SCIP_Bool mwpc = graph_pc_isPcMw(graph);
2844  const SCIP_Bool rpcmw = graph_pc_isRootedPcMw(graph);
2845 
2846  assert(!Is_term(graph->term[vertex]));
2847  assert(heurdata);
2848 
2849  for( int oedge = graph->outbeg[vertex]; oedge != EAT_LAST; oedge = graph->oeat[oedge] )
2850  {
2851  const int head = graph->head[oedge];
2852 
2853  if( !solNodes[head] )
2854  continue;
2855 
2856  /* skip dummy terminals */
2857  if( mwpc )
2858  {
2859  if( Is_term(graph->term[head]) && !graph_pc_knotIsFixedTerm(graph, head) )
2860  continue;
2861 
2862  /* note: for unrooted problems the (artificial) source is also a fixed terminal */
2863  if( !rpcmw && head == graph->source )
2864  continue;
2865  }
2866 
2867  assert(SCIPisLT(scip, graph->cost[oedge], FARAWAY));
2868 
2869  insertcands[insertcount++] = oedge;
2870  }
2871 
2872  if( insertcount >= 3 )
2873  SCIPrandomPermuteIntArray(heurdata->randnumgen, insertcands, 0, insertcount);
2874 
2875  *ncands = insertcount;
2876 }
2877 
2878 
2879 /** perform local vertex insertion heuristic on given Steiner tree */
2880 static
2882  SCIP* scip, /**< SCIP data structure */
2883  SCIP_HEURDATA* heurdata, /**< heuristic data or NULL */
2884  const GRAPH* graph, /**< graph data structure */
2885  STP_Bool* solNodes, /**< Steiner tree nodes */
2886  LCNODE* linkcutNodes, /**< Steiner tree nodes */
2887  int* solEdges /**< array indicating whether an arc is part of the solution (CONNECTED/UNKNOWN) */
2888  )
2889 {
2890  int* chainStarts;
2891  int* chainEnds;
2892  int* insertCands = NULL;
2893  int* addedEdges = NULL;
2894  int* cutedgesStart = NULL;
2895  int* cutedgesEnd = NULL;
2896  int* solDegree = NULL;
2897  int* vertices = NULL;
2898  int* blockedList = NULL;
2899  SCIP_Bool* nodeIsBlocked = NULL;
2900 
2901  int newnverts = 0;
2902  const int nnodes = graph->knots;
2903  const int nedges = graph->edges;
2904  const int root = graph->source;
2905  const int stp_type = graph->stp_type;
2906  const SCIP_Bool mw = (stp_type == STP_MWCSP || stp_type == STP_RMWCSP);
2907  const SCIP_Bool pc = graph_pc_isPc(graph);
2908  const SCIP_Bool mwpc = graph_pc_isPcMw(graph);
2909  SCIP_Bool solimproved = TRUE;
2910  SCIP_Real* edgecosts;
2911 
2912 
2913 #ifndef NDEBUG
2914  const SCIP_Real initialobj = solstp_getObjBounded(graph, solEdges, 0.0, nedges);
2915  SCIP_Real diff_total = 0.0;
2916 #endif
2917 
2918  if( stp_type != STP_SPG && stp_type != STP_RSMT && stp_type != STP_OARSMT && stp_type != STP_GSTP && !mwpc )
2919  {
2920  SCIPdebugMessage("vertex inclusion does not work for current problem type \n");
2921  return SCIP_OKAY;
2922  }
2923 
2924  if( pc )
2925  {
2926  SCIP_CALL( SCIPallocBufferArray(scip, &edgecosts, nedges) );
2927 
2928  graph_pc_getOrgCosts(scip, graph, edgecosts);
2929  }
2930  else
2931  {
2932  edgecosts = graph->cost;
2933  }
2934 
2935  SCIP_CALL( SCIPallocBufferArray(scip, &chainStarts, nnodes) );
2936  SCIP_CALL( SCIPallocBufferArray(scip, &chainEnds, nnodes) );
2937  SCIP_CALL( SCIPallocBufferArray(scip, &blockedList, nnodes) );
2938  SCIP_CALL( SCIPallocBufferArray(scip, &nodeIsBlocked, nnodes) );
2939  SCIP_CALL( SCIPallocBufferArray(scip, &vertices, nnodes) );
2940  SCIP_CALL( SCIPallocBufferArray(scip, &insertCands, nnodes) );
2941  SCIP_CALL( SCIPallocBufferArray(scip, &addedEdges, nnodes) );
2942  SCIP_CALL( SCIPallocBufferArray(scip, &cutedgesStart, nnodes) );
2943  SCIP_CALL( SCIPallocBufferArray(scip, &cutedgesEnd, nnodes) );
2944  SCIP_CALL( SCIPallocBufferArray(scip, &solDegree, nnodes) );
2945 
2946  insertionInit(heurdata, graph, solEdges, solDegree, nodeIsBlocked, vertices);
2947 
2948  assert(solDegIsValid(scip, graph, solDegree, linkcutNodes));
2949  assert(solNodeIsValid(scip, graph, solNodes, linkcutNodes));
2950  assert(linkcutNodes[graph->source].edge == -1);
2951 
2952  /* main loop */
2953  while( solimproved )
2954  {
2955  INSERT insertData = { .chainStarts = chainStarts, .chainEnds = chainEnds, .edgecosts = edgecosts,
2956  .solDegreeNonTerm = solDegree, .solNodes = solNodes, .nodeIsBlocked = nodeIsBlocked,
2957  .cutedgesStart = cutedgesStart, .cutedgesEnd = cutedgesEnd, .addedEdges = addedEdges,
2958  .blockedList = blockedList, .blockedListSize = 0, .nInsertions = 0, .insertionVertex = -1 };
2959 
2960  solimproved = FALSE;
2961 
2962  for( int i = 0; i < nnodes; i++ )
2963  {
2964  SCIP_Real diff;
2965  int ninsertcands = 0;
2966  const int v = vertices[i];
2967 
2968  assert(v >= 0 && v < graph->knots);
2969 
2970  if( solNodes[v] || graph->grad[v] <= 1 )
2971  continue;
2972 
2973  assert(!Is_term(graph->term[v]));
2974 
2975  insertionGetCandidateEdges(scip, heurdata, graph, solNodes, v, insertCands, &ninsertcands);
2976 
2977  /* if there are less than two edges connecting node i and the current tree, continue */
2978  if( ninsertcands <= 1 )
2979  continue;
2980 
2981  insertionInitInsert(scip, graph, v, insertCands[0], linkcutNodes, &insertData, &diff);
2982 
2983  /* try to add additional edges between new vertex and tree */
2984  for( int k = 1; k < ninsertcands; k++ )
2985  {
2986  int chainfirst = -1;
2987  int chainlast = -1;
2988  const int insertEdgeCurr = insertCands[k];
2989  const int insertHeadCurr = graph->head[insertEdgeCurr];
2990 
2991  if( nodeIsBlocked[insertHeadCurr] )
2992  {
2993  assert(k > 1);
2994 
2995  continue;
2996  }
2997 
2998  SCIPdebugMessage("insertHeadCurr=%d \n", insertHeadCurr);
2999 
3000  SCIPlinkcuttreeEvert(linkcutNodes, v);
3001 
3002  if( mw )
3003  {
3004  const SCIP_Real chainweight = SCIPlinkcuttreeFindMinChainMw(scip, linkcutNodes, graph->prize,
3005  graph->head, solDegree, nodeIsBlocked, insertHeadCurr,
3006  &chainfirst, &chainlast);
3007 
3008  if( SCIPisLT(scip, chainweight, 0.0) )
3009  {
3010  diff += chainweight;
3011  insertionReplaceChain(graph, insertEdgeCurr, linkcutNodes, &insertData, v, insertHeadCurr, chainfirst, chainlast);
3012  }
3013  }
3014  else
3015  {
3016  const SCIP_Real chainweight = SCIPlinkcuttreeFindMaxChain(scip, linkcutNodes, edgecosts, pc ? graph->prize : NULL,
3017  graph->head, solDegree, nodeIsBlocked, insertHeadCurr,
3018  &chainfirst, &chainlast);
3019 
3020  SCIPdebugMessage("chainweight=%f edgecost=%f \n", chainweight, edgecosts[insertEdgeCurr]);
3021 
3022  /* note: comparison needs to be strict to avoid (redundant) removal of current edge */
3023  if( SCIPisGT(scip, chainweight, edgecosts[insertEdgeCurr]) )
3024  {
3025  diff += edgecosts[insertEdgeCurr];
3026  diff -= chainweight;
3027  insertionReplaceChain(graph, insertEdgeCurr, linkcutNodes, &insertData, v, insertHeadCurr, chainfirst, chainlast);
3028  }
3029  }
3030  }
3031 
3032  /* is the new tree better? */
3033  if( SCIPisNegative(scip, diff) )
3034  {
3035  SCIPlinkcuttreeEvert(linkcutNodes, root);
3036  solimproved = TRUE;
3037  newnverts++;
3038  assert(solNodes[v]);
3039 
3040  insertionFinalizeReplacement(graph, v, linkcutNodes, &insertData);
3041 
3042  SCIPdebugMessage("Inclusion: ADDED VERTEX %d \n", v);
3043 #ifndef NDEBUG
3044  diff_total -= diff;
3045 #endif
3046  }
3047  else
3048  {
3049  /* restore the old tree */
3050  insertionRestoreTree(graph, insertCands, v, linkcutNodes, &insertData);
3051  }
3052 
3053  assert(solDegIsValid(scip, graph, solDegree, linkcutNodes));
3054  assert(solNodeIsValid(scip, graph, solNodes, linkcutNodes));
3055 
3056  } /* for i 0:nnodes */
3057  }
3058 
3059  /* free buffer memory */
3060  SCIPfreeBufferArray(scip, &solDegree);
3061  SCIPfreeBufferArray(scip, &cutedgesEnd);
3062  SCIPfreeBufferArray(scip, &cutedgesStart);
3063  SCIPfreeBufferArray(scip, &addedEdges);
3064  SCIPfreeBufferArray(scip, &insertCands);
3065  SCIPfreeBufferArray(scip, &vertices);
3066  SCIPfreeBufferArray(scip, &blockedList);
3067  SCIPfreeBufferArray(scip, &nodeIsBlocked);
3068  SCIPfreeBufferArray(scip, &chainEnds);
3069  SCIPfreeBufferArray(scip, &chainStarts);
3070 
3071  if( pc )
3072  SCIPfreeBufferArray(scip, &edgecosts);
3073 
3074  for( int e = 0; e < nedges; e++ )
3075  solEdges[e] = UNKNOWN;
3076 
3077  if( newnverts > 0 )
3078  {
3079  SCIP_CALL( solstp_prune(scip, graph, solEdges, solNodes) );
3080 
3081  for( int i = 0; i < nnodes; i++ )
3082  SCIPlinkcuttreeInitNode(&linkcutNodes[i]);
3083 
3084  /* create a link-cut tree representing the current Steiner tree */
3085  for( int e = 0; e < nedges; e++ )
3086  {
3087  if( solEdges[e] == CONNECT )
3088  {
3089  assert(solNodes[graph->tail[e]] && solNodes[graph->head[e]]);
3090 
3091  SCIPlinkcuttreeLink(linkcutNodes, graph->head[e], graph->tail[e], flipedge(e));
3092  }
3093  }
3094  SCIPlinkcuttreeEvert(linkcutNodes, root);
3095  }
3096  else
3097  {
3098  SCIPlinkcuttreeEvert(linkcutNodes, root);
3099  for( int i = 0; i < nnodes; i++ )
3100  {
3101  if( solNodes[i] && linkcutNodes[i].edge != -1 )
3102  solEdges[flipedge(linkcutNodes[i].edge)] = CONNECT;
3103  }
3104  }
3105 
3106 #ifndef NDEBUG
3107  {
3108  const SCIP_Real newobj = solstp_getObjBounded(graph, solEdges, 0.0, nedges);
3109  SCIPdebugMessage("vertex inclusion obj before/after: %f/%f \n", initialobj, newobj);
3110  assert(SCIPisLE(scip, newobj + diff_total, initialobj));
3111  assert(SCIPisGE(scip, diff_total, 0.0));
3112  }
3113 #endif
3114 
3115  assert(solstp_isValid(scip, graph, solEdges));
3116 
3117  return SCIP_OKAY;
3118 }
3119 
3120 /** perform local vertex insertion heuristic on given Steiner tree */
3121 static
3123  SCIP* scip, /**< SCIP data structure */
3124  GRAPH* graph, /**< graph data structure */
3125  SCIP_Bool useFast, /**< fast variant? */
3126  STP_Bool* solNodes, /**< Steiner tree nodes */
3127  LCNODE* linkcutNodes, /**< Steiner tree nodes */
3128  int* solEdges, /**< array indicating whether an arc is part of the solution (CONNECTED/UNKNOWN) */
3129  SCIP_Bool* success /**< solution improved? */
3130  )
3131 {
3132  UF uf; /* union-find */
3133  STP_Vectype(int)* blists_start = NULL;
3134  PATH* vnoipath = NULL;
3135  STP_Vectype(int)* lvledges_start = NULL; /* horizontal edges */
3136  PHNODE** boundpaths = NULL;
3137  PHNODE* pheapelements = NULL;
3138  SCIP_Real* memvdist = NULL;
3139  SCIP_Real* edgecostBiased_pc = NULL;
3140  SCIP_Real* edgecostOrg_pc = NULL;
3141  SCIP_Real* prize_pc = NULL;
3142  int* vnoibase = NULL;
3143  int* kpedges = NULL;
3144  int* kpnodes = NULL;
3145  int* dfstree = NULL;
3146  int* newedges = NULL;
3147  int* memvbase = NULL;
3148  int* pheapsize = NULL;
3149  int* boundedges = NULL;
3150  int* meminedges = NULL;
3151  int* supernodes = NULL;
3152  int* prizemarklist = NULL;
3153  STP_Bool* pinned = NULL;
3154  STP_Bool* scanned = NULL;
3155  STP_Bool* supernodesmark = NULL;
3156  STP_Bool* prizemark = NULL;
3157  const int root = graph->source;
3158  const int nnodes = graph->knots;
3159  const int nedges = graph->edges;
3160  const int maxnrestarts = (useFast ? LOCAL_MAXRESTARTS_FAST : LOCAL_MAXRESTARTS);
3161  const int solroot = pcmwGetSolRoot(graph, solEdges);
3162  const STP_Bool mwpc = graph_pc_isPcMw(graph);
3163  SCIP_Bool solimproved = FALSE;
3164 
3165 #ifndef NDEBUG
3166  const SCIP_Real initialobj = solstp_getObjBounded(graph, solEdges, 0.0, graph->edges);
3167  SCIP_Real objimprovement = 0.0;
3168 #endif
3169 
3170  *success = FALSE;
3171 
3172  assert(mwpc || solroot == -1);
3173 
3174  /* memory needed for both Key-Path Elimination and Exchange */
3175  SCIP_CALL( SCIPallocBufferArray(scip, &vnoipath, nnodes) );
3176  SCIP_CALL( SCIPallocBufferArray(scip, &vnoibase, nnodes) );
3177 
3178  /* only needed for Key-Path Elimination */
3179  SCIP_CALL( SCIPallocBufferArray(scip, &newedges, nedges) );
3180  SCIP_CALL( SCIPallocBufferArray(scip, &lvledges_start, nnodes) );
3181  SCIP_CALL( SCIPallocBufferArray(scip, &boundedges, nedges) );
3182 
3183  /* memory needed for both Key-Path Elimination and Exchange */
3184  if( mwpc )
3185  {
3186  if( graph_pc_isPc(graph) )
3187  {
3188  SCIP_CALL(SCIPallocBufferArray(scip, &edgecostOrg_pc, nedges));
3189 
3190  graph_pc_getOrgCosts(scip, graph, edgecostOrg_pc);
3191  }
3192 
3193  SCIP_CALL(SCIPallocBufferArray(scip, &edgecostBiased_pc, nedges));
3194  SCIP_CALL(SCIPallocBufferArray(scip, &prize_pc, nnodes));
3195  SCIP_CALL(SCIPallocBufferArray(scip, &prizemark, nnodes));
3196  SCIP_CALL(SCIPallocBufferArray(scip, &prizemarklist, nnodes));
3197 
3198  for( int k = 0; k < nnodes; k++ )
3199  prizemark[k] = FALSE;
3200 
3201  if( graph_pc_isPc(graph) )
3202  {
3203  graph_pc_getBiased(graph, edgecostBiased_pc, prize_pc);
3204  }
3205  else
3206  {
3207  BMScopyMemoryArray(edgecostBiased_pc, graph->cost, graph->edges);
3208  BMScopyMemoryArray(prize_pc, graph->prize, graph->knots);
3209  }
3210  }
3211 
3212  SCIP_CALL( SCIPallocBufferArray(scip, &scanned, nnodes) );
3213  SCIP_CALL( SCIPallocBufferArray(scip, &pheapsize, nnodes) );
3214  SCIP_CALL( SCIPallocBufferArray(scip, &blists_start, nnodes) );
3215  SCIP_CALL( SCIPallocBufferArray(scip, &memvbase, nnodes) );
3216  SCIP_CALL( SCIPallocBufferArray(scip, &memvdist, nnodes) );
3217  SCIP_CALL( SCIPallocBufferArray(scip, &meminedges, nnodes) );
3218  SCIP_CALL( SCIPallocBufferArray(scip, &boundpaths, nnodes) );
3219  SCIP_CALL( SCIPallocBufferArray(scip, &pinned, nnodes) );
3220  SCIP_CALL( SCIPallocBufferArray(scip, &dfstree, nnodes) );
3221  SCIP_CALL( SCIPallocBufferArray(scip, &supernodesmark, nnodes) );
3222  SCIP_CALL( SCIPallocBufferArray(scip, &supernodes, nnodes) );
3223  SCIP_CALL( SCIPallocBufferArray(scip, &kpnodes, nnodes) );
3224  SCIP_CALL( SCIPallocBufferArray(scip, &kpedges, nnodes) );
3225  SCIP_CALL( SCIPallocMemoryArray(scip, &pheapelements, nedges) );
3226 
3227  for( int k = 0; k < nnodes; k++ )
3228  graph->mark[k] = (graph->grad[k] > 0);
3229 
3230  graph->mark[root] = TRUE;
3231 
3232  SCIP_CALL( SCIPStpunionfindInit(scip, &uf, nnodes) );
3233 
3234  BMSclearMemoryArray(blists_start, nnodes);
3235  BMSclearMemoryArray(lvledges_start, nnodes);
3236 
3237  /* main loop */
3238  for( int nruns = 0, localmoves = 1; nruns < maxnrestarts && localmoves > 0; nruns++ )
3239  {
3240  VNOILOC vnoiData = { .vnoi_path = vnoipath, .vnoi_base = vnoibase, .memvdist = memvdist, .memvbase = memvbase,
3241  .meminedges = meminedges, .vnoi_nodestate = graph->path_state, .nmems = 0, .nkpnodes = -1 };
3242  KPATHS keypathsData = { .kpnodes = kpnodes, .kpedges = kpedges, .kpcost = 0.0, .nkpnodes = 0, .nkpedges = 0,
3243  .kptailnode = -1 };
3244  CONN connectivityData = { .blists_start = blists_start, .pheap_boundpaths = boundpaths, .lvledges_start = lvledges_start,
3245  .pheap_sizes = pheapsize, .uf = &uf, .boundaryedges = boundedges, .nboundaryedges = 0,
3246  .pheap_elements = pheapelements };
3247  SOLTREE soltreeData = { .solNodes = solNodes, .linkcutNodes = linkcutNodes, .solEdges = solEdges, .nodeIsPinned = pinned,
3248  .nodeIsScanned = scanned, .newedges = newedges };
3249  SGRAPH supergraphData = { .supernodes = supernodes, .nodeIsLowSupernode = supernodesmark,
3250  .mst = NULL, .mstcost = 0.0, .nsupernodes = 0 };
3251  PCMW pcmwData = { .prize_biased = prize_pc, .edgecost_biased = edgecostBiased_pc, .edgecost_org = edgecostOrg_pc,
3252  .prizemark = prizemark, .prizemarklist = prizemarklist, .nprizemarks = 0, .isActive = graph_pc_isPcMw(graph),
3253  .solroot = solroot };
3254  int nstnodes = 0;
3255 
3256  localmoves = 0;
3257 
3258  /* find a DFS order of the ST nodes */
3259  dfsorder(graph, solEdges, root, &nstnodes, dfstree);
3260 
3261  /* initialize data structures */
3262  for( int k = 0; k < nnodes; k++ )
3263  {
3264  pinned[k] = FALSE;
3265  scanned[k] = FALSE;
3266  supernodesmark[k] = FALSE;
3267  }
3268 
3269  for( int e = 0; e < nedges; e++ )
3270  pheapelements[e].element = -1;
3271 
3272  for( int e = 0; e < nedges; e++ )
3273  newedges[e] = UNKNOWN;
3274 
3275  if( mwpc )
3276  {
3277  assert(graph->extended);
3278  SCIP_CALL( pcmwUpdate(scip, graph, &soltreeData, &pcmwData) );
3279 
3280  /* compute a Voronoi diagram with the Steiner tree nodes as bases */
3281  graph_voronoi(graph, pcmwData.edgecost_biased, NULL, soltreeData.solNodes, vnoiData.vnoi_base, vnoiData.vnoi_path);
3282  }
3283  else
3284  {
3285  graph_voronoi(graph, graph->cost, NULL, soltreeData.solNodes, vnoiData.vnoi_base, vnoiData.vnoi_path);
3286  }
3287 
3288  SCIP_CALL( connectivityDataInit(scip, graph, &vnoiData, &soltreeData, &pcmwData, &connectivityData) );
3289 
3290  /* henceforth, the union-find structure will be used on the Steiner tree */
3291  assert(uf.nElements == nnodes);
3292  SCIPStpunionfindClear(scip, &uf);
3293 
3294  /* main loop visiting all nodes of the current Steiner tree in post-order */
3295  for( int dfstree_pos = 0; dfstree[dfstree_pos] != root; dfstree_pos++ )
3296  {
3297  /* current crucial node */
3298  const int crucnode = dfstree[dfstree_pos];
3299  int nheapelems = -1;
3300 
3301  assert(!scanned[crucnode]);
3302  scanned[crucnode] = TRUE;
3303 
3304  SCIPdebugMessage("iteration %d (crucial node: %d) \n", dfstree_pos, crucnode);
3305 
3306  if( solroot == crucnode )
3307  {
3308  assert(mwpc && !graph_pc_isRootedPcMw(graph));
3309  continue;
3310  }
3311 
3312  /* has the node been temporarily removed from the Steiner tree? */
3313  if( !graph->mark[crucnode] )
3314  continue;
3315 
3316  /* key vertex elimination: */
3317  /* is node 'crucnode' a removable crucial node? (i.e. not pinned or a terminal) */
3318  if( !pinned[crucnode] && !Is_term(graph->term[crucnode]) && nodeIsCrucial(graph, solEdges, crucnode) )
3319  {
3320  SCIP_Bool allgood;
3321 
3322  getKeyPathsStar(crucnode, graph, &connectivityData, &soltreeData, &pcmwData, &keypathsData, &supergraphData, &allgood);
3323 
3324  /* todo: check this! */
3325  if( !allgood )
3326  {
3327  *success = FALSE;
3328  localmoves = 0;
3329  SCIPdebugMessage("terminate key vertex heuristic \n");
3330  goto TERMINATE;
3331  }
3332 
3333  assert(keypathsData.nkpnodes != 0); /* if there are no key-path nodes, something has gone wrong */
3334 
3335  /* reset all nodes (referred to as 'C' henceforth) whose bases are internal nodes of the current key-paths */
3336  vnoiDataReset(&connectivityData, &keypathsData, graph->mark, &vnoiData);
3337 
3338  SCIP_CALL( connectivityDataKeyElimUpdate(scip, graph, &vnoiData, &supergraphData, crucnode, &connectivityData) );
3339 
3340  /* try to connect the nodes of C (directly) to COMP(C), as a preprocessing for graph_voronoiRepairMult */
3341  vnoiDataRepairPreprocess(scip, graph, &keypathsData, &connectivityData, &pcmwData, &vnoiData, &nheapelems);
3342 
3343  graph_voronoiRepairMult(scip, graph, (pcmwData.isActive? pcmwData.edgecost_biased : graph->cost),
3344  supernodesmark, &nheapelems, vnoibase, connectivityData.boundaryedges, &(connectivityData.nboundaryedges), &uf, vnoipath);
3345 
3346  SCIP_CALL( supergraphComputeMst(scip, graph, &connectivityData, &vnoiData, &keypathsData,
3347  crucnode, &soltreeData, &pcmwData, &supergraphData) );
3348 
3349  assert(crucnode == dfstree[dfstree_pos]);
3350 
3351  SCIPdebugMessage("kpcost=%f mstcost=%f \n", keypathsData.kpcost, supergraphData.mstcost);
3352 
3353  /* improving solution found? */
3354  if( SCIPisLT(scip, supergraphData.mstcost, keypathsData.kpcost) )
3355  {
3356  localmoves++;
3357  solimproved = TRUE;
3358 
3359  SCIPdebugMessage("found improving solution in KEY VERTEX ELIMINATION (round: %d) \n \n ", nruns);
3360 
3361  SCIP_CALL( soltreeElimKeyPathsStar(scip, graph, &connectivityData, &vnoiData, &keypathsData, &supergraphData,
3362  dfstree, scanned, dfstree_pos, &soltreeData) );
3363 
3364 #ifndef NDEBUG
3365  assert((keypathsData.kpcost - supergraphData.mstcost) >= 0.0);
3366  objimprovement += (keypathsData.kpcost - supergraphData.mstcost);
3367 #endif
3368  }
3369  else /* no improving solution has been found during the move */
3370  {
3371  /* meld the heap pertaining to 'crucnode' and all heaps pertaining to descendant key-paths of node 'crucnode' */
3372  for( int k = 0; k < keypathsData.rootpathstart; k++ )
3373  {
3374  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[kpnodes[k]], &pheapsize[crucnode], &pheapsize[kpnodes[k]]);
3375  }
3376  for( int k = 0; k < supergraphData.nsupernodes - 1; k++ )
3377  {
3378  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[supernodes[k]], &pheapsize[crucnode], &pheapsize[supernodes[k]]);
3379  SCIPStpunionfindUnion(&uf, crucnode, supernodes[k], FALSE);
3380  }
3381  }
3382 
3383  supergraphDataRestore(graph, &supergraphData);
3384 
3385  /* restore the original Voronoi diagram */
3386  vnoiDataRestore(&connectivityData, &keypathsData, &vnoiData);
3387  } /* end of key vertex elimination */
3388 
3389  /* Key-Path Exchange:
3390  * If the crucnode has just been eliminated, skip Key-Path Exchange */
3391  if( graph->mark[crucnode] )
3392  {
3393  SCIP_Real edgecost = FARAWAY;
3394  int e = UNKNOWN;
3395  int oldedge = UNKNOWN;
3396  int newedge = UNKNOWN;
3397 
3398  assert(graph->mark[crucnode]);
3399 
3400  /* is crucnode not a crucial node and not a pinned vertex? */
3401  if( (!nodeIsCrucial(graph, solEdges, crucnode) && !pinned[crucnode]) )
3402  continue;
3403 
3404  /* gets key path from crucnode towards tree root */
3405  getKeyPathUpper(scip, crucnode, graph, &soltreeData, &pcmwData, &connectivityData, &keypathsData);
3406 
3407 #ifndef NDEBUG
3408  for( int k = 0; k < nnodes; k++ )
3409  assert(graph->path_state[k] == CONNECT || !graph->mark[k]);
3410 #endif
3411 
3412  /* reset all nodes (henceforth referred to as 'C') whose bases are internal nodes of the current keypath */
3413  vnoiDataReset(&connectivityData, &keypathsData, graph->mark, &vnoiData);
3414 
3415  while( boundpaths[crucnode] != NULL )
3416  {
3417  int base;
3418  int node;
3419 
3420  SCIP_CALL( SCIPpairheapDeletemin(scip, &e, &edgecost, &boundpaths[crucnode], &(pheapsize[crucnode])) );
3421 
3422  assert(e != UNKNOWN);
3423  base = vnoibase[graph->head[e]];
3424 
3425  assert(graph->mark[vnoibase[graph->tail[e]]]);
3426  node = (base == UNKNOWN || !graph->mark[base] )? UNKNOWN : SCIPStpunionfindFind(&uf, base);
3427 
3428  /* does the boundary-path end in the root component? */
3429  if( node != UNKNOWN && node != crucnode && graph->mark[base] )
3430  {
3431  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[crucnode],
3432  pheapelements,
3433  e, edgecost, &(pheapsize[crucnode])) );
3434  break;
3435  }
3436  }
3437 
3438  if( boundpaths[crucnode] == NULL )
3439  oldedge = UNKNOWN;
3440  else
3441  oldedge = e;
3442 
3443  /* try to connect the nodes of C (directly) to COMP(C), as a preprocessing for Voronoi-repair */
3444  vnoiDataRepairPreprocess(scip, graph, &keypathsData, &connectivityData, &pcmwData, &vnoiData, &nheapelems);
3445 
3446  newedge = UNKNOWN;
3447 
3448  /* if there is no key path, nothing has to be repaired */
3449  if( keypathsData.nkpnodes > 0 )
3450  {
3451  graph_voronoiRepair(scip, graph, pcmwData.isActive? pcmwData.edgecost_biased : graph->cost,
3452  pcmwData.edgecost_org, &nheapelems, vnoibase, vnoipath, &newedge, crucnode, &uf);
3453  }
3454  else
3455  {
3456  newedge = linkcutNodes[crucnode].edge;
3457  }
3458 
3459  edgecost = getKeyPathReplaceCost(scip, graph, &vnoiData, &soltreeData, edgecost, oldedge, &pcmwData, &newedge);
3460 
3461  if( SCIPisLT(scip, edgecost, keypathsData.kpcost) )
3462  {
3463  localmoves++;
3464  solimproved = TRUE;
3465 
3466  SCIPdebugMessage( "ADDING NEW KEY PATH (%f )\n\n", edgecost - keypathsData.kpcost );
3467 #ifndef NDEBUG
3468  assert((keypathsData.kpcost - edgecost) >= 0.0);
3469  objimprovement += (keypathsData.kpcost - edgecost);
3470  assert(crucnode == dfstree[dfstree_pos]);
3471 #endif
3472 
3473  SCIP_CALL( soltreeExchangeKeyPath(scip, graph, &connectivityData, &vnoiData, &keypathsData,
3474  dfstree, scanned, dfstree_pos, newedge, &soltreeData) );
3475  }
3476 
3477  /* restore the original Voronoi diagram */
3478  vnoiDataRestore(&connectivityData, &keypathsData, &vnoiData);
3479  }
3480  }
3481 
3482 
3483  /**********************************************************/
3484 
3485  TERMINATE:
3486 
3487  assert(uf.nElements == nnodes);
3488  SCIPStpunionfindClear(scip, &uf);
3489 
3490  /* free data structures */
3491 
3492  for( int k = nnodes - 1; k >= 0; k-- )
3493  {
3494  if( boundpaths[k] )
3495  SCIPpairheapFree(scip, &boundpaths[k]);
3496  }
3497 
3498  /* has there been a move during this run? */
3499  if( localmoves > 0 )
3500  {
3501  for( int i = 0; i < nnodes; i++ )
3502  {
3503  solNodes[i] = FALSE;
3504  graph->mark[i] = (graph->grad[i] > 0);
3505  SCIPlinkcuttreeInitNode(&linkcutNodes[i]);
3506  }
3507 
3508  graph->mark[root] = TRUE;
3509 
3510  /* create a link-cut tree representing the current Steiner tree */
3511  for( int e = 0; e < nedges; e++ )
3512  {
3513  assert(graph->head[e] == graph->tail[flipedge(e)]);
3514 
3515  /* if edge e is in the tree, so are its incident vertices */
3516  if( solEdges[e] != -1 )
3517  {
3518  assert(CONNECT == solEdges[e]);
3519 
3520  solNodes[graph->tail[e]] = TRUE;
3521  solNodes[graph->head[e]] = TRUE;
3522  SCIPlinkcuttreeLink(linkcutNodes, graph->head[e], graph->tail[e], flipedge(e));
3523  }
3524  }
3525  assert( linkcutNodes[root].edge == -1 );
3526  linkcutNodes[root].edge = -1;
3527  }
3528  } /* main loop */
3529 
3530  /* free data structures */
3531  SCIPStpunionfindFreeMembers(scip, &uf);
3532  SCIPfreeMemoryArray(scip, &pheapelements);
3533 
3534  for( int k = nnodes - 1; k >= 0; k-- )
3535  {
3536  StpVecFree(scip, lvledges_start[k]);
3537  StpVecFree(scip, blists_start[k]);
3538  }
3539 
3540  SCIPfreeBufferArray(scip, &kpedges);
3541  SCIPfreeBufferArray(scip, &kpnodes);
3542  SCIPfreeBufferArray(scip, &supernodes);
3543  SCIPfreeBufferArray(scip, &supernodesmark);
3544  SCIPfreeBufferArray(scip, &dfstree);
3545  SCIPfreeBufferArray(scip, &pinned);
3546  SCIPfreeBufferArray(scip, &boundpaths);
3547  SCIPfreeBufferArray(scip, &meminedges);
3548  SCIPfreeBufferArray(scip, &memvdist);
3549  SCIPfreeBufferArray(scip, &memvbase);
3550  SCIPfreeBufferArray(scip, &blists_start);
3551  SCIPfreeBufferArray(scip, &pheapsize);
3552  SCIPfreeBufferArray(scip, &scanned);
3553  SCIPfreeBufferArrayNull(scip, &prizemarklist);
3554  SCIPfreeBufferArrayNull(scip, &prizemark);
3555  SCIPfreeBufferArrayNull(scip, &prize_pc);
3556  SCIPfreeBufferArrayNull(scip, &edgecostBiased_pc);
3557  SCIPfreeBufferArrayNull(scip, &edgecostOrg_pc);
3558  SCIPfreeBufferArray(scip, &boundedges);
3559  SCIPfreeBufferArray(scip, &lvledges_start);
3560  SCIPfreeBufferArray(scip, &newedges);
3561  SCIPfreeBufferArray(scip, &vnoibase);
3562  SCIPfreeBufferArray(scip, &vnoipath);
3563  /******/
3564 
3565  if( solimproved )
3566  {
3567  SCIP_CALL( solstp_pruneFromEdges(scip, graph, solEdges) );
3568  *success = TRUE;
3569  }
3570 
3571 #ifndef NDEBUG
3572  {
3573  const SCIP_Real newobj = solstp_getObjBounded(graph, solEdges, 0.0, nedges);
3574  SCIPdebugMessage("key vertex heuristic obj before/after: %f/%f (improvement=%f)\n", initialobj, newobj, objimprovement);
3575  assert(SCIPisLE(scip, newobj + objimprovement, initialobj));
3576  }
3577 #endif
3578 
3579  return SCIP_OKAY;
3580 }
3581 
3582 
3583 /*
3584  * Callback methods of primal heuristic
3585  */
3586 
3587 /** copy method for primal heuristic plugins (called when SCIP copies plugins) */
3588 static
3589 SCIP_DECL_HEURCOPY(heurCopyLocal)
3590 { /*lint --e{715}*/
3591  assert(scip != NULL);
3592  assert(heur != NULL);
3593  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
3594 
3595  /* call inclusion method of primal heuristic */
3597 
3598  return SCIP_OKAY;
3599 }
3600 
3601 /** destructor of primal heuristic to free user data (called when SCIP is exiting) */
3602 static
3603 SCIP_DECL_HEURFREE(heurFreeLocal)
3604 { /*lint --e{715}*/
3605  SCIP_HEURDATA* heurdata;
3606 
3607  assert(heur != NULL);
3608  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
3609  assert(scip != NULL);
3610 
3611  /* free heuristic data */
3612  heurdata = SCIPheurGetData(heur);
3613 
3614  assert(heurdata != NULL);
3615 
3616  SCIPfreeRandom(scip, &heurdata->randnumgen);
3617  SCIPfreeMemory(scip, &heurdata);
3618  SCIPheurSetData(heur, NULL);
3619 
3620  return SCIP_OKAY;
3621 }
3622 
3623 /** solving process initialization method of primal heuristic (called when branch and bound process is about to begin) */
3624 static
3625 SCIP_DECL_HEURINITSOL(heurInitsolLocal)
3626 { /*lint --e{715}*/
3627  SCIP_HEURDATA* heurdata;
3628 
3629  assert(heur != NULL);
3630  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
3631  assert(scip != NULL);
3632 
3633  /* free heuristic data */
3634  heurdata = SCIPheurGetData(heur);
3635 
3636  heurdata->nfails = 1;
3637  heurdata->nbestsols = -1;
3638  heurdata->nbestsols_min = -1;
3639 
3640  SCIP_CALL( SCIPallocMemoryArray(scip, &(heurdata->lastsolindices), heurdata->maxnsols) );
3641 
3642  for( int i = 0; i < heurdata->maxnsols; i++ )
3643  heurdata->lastsolindices[i] = -1;
3644 
3645  return SCIP_OKAY;
3646 }
3647 
3648 
3649 /** solving process deinitialization method of primal heuristic (called before branch and bound process data is freed) */
3650 static
3651 SCIP_DECL_HEUREXITSOL(heurExitsolLocal)
3652 { /*lint --e{715}*/
3653  SCIP_HEURDATA* heurdata;
3654 
3655  assert(heur != NULL);
3656  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
3657  assert(scip != NULL);
3658 
3659  /* free heuristic data */
3660  heurdata = SCIPheurGetData(heur);
3661  assert(heurdata != NULL);
3662  assert(heurdata->lastsolindices != NULL);
3663  SCIPfreeMemoryArray(scip, &(heurdata->lastsolindices));
3664 
3665  return SCIP_OKAY;
3666 }
3667 
3668 
3669 /** execution method of primal heuristic */
3670 static
3671 SCIP_DECL_HEUREXEC(heurExecLocal)
3672 { /*lint --e{715}*/
3673  SCIP_HEURDATA* heurdata;
3674  SCIP_PROBDATA* probdata;
3675  GRAPH* graph; /* graph structure */
3676  SCIP_SOL* initialsol; /* initial solution */
3677  SCIP_SOL** sols; /* solutions */
3678  int v;
3679  int nActiveSols;
3680  int nsols; /* number of all solutions found so far */
3681  int nedges;
3682  int* results;
3683  int* lastsolindices;
3684  SCIP_Real initialsol_obj;
3685 
3686  assert(heur != NULL);
3687  assert(scip != NULL);
3688  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
3689  assert(result != NULL);
3690 
3691  /* get heuristic data */
3692  heurdata = SCIPheurGetData(heur);
3693  assert(heurdata != NULL);
3694  lastsolindices = heurdata->lastsolindices;
3695  assert(lastsolindices != NULL);
3696 
3697  probdata = SCIPgetProbData(scip);
3698  assert(probdata != NULL);
3699 
3700  graph = SCIPprobdataGetGraph(probdata);
3701  assert(graph != NULL);
3702 
3703  *result = SCIP_DIDNOTRUN;
3704 
3705  /* the local heuristics may not work correctly for several problem variants */
3706  if( !probtypeIsValidForLocal(graph) )
3707  return SCIP_OKAY;
3708 
3709  /* don't run local in a Subscip */
3710  if( SCIPgetSubscipDepth(scip) > 0 )
3711  return SCIP_OKAY;
3712 
3713  /* no solution available? */
3714  if( SCIPgetBestSol(scip) == NULL )
3715  return SCIP_OKAY;
3716 
3717  sols = SCIPgetSols(scip);
3718  nsols = SCIPgetNSols(scip);
3719  nedges = graph->edges;
3720 
3721  assert(heurdata->maxnsols >= 0);
3722 
3723  nActiveSols = MIN(heurdata->maxnsols, nsols);
3724 
3725  /* only process each solution once */
3726  for( v = 0; v < nActiveSols; v++ )
3727  {
3728  if( SCIPsolGetIndex(sols[v]) != lastsolindices[v] )
3729  {
3730  /* shift all solution indices right of the new solution index */
3731  for( int i = nActiveSols - 1; i >= v + 1; i-- )
3732  lastsolindices[i] = lastsolindices[i - 1];
3733  break;
3734  }
3735  }
3736 
3737  /* no new solution available? */
3738  if( v == nActiveSols )
3739  return SCIP_OKAY;
3740 
3741  if( heurdata->nbestsols == -1 )
3742  initSolNumberBounds(scip, heurdata);
3743 
3744  initialsol = sols[v];
3745  lastsolindices[v] = SCIPsolGetIndex(initialsol);
3746 
3747  /* solution not good enough? */
3748  if( (v > heurdata->nbestsols && !(heurdata->maxfreq)) )
3749  return SCIP_OKAY;
3750 
3751  /* has the new solution been found by this very heuristic
3752  * and not among the elite solutions? (note that there can still be an improvement because heuristic is aborted early) */
3753  if( SCIPsolGetHeur(initialsol) == heur && v > DEFAULT_NELITESOLS )
3754  return SCIP_OKAY;
3755 
3756  *result = SCIP_DIDNOTFIND;
3757 
3758  assert(SCIPprobdataGetVars(scip) != NULL);
3759  assert(SCIPprobdataGetNVars(scip) == nedges);
3760 
3761  /* allocate memory */
3762  SCIP_CALL( SCIPallocBufferArray(scip, &results, nedges) );
3763 
3764  /* write solution into 'results' array */
3765  solGetStpSol(scip, graph, initialsol, results);
3766 
3767  if( !solstp_isValid(scip, graph, results) )
3768  {
3769  SCIPfreeBufferArray(scip, &results);
3770 
3771  return SCIP_OKAY;
3772  }
3773 
3774  /* initial pruning necessary? */
3775  if( solNeedsPruning(initialsol) )
3776  {
3777  SCIP_CALL( solPrune(scip, graph, results) );
3778  }
3779 
3780  initialsol_obj = solstp_getObjBounded(graph, results, 0.0, nedges);
3781 
3782  /* execute local heuristics */
3783  SCIP_CALL( SCIPStpHeurLocalRun(scip, graph, results) );
3784 
3785 #if 0
3786  if( graph_pc_isPcMw(graph) )
3787  SCIP_CALL( SCIPStpHeurLocalExtendPcMwImp(scip, graph, results) );
3788 #endif
3789 
3790 
3791  /* finally, try to add the solution to the SCIP pool */
3792  SCIP_CALL( solAddTry(scip, sols, heur, graph, initialsol_obj, initialsol, results, result) );
3793 
3794  SCIPfreeBufferArray(scip, &results);
3795 
3796  return SCIP_OKAY;
3797 }
3798 
3799 
3800 /** perform local heuristics on a given Steiner tree */
3801 static
3803  SCIP* scip, /**< SCIP data structure */
3804  GRAPH* graph, /**< graph data structure */
3805  SCIP_Bool useFast, /**< fast variant? */
3806  int* solEdges /**< array indicating whether an arc is part of the solution: CONNECTED/UNKNOWN (in/out) */
3807  )
3808 {
3809  SCIP_HEUR* heur;
3810  SCIP_HEURDATA* heurdata;
3811  LCNODE* linkcutNodes;
3812  const int root = graph->source;
3813  const int nnodes = graph->knots;
3814  STP_Bool* solNodes;
3815  const STP_Bool mwpc = graph_pc_isPcMw(graph);
3816  SCIP_Bool success = FALSE;
3817 #ifndef NDEBUG
3818  const SCIP_Real initialobj = solstp_getObjBounded(graph, solEdges, 0.0, graph->edges);
3819 #endif
3820 
3821  assert(graph && solEdges);
3822  assert(graph_valid(scip, graph));
3823 
3824  if( graph->grad[root] == 0 || graph->terms == 1 )
3825  return SCIP_OKAY;
3826 
3827  if( mwpc )
3828  {
3829  assert(graph->extended);
3830 
3831  if( solIsTrivialPcMw(graph, solEdges) )
3832  return SCIP_OKAY;
3833  }
3834 
3835  heur = SCIPfindHeur(scip, "local");
3836  assert(heur != NULL);
3837  heurdata = SCIPheurGetData(heur);
3838  assert(heurdata != NULL);
3839 
3840  SCIP_CALL( SCIPallocBufferArray(scip, &linkcutNodes, nnodes) );
3841  SCIP_CALL( SCIPallocBufferArray(scip, &solNodes, nnodes) );
3842 
3843  /* now call the actual heuristics */
3844 
3845  if( mwpc )
3846  {
3847  SCIP_CALL( SCIPStpHeurLocalExtendPcMw(scip, graph, graph->cost, solEdges, solNodes) );
3848  }
3849 
3850  markSolTreeNodes(scip, graph, solEdges, linkcutNodes, solNodes);
3851 
3852  SCIP_CALL( localVertexInsertion(scip, heurdata, graph, solNodes, linkcutNodes, solEdges) );
3853 
3854  SCIP_CALL( localKeyVertexHeuristics(scip, graph, useFast, solNodes, linkcutNodes, solEdges, &success) );
3855 
3856  if( success && !useFast )
3857  {
3858  markSolTreeNodes(scip, graph, solEdges, linkcutNodes, solNodes);
3859  SCIP_CALL( localVertexInsertion(scip, heurdata, graph, solNodes, linkcutNodes, solEdges) );
3860  }
3861 
3862  if( success && mwpc && !useFast )
3863  {
3864  SCIP_CALL( SCIPStpHeurLocalExtendPcMw(scip, graph, graph->cost, solEdges, solNodes) );
3865  }
3866 
3867 #ifndef NDEBUG
3868  {
3869  const SCIP_Real newobj = solstp_getObjBounded(graph, solEdges, 0.0, graph->edges);
3870  assert(SCIPisLE(scip, newobj, initialobj));
3871  assert(solstp_isValid(scip, graph, solEdges));
3872  }
3873 #endif
3874 
3875  SCIPfreeBufferArray(scip, &solNodes);
3876  SCIPfreeBufferArray(scip, &linkcutNodes);
3877 
3878  return SCIP_OKAY;
3879 }
3880 
3881 
3882 /** perform local heuristics on a given Steiner tree */
3884  SCIP* scip, /**< SCIP data structure */
3885  GRAPH* graph, /**< graph data structure */
3886  int* solEdges /**< array indicating whether an arc is part of the solution: CONNECTED/UNKNOWN (in/out) */
3887  )
3888 {
3889  SCIP_CALL( localRun(scip, graph, FALSE, solEdges) );
3890  return SCIP_OKAY;
3891 }
3892 
3893 
3894 /** perform local heuristics on a given Steiner tree */
3896  SCIP* scip, /**< SCIP data structure */
3897  GRAPH* graph, /**< graph data structure */
3898  int* solEdges /**< array indicating whether an arc is part of the solution: CONNECTED/UNKNOWN (in/out) */
3899  )
3900 {
3901  SCIP_CALL( localRun(scip, graph, TRUE, solEdges) );
3902  return SCIP_OKAY;
3903 }
3904 
3905 
3906 /** Implication based local heuristic for (R)PC and MW */
3908  SCIP* scip, /**< SCIP data structure */
3909  const GRAPH* graph, /**< graph data structure */
3910  int* result /**< array indicating whether an arc is part of the solution (CONNECTED/UNKNOWN) */
3911 )
3912 {
3913  const int* starts = SCIPStpGetPcImplStarts(scip);
3914  const int* verts = SCIPStpGetPcImplVerts(scip);
3915 
3916  assert(graph_pc_isPcMw(graph));
3917  assert(0 && "currently not supported");
3918 
3919  if( starts != NULL )
3920  {
3921  const int nnodes = graph->knots;
3922  STP_Bool* stvertex;
3923  int nfound = 0;
3924  int ptermcount = 0;
3925 
3926  assert(graph->extended);
3927  assert(verts != NULL);
3928 
3929  SCIPallocBufferArray(scip, &stvertex, nnodes);
3930 
3931  solstp_setVertexFromEdge(graph, result, stvertex);
3932 
3933  for( int i = 0; i < nnodes; i++ )
3934  {
3935  if( !(Is_term(graph->term[i]) || Is_pseudoTerm(graph->term[i])) )
3936  continue;
3937 
3938  assert(!graph_pc_knotIsFixedTerm(graph, i));
3939 
3940  ptermcount++;
3941 
3942  if( stvertex[i] )
3943  continue;
3944 
3945  for( int j = starts[ptermcount - 1]; j < starts[ptermcount]; j++ )
3946  {
3947  const int vert = verts[j];
3948  if( stvertex[vert] )
3949  {
3950  /* now connect the vertex */
3951 
3952  graph_knot_printInfo(graph, i);
3953  nfound++;
3954  break;
3955  }
3956  }
3957  }
3958 
3959  assert(ptermcount == graph_pc_nNonFixedTerms(graph));
3960 
3961  if( nfound > 0 )
3962  {
3963  printf("nfound: %d \n\n\n", nfound);
3964  }
3965  else
3966  printf("none %d \n", 0);
3967 
3968  SCIPfreeBufferArray(scip, &stvertex);
3969  }
3970 
3971  return SCIP_OKAY;
3972 }
3973 
3974 /** Greedy Extension local heuristic for (R)PC and MW */
3976  SCIP* scip, /**< SCIP data structure */
3977  GRAPH* graph, /**< graph data structure */
3978  const SCIP_Real* cost, /**< edge cost array */
3979  int* stedge, /**< initialized array to indicate whether an edge is part of the Steiner tree */
3980  STP_Bool* stvertex /**< uninitialized array to indicate whether a vertex is part of the Steiner tree */
3981  )
3982 {
3983  GNODE candidates[GREEDY_EXTENSIONS_MW];
3984  int candidatesup[GREEDY_EXTENSIONS_MW];
3985 
3986  PATH* path;
3987  PATH* orgpath;
3988  SCIP_PQUEUE* pqueue;
3989  SCIP_Real bestsolval;
3990 
3991  int nextensions;
3992  const int greedyextensions = (graph->stp_type == STP_MWCSP) ? GREEDY_EXTENSIONS_MW : GREEDY_EXTENSIONS;
3993  const int nedges = graph->edges;
3994  const int nnodes = graph->knots;
3995  const int root = graph->source;
3996  STP_Bool* stvertextmp;
3997  SCIP_Bool extensions = FALSE;
3998 
3999 #ifndef NDEBUG
4000  SCIP_Real objdiff = 0.0;
4001  const SCIP_Real initialobj = solstp_getObjBounded(graph, stedge, 0.0, nedges);
4002 #endif
4003 
4004  assert(GREEDY_EXTENSIONS_MW >= greedyextensions);
4005  assert(scip != NULL);
4006  assert(graph != NULL);
4007  assert(stedge != NULL);
4008  assert(cost != NULL);
4009  assert(stvertex != NULL);
4010  assert(graph->extended);
4011 
4012  graph_pc_2transcheck(scip, graph);
4013  SCIP_CALL( SCIPallocBufferArray(scip, &stvertextmp, nnodes) );
4014  SCIP_CALL( SCIPallocBufferArray(scip, &orgpath, nnodes) );
4015  SCIP_CALL( SCIPallocBufferArray(scip, &path, nnodes) );
4016 
4017  /* initialize solution vertex array with FALSE */
4018  BMSclearMemoryArray(stvertex, nnodes);
4019 
4020  stvertex[root] = TRUE;
4021 
4022  for( int j = 0; j < nnodes; j++ )
4023  path[j].edge = UNKNOWN;
4024 
4025  for( int e = 0; e < nedges; e++ )
4026  {
4027  if( stedge[e] == CONNECT )
4028  {
4029  path[graph->head[e]].edge = e;
4030  stvertex[graph->head[e]] = TRUE;
4031  }
4032  }
4033 
4034 #ifndef NDEBUG
4035  for( int e = 0; e < nedges; e++ )
4036  if( stedge[e] == CONNECT )
4037  assert(stvertex[graph->tail[e]]);
4038 #endif
4039 
4040  graph_path_st_pcmw_extend(scip, graph, cost, FALSE, path, stvertex, &extensions);
4041 
4042  BMScopyMemoryArray(orgpath, path, nnodes);
4043 
4044  /*** compute solution value and save greedyextensions many best unconnected nodes ***/
4045 
4046  SCIP_CALL( SCIPpqueueCreate(&pqueue, greedyextensions, 2.0, GNODECmpByDist, NULL) );
4047 
4048  assert(orgpath[root].edge == UNKNOWN);
4049 
4050  bestsolval = 0.0;
4051  nextensions = 0;
4052 
4053  for( int i = 0; i < nnodes; i++ )
4054  {
4055  if( graph->grad[i] == 0 || root == i )
4056  continue;
4057 
4058  if( Is_term(graph->term[i]) && !graph_pc_knotIsFixedTerm(graph, i) )
4059  continue;
4060 
4061  if( stvertex[i] )
4062  {
4063  assert(orgpath[i].edge >= 0);
4064 
4065  bestsolval += graph->cost[orgpath[i].edge];
4066 
4067  if( Is_pseudoTerm(graph->term[i]) )
4068  {
4069  bestsolval -= graph->prize[i];
4070  }
4071  }
4072  else if( orgpath[i].edge != UNKNOWN && Is_pseudoTerm(graph->term[i]) )
4073  {
4074  SCIP_CALL( addToCandidates(scip, graph, path, i, greedyextensions, &nextensions, candidates, pqueue) );
4075  }
4076  }
4077 
4078  for( int restartcount = 0; restartcount < GREEDY_MAXRESTARTS && !graph_pc_isRootedPcMw(graph); restartcount++ )
4079  {
4080  int l = 0;
4081  SCIP_Bool extensionstmp = FALSE;
4082  int extcount = nextensions;
4083 
4084  /* write extension candidates into array, from max to min */
4085  while( SCIPpqueueNElems(pqueue) > 0 )
4086  {
4087  GNODE* min = (GNODE*) SCIPpqueueRemove(pqueue);
4088  assert(extcount > 0);
4089  candidatesup[--extcount] = min->number;
4090  }
4091  assert(extcount == 0);
4092 
4093  /* iteratively insert new subpaths and try to improve solution */
4094  for( ; l < nextensions; l++ )
4095  {
4096  const int extensioncand = candidatesup[l];
4097  if( !stvertex[extensioncand] )
4098  {
4099  SCIP_Real newsolval = 0.0;
4100  int k = extensioncand;
4101 
4102  BMScopyMemoryArray(stvertextmp, stvertex, nnodes);
4103  BMScopyMemoryArray(path, orgpath, nnodes);
4104 
4105  /* add new extension */
4106  while( !stvertextmp[k] )
4107  {
4108  stvertextmp[k] = TRUE;
4109  assert(orgpath[k].edge != UNKNOWN);
4110  k = graph->tail[orgpath[k].edge];
4111  assert(k != extensioncand);
4112  }
4113 
4114  graph_path_st_pcmw_extend(scip, graph, cost, TRUE, path, stvertextmp, &extensionstmp);
4115 
4116 
4117  for( int j = 0; j < nnodes; j++ )
4118  {
4119  if( graph->grad[j] == 0 || root == j )
4120  continue;
4121 
4122  if( Is_term(graph->term[j]) && !graph_pc_knotIsFixedTerm(graph, j) )
4123  continue;
4124 
4125  if( stvertextmp[j] )
4126  {
4127  assert(path[j].edge >= 0);
4128 
4129  newsolval += graph->cost[path[j].edge];
4130 
4131  if( Is_pseudoTerm(graph->term[j]) )
4132  {
4133  newsolval -= graph->prize[j];
4134  }
4135  }
4136  }
4137 
4138  /* new solution value better than old one? */
4139  if( SCIPisLT(scip, newsolval, bestsolval) )
4140  {
4141  SCIPdebugMessage("newsolval=%f bestsolval=%f \n", newsolval, bestsolval);
4142 
4143  extensions = TRUE;
4144  bestsolval = newsolval;
4145 
4146 #ifdef SCIP_DEBUG
4147  for( int k2 = 0; k2 < nnodes; k2++ )
4148  {
4149  if( !stvertex[k2] && stvertextmp[k2] )
4150  {
4151  printf("added \n");
4152  graph_knot_printInfo(graph, k2);
4153  graph_edge_printInfo(graph, path[k2].edge);
4154  }
4155  }
4156 #endif
4157 
4158 #ifndef NDEBUG
4159  objdiff += bestsolval - newsolval;
4160 #endif
4161 
4162  BMScopyMemoryArray(stvertex, stvertextmp, nnodes);
4163  BMScopyMemoryArray(orgpath, path, nnodes);
4164 
4165  /* save greedyextensions many best unconnected nodes */
4166  nextensions = 0;
4167 
4168  for( int j = 0; j < nnodes; j++ )
4169  if( !stvertex[j] && Is_pseudoTerm(graph->term[j]) && path[j].edge != UNKNOWN )
4170  SCIP_CALL( addToCandidates(scip, graph, path, j, greedyextensions, &nextensions, candidates, pqueue) );
4171 
4172  break;
4173  } /* if new solution value better than old one? */
4174  } /* if !stvertex[i] */
4175  } /* for l < nextension */
4176 
4177  /* no more extensions performed? */
4178  if( l == nextensions )
4179  break;
4180  } /* main loop */
4181 
4182  /* have vertices been added? */
4183  if( extensions )
4184  {
4185  assert(graph_pc_isPcMw(graph));
4186 
4187  SCIP_CALL( solstp_prune(scip, graph, stedge, stvertex) );
4188  }
4189 
4190  assert(solstp_isValid(scip, graph, stedge));
4191 
4192  SCIPpqueueFree(&pqueue);
4193  SCIPfreeBufferArray(scip, &path);
4194  SCIPfreeBufferArray(scip, &orgpath);
4195  SCIPfreeBufferArray(scip, &stvertextmp);
4196 
4197 #ifndef NDEBUG
4198  {
4199  const SCIP_Real newsolval = solstp_getObjBounded(graph, stedge, 0.0, nedges);
4200  SCIPdebugMessage("pcmw extend obj: initial=%f, new=%f \n", initialobj, newsolval);
4201  assert(SCIPisLE(scip, newsolval + objdiff, initialobj));
4202  }
4203 #endif
4204 
4205  return SCIP_OKAY;
4206 }
4207 
4208 /** Greedy Extension local heuristic for (R)PC and MW */
4210  SCIP* scip, /**< SCIP data structure */
4211  GRAPH* graph, /**< graph data structure */
4212  int* stedge, /**< initialized array to indicate whether an edge is part of the Steiner tree */
4213  STP_Bool* stvertex /**< uninitialized array to indicate whether a vertex is part of the Steiner tree */
4214  )
4215 {
4216  int candidates[GREEDY_EXTENSIONS];
4217  int ncandidates;
4218  DHEAP* dheap;
4219  STP_Bool* stvertextmp;
4220  SCIP_Real* dist;
4221  int* pred;
4222  const int nnodes = graph->knots;
4223  SCIP_Bool extensions = FALSE;
4224  int maxnode;
4225  const SCIP_Bool isexended = graph->extended;
4226 
4227 #ifndef NDEBUG
4228  const int nedges = graph->edges;
4229  const SCIP_Real initialobj = solstp_getObjBounded(graph, stedge, 0.0, nedges);
4230 #endif
4231 
4232  assert(scip && graph && stedge && stvertex);
4233 
4234  graph_pc_2orgcheck(scip, graph);
4235 
4236  solstp_setVertexFromEdge(graph, stedge, stvertex);
4237 
4238  /* compute candidates for extension */
4239 
4240  maxnode = -1;
4241  ncandidates = 0;
4242 
4243  for( int k = 0; k < nnodes; k++ )
4244  if( graph->mark[k] && !stvertex[k] && Is_term(graph->term[k]) && !graph_pc_termIsNonLeafTerm(graph, k) )
4245  {
4246  assert(graph->mark[k]);
4247 
4248  if( maxnode == -1 || graph->prize[k] > graph->prize[maxnode] )
4249  maxnode = k;
4250  }
4251 
4252  if( maxnode != -1 )
4253  {
4254  SCIP_RANDNUMGEN* randnumgen;
4255  int shift;
4256 
4257  SCIP_CALL( SCIPcreateRandom(scip, &randnumgen, 1, TRUE) );
4258 
4259  SCIP_CALL( SCIPallocBufferArray(scip, &dist, nnodes) );
4260  SCIP_CALL( SCIPallocBufferArray(scip, &pred, nnodes) );
4261  SCIP_CALL( SCIPallocBufferArray(scip, &stvertextmp, nnodes) );
4262 
4263  graph_heap_create(scip, nnodes, NULL, NULL, &dheap);
4264  graph_init_csr(scip, graph);
4265 
4266  shift = SCIPrandomGetInt(randnumgen, 0, nnodes - 1);
4267  ncandidates = 1;
4268  candidates[0] = maxnode;
4269 
4270  for( int k = 0; k < nnodes && ncandidates < GREEDY_EXTENSIONS; k++ )
4271  {
4272  const int node = (k + shift) % nnodes;
4273  if( graph->mark[k] && !stvertex[node] && Is_term(graph->term[node])
4274  && !graph_pc_termIsNonLeafTerm(graph, node) && node != maxnode )
4275  {
4276  assert(graph->mark[node]);
4277  candidates[ncandidates++] = node;
4278  }
4279  }
4280 
4281  SCIPfreeRandom(scip, &randnumgen);
4282  }
4283 
4284  /* main loop */
4285  for( int k = 0; k < ncandidates; k++ )
4286  {
4287  const int cand = candidates[k];
4288  SCIP_Bool success = FALSE;
4289 
4290  if( stvertex[cand] )
4291  {
4292  assert(k > 0);
4293  continue;
4294  }
4295 
4296  graph_path_st_pcmw_extendOut(scip, graph, cand, stvertex, dist, pred, stvertextmp, dheap, &success);
4297 
4298  if( success )
4299  extensions = TRUE;
4300  }
4301 
4302  /* have vertices been added? */
4303  if( extensions )
4304  {
4305  graph_pc_2trans(scip, graph);
4306 
4307  assert(graph_pc_isPcMw(graph));
4308 
4309  SCIP_CALL( solstp_prune(scip, graph, stedge, stvertex) );
4310  }
4311 
4312  if( maxnode != -1 )
4313  {
4314  graph_heap_free(scip, TRUE, TRUE, &dheap);
4315  graph_free_csr(scip, graph);
4316 
4317  SCIPfreeBufferArray(scip, &stvertextmp);
4318  SCIPfreeBufferArray(scip, &pred);
4319  SCIPfreeBufferArray(scip, &dist);
4320  }
4321 
4322 #ifndef NDEBUG
4323  assert(SCIPisLE(scip, solstp_getObjBounded(graph, stedge, 0.0, nedges), initialobj));
4324 #endif
4325 
4326  if( isexended && !graph->extended )
4327  graph_pc_2trans(scip, graph);
4328 
4329  if( !isexended && graph->extended )
4330  graph_pc_2org(scip, graph);
4331 
4332  return SCIP_OKAY;
4333 }
4334 
4335 
4336 /** creates the local primal heuristic and includes it in SCIP */
4338  SCIP* scip /**< SCIP data structure */
4339  )
4340 {
4341  SCIP_HEURDATA* heurdata;
4342  SCIP_HEUR* heur;
4343 
4344  /* create Local primal heuristic data */
4345  SCIP_CALL( SCIPallocMemory(scip, &heurdata) );
4346 
4347  /* include primal heuristic */
4348  SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
4350  HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecLocal, heurdata) );
4351 
4352  assert(heur != NULL);
4353 
4354  /* set non-NULL pointers to callback methods */
4355  SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyLocal) );
4356  SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeLocal) );
4357  SCIP_CALL( SCIPsetHeurInitsol(scip, heur, heurInitsolLocal) );
4358  SCIP_CALL( SCIPsetHeurExitsol(scip, heur, heurExitsolLocal) );
4359 
4360  /* add local primal heuristic parameters */
4361  SCIP_CALL( SCIPaddBoolParam(scip, "stp/duringroot",
4362  "should the heuristic be called during the root node?",
4363  &heurdata->duringroot, FALSE, DEFAULT_DURINGROOT, NULL, NULL) );
4364 
4365  SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/"HEUR_NAME"/maxfreq",
4366  "should the heuristic be executed at maximum frequeny?",
4367  &heurdata->maxfreq, FALSE, DEFAULT_MAXFREQLOC, NULL, NULL) );
4368 
4369  SCIP_CALL( SCIPaddIntParam(scip, "heuristics/"HEUR_NAME"/maxnsols",
4370  "maximum number of best solutions to improve",
4371  &heurdata->maxnsols, FALSE, DEFAULT_MAXNBESTSOLS, DEFAULT_NBESTSOLS_HARD, 100, NULL, NULL) );
4372 
4373  SCIP_CALL( SCIPcreateRandom(scip, &heurdata->randnumgen, DEFAULT_RANDSEED, TRUE) );
4374 
4375  return SCIP_OKAY;
4376 }
enum SCIP_Result SCIP_RESULT
Definition: type_result.h:52
SCIP_Real SCIPlinkcuttreeFindMinChainMw(SCIP *scip, const LCNODE *tree, const SCIP_Real *nodeweight, const int *heads, const int *stdeg, const SCIP_Bool *nodeIsBlocked, int start, int *first, int *last)
Definition: misc_stp.c:381
SCIP_RETCODE SCIPsetHeurExitsol(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEUREXITSOL((*heurexitsol)))
Definition: scip_heur.c:233
#define HEUR_USESSUBSCIP
Definition: heur_local.c:52
#define STP_Vectype(type)
Definition: stpvector.h:44
void SCIPfreeRandom(SCIP *scip, SCIP_RANDNUMGEN **randnumgen)
SCIP_Real SCIPlinkcuttreeFindMaxChain(SCIP *scip, const LCNODE *tree, const SCIP_Real *edgecosts, const SCIP_Real *prizes, const int *heads, const int *nonTermDeg, const SCIP_Bool *nodeIsBlocked, int start, int *first, int *last)
Definition: misc_stp.c:445
SCIP_RETCODE graph_init_csr(SCIP *, GRAPH *)
Definition: graph_util.c:1581
static SCIP_RETCODE solPrune(SCIP *scip, const GRAPH *graph, int *results)
Definition: heur_local.c:220
int graph_get_nEdges(const GRAPH *)
Definition: graph_stats.c:548
int SCIPpqueueNElems(SCIP_PQUEUE *pqueue)
Definition: misc.c:1468
int nprizemarks
Definition: heur_local.c:162
#define HEUR_DESC
Definition: heur_local.c:44
static SCIP_Bool probtypeIsValidForLocal(const GRAPH *graph)
Definition: heur_local.c:401
static SCIP_RETCODE solAddTry(SCIP *scip, SCIP_SOL **sols, SCIP_HEUR *heur, const GRAPH *graph, SCIP_Real initialsol_obj, SCIP_SOL *initialsol, int *results, SCIP_RESULT *result)
Definition: heur_local.c:312
#define StpVecGetSize(vec)
Definition: stpvector.h:139
SCIP_RETCODE SCIPStpunionfindInit(SCIP *scip, UF *uf, int length)
Definition: misc_stp.c:817
struct pcmw_data PCMW
void graph_pc_getBiased(const GRAPH *, SCIP_Real *RESTRICT, SCIP_Real *RESTRICT)
#define StpVecClear(vec)
Definition: stpvector.h:134
struct connectivity_data CONN
int *RESTRICT head
Definition: graphdefs.h:224
int source
Definition: graphdefs.h:195
SCIP_Bool graph_pc_isRootedPcMw(const GRAPH *)
static void vnoiDataReset(const CONN *connectData, const KPATHS *keypathsData, const int *graphmark, VNOILOC *vnoiData)
Definition: heur_local.c:2293
#define Is_term(a)
Definition: graphdefs.h:102
static void vnoiDataRestore(const CONN *connectData, const KPATHS *keypathsData, VNOILOC *vnoiData)
Definition: heur_local.c:2255
static SCIP_Bool solDegIsValid(SCIP *scip, const GRAPH *graph, const int *solDegree, const LCNODE *linkcutNodes)
Definition: heur_local.c:2345
SCIP_Real mstcost
Definition: heur_local.c:149
Constraint handler for Steiner problems.
#define SCIPfreeMemoryArray(scip, ptr)
Definition: scip_mem.h:71
signed int edge
Definition: graphdefs.h:287
SCIP_Bool graph_valid(SCIP *, const GRAPH *)
Definition: graph_base.c:1480
void SCIPStpunionfindClear(SCIP *scip, UF *uf)
Definition: misc_stp.c:841
int terms
Definition: graphdefs.h:192
SCIP_VAR ** SCIPprobdataGetVars(SCIP *scip)
SCIP_Bool SCIPisPositive(SCIP *scip, SCIP_Real val)
#define SCIPallocMemoryArray(scip, ptr, num)
Definition: scip_mem.h:55
SCIP_Bool SCIPisGE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
int *const solDegreeNonTerm
Definition: heur_local.c:174
void graph_path_exec(SCIP *, const GRAPH *, int, int, const SCIP_Real *, PATH *)
Definition: graph_path.c:540
void SCIPpairheapMeldheaps(SCIP *scip, PHNODE **root1, PHNODE **root2, int *sizeroot1, int *sizeroot2)
Definition: misc_stp.c:708
#define GREEDY_EXTENSIONS
Definition: heur_local.c:69
includes methods for Steiner tree problem solutions
SCIP_RETCODE SCIPStpHeurLocalRunFast(SCIP *scip, GRAPH *graph, int *solEdges)
Definition: heur_local.c:3895
SCIP_RETCODE SCIPStpIncludeHeurLocal(SCIP *scip)
Definition: heur_local.c:4337
#define HEUR_PRIORITY
Definition: heur_local.c:46
void graph_knot_printInfo(const GRAPH *, int)
Definition: graph_stats.c:151
static SCIP_DECL_HEURCOPY(heurCopyLocal)
Definition: heur_local.c:3589
void graph_free(SCIP *, GRAPH **, SCIP_Bool)
Definition: graph_base.c:796
void graph_pc_2org(SCIP *, GRAPH *)
SCIP_RETCODE graph_path_init(SCIP *, GRAPH *)
Definition: graph_path.c:480
#define EAT_FREE
Definition: graphdefs.h:57
SCIP_SOL ** SCIPgetSols(SCIP *scip)
Definition: scip_sol.c:2254
static SCIP_Real getEdgeCostUnbiased(const GRAPH *graph, const PCMW *pcmwData, int edge)
Definition: heur_local.c:897
#define FALSE
Definition: def.h:87
int *RESTRICT inpbeg
Definition: graphdefs.h:202
static void supergraphDataRestore(const GRAPH *graph, SGRAPH *supergraphData)
Definition: heur_local.c:1782
int SCIPgetSubscipDepth(SCIP *scip)
Definition: scip_copy.c:2591
static void vnoiDataRepairPreprocess(SCIP *scip, const GRAPH *graph, const KPATHS *keypathsData, const CONN *connectData, const PCMW *pcmwData, VNOILOC *vnoiData, int *nheapelems)
Definition: heur_local.c:2188
int *RESTRICT path_state
Definition: graphdefs.h:256
SCIP_Bool isActive
Definition: heur_local.c:163
void graph_pathHeapAdd(const PATH *, int, int *RESTRICT, int *RESTRICT, int *)
static SCIP_Real getBoundaryPathCost(const GRAPH *graph, const VNOILOC *vnoiData, const PCMW *pcmwData, int boundaryedge)
Definition: heur_local.c:937
SCIP_Bool SCIPisNegative(SCIP *scip, SCIP_Real val)
Problem data for stp problem.
#define TRUE
Definition: def.h:86
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:54
SCIP_RETCODE SCIPStpHeurLocalExtendPcMwOut(SCIP *scip, GRAPH *graph, int *stedge, STP_Bool *stvertex)
Definition: heur_local.c:4209
static void markSolTreeNodes(SCIP *scip, const GRAPH *graph, const int *solEdges, LCNODE *linkcutNodes, STP_Bool *solNodes)
Definition: heur_local.c:455
static void initSolNumberBounds(SCIP *scip, SCIP_HEURDATA *heurdata)
Definition: heur_local.c:290
#define HEUR_MAXDEPTH
Definition: heur_local.c:49
void graph_voronoiRepairMult(SCIP *, const GRAPH *, const SCIP_Real *, const STP_Bool *, int *RESTRICT, int *RESTRICT, int *RESTRICT, int *RESTRICT, UF *RESTRICT, PATH *RESTRICT)
void * SCIPpqueueFirst(SCIP_PQUEUE *pqueue)
Definition: misc.c:1454
static void insertionDecrementSolDegree(const GRAPH *graph, int node, INSERT *insertData)
Definition: heur_local.c:2451
int SCIPprobdataGetNVars(SCIP *scip)
static void insertionGetCandidateEdges(SCIP *scip, SCIP_HEURDATA *heurdata, const GRAPH *graph, const STP_Bool *solNodes, int vertex, int *insertcands, int *ncands)
Definition: heur_local.c:2832
int *const cutedgesEnd
Definition: heur_local.c:177
struct SCIP_HeurData SCIP_HEURDATA
Definition: type_heur.h:67
int *const solEdges
Definition: heur_local.c:135
SCIP_RETCODE SCIPpairheapDeletemin(SCIP *scip, int *element, SCIP_Real *key, PHNODE **root, int *size)
Definition: misc_stp.c:671
static void insertionResetBlockedNodes(INSERT *insertData)
Definition: heur_local.c:2643
SCIP_RETCODE SCIPincludeHeurBasic(SCIP *scip, SCIP_HEUR **heur, const char *name, const char *desc, char dispchar, int priority, int freq, int freqofs, int maxdepth, SCIP_HEURTIMING timingmask, SCIP_Bool usessubscip, SCIP_DECL_HEUREXEC((*heurexec)), SCIP_HEURDATA *heurdata)
Definition: scip_heur.c:108
#define EAT_LAST
Definition: graphdefs.h:58
#define STP_SPG
Definition: graphdefs.h:42
void SCIPpqueueFree(SCIP_PQUEUE **pqueue)
Definition: misc.c:1263
#define SCIPdebugMessage
Definition: pub_message.h:87
int SCIPrandomGetInt(SCIP_RANDNUMGEN *randnumgen, int minrandval, int maxrandval)
Definition: misc.c:10003
int *const supernodes
Definition: heur_local.c:145
struct supergraph_data SGRAPH
#define FARAWAY
Definition: graphdefs.h:89
int *const addedEdges
Definition: heur_local.c:175
static SCIP_Bool pcmwDataIsClean(const GRAPH *graph, const PCMW *pcmwData)
Definition: heur_local.c:750
void SCIPpairheapFree(SCIP *scip, PHNODE **root)
Definition: misc_stp.c:736
SCIP_Bool SCIPisEQ(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
static void insertionFinalizeReplacement(const GRAPH *graph, int v_insert, LCNODE *linkcutNodes, INSERT *insertData)
Definition: heur_local.c:2592
#define SCIPfreeBufferArray(scip, ptr)
Definition: scip_mem.h:127
int number
Definition: misc_stp.h:55
void SCIPheurSetData(SCIP_HEUR *heur, SCIP_HEURDATA *heurdata)
Definition: heur.c:1362
void graph_voronoiRepair(SCIP *, const GRAPH *, const SCIP_Real *, const SCIP_Real *, int *, int *, PATH *, int *, int, UF *)
Definition: graph_vnoi.c:1100
SCIP_RETCODE SCIPaddIntParam(SCIP *scip, const char *name, const char *desc, int *valueptr, SCIP_Bool isadvanced, int defaultvalue, int minvalue, int maxvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip_param.c:74
header only, simple implementation of an STL like vector
GRAPH * SCIPprobdataGetGraph(SCIP_PROBDATA *probdata)
int *RESTRICT mark
Definition: graphdefs.h:198
SCIP_Bool graph_pc_termIsNonLeafTerm(const GRAPH *, int)
static const SCIP_Real * getEdgeCosts(const GRAPH *graph, const PCMW *pcmwData)
Definition: heur_local.c:1033
void SCIPlinkcuttreeEvert(LCNODE *RESTRICT tree, int root_new)
Definition: misc_stp.c:559
static SCIP_RETCODE localKeyVertexHeuristics(SCIP *scip, GRAPH *graph, SCIP_Bool useFast, STP_Bool *solNodes, LCNODE *linkcutNodes, int *solEdges, SCIP_Bool *success)
Definition: heur_local.c:3122
#define HEUR_FREQ
Definition: heur_local.c:47
static SCIP_RETCODE pcmwUpdate(SCIP *scip, GRAPH *graph, SOLTREE *soltreeData, PCMW *pcmwData)
Definition: heur_local.c:804
#define DEFAULT_DURINGROOT
Definition: heur_local.c:54
#define STP_RSMT
Definition: graphdefs.h:49
int *RESTRICT oeat
Definition: graphdefs.h:231
STP_Bool *const solNodes
Definition: heur_local.c:178
int *const newedges
Definition: heur_local.c:138
void SCIPlinkcuttreeInitNode(LCNODE *v)
Definition: misc_stp.c:347
SCIP_RETCODE SCIPsetHeurInitsol(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEURINITSOL((*heurinitsol)))
Definition: scip_heur.c:217
void graph_pc_2trans(SCIP *, GRAPH *)
const char * SCIPheurGetName(SCIP_HEUR *heur)
Definition: heur.c:1441
SCIP_Bool extended
Definition: graphdefs.h:267
SCIP_HEUR * SCIPfindHeur(SCIP *scip, const char *name)
Definition: scip_heur.c:249
int stp_type
Definition: graphdefs.h:264
SCIP_Bool *const nodeIsBlocked
Definition: heur_local.c:179
static SCIP_RETCODE getLowestCommonAncestors(SCIP *scip, const GRAPH *graph, const VNOILOC *vnoiData, const SOLTREE *soltreeData, CONN *connectData)
Definition: heur_local.c:600
SCIP_RETCODE SCIPsetHeurFree(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEURFREE((*heurfree)))
Definition: scip_heur.c:169
SCIP_Bool SCIPisLT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
static SCIP_DECL_HEURINITSOL(heurInitsolLocal)
Definition: heur_local.c:3625
void graph_pc_2transcheck(SCIP *, GRAPH *)
static void insertionRestoreTree(const GRAPH *graph, const int *insertCands, int lcvertex_insert, LCNODE *linkcutNodes, INSERT *insertData)
Definition: heur_local.c:2702
static STP_Bool nodeIsCrucial(const GRAPH *graph, const int *steineredges, int node)
Definition: heur_local.c:863
STP_Bool *const nodeIsScanned
Definition: heur_local.c:137
SCIP_Real SCIPprobdataGetOffset(SCIP *scip)
SCIP_Bool SCIPprobdataProbIsAdversarial(SCIP *scip)
static void getKeyPathsStar(int keyvertex, const GRAPH *graph, const CONN *connectData, const SOLTREE *soltreeData, const PCMW *pcmwData, KPATHS *keypathsData, SGRAPH *supergraphData, SCIP_Bool *success)
Definition: heur_local.c:1983
static SCIP_Real getKeyPathReplaceCost(SCIP *scip, const GRAPH *graph, const VNOILOC *vnoiData, const SOLTREE *soltreeData, SCIP_Real edgecost_old_in, int boundedge_old, PCMW *pcmwData, int *boundedge_new)
Definition: heur_local.c:1723
#define SCIPfreeBufferArrayNull(scip, ptr)
Definition: scip_mem.h:128
SCIP_Real dist
Definition: misc_stp.h:56
SCIP_Real * prize
Definition: graphdefs.h:210
SCIP_RETCODE SCIPprobdataAddNewSol(SCIP *scip, SCIP_Real *nval, SCIP_HEUR *heur, SCIP_Bool *success)
#define GREEDY_EXTENSIONS_MW
Definition: heur_local.c:68
void graph_edge_add(SCIP *, GRAPH *, int, int, double, double)
SCIP_Real dist
Definition: graphdefs.h:286
int *RESTRICT grad
Definition: graphdefs.h:201
SCIP_RETCODE SCIPpairheapInsert(SCIP *scip, PHNODE **root, PHNODE *pheapelems, int element, SCIP_Real key, int *size)
Definition: misc_stp.c:636
struct insertion_data INSERT
void graph_knot_add(GRAPH *, int)
Definition: graph_node.c:64
#define HEUR_TIMING
Definition: heur_local.c:50
static SCIP_RETCODE localVertexInsertion(SCIP *scip, SCIP_HEURDATA *heurdata, const GRAPH *graph, STP_Bool *solNodes, LCNODE *linkcutNodes, int *solEdges)
Definition: heur_local.c:2881
void graph_path_st_pcmw_extend(SCIP *, const GRAPH *, const SCIP_Real *, SCIP_Bool, PATH *, STP_Bool *, SCIP_Bool *)
Definition: graph_path.c:1705
#define DEFAULT_MAXFREQLOC
Definition: heur_local.c:55
#define NULL
Definition: lpi_spx1.cpp:155
SCIP_HEUR * SCIPsolGetHeur(SCIP_SOL *sol)
Definition: sol.c:2604
SCIP_RETCODE SCIPStpHeurLocalRun(SCIP *scip, GRAPH *graph, int *solEdges)
Definition: heur_local.c:3883
#define STP_RMWCSP
Definition: graphdefs.h:54
static SCIP_RETCODE soltreeExchangeKeyPath(SCIP *scip, GRAPH *graph, const CONN *connectData, const VNOILOC *vnoiData, const KPATHS *keypathsData, const int *dfstree, const STP_Bool *scanned, int dfstree_pos, int boundedge_new, SOLTREE *soltreeData)
Definition: heur_local.c:1376
int knots
Definition: graphdefs.h:190
#define SCIP_CALL(x)
Definition: def.h:384
void SCIPlinkcuttreeCutNode(LCNODE *v)
Definition: misc_stp.c:372
void SCIPStpunionfindUnion(UF *uf, int p, int q, SCIP_Bool compress)
Definition: misc_stp.c:912
static SCIP_Real pcmwGetNewEdgePrize(const GRAPH *graph, const STP_Bool *steinertree, int edge, PCMW *pcmwData)
Definition: heur_local.c:716
int * term2edge
Definition: graphdefs.h:208
#define MST_MODE
Definition: graphdefs.h:98
int SCIPStpunionfindFind(UF *uf, int element)
Definition: misc_stp.c:885
static void pcmwDataClean(const GRAPH *graph, PCMW *pcmwData)
Definition: heur_local.c:780
#define STP_PCSPG
Definition: graphdefs.h:44
#define STP_OARSMT
Definition: graphdefs.h:50
void * SCIPpqueueRemove(SCIP_PQUEUE *pqueue)
Definition: misc.c:1434
STP_Bool *const nodeIsPinned
Definition: heur_local.c:136
Improvement heuristic for Steiner problems.
SCIP_RETCODE SCIPStpValidateSol(SCIP *, const GRAPH *, const double *, SCIP_Bool, SCIP_Bool *)
Definition: validate.c:202
unsigned char STP_Bool
Definition: portab.h:34
static SCIP_DECL_HEURFREE(heurFreeLocal)
Definition: heur_local.c:3603
int GNODECmpByDist(void *first_arg, void *second_arg)
const SCIP_Real * edgecosts
Definition: heur_local.c:173
SCIP_RETCODE solstp_pruneFromEdges(SCIP *scip, const GRAPH *g, int *result)
Definition: solstp.c:1432
SCIP_RETCODE SCIPcreateRandom(SCIP *scip, SCIP_RANDNUMGEN **randnumgen, unsigned int initialseed, SCIP_Bool useglobalseed)
SCIP_Real solstp_getObjBounded(const GRAPH *g, const int *soledge, SCIP_Real offset, int nedges)
Definition: solstp.c:1833
static SCIP_Bool solNeedsPruning(SCIP_SOL *initialsol)
Definition: heur_local.c:194
#define SCIPallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:115
void graph_heap_free(SCIP *, SCIP_Bool, SCIP_Bool, DHEAP **)
Definition: graph_util.c:1034
static void solGetStpSol(SCIP *scip, const GRAPH *graph, SCIP_SOL *initialsol, int *results)
Definition: heur_local.c:258
#define SCIP_Bool
Definition: def.h:84
int *RESTRICT ieat
Definition: graphdefs.h:230
int *RESTRICT path_heap
Definition: graphdefs.h:255
#define DEFAULT_RANDSEED
Definition: heur_local.c:62
int *RESTRICT tail
Definition: graphdefs.h:223
#define LOCAL_MAXRESTARTS_FAST
Definition: heur_local.c:64
#define DEFAULT_MINNBESTSOLS
Definition: heur_local.c:60
static void insertData(SCIP *scip, STP_Bitset termsmark, STP_Bitset rootsmark, int64_t nsubsets, int node_pos, int split_pos, int *subrootpos, DPSTREE *dpstree)
Definition: dpterms_util.c:172
void SCIPrandomPermuteIntArray(SCIP_RANDNUMGEN *randnumgen, int *array, int begin, int end)
Definition: misc.c:10044
#define flipedge(edge)
Definition: graphdefs.h:84
static SCIP_RETCODE localRun(SCIP *scip, GRAPH *graph, SCIP_Bool useFast, int *solEdges)
Definition: heur_local.c:3802
static void insertionInit(SCIP_HEURDATA *heurdata, const GRAPH *graph, const int *solEdges, int *solDegreeNonTerm, SCIP_Bool *nodeIsBlocked, int *vertices)
Definition: heur_local.c:2500
static SCIP_DECL_HEUREXEC(heurExecLocal)
Definition: heur_local.c:3671
void SCIPlinkcuttreeLink(LCNODE *tree, int v, int w, int edge)
Definition: misc_stp.c:357
#define LOCAL_MAXRESTARTS
Definition: heur_local.c:63
int SCIPgetNSols(SCIP *scip)
Definition: scip_sol.c:2205
#define DEFAULT_NELITESOLS
Definition: heur_local.c:59
int *RESTRICT term
Definition: graphdefs.h:196
static void insertionIncrementSolDegree(const GRAPH *graph, int node, INSERT *insertData)
Definition: heur_local.c:2477
#define BMScopyMemoryArray(ptr, source, num)
Definition: memory.h:127
void graph_edge_printInfo(const GRAPH *, int)
Definition: graph_stats.c:133
static SCIP_DECL_HEUREXITSOL(heurExitsolLocal)
Definition: heur_local.c:3651
SCIP_Real SCIPgetSolOrigObj(SCIP *scip, SCIP_SOL *sol)
Definition: scip_sol.c:1435
SCIP_Bool graph_pc_knotIsNonLeafTerm(const GRAPH *, int)
void graph_path_st_pcmw_extendOut(SCIP *, const GRAPH *, int, STP_Bool *, SCIP_Real *, int *, STP_Bool *, DHEAP *, SCIP_Bool *)
Definition: graph_path.c:1050
static void soltreeUnmarkKpNodes(const GRAPH *graph, const KPATHS *keypathsData, SOLTREE *soltreeData)
Definition: heur_local.c:1328
const int * SCIPStpGetPcImplVerts(SCIP *scip)
Definition: cons_stp.c:1267
SCIP_RETCODE graph_trail_costAware(SCIP *, const GRAPH *, int, SCIP_Bool *)
Definition: graph_util.c:353
void graph_path_exit(SCIP *, GRAPH *)
Definition: graph_path.c:515
SCIP_Bool graph_pc_isPc(const GRAPH *)
static void getKeyPathUpper(SCIP *scip, int crucnode, const GRAPH *graph, const SOLTREE *soltreeData, const PCMW *pcmwData, CONN *connectData, KPATHS *keypathsData)
Definition: heur_local.c:1201
#define CONNECT
Definition: graphdefs.h:87
int *const cutedgesStart
Definition: heur_local.c:176
int * chainStarts
Definition: heur_local.c:171
struct Voronoi_data_structures VNOILOC
static SCIP_Real getBoundaryPathCostPrized(const GRAPH *graph, const VNOILOC *vnoiData, const SOLTREE *soltreeData, int boundaryedge, PCMW *pcmwData)
Definition: heur_local.c:974
int *const prizemarklist
Definition: heur_local.c:161
#define SCIPfreeMemory(scip, ptr)
Definition: scip_mem.h:69
SCIP_RETCODE SCIPpqueueCreate(SCIP_PQUEUE **pqueue, int initsize, SCIP_Real sizefac, SCIP_DECL_SORTPTRCOMP((*ptrcomp)), SCIP_DECL_PQUEUEELEMCHGPOS((*elemchgpos)))
Definition: misc.c:1236
static void insertionBlockChain(const GRAPH *graph, const LCNODE *lctree, int chainfirst_index, int chainlast_index, INSERT *insertData)
Definition: heur_local.c:2669
STP_Bool *const prizemark
Definition: heur_local.c:160
#define Is_pseudoTerm(a)
Definition: graphdefs.h:103
static SCIP_Bool solIsTrivialPcMw(const GRAPH *graph, const int *solEdges)
Definition: heur_local.c:499
static SCIP_RETCODE connectivityDataInit(SCIP *scip, const GRAPH *graph, const VNOILOC *vnoiData, const SOLTREE *soltreeData, const PCMW *pcmwData, CONN *connectData)
Definition: heur_local.c:1122
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip_sol.c:2304
#define DEFAULT_NBESTSOLS_HARD
Definition: heur_local.c:58
SCIP_Bool SCIPisGT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
static SCIP_Real getNewPrizeNode(const GRAPH *graph, const STP_Bool *steinertree, const SCIP_Real *prize_biased, const int *graphmark, int node, STP_Bool *prizemark, int *prizemarklist, int *prizemarkcount)
Definition: heur_local.c:659
struct SCIP_ProbData SCIP_PROBDATA
Definition: type_prob.h:44
int solroot
Definition: heur_local.c:164
SCIP_RETCODE SCIPpqueueInsert(SCIP_PQUEUE *pqueue, void *elem)
Definition: misc.c:1335
SCIP_PROBDATA * SCIPgetProbData(SCIP *scip)
Definition: scip_prob.c:962
void graph_pc_2orgcheck(SCIP *, GRAPH *)
static void insertionReplaceChain(const GRAPH *graph, int newedge, LCNODE *lctree, INSERT *insertData, int v_lc, int headCurr_lc, int chainfirst_index, int chainlast_index)
Definition: heur_local.c:2778
SCIP_Bool graph_pc_knotIsFixedTerm(const GRAPH *, int)
SCIP_Real * cost
Definition: graphdefs.h:221
static void dfsorder(const GRAPH *graph, const int *edges, int node, int *counter, int *dfst)
Definition: heur_local.c:633
#define DEFAULT_MINNBESTSOLS_HARD
Definition: heur_local.c:61
static int pcmwGetSolRoot(const GRAPH *graph, const int *solEdges)
Definition: heur_local.c:686
SCIP_RETCODE graph_heap_create(SCIP *, int, int *, DENTRY *, DHEAP **)
Definition: graph_util.c:992
struct solution_tree_data SOLTREE
SCIP_Real *const edgecost_org
Definition: heur_local.c:159
SCIP_Bool graph_pc_isPcMw(const GRAPH *)
#define GREEDY_MAXRESTARTS
Definition: heur_local.c:67
SCIP_Real * memvdist
Definition: heur_local.c:94
struct keypaths_data_structures KPATHS
#define HEUR_NAME
Definition: heur_local.c:43
const int * SCIPStpGetPcImplStarts(SCIP *scip)
Definition: cons_stp.c:1229
#define SCIP_Real
Definition: def.h:177
SCIP_RETCODE SCIPStpHeurLocalExtendPcMwImp(SCIP *scip, const GRAPH *graph, int *result)
Definition: heur_local.c:3907
#define StpVecFree(scip, vec)
Definition: stpvector.h:149
SCIP_RETCODE SCIPpairheapBuffarr(SCIP *scip, const PHNODE *root, int size, int **elements)
Definition: misc_stp.c:761
void graph_free_csr(SCIP *, GRAPH *)
Definition: graph_util.c:1623
void graph_voronoi(const GRAPH *, const SCIP_Real *, const SCIP_Real *, const STP_Bool *, int *, PATH *)
Definition: graph_vnoi.c:333
SCIP_Real * SCIPprobdataGetXval(SCIP *scip, SCIP_SOL *sol)
int *RESTRICT outbeg
Definition: graphdefs.h:204
#define DEFAULT_NBESTSOLS
Definition: heur_local.c:57
shortest paths based primal heuristics for Steiner problems
STP_Bool *const nodeIsLowSupernode
Definition: heur_local.c:146
int edges
Definition: graphdefs.h:219
SCIP_Real *const prize_biased
Definition: heur_local.c:157
static SCIP_Bool solNodeIsValid(SCIP *scip, const GRAPH *graph, const STP_Bool *solNodes, const LCNODE *linkcutNodes)
Definition: heur_local.c:2402
void graph_pc_getOrgCosts(SCIP *, const GRAPH *, SCIP_Real *)
#define STP_GSTP
Definition: graphdefs.h:53
void solstp_setVertexFromEdge(const GRAPH *g, const int *result, STP_Bool *solnode)
Definition: solstp.c:2078
int graph_pc_nNonFixedTerms(const GRAPH *)
SCIP_RETCODE SCIPStpHeurLocalExtendPcMw(SCIP *scip, GRAPH *graph, const SCIP_Real *cost, int *stedge, STP_Bool *stvertex)
Definition: heur_local.c:3975
#define SCIPallocMemory(scip, ptr)
Definition: scip_mem.h:51
SCIP_Bool graph_knot_isInRange(const GRAPH *, int)
Definition: graph_node.c:52
static void insertionInitInsert(SCIP *scip, const GRAPH *graph, int v_insert, int initialEdge, LCNODE *linkcutNodes, INSERT *insertData, SCIP_Real *diff)
Definition: heur_local.c:2539
SCIP_Bool SCIPisLE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
#define UNKNOWN
Definition: sepa_mcf.c:4095
static SCIP_RETCODE lca(SCIP *scip, const GRAPH *graph, int u, UF *uf, STP_Bool *nodesmark, const int *steineredges, STP_Vectype(int) *lcalists, PHNODE **boundpaths, const int *heapsize, const int *vbase)
Definition: heur_local.c:543
SCIP_RETCODE SCIPsetHeurCopy(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEURCOPY((*heurcopy)))
Definition: scip_heur.c:153
#define nnodes
Definition: gastrans.c:65
#define StpVecPushBack(scip, vec, value)
Definition: stpvector.h:163
SCIP_RETCODE graph_init(SCIP *, GRAPH **, int, int, int)
Definition: graph_base.c:607
#define BMSclearMemoryArray(ptr, num)
Definition: memory.h:123
SCIP_Bool SCIPStpunionfindIsClear(SCIP *scip, const UF *uf)
Definition: misc_stp.c:861
SCIP_Bool solstp_isValid(SCIP *scip, const GRAPH *graph, const int *result)
Definition: solstp.c:1650
STP_Bool *const solNodes
Definition: heur_local.c:133
#define SCIP_CALL_ABORT(x)
Definition: def.h:363
SCIP_HEURDATA * SCIPheurGetData(SCIP_HEUR *heur)
Definition: heur.c:1352
#define HEUR_FREQOFS
Definition: heur_local.c:48
SCIP_Real *const edgecost_biased
Definition: heur_local.c:158
SCIP_RETCODE solstp_prune(SCIP *scip, const GRAPH *g, int *result, STP_Bool *connected)
Definition: solstp.c:1366
LCNODE *const linkcutNodes
Definition: heur_local.c:134
void SCIPStpunionfindFreeMembers(SCIP *scip, UF *uf)
Definition: misc_stp.c:955
static SCIP_RETCODE supergraphComputeMst(SCIP *scip, const GRAPH *graph, const CONN *connectData, const VNOILOC *vnoiData, const KPATHS *keypathsData, int crucnode, SOLTREE *soltreeData, PCMW *pcmwData, SGRAPH *supergraphData)
Definition: heur_local.c:1804
#define HEUR_DISPCHAR
Definition: heur_local.c:45
#define STP_RPCSPG
Definition: graphdefs.h:45
int *const blockedList
Definition: heur_local.c:180
#define DEFAULT_MAXNBESTSOLS
Definition: heur_local.c:56
#define STP_MWCSP
Definition: graphdefs.h:51
static SCIP_RETCODE addToCandidates(SCIP *scip, const GRAPH *graph, const PATH *path, int i, int greedyextensions, int *nextensions, GNODE *candidates, SCIP_PQUEUE *pqueue)
Definition: heur_local.c:415
int SCIPsolGetIndex(SCIP_SOL *sol)
Definition: sol.c:2635
SCIP_RETCODE SCIPaddBoolParam(SCIP *scip, const char *name, const char *desc, SCIP_Bool *valueptr, SCIP_Bool isadvanced, SCIP_Bool defaultvalue, SCIP_DECL_PARAMCHGD((*paramchgd)), SCIP_PARAMDATA *paramdata)
Definition: scip_param.c:48
static SCIP_RETCODE soltreeElimKeyPathsStar(SCIP *scip, const GRAPH *graph, const CONN *connectData, const VNOILOC *vnoiData, const KPATHS *keypathsData, const SGRAPH *supergraphData, const int *dfstree, const STP_Bool *scanned, int dfstree_pos, SOLTREE *soltreeData)
Definition: heur_local.c:1525
static void soltreeMarkKpNodes(const GRAPH *graph, const KPATHS *keypathsData, SOLTREE *soltreeData)
Definition: heur_local.c:1352
static SCIP_RETCODE connectivityDataKeyElimUpdate(SCIP *scip, const GRAPH *graph, const VNOILOC *vnoiData, const SGRAPH *supergraphData, int crucnode, CONN *connectData)
Definition: heur_local.c:1046