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-2019 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 scip.zib.de. */
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 
29 #include <assert.h>
30 #include <string.h>
31 
32 #include "heur_local.h"
33 #include "heur_tm.h"
34 #include "probdata_stp.h"
35 
36 
37 /* @note if heuristic is running in root node timing is changed there to (SCIP_HEURTIMING_DURINGLPLOOP |
38  * SCIP_HEURTIMING_BEFORENODE), see SCIP_DECL_HEURINITSOL callback
39  */
40 
41 #define HEUR_NAME "local"
42 #define HEUR_DESC "improvement heuristic for STP"
43 #define HEUR_DISPCHAR '-'
44 #define HEUR_PRIORITY 100
45 #define HEUR_FREQ 1
46 #define HEUR_FREQOFS 0
47 #define HEUR_MAXDEPTH -1
48 #define HEUR_TIMING (SCIP_HEURTIMING_BEFORENODE | SCIP_HEURTIMING_DURINGLPLOOP | SCIP_HEURTIMING_AFTERLPLOOP | SCIP_HEURTIMING_AFTERNODE)
49 
50 #define HEUR_USESSUBSCIP FALSE /**< does the heuristic use a secondary SCIP instance? */
51 
52 #define DEFAULT_DURINGROOT TRUE
53 #define DEFAULT_MAXFREQLOC FALSE
54 #define DEFAULT_MAXNBESTSOLS 25
55 #define DEFAULT_NBESTSOLS 10
56 #define DEFAULT_MINNBESTSOLS 6
57 #define LOCAL_MAXRESTARTS 5
58 
59 #define GREEDY_MAXRESTARTS 3 /**< Max number of restarts for greedy PC/MW heuristic if improving solution has been found. */
60 #define GREEDY_EXTENSIONS_MW 6 /**< Number of extensions for greedy MW heuristic. MUST BE HIGHER THAN GREEDY_EXTENSIONS */
61 #define GREEDY_EXTENSIONS 5 /**< Number of extensions for greedy PC heuristic. */
62 
63 
64 /*
65  * Data structures
66  */
67 
68 /** primal heuristic data */
69 struct SCIP_HeurData
70 {
71  int nfails; /**< number of fails */
72  int maxnsols; /**< maximal number of best solutions to improve */
73  int nbestsols; /**< number of best solutions to improve */
74  int* lastsolindices; /**< indices of a number of best solutions already tried */
75  SCIP_Bool maxfreq; /**< should the heuristic be called with maximum frequency? */
76  SCIP_Bool duringroot; /**< should the heuristic be called during the root node? */
77 };
78 
79 
80 /*
81  * Local methods
82  */
83 
84 /** recursive methode for a DFS ordering of graph 'g' */
85 static
86 void dfsorder(
87  const GRAPH* graph,
88  const int* edges,
89  const int* node,
90  int* counter,
91  int* dfst
92  )
93 {
94  int oedge = graph->outbeg[*node];
95 
96  while( oedge >= 0 )
97  {
98  if( edges[oedge] >= 0 )
99  {
100  dfsorder(graph, edges, &(graph->head[oedge]), counter, dfst);
101  }
102  oedge = graph->oeat[oedge];
103  }
104 
105  dfst[*counter] = *node;
106  (*counter)++;
107 }
108 
109 
110 /** computes lowest common ancestors for all pairs {vbase(v), vbase(u)} such that {u,w} is a boundary edge,
111  * first call should be with u := root */
112 static
114  SCIP* scip,
115  const GRAPH* graph,
116  int u,
117  UF* uf,
118  STP_Bool* nodesmark,
119  int* steineredges,
120  IDX** lcalists,
121  PHNODE** boundpaths,
122  int* heapsize,
123  int* vbase
124  )
125 {
126  int* uboundpaths; /* boundary-paths (each one represented by its boundary edge) having node 'u' as an endpoint */
127  int ancestor;
128  int v;
129  int i;
130  int oedge; /* outgoing edge */
131  IDX* curr;
132  uf->parent[u] = u;
133 
134  for( oedge = graph->outbeg[u]; oedge != EAT_LAST; oedge = graph->oeat[oedge] )
135  {
136  v = graph->head[oedge];
137  if( steineredges[oedge] == 0 )
138  {
139  SCIP_CALL( lca(scip, graph, v, uf, nodesmark, steineredges, lcalists, boundpaths, heapsize, vbase) );
140  SCIPStpunionfindUnion(uf, u, v, FALSE);
141  uf->parent[SCIPStpunionfindFind(uf, u)] = u;
142  }
143  }
144  nodesmark[u] = TRUE;
145 
146  /* iterate through all boundary-paths having one endpoint in the voronoi region of node 'u' */
147  SCIP_CALL( SCIPpairheapBuffarr(scip, boundpaths[u], heapsize[u], &uboundpaths) );
148  for( i = 0; i < heapsize[u]; i++ )
149  {
150  oedge = uboundpaths[i];
151  v = vbase[graph->head[oedge]];
152  if( nodesmark[v] )
153  {
154  ancestor = uf->parent[SCIPStpunionfindFind(uf, v)];
155 
156  /* if the ancestor of 'u' and 'v' is one of the two, the boundary-edge is already in boundpaths[u] */
157  if( ancestor != u && ancestor != v)
158  {
159  SCIP_CALL( SCIPallocBlockMemory(scip, &curr) );
160  curr->index = oedge;
161  curr->parent = lcalists[ancestor];
162  lcalists[ancestor] = curr;
163  }
164  }
165  }
166 
167  /* free the boundary-paths array */
168  SCIPfreeBufferArray(scip, &uboundpaths);
169 
170  return SCIP_OKAY;
171 }
172 
173 /** checks whether node is crucial, i.e. a terminal or a vertex with degree at least 3 (w.r.t. the steinertree) */
174 static
176  const GRAPH* graph,
177  int* steineredges,
178  int node
179  )
180 {
181  assert(graph != NULL);
182  assert(steineredges != NULL);
183 
184  if( graph->term[node] == -1 )
185  {
186  int counter = 0;
187  int e = graph->outbeg[node];
188  while( e >= 0 )
189  {
190  /* check if the adjacent node is in the ST */
191  if( steineredges[e] > -1 || steineredges[flipedge(e)] > -1 )
192  {
193  counter++;
194  }
195  e = graph->oeat[e];
196  }
197 
198  if( counter < 3 )
199  {
200  return FALSE;
201  }
202  }
203 
204  return TRUE;
205 }
206 
207 /** perform local heuristics on a given Steiner tree todo delete cost parameter */
209  SCIP* scip, /**< SCIP data structure */
210  GRAPH* graph, /**< graph data structure */
211  const SCIP_Real* cost, /**< arc cost array todo delete this parameter and use graph->cost */
212  int* best_result /**< array indicating whether an arc is part of the solution (CONNECTED/UNKNOWN) */
213  )
214 {
215  NODE* nodes;
216  PATH* vnoi;
217  int* vbase;
218  int* const graphmark = graph->mark;
219  int e;
220  int i;
221  int k;
222  const int root = graph->source;
223  const int nnodes = graph->knots;
224  const int nedges = graph->edges;
225  const int probtype = graph->stp_type;
226  int newnverts;
227  STP_Bool* steinertree;
228  const STP_Bool pc = ((probtype == STP_PCSPG) || (probtype == STP_RPCSPG));
229  const STP_Bool mw = (probtype == STP_MWCSP);
230  const STP_Bool mwpc = (pc || mw);
231 #ifndef NDEBUG
232  const SCIP_Real initialobj = graph_sol_getObj(graph->cost, best_result, 0.0, nedges);
233 #endif
234 
235  assert(graph != NULL);
236  assert(cost != NULL);
237  assert(best_result != NULL);
238  assert(graph_valid(graph));
239 
240  newnverts = 0;
241 
242  if( graph->grad[root] == 0 || graph->terms == 1 )
243  return SCIP_OKAY;
244 
245  /* for PC variants test whether solution is trivial */
246  if( mwpc )
247  {
248  for( e = graph->outbeg[root]; e != EAT_LAST; e = graph->oeat[e] )
249  if( !Is_term(graph->term[graph->head[e]]) && best_result[e] )
250  break;
251 
252  if( e == EAT_LAST )
253  {
254  SCIPdebugMessage("Local heuristic: return trivial \n");
255  return SCIP_OKAY;
256  }
257  }
258 
259  SCIP_CALL( SCIPallocBufferArray(scip, &vnoi, nnodes) );
260  SCIP_CALL( SCIPallocBufferArray(scip, &vbase, nnodes) );
261  SCIP_CALL( SCIPallocBufferArray(scip, &nodes, nnodes) );
262  SCIP_CALL( SCIPallocBufferArray(scip, &steinertree, nnodes) );
263 
264  if( mwpc )
265  {
266  SCIP_Bool dummy;
267  SCIP_CALL( SCIPStpHeurLocalExtendPcMw(scip, graph, cost, vnoi, best_result, steinertree, &dummy) );
268  }
269 
270  for( i = 0; i < nnodes; i++ )
271  {
272  steinertree[i] = FALSE;
273  SCIPlinkcuttreeInit(&nodes[i]);
274  }
275 
276  /* create a link-cut tree representing the current Steiner tree */
277  for( e = 0; e < nedges; e++ )
278  {
279  assert(graph->head[e] == graph->tail[flipedge(e)]);
280 
281  /* if edge e is in the tree, so are its incident vertices */
282  if( best_result[e] == CONNECT )
283  {
284  steinertree[graph->tail[e]] = TRUE;
285  steinertree[graph->head[e]] = TRUE;
286  SCIPlinkcuttreeLink(&nodes[graph->head[e]], &nodes[graph->tail[e]], flipedge(e));
287  }
288  }
289 
290  assert( nodes[root].edge == -1 );
291  nodes[root].edge = -1;
292 
293  /* VERTEX INSERTION */
294  if( probtype == STP_SPG || probtype == STP_RSMT || probtype == STP_OARSMT || probtype == STP_GSTP || (mwpc) )
295  {
296  int newnode;
297  int* insert;
298  int* adds;
299  int* cuts;
300  int* cuts2;
301  int* stdeg;
302 
303  SCIP_CALL( SCIPallocBufferArray(scip, &insert, nnodes) );
304  SCIP_CALL( SCIPallocBufferArray(scip, &adds, nnodes) );
305  SCIP_CALL( SCIPallocBufferArray(scip, &cuts, nnodes) );
306 
307  if( mw )
308  {
309  SCIP_CALL( SCIPallocBufferArray(scip, &cuts2, nnodes) );
310  SCIP_CALL( SCIPallocBufferArray(scip, &stdeg, nnodes) );
311 
312  BMSclearMemoryArray(stdeg, nnodes);
313 
314  for( e = 0; e < nedges; e++ )
315  if( best_result[e] == CONNECT )
316  {
317  stdeg[graph->tail[e]]++;
318  stdeg[graph->head[e]]++;
319  }
320  }
321  else
322  {
323  cuts2 = NULL;
324  stdeg = NULL;
325  }
326 
327  i = 0;
328  newnode = 0;
329 
330  for( ;; )
331  {
332  SCIP_Real diff;
333 
334  /* if vertex i is not in the current ST and has at least two adjacent nodes, it might be added */
335  if( !steinertree[i] && graph->grad[i] > 1 && (!mwpc || !Is_term(graph->term[i])) )
336  {
337  NODE* v;
338  int counter;
339  int lastnodeidx;
340  int insertcount = 0;
341 
342  /* if an outgoing edge of vertex i points to the current ST, SCIPlinkcuttreeLink the edge to a list */
343  for( int oedge = graph->outbeg[i]; oedge != EAT_LAST; oedge = graph->oeat[oedge])
344  if( steinertree[graph->head[oedge]] && (!mwpc || !Is_term(graph->term[graph->head[oedge]])) )
345  insert[insertcount++] = oedge;
346 
347  /* if there are less than two edges connecting node i and the current tree, continue */
348  if( insertcount <= 1 )
349  goto ENDOFLOOP;
350 
351  if( mw )
352  SCIPlinkcuttreeInit(&nodes[i]);
353 
354  /* the node to insert */
355  v = &nodes[i];
356 
357  SCIPlinkcuttreeLink(v, &nodes[graph->head[insert[0]]], insert[0]);
358 
359  lastnodeidx = graph->head[insert[0]];
360 
361  if( mw )
362  {
363  assert(!SCIPisPositive(scip, graph->prize[i]));
364 
365  diff = -1.0;
366  assert(stdeg != NULL);
367  stdeg[lastnodeidx]++;
368  }
369  else
370  diff = graph->cost[v->edge];
371 
372  counter = 0;
373 
374  /* try to add edges between new vertex and tree */
375  for( k = 1; k < insertcount; k++ )
376  {
377  NODE* firstnode;
378  int firstnodidx;
380 
381  /* next vertex in the current Steiner tree adjacent to vertex i resp. v (the one being scrutinized for possible insertion) */
382  firstnodidx = graph->head[insert[k]];
383  firstnode = &nodes[firstnodidx];
384 
385  if( mw )
386  {
387  NODE* chainfirst;
388  NODE* chainlast;
389  SCIP_Real minweight;
390 
391  assert(stdeg != NULL);
392 
393  minweight = SCIPlinkcuttreeFindMinChain(scip, graph->prize, graph->head, stdeg, firstnode, &chainfirst, &chainlast);
394 
395  if( SCIPisLT(scip, minweight, graph->prize[i]) )
396  {
397  assert(chainfirst != NULL && chainlast != NULL);
398  for( NODE* mynode = chainfirst; mynode != chainlast; mynode = mynode->parent )
399  {
400  int mynodeidx = graph->head[mynode->edge];
401  steinertree[mynodeidx] = FALSE;
402  stdeg[mynodeidx] = 0;
403  }
404 
405  SCIPlinkcuttreeCut(chainfirst);
406  SCIPlinkcuttreeCut(chainlast);
407 
408  SCIPlinkcuttreeLink(v, firstnode, insert[k]);
409  stdeg[graph->head[insert[k]]]++;
410 
411  diff = graph->prize[i] - minweight;
412  break;
413  }
414  }
415  else
416  {
417  /* if there is an edge with cost greater than that of the current edge... */
418  NODE* max = SCIPlinkcuttreeFindMax(scip, graph->cost, firstnode);
419  if( SCIPisGT(scip, graph->cost[max->edge], graph->cost[insert[k]]) )
420  {
421  diff += graph->cost[insert[k]];
422  diff -= graph->cost[max->edge];
423  cuts[counter] = max->edge;
424  SCIPlinkcuttreeCut(max);
425  SCIPlinkcuttreeLink(v, firstnode, insert[k]);
426  assert(v->edge == insert[k]);
427  adds[counter++] = v->edge;
428  }
429  }
430  }
431 
432  if( pc && Is_pterm(graph->term[i]) )
433  diff -= graph->prize[i];
434 
435  /* if the new tree is more expensive than the old one, restore the latter */
436  if( mw )
437  {
438  if( SCIPisLT(scip, diff, 0.0) )
439  {
440  assert(stdeg != NULL);
441 
443  stdeg[lastnodeidx]--;
444  SCIPlinkcuttreeCut(&nodes[graph->head[insert[0]]]);
445  }
446  else
447  {
448  steinertree[i] = TRUE;
449  newnverts++;
450  }
451  }
452  else
453  {
454  if( !SCIPisNegative(scip, diff) )
455  {
457  for (k = counter - 1; k >= 0; k--)
458  {
459  SCIPlinkcuttreeCut(&nodes[graph->head[adds[k]]]);
460  SCIPlinkcuttreeEvert(&nodes[graph->tail[cuts[k]]]);
461  SCIPlinkcuttreeLink(&nodes[graph->tail[cuts[k]]],
462  &nodes[graph->head[cuts[k]]], cuts[k]);
463  }
464 
465  /* finally, cut the edge added first (if it had been cut during the insertion process, it would have been restored above) */
467  SCIPlinkcuttreeCut(&nodes[graph->head[insert[0]]]);
468  }
469  else
470  {
471  SCIPlinkcuttreeEvert(&nodes[root]);
472  adds[counter] = insert[0];
473  newnode = i;
474  steinertree[i] = TRUE;
475  newnverts++;
476  SCIPdebugMessage("ADDED VERTEX \n");
477  }
478  }
479  }
480 
481  ENDOFLOOP:
482 
483  if( i < nnodes - 1 )
484  i++;
485  else
486  i = 0;
487 
488  if( newnode == i )
489  break;
490  }
491 
492  /* free buffer memory */
493  if( mw )
494  {
495  SCIPfreeBufferArray(scip, &stdeg);
496  SCIPfreeBufferArray(scip, &cuts2);
497  }
498  SCIPfreeBufferArray(scip, &cuts);
499  SCIPfreeBufferArray(scip, &adds);
500  SCIPfreeBufferArray(scip, &insert);
501 
502  for( e = 0; e < nedges; e++ )
503  best_result[e] = UNKNOWN;
504 
505  if( newnverts > 0 )
506  {
507  if( mwpc )
508  SCIP_CALL( SCIPStpHeurTMPrunePc(scip, graph, graph->cost, best_result, steinertree) );
509  else
510  SCIP_CALL( SCIPStpHeurTMPrune(scip, graph, graph->cost, 0, best_result, steinertree) );
511 
512  for( i = 0; i < nnodes; i++ )
513  SCIPlinkcuttreeInit(&nodes[i]);
514 
515  /* create a link-cut tree representing the current Steiner tree */
516  for( e = 0; e < nedges; e++ )
517  {
518  if( best_result[e] == CONNECT )
519  {
520  assert(steinertree[graph->tail[e]]);
521  assert(steinertree[graph->head[e]]);
522  SCIPlinkcuttreeLink(&nodes[graph->head[e]], &nodes[graph->tail[e]], flipedge(e));
523  }
524  }
525  SCIPlinkcuttreeEvert(&nodes[root]);
526  }
527  else
528  {
529  SCIPlinkcuttreeEvert(&nodes[root]);
530  for( i = 0; i < nnodes; i++ )
531  {
532  if( steinertree[i] && nodes[i].edge != -1 )
533  best_result[flipedge(nodes[i].edge)] = 0;
534  }
535  }
536 
537 #ifdef SCIP_DEBUG
538  {
539  SCIP_Real obj = 0.0;
540  for( e = 0; e < nedges; e++ )
541  obj += (best_result[e] > -1) ? graph->cost[e] : 0.0;
542  printf("ObjAfterVertexInsertion=%.12e\n", obj);
543 
544  }
545 #endif
546  }
547 
548  assert(graph_sol_valid(scip, graph, best_result));
549 
550  /* Key-Vertex Elimination & Key-Path Exchange */
551  if( !mw )
552  {
553  IDX* blists_curr;
554  IDX** blists_start; /* array [1,..,nnodes],
555  * if node i is in the current ST, blists_start[i] points to a linked list of all nodes having i as their base */
556  PATH* mst; /* minimal spanning tree structure */
557  GRAPH* supergraph;
558  IDX** lvledges_start; /* horizontal edges */
559  IDX* lvledges_curr;
560  PHNODE** boundpaths;
561  UF uf; /* union-find*/
562  SCIP_Real* memdist;
563  SCIP_Real kpcost;
564  SCIP_Real mstcost;
565  SCIP_Real edgecost;
566  SCIP_Real kpathcost;
567  int* state;
568  int* kpedges;
569  int* kpnodes;
570  int* dfstree;
571  int* newedges;
572  int* memvbase;
573  int* heapsize;
574  int* boundedges;
575  int* meminedges;
576  int* supernodes;
577  int* supernodesid;
578  int l;
579  int node;
580  int edge;
581  int count;
582  int nruns;
583  int newedge;
584  int oldedge;
585  int adjnode;
586  int nkpedges;
587  int nstnodes;
588  int nkpnodes;
589  int crucnode; /* current crucial node*/
590  int nresnodes;
591  int kptailnode; /* tail node of the current keypath*/
592  int localmoves = 2;
593  int newpathend = -1;
594  int nsupernodes;
595  int nboundedges;
596  int rootpathstart;
597  STP_Bool* pinned;
598  STP_Bool* scanned;
599  STP_Bool* nodesmark;
600 
601 #ifdef SCIP_DEBUG
602  SCIP_Real obj = 0.0;
603  for( e = 0; e < nedges; e++ )
604  obj += (best_result[e] > -1) ? graph->cost[e] : 0.0;
605  printf(" ObjBEFKEYVertexELimination=%.12e\n", obj);
606 #endif
607 
608  for( k = 0; k < nnodes; k++ )
609  graphmark[k] = (graph->grad[k] > 0);
610 
611  graphmark[root] = TRUE;
612 
613  /* allocate memory */
614 
615  /* only needed for Key-Path Elimination */
616  SCIP_CALL( SCIPallocBufferArray(scip, &newedges, nedges) );
617  SCIP_CALL( SCIPallocBufferArray(scip, &lvledges_start, nnodes) );
618  SCIP_CALL( SCIPallocBufferArray(scip, &boundedges, nedges) );
619  SCIP_CALL( SCIPallocBufferArray(scip, &supernodesid, nnodes) );
620 
621  /* only needed for Key-Path Exchange */
622 
623  /* memory needed for both Key-Path Elimination and Exchange */
624  SCIP_CALL( SCIPallocBufferArray(scip, &scanned, nnodes) );
625  SCIP_CALL( SCIPallocBufferArray(scip, &heapsize, nnodes) );
626  SCIP_CALL( SCIPallocBufferArray(scip, &blists_start, nnodes) );
627  SCIP_CALL( SCIPallocBufferArray(scip, &memvbase, nnodes) );
628  SCIP_CALL( SCIPallocBufferArray(scip, &memdist, nnodes) );
629  SCIP_CALL( SCIPallocBufferArray(scip, &meminedges, nnodes) );
630  SCIP_CALL( SCIPallocBufferArray(scip, &boundpaths, nnodes) );
631  SCIP_CALL( SCIPallocBufferArray(scip, &pinned, nnodes) );
632  SCIP_CALL( SCIPallocBufferArray(scip, &dfstree, nnodes) );
633  SCIP_CALL( SCIPallocBufferArray(scip, &nodesmark, nnodes) );
634 
635  /* initialize data structures */
636  SCIP_CALL( SCIPStpunionfindInit(scip, &uf, nnodes) );
637 
638  for( nruns = 0; nruns < LOCAL_MAXRESTARTS && localmoves > 0; nruns++ )
639  {
640  localmoves = 0;
641 
642  BMSclearMemoryArray(blists_start, nnodes);
643 
644  /* find a DFS order of the ST nodes */
645  nstnodes = 0;
646  dfsorder(graph, best_result, &(root), &nstnodes, dfstree);
647  assert(root == graph->source);
648 
649  /* compute a voronoi diagram with the ST nodes as bases */
650  graph_voronoi(scip, graph, graph->cost, graph->cost, steinertree, vbase, vnoi);
651 
652  state = graph->path_state;
653 
654  /* initialize data structures */
655  for( k = 0; k < nnodes; k++ )
656  {
657  assert(state[k] == CONNECT || !graphmark[k]);
658 
659  pinned[k] = FALSE;
660  scanned[k] = FALSE;
661  nodesmark[k] = FALSE;
662 
663  /* initialize pairing heaps */
664  heapsize[k] = 0;
665  boundpaths[k] = NULL;
666 
667  lvledges_start[k] = NULL;
668 
669  if( !graphmark[k] )
670  continue;
671 
672  /* link all nodes to their (respective) voronoi base */
673  SCIP_CALL( SCIPallocBlockMemory(scip, &blists_curr) );
674  blists_curr->index = k;
675  blists_curr->parent = blists_start[vbase[k]];
676  blists_start[vbase[k]] = blists_curr;
677  }
678 
679  SCIP_CALL( SCIPallocBufferArray(scip, &supernodes, nstnodes) );
680  SCIP_CALL( SCIPallocBufferArray(scip, &kpnodes, nstnodes) );
681  SCIP_CALL( SCIPallocBufferArray(scip, &kpedges, nstnodes) );
682 
683  if( mwpc )
684  {
685  for( e = graph->outbeg[root]; e != EAT_LAST; e = graph->oeat[e] )
686  {
687  k = graph->head[e];
688  if( Is_term(graph->term[k]) )
689  {
690  graphmark[k] = FALSE;
691  for( l = graph->outbeg[k]; l != EAT_LAST; l = graph->oeat[l] )
692  {
693  if( !Is_term(graph->term[graph->head[l]]) )
694  {
695  assert(graph->head[l] != root);
696  pinned[graph->head[l]] = TRUE;
697  }
698  }
699  }
700  }
701  if( probtype != STP_RPCSPG )
702  graphmark[root] = FALSE;
703  }
704 
705  /* for each node, store all of its outgoing boundary-edges in a (respective) heap*/
706  for( e = 0; e < nedges; e += 2 )
707  {
708  if( graph->oeat[e] == EAT_FREE )
709  continue;
710 
711  node = graph->tail[e];
712  adjnode = graph->head[e];
713  newedges[e] = UNKNOWN;
714  newedges[e + 1] = UNKNOWN;
715 
716  /* is edge 'e' a boundary-edge? */
717  if( vbase[node] != vbase[adjnode] && graphmark[node] && graphmark[adjnode] )
718  {
719  edgecost = vnoi[node].dist + graph->cost[e] + vnoi[adjnode].dist;
720 
721  /* add the boundary-edge 'e' and its reversed to the corresponding heaps */
722  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[vbase[node]], e, edgecost, &(heapsize[vbase[node]])) );
723  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[vbase[adjnode]], flipedge(e), edgecost, &(heapsize[vbase[adjnode]])) );
724  }
725  }
726 
727  /* find LCAs for all edges */
728  SCIP_CALL( lca(scip, graph, root, &uf, nodesmark, best_result, lvledges_start, boundpaths, heapsize, vbase) );
729 
730  /* henceforth, the union-find structure will be used on the ST */
731  SCIPStpunionfindClear(scip, &uf, nnodes);
732 
733  /* henceforth, nodesmark will be used to mark the current supervertices (except for the one representing the root-component) */
734  for( i = 0; dfstree[i] != root; i++ )
735  nodesmark[dfstree[i]] = FALSE;
736  nodesmark[dfstree[i]] = FALSE;
737 
738  for( k = 0; k < nnodes; k++ )
739  assert(!nodesmark[k]);
740 
741  /* main loop visiting all nodes of the current ST in post-order */
742  for( i = 0; dfstree[i] != root; i++ )
743  {
744  crucnode = dfstree[i];
745  scanned[crucnode] = TRUE;
746 
747  SCIPdebugMessage("iteration %d (crucial node: %d) \n", i, crucnode);
748 
749  /* has the node been temporarily removed from the ST? */
750  if( !graphmark[crucnode] )
751  continue;
752 
753  /* is node 'crucnode' a removable crucial node? (i.e. not pinned or a terminal) */
754  if( !pinned[crucnode] && !Is_term(graph->term[crucnode]) && nodeIsCrucial(graph, best_result, crucnode) )
755  {
756  for( k = 0; k < nnodes; k++ )
757  assert(state[k] == CONNECT || !graphmark[k]);
758 
759  /* find all (unique) key-paths starting in node 'crucnode' */
760  k = UNKNOWN;
761  kpcost = 0.0;
762  nkpnodes = 0;
763  nkpedges = 0;
764  nsupernodes = 0;
765  for( edge = graph->outbeg[crucnode]; edge != EAT_LAST; edge = graph->oeat[edge] )
766  {
767  /* check whether the outgoing edge is in the ST */
768  if( (best_result[edge] > -1 && steinertree[graph->head[edge]]) || (best_result[flipedge(edge)] > -1 && steinertree[graph->tail[edge]]) )
769  {
770  kpcost += graph->cost[edge];
771 
772  /* check whether the current edge leads to the ST root*/
773  if( best_result[flipedge(edge)] > -1 )
774  {
775  k = flipedge(edge);
776  kpedges[nkpedges++] = k;
777  assert( edge == nodes[crucnode].edge );
778  }
779  else
780  {
781  kpedges[nkpedges++] = edge;
782  adjnode = graph->head[edge];
783  e = edge;
784 
785  /* move along the key-path until its end (i.e. a crucial or pinned node) is reached */
786  while( !pinned[adjnode] && !nodeIsCrucial(graph, best_result, adjnode) && steinertree[adjnode] )
787  {
788  /* update the union-find data structure */
789  SCIPStpunionfindUnion(&uf, crucnode, adjnode, FALSE);
790 
791  kpnodes[nkpnodes++] = adjnode;
792 
793  for( e = graph->outbeg[adjnode]; e != EAT_LAST; e = graph->oeat[e] )
794  {
795  if( best_result[e] > -1 )
796  {
797  kpcost += graph->cost[e];
798  kpedges[nkpedges++] = e;
799  break;
800  }
801  }
802 
803  /* assert that each leaf of the ST is a terminal */
804 
805 
806  if( e == EAT_LAST )
807  {
808  localmoves = 0;
809 
810  goto TERMINATE;
811  }
812  assert(e != EAT_LAST);
813  adjnode = graph->head[e];
814  }
815  /* does the last node on the path belong to a removed component? */
816  if( !steinertree[adjnode] )
817  {
818  kpcost -= graph->cost[e];
819  nkpedges--;
820  adjnode = graph->tail[e];
821  if( adjnode != crucnode )
822  {
823  supernodes[nsupernodes++] = adjnode;
824  nodesmark[adjnode] = TRUE;
825  }
826  }
827  else
828  {
829  supernodes[nsupernodes++] = adjnode;
830  nodesmark[adjnode] = TRUE;
831  }
832  }
833  }
834  }
835 
836  /* traverse the key-path leading to the root-component */
837  rootpathstart = nkpnodes;
838  if( k != -1 )
839  {
840  /* begin with the edge starting in the root-component of node 'crucnode' */
841  e = k;
842  adjnode = graph->tail[e];
843  while( !pinned[adjnode] && !nodeIsCrucial(graph, best_result, adjnode) && steinertree[adjnode] )
844  {
845  /* update the union-find data structure */
846  kpnodes[nkpnodes++] = adjnode;
847 
848  for( e = graph->inpbeg[adjnode]; e != EAT_LAST; e = graph->ieat[e] )
849  {
850  if( best_result[e] > -1 )
851  {
852  assert(steinertree[graph->tail[e]]);
853  kpcost += graph->cost[e];
854  kpedges[nkpedges++] = e;
855  break;
856  }
857  }
858 
859  assert( e != EAT_LAST );
860  adjnode = graph->tail[e];
861  }
862  supernodes[nsupernodes++] = adjnode;
863  }
864 
865  /* the last of the key-path nodes to be stored is the current key-node */
866  kpnodes[nkpnodes++] = crucnode;
867 
868  /* number of reset nodes */
869  nresnodes = 0;
870 
871  /* reset all nodes (referred to as 'C' henceforth) whose bases are internal nodes of the current key-paths */
872  for( k = 0; k < nkpnodes; k++ )
873  {
874  /* reset all nodes having the current (internal) keypath node as their voronoi base */
875  blists_curr = blists_start[kpnodes[k]];
876  while( blists_curr != NULL )
877  {
878  node = blists_curr->index;
879  assert(graphmark[node]);
880 
881  /* store all relevant data */
882  memvbase[nresnodes] = vbase[node];
883  memdist[nresnodes] = vnoi[node].dist;
884  meminedges[nresnodes] = vnoi[node].edge;
885  nresnodes++;
886 
887  /* reset data */
888  vbase[node] = UNKNOWN;
889  vnoi[node].dist = FARAWAY;
890  vnoi[node].edge = UNKNOWN;
891  state[node] = UNKNOWN;
892  blists_curr = blists_curr->parent;
893  }
894  }
895 
896  /* add vertical boundary-paths between the child components and the root-component (wrt node 'crucnode') */
897  nboundedges = 0;
898  for( k = 0; k < nsupernodes - 1; k++ )
899  {
900  l = supernodes[k];
901  edge = UNKNOWN;
902  while( boundpaths[l] != NULL )
903  {
904  SCIP_CALL( SCIPpairheapDeletemin(scip, &edge, &edgecost, &boundpaths[l], &heapsize[l]) );
905 
906  node = (vbase[graph->head[edge]] == UNKNOWN)? UNKNOWN : SCIPStpunionfindFind(&uf, vbase[graph->head[edge]]);
907  assert( (vbase[graph->tail[edge]] == UNKNOWN)? UNKNOWN : SCIPStpunionfindFind(&uf, vbase[graph->tail[edge]]) == l );
908 
909  /* check whether edge 'edge' represents a boundary-path having an endpoint in the kth-component and in the root-component respectively */
910  if( node != UNKNOWN && !nodesmark[node] && graphmark[node] )
911  {
912  boundedges[nboundedges++] = edge;
913  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[l], edge, edgecost, &heapsize[l]) );
914  break;
915  }
916  }
917  }
918 
919  /* add horizontal boundary-paths (between the child-components) */
920  lvledges_curr = lvledges_start[crucnode];
921  while( lvledges_curr != NULL )
922  {
923  edge = lvledges_curr->index;
924  k = vbase[graph->tail[edge]];
925  l = vbase[graph->head[edge]];
926  node = (l == UNKNOWN)? UNKNOWN : SCIPStpunionfindFind(&uf, l);
927  adjnode = (k == UNKNOWN)? UNKNOWN : SCIPStpunionfindFind(&uf, k);
928 
929  /* check whether the current boundary-path connects two child components */
930  if( node != UNKNOWN && nodesmark[node] && adjnode != UNKNOWN && nodesmark[adjnode] )
931  {
932  assert(graphmark[node]);
933  assert(graphmark[adjnode]);
934  boundedges[nboundedges++] = edge;
935  }
936  lvledges_curr = lvledges_curr->parent;
937  }
938 
939  /* try to connect the nodes of C (directly) to COMP(C), as a preprocessing for graph_voronoiRepair */
940  count = 0;
941  for( k = 0; k < nkpnodes; k++ )
942  {
943  blists_curr = blists_start[kpnodes[k]];
944  assert( blists_curr != NULL );
945  while( blists_curr != NULL )
946  {
947  node = blists_curr->index;
948 
949  /* iterate through all outgoing edges of 'node' */
950  for( edge = graph->inpbeg[node]; edge != EAT_LAST; edge = graph->ieat[edge] )
951  {
952  adjnode = graph->tail[edge];
953 
954  /* check whether the adjacent node is not in C and allows a better voronoi assignment of the current node */
955  if( state[adjnode] == CONNECT && SCIPisGT(scip, vnoi[node].dist, vnoi[adjnode].dist + graph->cost[edge])
956  && graphmark[vbase[adjnode]] && graphmark[adjnode] )
957  {
958  vnoi[node].dist = vnoi[adjnode].dist + graph->cost[edge];
959  vbase[node] = vbase[adjnode];
960  vnoi[node].edge = edge;
961  }
962  }
963  if( vbase[node] != UNKNOWN )
964  {
965  heap_add(graph->path_heap, state, &count, node, vnoi);
966  }
967  blists_curr = blists_curr->parent;
968  }
969  }
970 
971  /* if there are no key-path nodes, something has gone wrong */
972  assert(nkpnodes != 0);
973 
974  graph_voronoiRepairMult(scip, graph, graph->cost, &count, vbase, boundedges, &nboundedges, nodesmark, &uf, vnoi);
975 
976  /* create a supergraph, having the endpoints of the key-paths incident to the current crucial node as (super-) vertices */
977  SCIP_CALL( graph_init(scip, &supergraph, nsupernodes, nboundedges * 2, 1) );
978  supergraph->stp_type = STP_SPG;
979 
980  /* add vertices to the supergraph */
981  for( k = 0; k < nsupernodes; k++ )
982  {
983  supernodesid[supernodes[k]] = k;
984  graph_knot_add(supergraph, graph->term[supernodes[k]]);
985  }
986 
987  /* the (super-) vertex representing the current root-component of the ST */
988  k = supernodes[nsupernodes - 1];
989 
990  /* add edges to the supergraph */
991  for( l = 0; l < nboundedges; l++ )
992  {
993  edge = boundedges[l];
994  node = SCIPStpunionfindFind(&uf, vbase[graph->tail[edge]]);
995  adjnode = SCIPStpunionfindFind(&uf, vbase[graph->head[edge]]);
996 
997  /* if node 'node' or 'adjnode' belongs to the root-component, take the (temporary) root-component identifier instead */
998  node = ((nodesmark[node])? node : k);
999  adjnode = ((nodesmark[adjnode])? adjnode : k);
1000 
1001  /* compute the cost of the boundary-path pertaining to the boundary-edge 'edge' */
1002  edgecost = vnoi[graph->tail[edge]].dist + graph->cost[edge] + vnoi[graph->head[edge]].dist;
1003  graph_edge_add(scip, supergraph, supernodesid[node], supernodesid[adjnode], edgecost, edgecost);
1004  }
1005 
1006  /* compute a MST on the supergraph */
1007  SCIP_CALL( SCIPallocBufferArray(scip, &mst, nsupernodes) );
1008  SCIP_CALL( graph_path_init(scip, supergraph) );
1009  graph_path_exec(scip, supergraph, MST_MODE, nsupernodes - 1, supergraph->cost, mst);
1010 
1011  /* compute the cost of the MST */
1012  mstcost = 0.0;
1013 
1014  /* compute the cost of the MST */
1015  for( l = 0; l < nsupernodes - 1; l++ )
1016  {
1017  /* compute the edge in the original graph corresponding to the current MST edge */
1018  if( mst[l].edge % 2 == 0 )
1019  edge = boundedges[mst[l].edge / 2 ];
1020  else
1021  edge = flipedge(boundedges[mst[l].edge / 2 ]);
1022 
1023  mstcost += graph->cost[edge];
1024  assert( newedges[edge] != crucnode && newedges[flipedge(edge)] != crucnode );
1025 
1026  /* mark the edge (in the original graph) as visited */
1027  newedges[edge] = crucnode;
1028 
1029  /* traverse along the boundary-path belonging to the boundary-edge 'edge' */
1030  for( node = graph->tail[edge]; node != vbase[node]; node = graph->tail[vnoi[node].edge] )
1031  {
1032  e = vnoi[node].edge;
1033 
1034  /* if edge 'e' and its reversed have not been visited yet */
1035  if( newedges[e] != crucnode && newedges[flipedge(e)] != crucnode )
1036  {
1037  newedges[e] = crucnode;
1038  mstcost += graph->cost[e];
1039  }
1040  }
1041  for( node = graph->head[edge]; node != vbase[node]; node = graph->tail[vnoi[node].edge] )
1042  {
1043  e = flipedge(vnoi[node].edge);
1044 
1045  /* if edge 'e' and its reversed have not been visited yet */
1046  if( newedges[vnoi[node].edge] != crucnode && newedges[e] != crucnode )
1047  {
1048  newedges[e] = crucnode;
1049  mstcost += graph->cost[e];
1050  }
1051  }
1052  }
1053 
1054  if( SCIPisLT(scip, mstcost, kpcost) )
1055  {
1056  localmoves++;
1057  SCIPdebugMessage("found improving solution in KEY VERTEX ELIMINATION (round: %d) \n ", nruns);
1058 
1059  /* unmark the original edges spanning the supergraph */
1060  for( e = 0; e < nkpedges; e++ )
1061  {
1062  assert(best_result[kpedges[e]] != -1);
1063  best_result[kpedges[e]] = -1;
1064  }
1065 
1066  /* mark all ST nodes except for those belonging to the root-component as forbidden */
1067  for( k = rootpathstart; k < nkpnodes; k++ )
1068  {
1069  graphmark[kpnodes[k]] = FALSE;
1070  steinertree[kpnodes[k]] = FALSE;
1071  }
1072 
1073  for( k = 0; k < i; k++ )
1074  {
1075  node = SCIPStpunionfindFind(&uf, dfstree[k]);
1076  if( nodesmark[node] || node == crucnode )
1077  {
1078  graphmark[dfstree[k]] = FALSE;
1079  steinertree[dfstree[k]] = FALSE;
1080  }
1081  }
1082 
1083  /* add the new edges reconnecting the (super-) components */
1084  for( l = 0; l < nsupernodes - 1; l++ )
1085  {
1086  if( mst[l].edge % 2 == 0 )
1087  edge = boundedges[mst[l].edge / 2 ];
1088  else
1089  edge = flipedge(boundedges[mst[l].edge / 2 ]);
1090 
1091  /* change the orientation within the target-component if necessary */
1092  if( !nodesmark[vbase[graph->head[edge]]] )
1093  {
1094  node = vbase[graph->head[edge]];
1095  k = SCIPStpunionfindFind(&uf, node);
1096  assert(nodesmark[k]);
1097  while( node != k )
1098  {
1099  /* the ST edge pointing towards the root */
1100  e = nodes[node].edge;
1101 
1102  assert(best_result[e] == -1 && best_result[flipedge(e)] != -1 );
1103  best_result[e] = CONNECT;
1104  best_result[flipedge(e)] = UNKNOWN;
1105  node = graph->head[e];
1106  }
1107  }
1108 
1109  /* is the vbase of the current boundary-edge tail in the root-component? */
1110  if( !nodesmark[SCIPStpunionfindFind(&uf, vbase[graph->tail[edge]])] )
1111  {
1112 
1113  best_result[edge] = CONNECT;
1114 
1115  for( node = graph->tail[edge], adjnode = graph->head[edge]; node != vbase[node]; adjnode = node, node = graph->tail[vnoi[node].edge] )
1116  {
1117  graphmark[node] = FALSE;
1118 
1119  if( best_result[flipedge(vnoi[node].edge)] == CONNECT )
1120  best_result[flipedge(vnoi[node].edge)] = UNKNOWN;
1121 
1122  best_result[vnoi[node].edge] = CONNECT;
1123  }
1124 
1125  assert(!nodesmark[node] && vbase[node] == node);
1126  assert( graphmark[node] == TRUE );
1127 
1128  /* is the pinned node its own component identifier? */
1129  if( !Is_term(graph->term[node]) && scanned[node] && !pinned[node] && SCIPStpunionfindFind(&uf, node) == node )
1130  {
1131  graphmark[graph->head[edge]] = FALSE;
1132  oldedge = edge;
1133 
1134  for( edge = graph->outbeg[node]; edge != EAT_LAST; edge = graph->oeat[edge] )
1135  {
1136  adjnode = graph->head[edge];
1137  /* check whether edge 'edge' leads to an ancestor of terminal 'node' */
1138  if( best_result[edge] == CONNECT && graphmark[adjnode] && steinertree[adjnode] && SCIPStpunionfindFind(&uf, adjnode) != node )
1139  {
1140 
1141  assert(scanned[adjnode]);
1142  /* meld the heaps */
1143  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[adjnode], &heapsize[node], &heapsize[adjnode]);
1144 
1145  /* update the union-find data structure */
1146  SCIPStpunionfindUnion(&uf, node, adjnode, FALSE);
1147 
1148  /* move along the key-path until its end (i.e. until a crucial node is reached) */
1149  while( !nodeIsCrucial(graph, best_result, adjnode) && !pinned[adjnode] )
1150  {
1151  for( e = graph->outbeg[adjnode]; e != EAT_LAST; e = graph->oeat[e] )
1152  {
1153  if( best_result[e] != -1 )
1154  break;
1155  }
1156 
1157  /* assert that each leaf of the ST is a terminal */
1158  assert( e != EAT_LAST );
1159  adjnode = graph->head[e];
1160  if( !steinertree[adjnode] )
1161  break;
1162  assert(scanned[adjnode]);
1163  assert(SCIPStpunionfindFind(&uf, adjnode) != node);
1164 
1165  /* update the union-find data structure */
1166  SCIPStpunionfindUnion(&uf, node, adjnode, FALSE);
1167 
1168  /* meld the heaps */
1169  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[adjnode], &heapsize[node], &heapsize[adjnode]);
1170  }
1171  }
1172  }
1173  edge = oldedge;
1174  }
1175 
1176  /* mark the start node (lying in the root-component of the ST) of the current boundary-path as pinned,
1177  * so that it may not be removed later on */
1178  pinned[node] = TRUE;
1179 
1180  for( node = graph->head[edge]; node != vbase[node]; node = graph->tail[vnoi[node].edge] )
1181  {
1182  graphmark[node] = FALSE;
1183  if( best_result[vnoi[node].edge] == CONNECT )
1184  best_result[vnoi[node].edge] = -1;
1185 
1186  best_result[flipedge(vnoi[node].edge)] = CONNECT;
1187 
1188  }
1189  }
1190  else
1191  {
1192 
1193  best_result[edge] = CONNECT;
1194 
1195  for( node = graph->tail[edge]; node != vbase[node]; node = graph->tail[vnoi[node].edge] )
1196  {
1197  graphmark[node] = FALSE;
1198  if( best_result[vnoi[node].edge] != CONNECT && best_result[flipedge(vnoi[node].edge)] != CONNECT )
1199  best_result[vnoi[node].edge] = CONNECT;
1200 
1201  }
1202 
1203  for( node = graph->head[edge]; node != vbase[node]; node = graph->tail[vnoi[node].edge] )
1204  {
1205  graphmark[node] = FALSE;
1206 
1207  best_result[flipedge(vnoi[node].edge)] = CONNECT;
1208  best_result[vnoi[node].edge] = UNKNOWN;
1209  }
1210  }
1211  }
1212 
1213  for( k = 0; k < nkpnodes; k++ )
1214  {
1215  assert(graphmark[kpnodes[k]] == FALSE);
1216  assert(steinertree[kpnodes[k]] == FALSE);
1217  }
1218  assert(!graphmark[crucnode]);
1219  }
1220  else
1221  {
1222  /* no improving solution has been found during the move */
1223 
1224  /* meld the heap pertaining to 'crucnode' and all heaps pertaining to descendant key-paths of node 'crucnode' */
1225  for( k = 0; k < rootpathstart; k++ )
1226  {
1227  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[kpnodes[k]], &heapsize[crucnode], &heapsize[kpnodes[k]]);
1228  }
1229  for( k = 0; k < nsupernodes - 1; k++ )
1230  {
1231  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[supernodes[k]], &heapsize[crucnode], &heapsize[supernodes[k]]);
1232 
1233  /* update the union-find data structure */
1234  SCIPStpunionfindUnion(&uf, crucnode, supernodes[k], FALSE);
1235  }
1236  }
1237 
1238  /* free the supergraph and the MST data structure */
1239  graph_path_exit(scip, supergraph);
1240  graph_free(scip, &supergraph, TRUE);
1241  SCIPfreeBufferArray(scip, &mst);
1242 
1243  /* unmark the descendant supervertices */
1244  for( k = 0; k < nsupernodes - 1; k++ )
1245  {
1246  nodesmark[supernodes[k]] = FALSE;
1247  }
1248 
1249  /* debug test; to be deleted later on */
1250  for( k = 0; k < nnodes; k++ )
1251  {
1252  assert( !nodesmark[k] );
1253  }
1254 
1255  /* restore the original voronoi diagram */
1256  l = 0;
1257  for( k = 0; k < nkpnodes; k++ )
1258  {
1259  /* restore data of all nodes having the current (internal) key-path node as their voronoi base */
1260  blists_curr = blists_start[kpnodes[k]];
1261  while( blists_curr != NULL )
1262  {
1263  node = blists_curr->index;
1264  vbase[node] = memvbase[l];
1265  vnoi[node].dist = memdist[l];
1266  vnoi[node].edge = meminedges[l];
1267  l++;
1268  blists_curr = blists_curr->parent;
1269  }
1270  }
1271 
1272  assert(l == nresnodes);
1273  }
1274 
1275  /** Key-Path Exchange */
1276  if( probtype != STP_MWCSP )
1277  {
1278  /* if the node has just been eliminated, skip Key-Path Exchange */
1279  if( !graphmark[crucnode] )
1280  continue;
1281 
1282  /* is crucnode a crucial or pinned vertex? */
1283  if( (!nodeIsCrucial(graph, best_result, crucnode) && !pinned[crucnode]) )
1284  continue;
1285 
1286  if( Is_term(graph->term[crucnode]) || pinned[crucnode] )
1287  {
1288  for( edge = graph->outbeg[crucnode]; edge != EAT_LAST; edge = graph->oeat[edge] )
1289  {
1290  adjnode = graph->head[edge];
1291  /* check whether edge 'edge' leads to an ancestor of terminal 'crucnode' */
1292  if( best_result[edge] == CONNECT && steinertree[adjnode] && graphmark[adjnode] )
1293  {
1294  assert( SCIPStpunionfindFind(&uf, adjnode) != crucnode);
1295  assert(scanned[adjnode]);
1296  /* meld the heaps */
1297  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[adjnode], &heapsize[crucnode], &heapsize[adjnode]);
1298 
1299  /* update the union-find data structure */
1300  SCIPStpunionfindUnion(&uf, crucnode, adjnode, FALSE);
1301 
1302  /* move along the key-path until its end (i.e. until a crucial node is reached) */
1303  while( !nodeIsCrucial(graph, best_result, adjnode) && !pinned[adjnode] )
1304  {
1305  for( e = graph->outbeg[adjnode]; e != EAT_LAST; e = graph->oeat[e] )
1306  {
1307  if( best_result[e] != -1 )
1308  break;
1309  }
1310 
1311  /* assert that each leaf of the ST is a terminal */
1312  assert( e != EAT_LAST );
1313  adjnode = graph->head[e];
1314  if( !steinertree[adjnode] || !graphmark[adjnode] )
1315  break;
1316  assert(scanned[adjnode]);
1317  assert(SCIPStpunionfindFind(&uf, adjnode) != crucnode);
1318 
1319  /* update the union-find data structure */
1320  SCIPStpunionfindUnion(&uf, crucnode, adjnode, FALSE);
1321 
1322  /* meld the heaps */
1323  SCIPpairheapMeldheaps(scip, &boundpaths[crucnode], &boundpaths[adjnode], &heapsize[crucnode], &heapsize[adjnode]);
1324  }
1325  }
1326  }
1327 
1328  }
1329 
1330  /* counts the internal nodes of the keypath */
1331  nkpnodes = 0;
1332 
1333  for( k = 0; k < nnodes; k++ )
1334  assert(state[k] == CONNECT || !graphmark[k]);
1335 
1336  /* find the (unique) key-path containing the parent of the current crucial node 'crucnode' */
1337  kptailnode = graph->head[nodes[crucnode].edge];
1338  kpathcost = graph->cost[nodes[crucnode].edge];
1339 
1340  while( !nodeIsCrucial(graph, best_result, kptailnode) && !pinned[kptailnode] )
1341  {
1342  kpathcost += graph->cost[nodes[kptailnode].edge];
1343 
1344  kpnodes[nkpnodes++] = kptailnode;
1345  kptailnode = graph->head[nodes[kptailnode].edge];
1346  }
1347 
1348  /* counts the reset nodes during voronoi repair */
1349  nresnodes = 0;
1350 
1351  /* reset all nodes (henceforth referred to as 'C') whose bases are internal nodes of the current keypath */
1352  for( k = 0; k < nkpnodes; k++ )
1353  {
1354  /* reset all nodes having the current (internal) keypath node as their voronoi base */
1355  blists_curr = blists_start[kpnodes[k]];
1356  while( blists_curr != NULL )
1357  {
1358  node = blists_curr->index;
1359  memvbase[nresnodes] = vbase[node];
1360  memdist[nresnodes] = vnoi[node].dist;
1361  meminedges[nresnodes] = vnoi[node].edge;
1362  nresnodes++;
1363  vbase[node] = UNKNOWN;
1364  vnoi[node].dist = FARAWAY;
1365  vnoi[node].edge = UNKNOWN;
1366  state[node] = UNKNOWN;
1367  blists_curr = blists_curr->parent;
1368  }
1369  }
1370 
1371  edgecost = UNKNOWN;
1372  e = UNKNOWN;
1373  while( boundpaths[crucnode] != NULL )
1374  {
1375  SCIP_CALL( SCIPpairheapDeletemin(scip, &e, &edgecost, &boundpaths[crucnode], &(heapsize[crucnode])) );
1376  assert( e != UNKNOWN );
1377  k = vbase[graph->tail[e]];
1378  l = vbase[graph->head[e]];
1379 
1380  assert(graphmark[k]);
1381  node = (l == UNKNOWN || !graphmark[l] )? UNKNOWN : SCIPStpunionfindFind(&uf, l);
1382 
1383  /* does the boundary-path end in the root component? */
1384  if( node != UNKNOWN && node != crucnode && graphmark[l] )
1385  {
1386  SCIP_CALL( SCIPpairheapInsert(scip, &boundpaths[crucnode], e, edgecost, &(heapsize[crucnode])) );
1387  break;
1388  }
1389  }
1390 
1391  if( boundpaths[crucnode] == NULL )
1392  {
1393  oldedge = UNKNOWN;
1394  }
1395  else
1396  {
1397  oldedge = e;
1398  }
1399 
1400  /* counts the nodes connected during the following 'preprocessing' */
1401  count = 0;
1402 
1403  /* try to connect the nodes of C (directly) to COMP(C), as a preprocessing for voronoi-repair */
1404  for( k = 0; k < nkpnodes; k++ )
1405  {
1406  blists_curr = blists_start[kpnodes[k]];
1407  assert( blists_curr != NULL );
1408  while( blists_curr != NULL )
1409  {
1410  node = blists_curr->index;
1411 
1412  /* iterate through all outgoing edges of 'node' */
1413  for( edge = graph->inpbeg[node]; edge != EAT_LAST; edge = graph->ieat[edge] )
1414  {
1415  adjnode = graph->tail[edge];
1416 
1417  /* check whether the adjacent node is not in C and allows a better voronoi assignment of the current node */
1418  if( state[adjnode] == CONNECT && SCIPisGT(scip, vnoi[node].dist, vnoi[adjnode].dist + graph->cost[edge])
1419  && graphmark[vbase[adjnode]] && graphmark[adjnode] )
1420  {
1421  vnoi[node].dist = vnoi[adjnode].dist + graph->cost[edge];
1422  vbase[node] = vbase[adjnode];
1423  vnoi[node].edge = edge;
1424  }
1425  }
1426  if( vbase[node] != UNKNOWN )
1427  {
1428  heap_add(graph->path_heap, state, &count, node, vnoi);
1429  }
1430  blists_curr = blists_curr->parent;
1431  }
1432  }
1433  if( nkpnodes > 0 )
1434  assert(count > 0);
1435  newedge = UNKNOWN;
1436 
1437  /* if there is no key path, nothing has to be repaired */
1438  if( nkpnodes > 0 )
1439  graph_voronoiRepair(scip, graph, graph->cost, &count, vbase, vnoi, &newedge, crucnode, &uf);
1440  else
1441  newedge = nodes[crucnode].edge;
1442 
1443  if( oldedge != UNKNOWN && newedge != UNKNOWN && SCIPisLT(scip, edgecost,
1444  vnoi[graph->tail[newedge]].dist + graph->cost[newedge] + vnoi[graph->head[newedge]].dist) )
1445  newedge = oldedge;
1446  if( oldedge != UNKNOWN && newedge == UNKNOWN )
1447  newedge = oldedge;
1448 
1449  assert( newedge != UNKNOWN );
1450  edgecost = vnoi[graph->tail[newedge]].dist + graph->cost[newedge] + vnoi[graph->head[newedge]].dist;
1451  if( SCIPisLT(scip, edgecost, kpathcost) )
1452  {
1453  node = SCIPStpunionfindFind(&uf, vbase[graph->head[newedge]]);
1454 
1455  SCIPdebugMessage( "ADDING NEW KEY PATH (%f )\n", edgecost - kpathcost );
1456 
1457  localmoves++;
1458 
1459  /* remove old keypath */
1460  assert( best_result[flipedge(nodes[crucnode].edge)] != UNKNOWN );
1461  best_result[flipedge(nodes[crucnode].edge)] = UNKNOWN;
1462  steinertree[crucnode] = FALSE;
1463  graphmark[crucnode] = FALSE;
1464 
1465  for( k = 0; k < nkpnodes; k++ )
1466  {
1467  assert( best_result[flipedge(nodes[kpnodes[k]].edge)] != UNKNOWN );
1468  best_result[flipedge(nodes[kpnodes[k]].edge)] = UNKNOWN;
1469  steinertree[kpnodes[k]] = FALSE;
1470  graphmark[kpnodes[k]] = FALSE;
1471  }
1472  assert(graphmark[kptailnode]);
1473 
1474  if( node == crucnode )
1475  newedge = flipedge(newedge);
1476 
1477  for( node = graph->tail[newedge]; node != vbase[node]; node = graph->tail[vnoi[node].edge] )
1478  {
1479  graphmark[node] = FALSE;
1480 
1481  best_result[flipedge(vnoi[node].edge)] = CONNECT;
1482  best_result[vnoi[node].edge] = UNKNOWN;
1483  }
1484 
1485  for( node = graph->head[newedge]; node != vbase[node]; node = graph->tail[vnoi[node].edge] )
1486  {
1487  graphmark[node] = FALSE;
1488 
1489  best_result[vnoi[node].edge] = CONNECT;
1490  }
1491 
1492  best_result[flipedge(newedge)] = CONNECT;
1493 
1494  newpathend = vbase[graph->tail[newedge]];
1495  assert(node == vbase[graph->head[newedge]] );
1496 
1497  /* flip all edges on the ST path between the endnode of the new key-path and the current crucial node */
1498  k = newpathend;
1499  assert(SCIPStpunionfindFind(&uf, newpathend) == crucnode);
1500 
1501  while( k != crucnode )
1502  {
1503  assert(graphmark[k]);
1504  assert( best_result[flipedge(nodes[k].edge)] != -1);
1505  best_result[flipedge(nodes[k].edge)] = UNKNOWN;
1506 
1507  best_result[nodes[k].edge] = CONNECT;
1508 
1509  k = graph->head[nodes[k].edge];
1510  }
1511 
1512  for( k = 0; k < i; k++ )
1513  {
1514  if( crucnode == SCIPStpunionfindFind(&uf, dfstree[k]) )
1515  {
1516  graphmark[dfstree[k]] = FALSE;
1517  steinertree[dfstree[k]] = FALSE;
1518  }
1519  }
1520 
1521  /* update union find */
1522  if( !Is_term(graph->term[node]) && scanned[node] && !pinned[node] && SCIPStpunionfindFind(&uf, node) == node )
1523  {
1524  for( edge = graph->outbeg[node]; edge != EAT_LAST; edge = graph->oeat[edge] )
1525  {
1526  adjnode = graph->head[edge];
1527  /* check whether edge 'edge' leads to an ancestor of terminal 'node' */
1528  if( best_result[edge] == CONNECT && steinertree[adjnode] && graphmark[adjnode] && SCIPStpunionfindFind(&uf, adjnode) != node )
1529  {
1530  assert(scanned[adjnode]);
1531  /* meld the heaps */
1532  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[adjnode], &heapsize[node], &heapsize[adjnode]);
1533 
1534  /* update the union-find data structure */
1535  SCIPStpunionfindUnion(&uf, node, adjnode, FALSE);
1536 
1537  /* move along the key-path until its end (i.e. until a crucial node is reached) */
1538  while( !nodeIsCrucial(graph, best_result, adjnode) && !pinned[adjnode] )
1539  {
1540  for( e = graph->outbeg[adjnode]; e != EAT_LAST; e = graph->oeat[e] )
1541  {
1542  if( best_result[e] != -1 )
1543  break;
1544  }
1545 
1546  /* assert that each leaf of the ST is a terminal */
1547  assert( e != EAT_LAST );
1548  adjnode = graph->head[e];
1549  if( !steinertree[adjnode] )
1550  break;
1551  assert(scanned[adjnode]);
1552  assert(SCIPStpunionfindFind(&uf, adjnode) != node);
1553 
1554  /* update the union-find data structure */
1555  SCIPStpunionfindUnion(&uf, node, adjnode, FALSE);
1556 
1557  /* meld the heaps */
1558  SCIPpairheapMeldheaps(scip, &boundpaths[node], &boundpaths[adjnode], &heapsize[node], &heapsize[adjnode]);
1559  }
1560  }
1561  }
1562 
1563  }
1564  pinned[node] = TRUE;
1565  }
1566 
1567  /* restore the original voronoi digram */
1568  l = 0;
1569  for( k = 0; k < nkpnodes; k++ )
1570  {
1571  /* reset all nodes having the current (internal) keypath node as their voronoi base */
1572  blists_curr = blists_start[kpnodes[k]];
1573  while( blists_curr != NULL )
1574  {
1575  node = blists_curr->index;
1576  vbase[node] = memvbase[l];
1577  vnoi[node].dist = memdist[l];
1578  vnoi[node].edge = meminedges[l];
1579  l++;
1580  blists_curr = blists_curr->parent;
1581  }
1582  }
1583  assert(l == nresnodes);
1584  }
1585  }
1586 
1587 
1588  /**********************************************************/
1589 
1590  TERMINATE:
1591 
1592  SCIPStpunionfindClear(scip, &uf, nnodes);
1593 
1594  /* free data structures */
1595  SCIPfreeBufferArray(scip, &kpedges);
1596  SCIPfreeBufferArray(scip, &kpnodes);
1597  SCIPfreeBufferArray(scip, &supernodes);
1598 
1599  for( k = nnodes - 1; k >= 0; k-- )
1600  {
1601  if( boundpaths[k] != NULL )
1602  SCIPpairheapFree(scip, &boundpaths[k]);
1603 
1604  blists_curr = blists_start[k];
1605  lvledges_curr = lvledges_start[k];
1606  while( lvledges_curr != NULL )
1607  {
1608  lvledges_start[k] = lvledges_curr->parent;
1609  SCIPfreeBlockMemory(scip, &lvledges_curr);
1610  lvledges_curr = lvledges_start[k];
1611  }
1612 
1613  while( blists_curr != NULL )
1614  {
1615  blists_start[k] = blists_curr->parent;
1616  SCIPfreeBlockMemory(scip, &blists_curr);
1617  blists_curr = blists_start[k];
1618  }
1619  }
1620 
1621  /* has there been a move during this run? */
1622  if( localmoves > 0 )
1623  {
1624  for( i = 0; i < nnodes; i++ )
1625  {
1626  steinertree[i] = FALSE;
1627  graphmark[i] = (graph->grad[i] > 0);
1628  SCIPlinkcuttreeInit(&nodes[i]);
1629  }
1630 
1631  graphmark[root] = TRUE;
1632 
1633  /* create a link-cut tree representing the current Steiner tree */
1634  for( e = 0; e < nedges; e++ )
1635  {
1636  assert(graph->head[e] == graph->tail[flipedge(e)]);
1637 
1638  /* if edge e is in the tree, so are its incident vertices */
1639  if( best_result[e] != -1 )
1640  {
1641  steinertree[graph->tail[e]] = TRUE;
1642  steinertree[graph->head[e]] = TRUE;
1643  SCIPlinkcuttreeLink(&nodes[graph->head[e]], &nodes[graph->tail[e]], flipedge(e));
1644  }
1645  }
1646  assert( nodes[root].edge == -1 );
1647  nodes[root].edge = -1;
1648  }
1649  }
1650 
1651  /* free data structures */
1652  SCIPStpunionfindFree(scip, &uf);
1653  SCIPfreeBufferArray(scip, &nodesmark);
1654  SCIPfreeBufferArray(scip, &dfstree);
1655  SCIPfreeBufferArray(scip, &pinned);
1656  SCIPfreeBufferArray(scip, &boundpaths);
1657  SCIPfreeBufferArray(scip, &meminedges);
1658  SCIPfreeBufferArray(scip, &memdist);
1659  SCIPfreeBufferArray(scip, &memvbase);
1660  SCIPfreeBufferArray(scip, &blists_start);
1661  SCIPfreeBufferArray(scip, &heapsize);
1662  SCIPfreeBufferArray(scip, &scanned);
1663  SCIPfreeBufferArray(scip, &supernodesid);
1664  SCIPfreeBufferArray(scip, &boundedges);
1665  SCIPfreeBufferArray(scip, &lvledges_start);
1666  SCIPfreeBufferArray(scip, &newedges);
1667  /******/
1668  }
1669 
1670 #ifdef SCIP_DEBUG
1671  {
1672  SCIP_Real obj = 0.0;
1673  for( e = 0; e < nedges; e++ )
1674  obj += (best_result[e] > -1) ? graph->cost[e] : 0.0;
1675 
1676  printf(" ObjAfterHeurLocal=%.12e\n", obj);
1677  }
1678 #endif
1679 
1680 #ifndef NDEBUG
1681  assert(SCIPisLE(scip, graph_sol_getObj(graph->cost, best_result, 0.0, nedges), initialobj));
1682 #endif
1683 
1684  SCIPfreeBufferArray(scip, &steinertree);
1685  SCIPfreeBufferArray(scip, &nodes);
1686  SCIPfreeBufferArray(scip, &vbase);
1687  SCIPfreeBufferArray(scip, &vnoi);
1688 
1689  return SCIP_OKAY;
1690 }
1691 
1692 /** Greedy Extension local heuristic for (R)PC and MW */
1694  SCIP* scip, /**< SCIP data structure */
1695  GRAPH* graph, /**< graph data structure */
1696  const SCIP_Real* cost, /**< edge cost array*/
1697  PATH* path, /**< shortest data structure array */
1698  int* stedge, /**< initialized array to indicate whether an edge is part of the Steiner tree */
1699  STP_Bool* stvertex, /**< uninitialized array to indicate whether a vertex is part of the Steiner tree */
1700  SCIP_Bool* extensions /**< pointer to store whether extensions have been made */
1701  )
1702 {
1704  int candidatesup[MAX(GREEDY_EXTENSIONS, GREEDY_EXTENSIONS_MW)];
1705 
1706  PATH* orgpath;
1707  SCIP_PQUEUE* pqueue;
1708  SCIP_Real bestsolval;
1709  int i;
1710  int root;
1711  int nedges;
1712  int nnodes;
1713  int nextensions;
1714  int greedyextensions;
1715  STP_Bool* stvertextmp;
1716 
1717  assert(scip != NULL);
1718  assert(path != NULL);
1719  assert(graph != NULL);
1720  assert(stedge != NULL);
1721  assert(cost != NULL);
1722  assert(stvertex != NULL);
1723 
1724  root = graph->source;
1725  nnodes = graph->knots;
1726  nedges = graph->edges;
1727 
1728  graph_pc_2transcheck(graph);
1729  SCIP_CALL( SCIPallocBufferArray(scip, &stvertextmp, nnodes) );
1730  SCIP_CALL( SCIPallocBufferArray(scip, &orgpath, nnodes) );
1731 
1732  /* initialize solution vertex array with FALSE */
1733  BMSclearMemoryArray(stvertex, nnodes);
1734 
1735  stvertex[graph->source] = TRUE;
1736  path[graph->source].edge = UNKNOWN;
1737 
1738  for( int e = 0; e < nedges; e++ )
1739  if( stedge[e] == CONNECT )
1740  {
1741  i = graph->head[e];
1742  path[i].edge = e;
1743  stvertex[i] = TRUE;
1744  }
1745 
1746  graph_path_st_pcmw_extend(scip, graph, cost, path, stvertex, extensions);
1747 
1748  BMScopyMemoryArray(orgpath, path, nnodes);
1749 
1750  if( graph->stp_type == STP_MWCSP )
1751  greedyextensions = GREEDY_EXTENSIONS_MW;
1752  else
1753  greedyextensions = GREEDY_EXTENSIONS;
1754 
1755  /*** compute solution value and save greedyextensions many best unconnected nodes ***/
1756 
1757  SCIP_CALL( SCIPpqueueCreate(&pqueue, greedyextensions, 2.0, GNODECmpByDist) );
1758 
1759  bestsolval = 0.0;
1760  nextensions = 0;
1761  for( i = 0; i < nnodes; i++ )
1762  {
1763  if( graph->grad[i] == 0 )
1764  continue;
1765 
1766  if( Is_term(graph->term[i]) )
1767  {
1768  if( i != root )
1769  {
1770  int e;
1771  SCIP_Bool connect = FALSE;
1772  SCIP_Real ecost = -1.0;
1773  for( e = graph->inpbeg[i]; e != EAT_LAST; e = graph->ieat[e] )
1774  {
1775  if( root == graph->tail[e] )
1776  ecost = graph->cost[e];
1777  else if( stvertex[graph->tail[e]] )
1778  connect = TRUE;
1779  }
1780 
1781  if( !connect )
1782  bestsolval += ecost;
1783  }
1784  }
1785  else if( stvertex[i] )
1786  {
1787  bestsolval += graph->cost[orgpath[i].edge];
1788  }
1789  else if( orgpath[i].edge != UNKNOWN && Is_pterm(graph->term[i]) )
1790  {
1791  if( nextensions < greedyextensions )
1792  {
1793  candidates[nextensions].dist = graph->prize[i] - path[i].dist;
1794  candidates[nextensions].number = i;
1795 
1796  SCIP_CALL( SCIPpqueueInsert(pqueue, &(candidates[nextensions++])) );
1797  }
1798  else
1799  {
1800  /* get candidate vertex of minimum value */
1801  GNODE* min = (GNODE*) SCIPpqueueFirst(pqueue);
1802  if( SCIPisLT(scip, min->dist, graph->prize[i] - path[i].dist) )
1803  {
1804  min = (GNODE*) SCIPpqueueRemove(pqueue);
1805  min->dist = graph->prize[i] - path[i].dist;
1806  min->number = i;
1807  SCIP_CALL( SCIPpqueueInsert(pqueue, min) );
1808  }
1809  }
1810  }
1811  }
1812 
1813  for( int restartcount = 0; restartcount < GREEDY_MAXRESTARTS; restartcount++ )
1814  {
1815  int l = 0;
1816  SCIP_Bool extensionstmp = FALSE;
1817 
1818  i = nextensions;
1819 
1820  /* write extension candidates into array, from max to min */
1821  while( SCIPpqueueNElems(pqueue) > 0 )
1822  {
1823  GNODE* min = (GNODE*) SCIPpqueueRemove(pqueue);
1824  assert(i > 0);
1825  candidatesup[--i] = min->number;
1826  }
1827  assert(i == 0);
1828 
1829  /* iteratively insert new subpaths and try to improve solution */
1830  for( ; l < nextensions; l++ )
1831  {
1832  i = candidatesup[l];
1833  if( !stvertex[i] )
1834  {
1835  SCIP_Real newsolval = 0.0;
1836  int k = i;
1837 
1838  BMScopyMemoryArray(stvertextmp, stvertex, nnodes);
1839  BMScopyMemoryArray(path, orgpath, nnodes);
1840 
1841  /* add new extension */
1842  while( !stvertextmp[k] )
1843  {
1844  stvertextmp[k] = TRUE;
1845  assert(orgpath[k].edge != UNKNOWN);
1846  k = graph->tail[orgpath[k].edge];
1847  }
1848 
1849  graph_path_st_pcmw_extend(scip, graph, cost, path, stvertextmp, &extensionstmp);
1850 
1851  for( int j = 0; j < nnodes; j++ )
1852  {
1853  if( graph->grad[j] == 0 )
1854  continue;
1855 
1856  if( Is_term(graph->term[j]) )
1857  {
1858  if( j != root )
1859  {
1860  int e;
1861  SCIP_Bool connect = FALSE;
1862  SCIP_Real ecost = -1.0;
1863  for( e = graph->inpbeg[j]; e != EAT_LAST; e = graph->ieat[e] )
1864  {
1865  if( root == graph->tail[e] )
1866  ecost = graph->cost[e];
1867  else if( stvertextmp[graph->tail[e]] )
1868  connect = TRUE;
1869  }
1870 
1871  if( !connect )
1872  newsolval += ecost;
1873  }
1874  }
1875  else if( stvertextmp[j] )
1876  {
1877  newsolval += graph->cost[path[j].edge];
1878  }
1879  }
1880 
1881  /* new solution value better than old one? */
1882  if( SCIPisLT(scip, newsolval, bestsolval) )
1883  {
1884  *extensions = TRUE;
1885 
1886  bestsolval = newsolval;
1887  BMScopyMemoryArray(stvertex, stvertextmp, nnodes);
1888 
1889  /* save greedyextensions many best unconnected nodes */
1890  nextensions = 0;
1891 
1892  for( int j = 0; j < nnodes; j++ )
1893  {
1894  orgpath[j].edge = path[j].edge;
1895  orgpath[j].dist = path[j].dist;
1896  if( !stvertex[j] && Is_pterm(graph->term[j]) && path[j].edge != UNKNOWN )
1897  {
1898  if( nextensions < greedyextensions )
1899  {
1900  candidates[nextensions].dist = graph->prize[j] - path[j].dist;
1901  candidates[nextensions].number = j;
1902 
1903  SCIP_CALL( SCIPpqueueInsert(pqueue, &(candidates[nextensions++])) );
1904  }
1905  else
1906  {
1907  /* get candidate vertex of minimum value */
1908  GNODE* min = (GNODE*) SCIPpqueueFirst(pqueue);
1909  if( SCIPisLT(scip, min->dist, graph->prize[j] - path[j].dist) )
1910  {
1911  min = (GNODE*) SCIPpqueueRemove(pqueue);
1912  min->dist = graph->prize[j] - path[j].dist;
1913  min->number = j;
1914  SCIP_CALL( SCIPpqueueInsert(pqueue, min) );
1915  }
1916  }
1917  }
1918  }
1919 
1920  break;
1921  } /* if new solution value better than old one? */
1922  } /* if !stvertex[i] */
1923  } /* for l < nextension */
1924  /* no more extensions performed? */
1925  if( l == nextensions )
1926  break;
1927  } /* main loop */
1928 
1929  /* have vertices been added? */
1930  if( *extensions )
1931  {
1932  SCIPdebugMessage("SCIPStpHeurLocalExtendPcMw found extensions \n");
1933  for( int e = 0; e < nedges; e++ )
1934  stedge[e] = UNKNOWN;
1935  SCIP_CALL( SCIPStpHeurTMPrunePc(scip, graph, cost, stedge, stvertex) );
1936  }
1937 
1938  SCIPpqueueFree(&pqueue);
1939  SCIPfreeBufferArray(scip, &orgpath);
1940  SCIPfreeBufferArray(scip, &stvertextmp);
1941 
1942 #ifdef SCIP_DEBUG
1943  {
1944  SCIP_Real t = 0.0;
1945  for( int e = 0; e < nedges; e++ )
1946  if( stedge[e] == CONNECT )
1947  t += graph->cost[e];
1948 
1949  SCIPdebugMessage("SCIPStpHeurLocalExtendPcMw: exit real cost %f \n", t);
1950  }
1951 #endif
1952 
1953  return SCIP_OKAY;
1954 }
1955 
1956 
1957 
1958 /*
1959  * Callback methods of primal heuristic
1960  */
1961 
1962 /** copy method for primal heuristic plugins (called when SCIP copies plugins) */
1963 static
1964 SCIP_DECL_HEURCOPY(heurCopyLocal)
1965 { /*lint --e{715}*/
1966  assert(scip != NULL);
1967  assert(heur != NULL);
1968  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
1969 
1970  /* call inclusion method of primal heuristic */
1972 
1973  return SCIP_OKAY;
1974 }
1975 
1976 /** destructor of primal heuristic to free user data (called when SCIP is exiting) */
1977 static
1978 SCIP_DECL_HEURFREE(heurFreeLocal)
1979 { /*lint --e{715}*/
1980  SCIP_HEURDATA* heurdata;
1981 
1982  assert(heur != NULL);
1983  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
1984  assert(scip != NULL);
1985 
1986  /* free heuristic data */
1987  heurdata = SCIPheurGetData(heur);
1988  assert(heurdata != NULL);
1989  SCIPfreeMemory(scip, &heurdata);
1990  SCIPheurSetData(heur, NULL);
1991 
1992  return SCIP_OKAY;
1993 }
1994 
1995 /** solving process initialization method of primal heuristic (called when branch and bound process is about to begin) */
1996 static
1997 SCIP_DECL_HEURINITSOL(heurInitsolLocal)
1998 { /*lint --e{715}*/
1999  SCIP_HEURDATA* heurdata;
2000 
2001  assert(heur != NULL);
2002  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
2003  assert(scip != NULL);
2004 
2005  /* free heuristic data */
2006  heurdata = SCIPheurGetData(heur);
2007 
2008  heurdata->nfails = 1;
2009  heurdata->nbestsols = DEFAULT_NBESTSOLS;
2010 
2011  SCIP_CALL( SCIPallocMemoryArray(scip, &(heurdata->lastsolindices), heurdata->maxnsols) );
2012 
2013  for( int i = 0; i < heurdata->maxnsols; i++ )
2014  heurdata->lastsolindices[i] = -1;
2015 
2016  return SCIP_OKAY;
2017 }
2018 
2019 
2020 /** solving process deinitialization method of primal heuristic (called before branch and bound process data is freed) */
2021 static
2022 SCIP_DECL_HEUREXITSOL(heurExitsolLocal)
2023 { /*lint --e{715}*/
2024  SCIP_HEURDATA* heurdata;
2025 
2026  assert(heur != NULL);
2027  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
2028  assert(scip != NULL);
2029 
2030  /* free heuristic data */
2031  heurdata = SCIPheurGetData(heur);
2032  assert(heurdata != NULL);
2033  assert(heurdata->lastsolindices != NULL);
2034  SCIPfreeMemoryArray(scip, &(heurdata->lastsolindices));
2035 
2036  return SCIP_OKAY;
2037 }
2038 
2039 
2040 
2041 /** execution method of primal heuristic */
2042 static
2043 SCIP_DECL_HEUREXEC(heurExecLocal)
2044 { /*lint --e{715}*/
2045  SCIP_HEURDATA* heurdata;
2046  SCIP_PROBDATA* probdata;
2047  GRAPH* graph; /* graph structure */
2048  SCIP_SOL* newsol; /* new solution */
2049  SCIP_SOL* impsol; /* new improved solution */
2050  SCIP_SOL** sols; /* solutions */
2051  SCIP_VAR** vars; /* SCIP variables */
2052  SCIP_Real pobj;
2053  SCIP_Real* nval;
2054  SCIP_Real* xval;
2055  int i;
2056  int e;
2057  int v;
2058  int min;
2059  int root;
2060  int nvars;
2061  int nsols; /* number of all solutions found so far */
2062  int nedges;
2063  int* results;
2064  int* lastsolindices;
2065  SCIP_Bool feasible;
2066 
2067  assert(heur != NULL);
2068  assert(scip != NULL);
2069  assert(strcmp(SCIPheurGetName(heur), HEUR_NAME) == 0);
2070  assert(result != NULL);
2071 
2072  /* get heuristic data */
2073  heurdata = SCIPheurGetData(heur);
2074  assert(heurdata != NULL);
2075  lastsolindices = heurdata->lastsolindices;
2076  assert(lastsolindices != NULL);
2077 
2078  probdata = SCIPgetProbData(scip);
2079  assert(probdata != NULL);
2080 
2081  graph = SCIPprobdataGetGraph(probdata);
2082  assert(graph != NULL);
2083 
2084  *result = SCIP_DIDNOTRUN;
2085 
2086  /* the local heuristics may not work correctly for several problem variants*/
2087  if( graph->stp_type != STP_SPG && graph->stp_type != STP_RSMT && graph->stp_type != STP_OARSMT &&
2088  graph->stp_type != STP_PCSPG && graph->stp_type != STP_RPCSPG && graph->stp_type != STP_GSTP
2089  && graph->stp_type != STP_MWCSP )
2090  return SCIP_OKAY;
2091 
2092  /* don't run local in a Subscip */
2093  if( SCIPgetSubscipDepth(scip) > 0 )
2094  return SCIP_OKAY;
2095 
2096  /* no solution available? */
2097  if( SCIPgetBestSol(scip) == NULL )
2098  return SCIP_OKAY;
2099 
2100  root = graph->source;
2101  sols = SCIPgetSols(scip);
2102  nsols = SCIPgetNSols(scip);
2103  nedges = graph->edges;
2104 
2105  assert(heurdata->maxnsols >= 0);
2106 
2107  min = MIN(heurdata->maxnsols, nsols);
2108 
2109  /* only process each solution once */
2110  for( v = 0; v < min; v++ )
2111  {
2112  if( SCIPsolGetIndex(sols[v]) != lastsolindices[v] )
2113  {
2114  /* shift all solution indices right of the new solution index */
2115  for( i = min - 1; i >= v + 1; i-- )
2116  lastsolindices[i] = lastsolindices[i - 1];
2117  break;
2118  }
2119  }
2120 
2121  /* no new solution available? */
2122  if( v == min )
2123  return SCIP_OKAY;
2124 
2125  newsol = sols[v];
2126  lastsolindices[v] = SCIPsolGetIndex(newsol);
2127 
2128  /* solution not good enough? */
2129  if( (v > heurdata->nbestsols && !(heurdata->maxfreq)) && graph->stp_type != STP_MWCSP )
2130  return SCIP_OKAY;
2131 
2132  /* has the new solution been found by this very heuristic? */
2133  if( SCIPsolGetHeur(newsol) == heur )
2134  return SCIP_OKAY;
2135 
2136  *result = SCIP_DIDNOTFIND;
2137 
2138  vars = SCIPprobdataGetVars(scip);
2139  nvars = SCIPprobdataGetNVars(scip);
2140  xval = SCIPprobdataGetXval(scip, newsol);
2141 
2142  if( vars == NULL )
2143  return SCIP_OKAY;
2144 
2145  assert(vars != NULL);
2146  assert(xval != NULL);
2147 
2148  /* for PC variants: test whether solution is trivial */
2149  if( graph->stp_type == STP_PCSPG || graph->stp_type == STP_RPCSPG || graph->stp_type == STP_MWCSP )
2150  {
2151  for( e = graph->outbeg[root]; e != EAT_LAST; e = graph->oeat[e] )
2152  if( !Is_term(graph->term[graph->head[e]]) && SCIPisEQ(scip, xval[e], 1.0) )
2153  break;
2154  if( e == EAT_LAST )
2155  return SCIP_OKAY;
2156  }
2157 
2158  /* allocate memory */
2159  SCIP_CALL( SCIPallocBufferArray(scip, &results, nedges) );
2160  SCIP_CALL( SCIPallocBufferArray(scip, &nval, nvars) );
2161 
2162  /* set solution array */
2163  for( e = 0; e < nedges; e++ )
2164  {
2165  if( SCIPisEQ(scip, xval[e], 1.0) )
2166  results[e] = CONNECT;
2167  else
2168  results[e] = UNKNOWN;
2169  }
2170 
2171  if( !graph_sol_valid(scip, graph, results) )
2172  {
2173  SCIPfreeBufferArray(scip, &nval);
2174  SCIPfreeBufferArray(scip, &results);
2175  return SCIP_OKAY;
2176  }
2177 
2178  /* pruning necessary? */
2179  if( SCIPsolGetHeur(newsol) == NULL ||
2180  !(strcmp(SCIPheurGetName(SCIPsolGetHeur(newsol)), "rec") == 0 ||
2181  strcmp(SCIPheurGetName(SCIPsolGetHeur(newsol)), "TM") == 0) )
2182  {
2183  const int nnodes = graph->knots;
2184  STP_Bool* steinertree;
2185  SCIP_CALL( SCIPallocBufferArray(scip, &steinertree, nnodes) );
2186  assert(graph_sol_valid(scip, graph, results));
2187 
2188  for( v = nnodes - 1; v >= 0; --v )
2189  steinertree[v] = FALSE;
2190 
2191  for( e = nedges - 1; e >= 0; --e )
2192  {
2193  if( results[e] == CONNECT )
2194  {
2195  steinertree[graph->tail[e]] = TRUE;
2196  steinertree[graph->head[e]] = TRUE;
2197  }
2198  results[e] = UNKNOWN;
2199  }
2200 
2201  if( graph->stp_type == STP_PCSPG || graph->stp_type == STP_RPCSPG || graph->stp_type == STP_MWCSP )
2202  SCIP_CALL( SCIPStpHeurTMPrunePc(scip, graph, graph->cost, results, steinertree) );
2203  else
2204  SCIP_CALL( SCIPStpHeurTMPrune(scip, graph, graph->cost, 0, results, steinertree) );
2205 
2206  SCIPfreeBufferArray(scip, &steinertree);
2207  }
2208 
2209  /* execute local heuristics */
2210  SCIP_CALL( SCIPStpHeurLocalRun(scip, graph, graph->cost, results) );
2211 
2212  /* can we connect the network */
2213  for( v = 0; v < nvars; v++ )
2214  nval[v] = (results[v % nedges] == (v / nedges)) ? 1.0 : 0.0;
2215 
2216  SCIP_CALL( SCIPStpValidateSol(scip, graph, nval, &feasible) );
2217 
2218  /* solution feasible? */
2219  if( feasible )
2220  {
2221  assert(nedges == nvars);
2222 
2223  pobj = 0.0;
2224 
2225  for( v = 0; v < nedges; v++ )
2226  pobj += graph->cost[v] * nval[v];
2227 
2228  /* has solution been improved? */
2229  if( SCIPisGT(scip, SCIPgetSolOrigObj(scip, newsol) - SCIPprobdataGetOffset(scip), pobj) )
2230  {
2231  SCIP_SOL* bestsol;
2232  SCIP_Bool success;
2233 
2234  bestsol = sols[0];
2235  impsol = NULL;
2236  SCIP_CALL( SCIPprobdataAddNewSol(scip, nval, impsol, heur, &success) );
2237 
2238  if( success )
2239  {
2240  *result = SCIP_FOUNDSOL;
2241 
2242  if( heurdata->nbestsols < heurdata->maxnsols && SCIPisGT(scip, SCIPgetSolOrigObj(scip, bestsol) - SCIPprobdataGetOffset(scip), pobj) )
2243  {
2244  heurdata->nfails = 0;
2245  heurdata->nbestsols++;
2246  }
2247  SCIPdebugMessage("success in local: old: %f new: %f \n", (SCIPgetSolOrigObj(scip, bestsol) - SCIPprobdataGetOffset(scip)), pobj);
2248  }
2249  }
2250  }
2251 
2252  if( *result != SCIP_FOUNDSOL )
2253  {
2254  heurdata->nfails++;
2255  if( heurdata->nbestsols > DEFAULT_MINNBESTSOLS && heurdata->nfails > 1 && graph->stp_type != STP_MWCSP )
2256  heurdata->nbestsols--;
2257 
2258  SCIPdebugMessage("fail! %d \n", heurdata->nbestsols);
2259  }
2260 
2261  SCIPfreeBufferArray(scip, &nval);
2262  SCIPfreeBufferArray(scip, &results);
2263 
2264  return SCIP_OKAY;
2265 }
2266 
2267 /*
2268  * primal heuristic specific interface methods
2269  */
2270 
2271 
2272 /** creates the local primal heuristic and includes it in SCIP */
2274  SCIP* scip /**< SCIP data structure */
2275  )
2276 {
2277  SCIP_HEURDATA* heurdata;
2278  SCIP_HEUR* heur;
2279 
2280  /* create Local primal heuristic data */
2281  SCIP_CALL( SCIPallocMemory(scip, &heurdata) );
2282 
2283  /* include primal heuristic */
2284  SCIP_CALL( SCIPincludeHeurBasic(scip, &heur,
2286  HEUR_MAXDEPTH, HEUR_TIMING, HEUR_USESSUBSCIP, heurExecLocal, heurdata) );
2287 
2288  assert(heur != NULL);
2289 
2290  /* set non-NULL pointers to callback methods */
2291  SCIP_CALL( SCIPsetHeurCopy(scip, heur, heurCopyLocal) );
2292  SCIP_CALL( SCIPsetHeurFree(scip, heur, heurFreeLocal) );
2293  SCIP_CALL( SCIPsetHeurInitsol(scip, heur, heurInitsolLocal) );
2294  SCIP_CALL( SCIPsetHeurExitsol(scip, heur, heurExitsolLocal) );
2295 
2296  /* add local primal heuristic parameters */
2297  SCIP_CALL( SCIPaddBoolParam(scip, "stp/duringroot",
2298  "should the heuristic be called during the root node?",
2299  &heurdata->duringroot, FALSE, DEFAULT_DURINGROOT, NULL, NULL) );
2300 
2301  SCIP_CALL( SCIPaddBoolParam(scip, "heuristics/"HEUR_NAME"/maxfreq",
2302  "should the heuristic be executed at maximum frequeny?",
2303  &heurdata->maxfreq, FALSE, DEFAULT_MAXFREQLOC, NULL, NULL) );
2304 
2305  SCIP_CALL( SCIPaddIntParam(scip, "heuristics/"HEUR_NAME"/maxnsols",
2306  "maximum number of best solutions to improve",
2307  &heurdata->maxnsols, FALSE, DEFAULT_MAXNBESTSOLS, 1, 50, NULL, NULL) );
2308 
2309  return SCIP_OKAY;
2310 }
SCIP_RETCODE graph_path_init(SCIP *, GRAPH *)
Definition: grphpath.c:444
#define HEUR_USESSUBSCIP
Definition: heur_local.c:50
SCIP_Bool SCIPisEQ(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
void * SCIPpqueueRemove(SCIP_PQUEUE *pqueue)
Definition: misc.c:1307
#define HEUR_DESC
Definition: heur_local.c:42
SCIP_Bool SCIPisPositive(SCIP *scip, SCIP_Real val)
#define NULL
Definition: def.h:253
void graph_voronoi(SCIP *scip, const GRAPH *, SCIP_Real *, SCIP_Real *, STP_Bool *, int *, PATH *)
Definition: grphpath.c:1751
SCIP_RETCODE SCIPStpunionfindInit(SCIP *scip, UF *uf, int length)
Definition: misc_stp.c:689
void SCIPpqueueFree(SCIP_PQUEUE **pqueue)
Definition: misc.c:1259
SCIP_Real SCIPlinkcuttreeFindMinChain(SCIP *scip, const SCIP_Real *nodeweight, const int *head, const int *stdeg, const NODE *start, NODE **first, NODE **last)
Definition: misc_stp.c:258
const char * SCIPheurGetName(SCIP_HEUR *heur)
Definition: heur.c:1254
int *RESTRICT head
Definition: grph.h:96
Definition: grph.h:57
int source
Definition: grph.h:67
void * SCIPpqueueFirst(SCIP_PQUEUE *pqueue)
Definition: misc.c:1348
SCIP_RETCODE SCIPsetHeurExitsol(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEUREXITSOL((*heurexitsol)))
Definition: scip_heur.c:232
#define STP_GSTP
Definition: grph.h:49
SCIP_HEURDATA * SCIPheurGetData(SCIP_HEUR *heur)
Definition: heur.c:1165
#define SCIPfreeMemoryArray(scip, ptr)
Definition: scip_mem.h:69
signed int edge
Definition: grph.h:146
#define MST_MODE
Definition: grph.h:162
SCIP_Bool graph_valid(const GRAPH *)
Definition: grphbase.c:4324
int terms
Definition: grph.h:64
SCIP_VAR ** SCIPprobdataGetVars(SCIP *scip)
struct ST_Node * parent
Definition: misc_stp.h:60
#define SCIPallocMemoryArray(scip, ptr, num)
Definition: scip_mem.h:53
void SCIPpairheapMeldheaps(SCIP *scip, PHNODE **root1, PHNODE **root2, int *sizeroot1, int *sizeroot2)
Definition: misc_stp.c:599
#define GREEDY_EXTENSIONS
Definition: heur_local.c:61
#define EAT_LAST
Definition: grph.h:31
SCIP_RETCODE SCIPStpIncludeHeurLocal(SCIP *scip)
Definition: heur_local.c:2273
#define HEUR_PRIORITY
Definition: heur_local.c:44
static SCIP_DECL_HEURCOPY(heurCopyLocal)
Definition: heur_local.c:1964
#define FALSE
Definition: def.h:73
int *RESTRICT inpbeg
Definition: grph.h:74
SCIP_EXPORT SCIP_HEUR * SCIPsolGetHeur(SCIP_SOL *sol)
Definition: sol.c:2553
int *RESTRICT path_state
Definition: grph.h:119
void graph_free(SCIP *, GRAPH **, SCIP_Bool)
Definition: grphbase.c:3674
Problem data for stp problem.
#define TRUE
Definition: def.h:72
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:53
void graph_path_exit(SCIP *, GRAPH *)
Definition: grphpath.c:466
int SCIPpqueueNElems(SCIP_PQUEUE *pqueue)
Definition: misc.c:1362
#define HEUR_MAXDEPTH
Definition: heur_local.c:47
SCIP_PROBDATA * SCIPgetProbData(SCIP *scip)
Definition: scip_prob.c:963
int SCIPprobdataGetNVars(SCIP *scip)
struct SCIP_HeurData SCIP_HEURDATA
Definition: type_heur.h:51
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:47
SCIP_RETCODE SCIPpairheapDeletemin(SCIP *scip, int *element, SCIP_Real *key, PHNODE **root, int *size)
Definition: misc_stp.c:562
#define SCIPfreeBlockMemory(scip, ptr)
Definition: scip_mem.h:95
#define STP_PCSPG
Definition: grph.h:40
#define SCIPdebugMessage
Definition: pub_message.h:77
void graph_voronoiRepair(SCIP *, const GRAPH *, SCIP_Real *, int *, int *, PATH *, int *, int, UF *)
Definition: grphpath.c:2979
void SCIPpairheapFree(SCIP *scip, PHNODE **root)
Definition: misc_stp.c:627
#define SCIPfreeBufferArray(scip, ptr)
Definition: scip_mem.h:123
int number
Definition: misc_stp.h:43
#define SCIPallocBlockMemory(scip, ptr)
Definition: scip_mem.h:78
NODE * SCIPlinkcuttreeFindMax(SCIP *scip, const SCIP_Real *cost, NODE *v)
Definition: misc_stp.c:327
GRAPH * SCIPprobdataGetGraph(SCIP_PROBDATA *probdata)
int *RESTRICT mark
Definition: grph.h:70
void graph_voronoiRepairMult(SCIP *, const GRAPH *, SCIP_Real *, int *, int *, int *, int *, STP_Bool *, UF *, PATH *)
Definition: grphpath.c:3057
int SCIPgetNSols(SCIP *scip)
Definition: scip_sol.c:2205
SCIP_Bool SCIPisLE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
#define HEUR_FREQ
Definition: heur_local.c:45
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:107
#define DEFAULT_DURINGROOT
Definition: heur_local.c:52
int *RESTRICT oeat
Definition: grph.h:103
#define CONNECT
Definition: grph.h:154
int stp_type
Definition: grph.h:127
void graph_path_st_pcmw_extend(SCIP *, const GRAPH *, const SCIP_Real *, PATH *, STP_Bool *, SCIP_Bool *)
Definition: grphpath.c:1410
void heap_add(int *, int *, int *, int, PATH *)
Definition: grphpath.c:244
#define Is_pterm(a)
Definition: grph.h:169
static SCIP_DECL_HEURINITSOL(heurInitsolLocal)
Definition: heur_local.c:1997
unsigned char STP_Bool
Definition: grph.h:52
SCIP_RETCODE SCIPprobdataAddNewSol(SCIP *scip, SCIP_Real *nval, SCIP_SOL *sol, SCIP_HEUR *heur, SCIP_Bool *success)
SCIP_Real SCIPprobdataGetOffset(SCIP *scip)
SCIP_Real dist
Definition: misc_stp.h:44
SCIP_Real * prize
Definition: grph.h:82
int edge
Definition: misc_stp.h:59
#define GREEDY_EXTENSIONS_MW
Definition: heur_local.c:60
SCIP_Real dist
Definition: grph.h:145
int *RESTRICT grad
Definition: grph.h:73
#define HEUR_TIMING
Definition: heur_local.c:48
#define DEFAULT_MAXFREQLOC
Definition: heur_local.c:53
SCIP_Bool graph_sol_valid(SCIP *, const GRAPH *, const int *)
Definition: grphbase.c:3066
void SCIPheurSetData(SCIP_HEUR *heur, SCIP_HEURDATA *heurdata)
Definition: heur.c:1175
void SCIPStpunionfindClear(SCIP *scip, UF *uf, int length)
Definition: misc_stp.c:708
void graph_path_exec(SCIP *, const GRAPH *, const int, int, const SCIP_Real *, PATH *)
Definition: grphpath.c:491
int knots
Definition: grph.h:62
#define SCIP_CALL(x)
Definition: def.h:365
void SCIPStpunionfindUnion(UF *uf, int p, int q, SCIP_Bool compress)
Definition: misc_stp.c:752
int SCIPStpunionfindFind(UF *uf, int element)
Definition: misc_stp.c:728
SCIP_RETCODE SCIPStpHeurLocalExtendPcMw(SCIP *scip, GRAPH *graph, const SCIP_Real *cost, PATH *path, int *stedge, STP_Bool *stvertex, SCIP_Bool *extensions)
Definition: heur_local.c:1693
void SCIPlinkcuttreeInit(NODE *v)
Definition: misc_stp.c:227
Improvement heuristic for Steiner problems.
void SCIPlinkcuttreeEvert(NODE *v)
Definition: misc_stp.c:352
static SCIP_DECL_HEURFREE(heurFreeLocal)
Definition: heur_local.c:1978
#define FARAWAY
Definition: grph.h:156
SCIP_RETCODE SCIPpqueueCreate(SCIP_PQUEUE **pqueue, int initsize, SCIP_Real sizefac, SCIP_DECL_SORTPTRCOMP((*ptrcomp)))
Definition: misc.c:1234
int GNODECmpByDist(void *first_arg, void *second_arg)
#define STP_SPG
Definition: grph.h:38
#define SCIPallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:111
void SCIPlinkcuttreeLink(NODE *v, NODE *w, int edge)
Definition: misc_stp.c:236
#define SCIP_Bool
Definition: def.h:70
int *RESTRICT ieat
Definition: grph.h:102
int *RESTRICT path_heap
Definition: grph.h:118
#define STP_MWCSP
Definition: grph.h:47
int *RESTRICT tail
Definition: grph.h:95
#define DEFAULT_MINNBESTSOLS
Definition: heur_local.c:56
#define MIN(x, y)
Definition: def.h:223
SCIP_EXPORT int SCIPsolGetIndex(SCIP_SOL *sol)
Definition: sol.c:2584
SCIP_Bool SCIPisNegative(SCIP *scip, SCIP_Real val)
static SCIP_DECL_HEUREXEC(heurExecLocal)
Definition: heur_local.c:2043
SCIP_RETCODE SCIPStpValidateSol(SCIP *, const GRAPH *, const double *, SCIP_Bool *)
Definition: validate.c:136
int *RESTRICT term
Definition: grph.h:68
SCIP_Real graph_sol_getObj(const SCIP_Real *, const int *, SCIP_Real, int)
Definition: grphbase.c:3196
#define BMScopyMemoryArray(ptr, source, num)
Definition: memory.h:124
static SCIP_DECL_HEUREXITSOL(heurExitsolLocal)
Definition: heur_local.c:2022
#define SCIPfreeMemory(scip, ptr)
Definition: scip_mem.h:67
SCIP_RETCODE SCIPpairheapInsert(SCIP *scip, PHNODE **root, int element, SCIP_Real key, int *size)
Definition: misc_stp.c:528
static STP_Bool nodeIsCrucial(const GRAPH *graph, int *steineredges, int node)
Definition: heur_local.c:175
#define Is_term(a)
Definition: grph.h:168
#define MAX(x, y)
Definition: def.h:222
#define EAT_FREE
Definition: grph.h:30
struct SCIP_ProbData SCIP_PROBDATA
Definition: type_prob.h:44
SCIP_SOL ** SCIPgetSols(SCIP *scip)
Definition: scip_sol.c:2254
SCIP_Real * cost
Definition: grph.h:94
SCIP_RETCODE SCIPsetHeurCopy(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEURCOPY((*heurcopy)))
Definition: scip_heur.c:152
void SCIPlinkcuttreeCut(NODE *v)
Definition: misc_stp.c:249
#define GREEDY_MAXRESTARTS
Definition: heur_local.c:59
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:73
#define HEUR_NAME
Definition: heur_local.c:41
void SCIPStpunionfindFree(SCIP *scip, UF *uf)
Definition: misc_stp.c:795
static void dfsorder(const GRAPH *graph, const int *edges, const int *node, int *counter, int *dfst)
Definition: heur_local.c:86
#define SCIP_Real
Definition: def.h:164
SCIP_Bool SCIPisLT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Bool SCIPisGT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Real * SCIPprobdataGetXval(SCIP *scip, SCIP_SOL *sol)
int *RESTRICT outbeg
Definition: grph.h:76
SCIP_RETCODE SCIPStpHeurTMPrunePc(SCIP *scip, const GRAPH *g, const SCIP_Real *cost, int *result, STP_Bool *connected)
Definition: heur_tm.c:167
#define DEFAULT_NBESTSOLS
Definition: heur_local.c:55
SCIP_RETCODE SCIPStpHeurLocalRun(SCIP *scip, GRAPH *graph, const SCIP_Real *cost, int *best_result)
Definition: heur_local.c:208
shortest paths based primal heuristics for Steiner problems
SCIP_RETCODE SCIPsetHeurInitsol(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEURINITSOL((*heurinitsol)))
Definition: scip_heur.c:216
int edges
Definition: grph.h:92
#define flipedge(edge)
Definition: grph.h:150
void graph_knot_add(GRAPH *, int)
Definition: grphbase.c:2196
SCIP_RETCODE SCIPsetHeurFree(SCIP *scip, SCIP_HEUR *heur, SCIP_DECL_HEURFREE((*heurfree)))
Definition: scip_heur.c:168
#define SCIPallocMemory(scip, ptr)
Definition: scip_mem.h:51
#define UNKNOWN
Definition: sepa_mcf.c:4081
#define STP_RSMT
Definition: grph.h:45
#define nnodes
Definition: gastrans.c:65
#define STP_OARSMT
Definition: grph.h:46
SCIP_RETCODE SCIPpairheapBuffarr(SCIP *scip, PHNODE *root, int size, int **elements)
Definition: misc_stp.c:670
#define BMSclearMemoryArray(ptr, num)
Definition: memory.h:120
static SCIP_RETCODE lca(SCIP *scip, const GRAPH *graph, int u, UF *uf, STP_Bool *nodesmark, int *steineredges, IDX **lcalists, PHNODE **boundpaths, int *heapsize, int *vbase)
Definition: heur_local.c:113
struct Int_List_Node * parent
Definition: misc_stp.h:76
#define HEUR_FREQOFS
Definition: heur_local.c:46
#define HEUR_DISPCHAR
Definition: heur_local.c:43
SCIP_RETCODE SCIPStpHeurTMPrune(SCIP *scip, const GRAPH *g, const SCIP_Real *cost, int layer, int *result, STP_Bool *connected)
Definition: heur_tm.c:555
void graph_edge_add(SCIP *, GRAPH *, int, int, double, double)
#define STP_RPCSPG
Definition: grph.h:41
#define DEFAULT_MAXNBESTSOLS
Definition: heur_local.c:54
SCIP_SOL * SCIPgetBestSol(SCIP *scip)
Definition: scip_sol.c:2304
int SCIPgetSubscipDepth(SCIP *scip)
Definition: scip_copy.c:2289
SCIP_RETCODE graph_init(SCIP *, GRAPH **, int, int, int)
Definition: grphbase.c:3491
void graph_pc_2transcheck(GRAPH *)
Definition: grphbase.c:1041
SCIP_Real SCIPgetSolOrigObj(SCIP *scip, SCIP_SOL *sol)
Definition: scip_sol.c:1435
SCIP_RETCODE SCIPpqueueInsert(SCIP_PQUEUE *pqueue, void *elem)
Definition: misc.c:1280