Scippy

SCIP

Solving Constraint Integer Programs

dualascent.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-2022 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 dualascent.c
17  * @brief Dual-ascent dual heuristic for Steiner problems
18  * @author Daniel Rehfeldt
19  *
20  * This file includes dual-ascent for classic Steiner tree and some variants.
21  *
22  */
23 
24 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
25 
26 #include <assert.h>
27 #include "scip/cons_linear.h"
28 #include "scip/cons_logicor.h"
29 #include "scip/cons_setppc.h"
30 #include "dualascent.h"
31 #include "probdata_stp.h"
32 #include "graph.h"
33 #include "heur_ascendprune.h"
34 #include "scip/scip.h"
35 #include "scip/misc.h"
36 
37 
38 
39 #define DEFAULT_DAMAXDEVIATION 0.25 /**< max deviation for dual ascent */
40 #define DA_MAXDEVIATION_LOWER 0.01 /**< lower bound for max deviation for dual ascent */
41 #define DA_MAXDEVIATION_UPPER 0.9 /**< upper bound for max deviation for dual ascent */
42 #define DA_EPS (5e-7)
43 #define DAPATHS_HITLIMIT_PCMW 20
44 #define DAPATHS_HITLIMIT 5
45 
46 /* do depth-first search */
47 #define DFS
48 
49 
50 #ifdef BITFIELDSARRAY
51 #define ARRLENGTH 32
52 #define SetBit(Arr, pos) ( Arr[(pos/ARRLENGTH)] |= (1 << (pos%ARRLENGTH)) )
53 #define CleanBit(Arr, pos) ( Arr[(pos/ARRLENGTH)] &= ~(1 << (pos%ARRLENGTH)) )
54 #define BitTrue(Arr, pos) ( Arr[(pos/ARRLENGTH)] & (1 << (pos%ARRLENGTH)) )
55 #endif
56 
58 
59 
60 /** internal data for path based dual-ascent */
61 typedef struct dual_ascent_paths
62 {
63  DIJK* dijklimited; /**< Dijkstra data */
64  int* startnodes; /**< start nodes */
65  int* nodes_hits; /**< counts how often a node has been hit */
66  SCIP_Bool* nodes_abort; /**< nodes to abort at */
67  SCIP_Real maxdist; /**< current max distance */
68  int nstartnodes; /**< number of */
69  int centernode; /**< node for PC/MW */
70  SCIP_Bool isUnrootedPcMw; /**< un-rooted PC or MW? */
71 } DAPATHS;
72 
73 
74 /**@name Local methods
75  *
76  * @{
77  */
78 
79 /** creates empty constraint */
80 static
82  SCIP* scip, /**< SCIP data structure */
83  enum DACONS_Type constype, /**< constraint type to be used */
84  SCIP_Bool consUseInital, /**< use dual-ascent cuts as initial constraints? */
85  SCIP_CONS** cons /**< to be initialized */
86  )
87 {
88  // todo PcMw did not propagate before, keep it?
89 
90  if( constype == dacons_linear )
91  {
92  SCIP_CALL( SCIPcreateConsLinear(scip, cons, "da", 0, NULL, NULL, 1.0, SCIPinfinity(scip),
93  consUseInital, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
94  }
95  else if( constype == dacons_logicor )
96  {
97  SCIP_CALL( SCIPcreateConsLogicor(scip, cons, "da", 0, NULL,
98  consUseInital, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
99  }
100  else
101  {
102  SCIP_CALL( SCIPcreateConsSetcover(scip, cons, "da", 0, NULL,
103  consUseInital, TRUE, TRUE, TRUE, TRUE, FALSE, FALSE, TRUE, TRUE, FALSE) );
104  }
105 
106  return SCIP_OKAY;
107 }
108 
109 
110 /** gets parameter */
111 static
113  SCIP* scip, /**< SCIP data structure */
114  enum DACONS_Type* constype, /**< pointer: constraint type to be used (OUT) */
115  SCIP_Bool* consUseInital /**< pointer: use dual-ascent cuts as initial constraints? (OUT) */
116  )
117 {
118  int type;
119  SCIP_CALL( SCIPgetBoolParam(scip, "stp/usedacutsinitial", consUseInital) );
120  SCIP_CALL( SCIPgetIntParam(scip, "stp/dacutstype", &type) );
121 
122  if( type == 0 )
123  {
124  *constype = dacons_linear;
125  }
126  else if( type == 1 )
127  {
128  *constype = dacons_logicor;
129  }
130  else
131  {
132  assert(type == 2);
133  *constype = dacons_setppc;
134  }
135 
136  return SCIP_OKAY;
137 
138 }
139 
140 
141 /** sorts according to distance in solution */
142 static
144  SCIP* scip, /**< SCIP data structure */
145  const GRAPH* graph, /**< graph */
146  const int* result, /**< solution array */
147  DAPATHS* dapaths /**< to be initialized */
148  )
149 {
150  SCIP_Real* RESTRICT starts_prio;
151  SCIP_Real* RESTRICT nodes_dist;
152  SCIP_Bool* RESTRICT nodes_visisted;
153  int* RESTRICT nodes_startid;
154  int* RESTRICT stackarr;
155  int* const startnodes = dapaths->startnodes;
156  const int nstartnodes = dapaths->nstartnodes;
157  const int nnodes = graph_get_nNodes(graph);
158  int stacksize;
159  int startcount = 0;
160 
161  assert(result);
162  assert(nstartnodes > 0);
163 
164  SCIP_CALL( SCIPallocBufferArray(scip, &starts_prio, nstartnodes) );
165  SCIP_CALL( SCIPallocBufferArray(scip, &nodes_dist, nnodes) );
166  SCIP_CALL( SCIPallocBufferArray(scip, &nodes_visisted, nnodes) );
167  SCIP_CALL( SCIPallocBufferArray(scip, &nodes_startid, nnodes) );
168  SCIP_CALL( SCIPallocBufferArray(scip, &stackarr, nnodes) );
169 
170  BMSclearMemoryArray(nodes_visisted, nnodes);
171  nodes_visisted[graph->source] = TRUE;
172  nodes_dist[graph->source] = 0.0;
173 
174 #ifndef NDEBUG
175  for( int i = 0; i < nstartnodes; i++ )
176  starts_prio[i] = FARAWAY;
177 #endif
178 
179  for( int i = 0; i < nnodes; i++ )
180  {
181  if( Is_term(graph->term[i]) && i != graph->source )
182  {
183  nodes_startid[i] = startcount++;
184  }
185  else
186  {
187  nodes_startid[i] = -1;
188  }
189  }
190  assert(startcount == nstartnodes);
191 
192  stacksize = 0;
193  stackarr[stacksize++] = graph->source;
194 
195  /* DFS loop */
196  while( stacksize != 0 )
197  {
198  const int node = stackarr[--stacksize];
199 
200  for( int a = graph->outbeg[node]; a != EAT_LAST; a = graph->oeat[a] )
201  {
202  if( result[a] == CONNECT )
203  {
204  const int head = graph->head[a];
205 
206  if( !nodes_visisted[head] )
207  {
208  nodes_dist[head] = nodes_dist[node] + graph->cost[a];
209  nodes_visisted[head] = TRUE;
210  stackarr[stacksize++] = head;
211 
212  if( Is_term(graph->term[head]) )
213  {
214  const int start = nodes_startid[head];
215  assert(start >= 0 && start < nstartnodes);
216 
217  starts_prio[start] = -nodes_dist[head];
218  }
219  }
220  }
221  }
222  }
223 
224 #ifndef NDEBUG
225  for( int i = 0; i < nstartnodes; i++ )
226  {
227  assert(!EQ(starts_prio[i], FARAWAY));
228  }
229 #endif
230 
231 #ifdef SCIP_DEBUG
232  for( int i = 0; i < nstartnodes; i++ )
233  {
234  graph_knot_printInfo(graph, startnodes[i]);
235  printf("...%f \n", starts_prio[i]);
236  }
237 #endif
238 
239  SCIPsortDownRealInt(starts_prio, startnodes, nstartnodes);
240 
241 #ifdef SCIP_DEBUG
242  printf("after \n");
243  for( int i = 0; i < nstartnodes; i++ )
244  {
245  graph_knot_printInfo(graph, startnodes[i]);
246  printf("...%f \n", starts_prio[i]);
247  }
248 #endif
249 
250  SCIPfreeBufferArray(scip, &stackarr);
251  SCIPfreeBufferArray(scip, &nodes_startid);
252  SCIPfreeBufferArray(scip, &nodes_visisted);
253  SCIPfreeBufferArray(scip, &nodes_dist);
254  SCIPfreeBufferArray(scip, &starts_prio);
255 
256  return SCIP_OKAY;
257 }
258 
259 
260 /** sets shortest path parameters: start node and abort nodes */
261 static
263  const GRAPH* graph, /**< graph */
264  DAPATHS* dapaths /**< to be initialized */
265  )
266 {
267  SCIP_Bool* RESTRICT nodes_abort = dapaths->nodes_abort;
268  int* RESTRICT startnodes = dapaths->startnodes;
269  const int nnodes = graph_get_nNodes(graph);
270  const int root = graph->source;
271  int centernode = UNKNOWN;
272  int nstartnodes = 0;
273 
274  BMSclearMemoryArray(nodes_abort, nnodes);
275  nodes_abort[root] = TRUE;
276 
277  for( int i = 0; i < nnodes; i++ )
278  {
279  if( !Is_term(graph->term[i]) )
280  continue;
281 
282  if( i == root )
283  continue;
284 
285  if( dapaths->isUnrootedPcMw )
286  {
287  assert(graph->grad[i] == 2);
288 
289  if( nstartnodes == 0 )
290  {
291  for( int e = graph->inpbeg[i]; e != EAT_LAST; e = graph->ieat[e] )
292  {
293  if( GT(graph->cost[e], 0.0) )
294  {
295  centernode = graph->tail[e];
296  assert(graph->grad[centernode] == 2 * (graph->terms - 1));
297  nodes_abort[centernode] = TRUE;
298  }
299  }
300  }
301  }
302 
303  startnodes[nstartnodes++] = i;
304  }
305 
306  assert(nstartnodes > 0);
307  assert(centernode >= 0 || !dapaths->isUnrootedPcMw);
308 
309  dapaths->nstartnodes = nstartnodes;
310  dapaths->maxdist = -FARAWAY;
311  dapaths->centernode = centernode;
312  assert(graph->outbeg[root] >= 0);
313 }
314 
315 
316 /** initializes */
317 static
319  SCIP* scip, /**< SCIP data structure */
320  const GRAPH* graph, /**< graph */
321  DAPATHS* dapaths /**< to be initialized */
322  )
323 {
324  const int nnodes = graph_get_nNodes(graph);
325  SCIP_CALL( graph_dijkLimited_init(scip, graph, &(dapaths->dijklimited)) );
326 
327  SCIP_CALL( SCIPallocBufferArray(scip, &(dapaths->nodes_hits), nnodes) );
328  SCIP_CALL( SCIPallocBufferArray(scip, &(dapaths->nodes_abort), nnodes) );
329  SCIP_CALL( SCIPallocBufferArray(scip, &(dapaths->startnodes), nnodes) );
330 
331  BMSclearMemoryArray(dapaths->nodes_hits, nnodes);
332 
333  dapathsSetRunParams(graph, dapaths);
334 
335  return SCIP_OKAY;
336 }
337 
338 
339 /** initializes reduced costs */
340 static
342  const GRAPH* graph, /**< graph; possibly transformed SAP graph */
343  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs */
344  SCIP_Real* objval /**< pointer to store (dual) objective value */
345 )
346 {
347  const int nedges = graph_get_nEdges(graph);
348 
349  BMScopyMemoryArray(redcost, graph->cost, nedges);
350  *objval = 0.0;
351 }
352 
353 
354 /** updates reduced costs and hit count */
355 static
357  const GRAPH* g, /**< graph */
358  const DAPATHS* dapaths, /**< DA paths data */
359  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs */
360  SCIP_Real* objval /**< pointer to store (dual) objective value */
361 )
362 {
363  const DIJK* dijklimited = dapaths->dijklimited;
364  const SCIP_Real* const node_distance = dijklimited->node_distance;
365  const int* const visitlist = dijklimited->visitlist;
366  const SCIP_Real maxdist = dapaths->maxdist;
367  const int nvisits = dijklimited->nvisits;
368  int* RESTRICT nodes_hits = dapaths->nodes_hits;
369  SCIP_Bool* RESTRICT nodes_abort = dapaths->nodes_abort;
370  const int hitlimit = dapaths->isUnrootedPcMw ? DAPATHS_HITLIMIT_PCMW : DAPATHS_HITLIMIT;
371 
372  assert(GE(maxdist, 0.0));
373  assert(GE(*objval, 0.0));
374  assert(nvisits >= 0);
375 
376  // go over all visited...in and out separated.. todo maybe use visited edge list if slow
377  for( int k = 0; k < nvisits; k++ )
378  {
379  const int node_i = visitlist[k];
380  const SCIP_Real dist_i = MIN(node_distance[node_i], maxdist);
381 
382  assert(graph_knot_isInRange(g, node_i));
383 
384  if( nodes_hits[node_i] > hitlimit )
385  {
386  nodes_abort[node_i] = TRUE;
387  }
388  else
389  {
390  nodes_hits[node_i]++;
391  }
392 
393  for( int e = g->outbeg[node_i]; e >= 0; e = g->oeat[e] )
394  {
395  const int node_j = g->head[e];
396  const SCIP_Real dist_j = node_distance[node_j];
397 
398  if( LT(dist_j, maxdist) )
399  {
400  const SCIP_Real offset = MAX(0.0, dist_i - dist_j);
401  assert(LE_FEAS_EPS(offset, redcost[e], EPSILON));
402 
403  redcost[e] -= offset;
404 
405  if( redcost[e] < 0.0 )
406  redcost[e] = 0.0;
407  }
408  }
409  }
410 
411  *objval += maxdist;
412 }
413 
414 
415 /** runs */
416 static
418  const GRAPH* graph, /**< graph; possibly transformed SAP graph */
419  DAPATHS* dapaths, /**< to be initialized */
420  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs */
421  SCIP_Real* objval /**< objective */
422  )
423 {
424  const int* startnodes = dapaths->startnodes;
425  SCIP_Bool* nodes_abort = dapaths->nodes_abort;
426  DIJK* dijklimited = dapaths->dijklimited;
427  const int nstarts = dapaths->nstartnodes;
428 
429  dapathsInitRedCosts(graph, redcost, objval);
430 
431  for( int i = 0; i < nstarts; i++ )
432  {
433  const int start = startnodes[i];
434  assert(graph_knot_isInRange(graph, start) && Is_term(graph->term[start]));
435  assert(nodes_abort[graph->source]);
436 
437  graph_pathInLimitedExec(graph, redcost, nodes_abort, start, dijklimited, &(dapaths->maxdist));
438  dapathsUpdate(graph, dapaths, redcost, objval);
439 
440  graph_dijkLimited_reset(graph, dijklimited);
441  }
442 
443  /* for PC/MW we make sure that we reach the root */
444  if( dapaths->isUnrootedPcMw )
445  {
446  assert(nodes_abort[dapaths->centernode]);
447 
448  BMSclearMemoryArray(dapaths->nodes_abort, graph->knots);
449  nodes_abort[graph->source] = TRUE;
450  graph_pathInLimitedExec(graph, redcost, nodes_abort, dapaths->centernode, dijklimited, &(dapaths->maxdist));
451  dapathsUpdate(graph, dapaths, redcost, objval);
452  }
453 }
454 
455 
456 /** frees */
457 static
459  SCIP* scip, /**< SCIP data structure */
460  DAPATHS* dapaths /**< to be initialized */
461  )
462 {
463  SCIPfreeBufferArray(scip, &(dapaths->startnodes));
464  SCIPfreeBufferArray(scip, &(dapaths->nodes_abort));
465  SCIPfreeBufferArray(scip, &(dapaths->nodes_hits));
466  graph_dijkLimited_free(scip, &(dapaths->dijklimited));
467 }
468 
469 
470 /** returns whether node realtail is active or leads to active node other than dfsbase */
471 static inline
473  const int* active, /**< active nodes array */
474  int realtail, /**< vertex to start from */
475  int dfsbase /**< DFS source vertex */
476  )
477 {
478  int curr;
479 
480  for( curr = active[realtail]; curr != 0 && curr != dfsbase + 1; curr = active[curr - 1] )
481  {
482  assert(curr >= 0);
483  }
484 
485  return (curr == 0);
486 }
487 
488 
489 /** initializes */
490 static
492  const GRAPH* g, /**< graph data structure */
493  const int* edgemap, /**< CSR ancestor edge array */
494  int ncsredges, /**< number of CSR edges */
495  SCIP_Real* RESTRICT rescap, /**< residual capacity */
496  SCIP_Real* dualobj /**< dual objective */
497 )
498 {
499  for( int i = 0; i < ncsredges; i++ )
500  rescap[i] = g->cost[edgemap[i]];
501 
502  *dualobj = 0.0;
503 }
504 
505 
506 /** updates */
507 static
509  SCIP* scip, /**< SCIP */
510  const GRAPH* g, /**< graph data structure */
511  const int* edgemap, /**< CSR ancestor edge array */
512  int ncsredges, /**< number of CSR edges */
513  SCIP_Real* rescap /**< residual capacity */
514 )
515 {
516  const int nedges = graph_get_nEdges(g);
517  SCIP_Real* rescaps_org;
518 
519  SCIP_CALL( SCIPallocBufferArray(scip, &rescaps_org, nedges) );
520  BMScopyMemoryArray(rescaps_org, rescap, nedges);
521 
522  for( int i = 0; i < ncsredges; i++ )
523  {
524  const int edgeorg = edgemap[i];
525  assert(graph_edge_isInRange(g, edgeorg));
526 
527  rescap[i] = rescaps_org[edgeorg];
528  }
529 
530  SCIPfreeBufferArray(scip, &rescaps_org);
531 
532  return SCIP_OKAY;
533 }
534 
535 
536 /** initializes */
537 static
539  SCIP* scip, /**< SCIP */
540  const GRAPH* g, /**< graph data structure */
541  int root, /**< the root */
542  SCIP_Bool is_pseudoroot, /**< is the root a pseudo root? */
543  int* gmark, /**< array for marking nodes */
544  int* RESTRICT active, /**< active vertices mark */
545  SCIP_PQUEUE* pqueue, /**< priority queue */
546  GNODE* gnodearr, /**< array containing terminal nodes*/
547  int* augmentingcomponent /**< augmenting component */
548 )
549 {
550  const int nnodes = g->knots;
551 
552  *augmentingcomponent = -1;
553 
554  /* mark terminals as active, add all except root to pqueue */
555  for( int i = 0, termcount = 0; i < nnodes; i++ )
556  {
557  if( !Is_term(g->term[i]) )
558  {
559  active[i] = -1;
560  continue;
561  }
562 
563  active[i] = 0;
564  assert(g->grad[i] > 0);
565 
566  if( i != root )
567  {
568  assert(termcount < g->terms - 1);
569  assert(gnodearr);
570 
571  gnodearr[termcount].number = i;
572  gnodearr[termcount].dist = g->grad[i];
573 
574  /* for variants with dummy terminals */
575  if( g->grad[i] == 2 )
576  {
577  int a;
578 
579  for( a = g->inpbeg[i]; a != EAT_LAST; a = g->ieat[a] )
580  if( SCIPisZero(scip, g->cost[a]) )
581  break;
582 
583  if( a != EAT_LAST )
584  {
585  const int tail = g->tail[a];
586  gnodearr[termcount].dist += g->grad[tail] - 1;
587 
588  if( is_pseudoroot )
589  {
590  for( a = g->inpbeg[tail]; a != EAT_LAST; a = g->ieat[a] )
591  {
592  if( SCIPisZero(scip, g->cost[a]) )
593  {
594  gnodearr[termcount].dist += g->grad[g->tail[a]] - 1;
595  }
596  }
597  }
598  }
599 
600  assert(gnodearr[termcount].dist > 0);
601  }
602 
603  SCIP_CALL(SCIPpqueueInsert(pqueue, &(gnodearr[termcount])));
604 
605  if( *augmentingcomponent == -1 )
606  *augmentingcomponent = i;
607 
608  termcount++;
609  }
610  }
611 
612  assert(*augmentingcomponent >= 0);
613 
614  for( int i = 0; i < nnodes + 1; i++ )
615  gmark[i] = FALSE;
616 
617  return SCIP_OKAY;
618 }
619 
620 
621 /** dual ascent heuristic */
622 static
624  SCIP* scip, /**< SCIP data structure */
625  const GRAPH* g, /**< graph data structure */
626  const int* result, /**< solution array or NULL */
627  const DAPARAMS* daparams, /**< parameter */
628  SCIP_Bool updateRescaps, /**< update? */
629  SCIP_Real* RESTRICT rescap, /**< residual capacities aka reduced costs */
630  SCIP_Real* objval /**< pointer to store objective value */
631 )
632 {
633  SCIP_CONSHDLR* conshdlr = NULL;
634  SCIP_PQUEUE* pqueue;
635  SCIP_VAR** vars;
636  GNODE* gnodearr;
637  int* RESTRICT edgearr;
638  int* RESTRICT tailarr;
639  int* RESTRICT start;
640  int* RESTRICT stackarr;
641  int* RESTRICT cutverts;
642  int* RESTRICT unsatarcs;
643  int* RESTRICT unsattails;
644  int* RESTRICT gmark;
645  int* RESTRICT active;
646  SCIP_Real dualobj;
647  SCIP_Real currscore;
648  const SCIP_Real maxdeviation = (daparams->damaxdeviation > 0.0) ? daparams->damaxdeviation : DEFAULT_DAMAXDEVIATION;
649  const int nnodes = g->knots;
650  const int nterms = g->terms;
651  const int nedges = g->edges;
652  int ncsredges;
653  int norgcutverts;
654  int stacklength;
655  int augmentingcomponent;
656  const SCIP_Bool addconss = (SCIPgetStage(scip) < SCIP_STAGE_INITSOLVE);
657  int root = daparams->root;
658  const SCIP_Bool addcuts = daparams->addcuts;
659  const SCIP_Bool is_pseudoroot = daparams->is_pseudoroot;
660  const SCIP_Bool withInfinityArcs = graph_typeIsDirected(g) || g->stp_type == STP_DCSTP;
661  SCIP_Bool consUseInital = TRUE;
662  enum DACONS_Type constype = dacons_logicor;
663 
664  assert(rescap);
665  assert(addconss || !addcuts); /* should currently not be activated */
666  assert(maxdeviation >= DA_MAXDEVIATION_LOWER && maxdeviation <= DA_MAXDEVIATION_UPPER);
667  assert(daparams->damaxdeviation == -1.0 || daparams->damaxdeviation > 0.0);
668 
669  /* if specified root is not a terminal, take default root */
670  if( root < 0 || !Is_term(g->term[root]) )
671  root = g->source;
672 
673 
674  if( addcuts )
675  {
676  SCIP_CALL( daconsGetParams(scip, &constype, &consUseInital) );
677 
678  vars = SCIPprobdataGetVars(scip);
679  assert(vars != NULL);
680 
681  if( !addconss )
682  {
683  conshdlr = SCIPfindConshdlr(scip, "stp");
684  assert(conshdlr != NULL);
685  }
686  }
687  else
688  {
689  vars = NULL;
690  }
691 
692 #ifdef BITFIELDSARRAY
693  u_int32_t* bitarr;
694  SCIP_CALL( SCIPallocBufferArray(scip, &bitarr, nedges / ARRLENGTH + 1) );
695 #endif
696 
697  stacklength = 0;
698 
699  if( nterms > 1 )
700  SCIP_CALL( SCIPallocMemoryArray(scip, &gnodearr, nterms - 1) );
701  else
702  gnodearr = NULL;
703 
704  SCIP_CALL( SCIPpqueueCreate(&pqueue, nterms, 2.0, GNODECmpByDist, NULL) );
705 
706  SCIP_CALL( SCIPallocMemoryArray(scip, &active, nnodes) );
707  SCIP_CALL( SCIPallocMemoryArray(scip, &edgearr, nedges) );
708  SCIP_CALL( SCIPallocMemoryArray(scip, &tailarr, nedges) );
709  SCIP_CALL( SCIPallocMemoryArray(scip, &start, nnodes + 1) );
710  SCIP_CALL( SCIPallocMemoryArray(scip, &gmark, nnodes + 1) );
711  SCIP_CALL( SCIPallocMemoryArray(scip, &stackarr, nnodes) );
712 
713  /* fill auxiliary adjacent vertex/edges arrays */
714  graph_getCsr(g, edgearr, tailarr, start, &ncsredges);
715 
716  if( updateRescaps )
717  {
718  dualobj = *objval;
719  SCIP_CALL( daUpdateRescaps(scip, g, edgearr, ncsredges, rescap) );
720  }
721  else
722  {
723  daInitRescaps(g, edgearr, ncsredges, rescap, &dualobj);
724  }
725 
726  SCIP_CALL( SCIPallocBufferArray(scip, &unsattails, nedges) );
727  SCIP_CALL( SCIPallocBufferArray(scip, &cutverts, nnodes) );
728  SCIP_CALL( SCIPallocBufferArray(scip, &unsatarcs, nedges) );
729 
730  SCIP_CALL( daInit(scip, g, root, is_pseudoroot, gmark, active, pqueue,
731  gnodearr, &augmentingcomponent) );
732 
733  /* mark whether an arc is satisfied (has capacity 0) */
734  for( int i = 0; i < ncsredges; i++ )
735  {
736 #ifdef BITFIELDSARRAY
737  if( SCIPisZero(scip, rescap[i]) )
738  SetBit(bitarr, i);
739  else
740  CleanBit(bitarr, i);
741 #else
742  if( EQ(rescap[i], 0.0) )
743  {
744  rescap[i] = 0.0;
745 
746  if( active[tailarr[i] - 1] == 0 )
747  tailarr[i] = 0;
748  else
749  tailarr[i] *= -1;
750  }
751 #endif
752  }
753 
754  norgcutverts = 0;
755 
756  /* (main) dual ascent loop */
757  while( SCIPpqueueNElems(pqueue) > 0 && !SCIPisStopped(scip) )
758  {
759  /* get active vertex of minimum score */
760  GNODE* const gnodeact = (GNODE*) SCIPpqueueRemove(pqueue);
761  const SCIP_Real prio1 = gnodeact->dist;
762  const SCIP_Real prio2 = (SCIPpqueueNElems(pqueue) > 0) ? ((GNODE*) SCIPpqueueFirst(pqueue))->dist : FARAWAY;
763  const int v = gnodeact->number;
764  SCIP_Real degsum = g->grad[v];
765  int ncutverts = 0;
766  int nunsatarcs = 0;
767 
768  SCIP_Bool firstrun = TRUE;
769 
770  SCIPdebugMessage("DA: START WITH v %d prio1 %f prio2 %f \n", v, prio1, prio2);
771 
772  /* perform augmentation as long as priority of root component does not exceed max deviation */
773  for( ; ; )
774  {
775  assert(stacklength == 0);
776 
777  /* 1. step: BFS from v (or connected component) on saturated, incoming arcs */
778 
779  if( firstrun )
780  {
781  firstrun = FALSE;
782  gmark[v + 1] = TRUE;
783  cutverts[ncutverts++] = v;
784  assert(stacklength < nnodes);
785  stackarr[stacklength++] = v;
786  }
787  /* not in first processing of root component: */
788  else
789  {
790  for( int i = norgcutverts; i < ncutverts; i++ )
791  {
792  const int s = cutverts[i];
793 
794  assert(gmark[s + 1]);
795  assert(active[s] != 0);
796  assert(stacklength < nnodes);
797 
798  stackarr[stacklength++] = s;
799  }
800  }
801 #ifdef DFS
802  while( stacklength )
803  {
804  const int node = stackarr[--stacklength];
805 #else
806  for( int n = 0; n < stacklength; n++ )
807  {
808  int end;
809 
810  assert(n < nnodes);
811  node = stackarr[n];
812 #endif
813 
814  /* traverse incoming arcs */
815  for( int i = start[node], end = start[node + 1]; i != end; i++ )
816  {
817  int tail = tailarr[i];
818 
819  /* zero reduced-cost arc? */
820  if( tail <= 0 )
821  {
822  tail *= -1;
823  if( !gmark[tail] )
824  {
825  /* if an active vertex has been hit (other than v), break */
826  if( 0 == tail )
827  {
828  const int realtail = g->tail[edgearr[i]];
829 
830  /* v should not be processed */
831  if( realtail == v )
832  continue;
833 
834  /* is realtail active or does realtail lead to an active vertex other than v? */
835  if( daNodeIsActive(active, realtail, v) )
836  {
837  active[v] = realtail + 1;
838  stacklength = 0;
839  goto ENDOFLOOP;
840  }
841 
842  tail = realtail + 1;
843 
844  /* have we processed tail already? */
845  if( gmark[tail] )
846  continue;
847  }
848 
849  assert(tail > 0);
850 
851  gmark[tail] = TRUE;
852  tail--;
853  cutverts[ncutverts++] = tail;
854  degsum += g->grad[tail];
855 
856  assert(stacklength < nnodes);
857  stackarr[stacklength++] = tail;
858  } /* marked */
859  } /* zero reduced-cost arc */
860  else if( !gmark[tail] )
861  {
862  unsattails[nunsatarcs] = tail;
863  unsatarcs[nunsatarcs++] = i;
864  }
865  }
866  }
867 #ifndef DFS
868  stacklength = 0;
869 #endif
870  currscore = degsum - (ncutverts - 1);
871 
872  /* guiding solution provided? */
873  if( result != NULL )
874  {
875  int nsolarcs = 0;
876  for( int i = 0; i < nunsatarcs; i++ )
877  {
878  const int a = unsatarcs[i];
879 
880  assert(tailarr[a] > 0);
881 
882  if( !(gmark[tailarr[a]]) )
883  {
884  if( result[edgearr[a]] == CONNECT )
885  nsolarcs++;
886  }
887  }
888 
889  assert(nsolarcs > 0);
890  assert(currscore <= nedges);
891 
892  if( nsolarcs > 1 )
893  currscore += (SCIP_Real) ((nsolarcs - 1) * (g->knots * 2.0));
894  }
895  else
896  {
897  assert(SCIPisGE(scip, currscore, prio1));
898  }
899 
900  SCIPdebugMessage("DA: deviation %f \n", (currscore - prio1) / prio1);
901  SCIPdebugMessage("DA: currscore %f prio1 %f prio2 %f \n", currscore, prio1, prio2);
902 
903  /* augmentation criteria met? */
904  if( ((currscore - prio1) / prio1) <= maxdeviation || currscore <= prio2 )
905  {
906  SCIP_CONS* cons = NULL;
907  SCIP_ROW* row = NULL;
908 
909  int shift = 0;
910  SCIP_Real min = FARAWAY;
911  SCIP_Bool isactive = FALSE;
912 
913  /* 2. step: get minimum residual capacity among cut-arcs */
914 
915  /* adjust array of unsatisfied arcs */
916 
917  for( int i = 0; i < nunsatarcs; i++ )
918  {
919  const int tail = unsattails[i];
920 
921  if( gmark[tail] )
922  {
923  shift++;
924  }
925  else
926  {
927  const int a = unsatarcs[i];
928 
929  assert(tailarr[a] > 0);
930  assert(rescap[a] > 0);
931 
932  if( rescap[a] < min )
933  min = rescap[a];
934  if( shift )
935  {
936  unsattails[i - shift] = tail;
937  unsatarcs[i - shift] = a;
938  }
939  }
940  }
941 
942  assert(SCIPisLT(scip, min, FARAWAY));
943  nunsatarcs -= shift;
944 
945  norgcutverts = ncutverts;
946 
947  /* 3. step: perform augmentation */
948 
949  /* create constraints/cuts ? */
950  if( addcuts )
951  {
952  if( addconss )
953  {
954  SCIP_CALL( daconsCreateEmpty(scip, constype, consUseInital, &cons) );
955  }
956  else
957  {
958  SCIP_CALL( SCIPcreateEmptyRowConshdlr(scip, &row, conshdlr, "da", 1.0,
959  SCIPinfinity(scip), FALSE, FALSE, TRUE) );
960 
961  SCIP_CALL( SCIPcacheRowExtensions(scip, row) );
962  }
963  }
964 
965  shift = 0;
966 
967  /* update (dual) objective */
968  dualobj += min;
969 
970  for( int i = 0; i < nunsatarcs; i++ )
971  {
972  const int a = unsatarcs[i];
973  assert(a >= 0);
974 
975  if( addcuts )
976  {
977  assert(vars != NULL);
978 
979  if( !withInfinityArcs || LT(g->cost[edgearr[a]], FARAWAY) )
980  {
981  assert(addconss);
982  if( constype == dacons_linear )
983  SCIP_CALL( SCIPaddCoefLinear(scip, cons, vars[edgearr[a]], 1.0) );
984  else if( constype == dacons_logicor )
985  SCIP_CALL( SCIPaddCoefLogicor(scip, cons, vars[edgearr[a]]) );
986  else
987  SCIP_CALL( SCIPaddCoefSetppc(scip, cons, vars[edgearr[a]]) );
988 
989 #ifdef SCIP_DISABLED
990  if( addconss )
991  SCIP_CALL( SCIPaddCoefLogicor(scip, cons, vars[edgearr[a]]) );//, 1.0) );
992  else
993  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[edgearr[a]], 1.0) );
994 #endif
995  }
996  }
997  rescap[a] -= min;
998 
999  assert(SCIPisGE(scip, rescap[a], 0.0));
1000 
1001  if( rescap[a] <= DA_EPS )
1002  {
1003  int tail = unsattails[i];
1004 
1005  rescap[a] = 0.0;
1006 
1007  assert(tail > 0);
1008  assert(tailarr[a] > 0);
1009 
1010  tailarr[a] *= -1;
1011 
1012  if( active[tail - 1] >= 0 && daNodeIsActive(active, tail - 1, v) )
1013  {
1014  assert(tail - 1 != v);
1015  tailarr[a] = 0;
1016  if( !isactive )
1017  {
1018  isactive = TRUE;
1019  active[v] = tail;
1020  }
1021  }
1022 
1023 
1024  if( !(gmark[tail]) )
1025  {
1026  assert(tail != 0);
1027 
1028  gmark[tail] = TRUE;
1029  tail--;
1030  degsum += g->grad[tail];
1031  cutverts[ncutverts++] = tail;
1032  }
1033 
1034  shift++;
1035  }
1036  else if( shift )
1037  {
1038  unsattails[i - shift] = unsattails[i];
1039  unsatarcs[i - shift] = a;
1040  }
1041  }
1042 
1043  if( addcuts )
1044  {
1045  if( addconss )
1046  {
1047  SCIP_CALL( SCIPaddCons(scip, cons) );
1048  SCIP_CALL( SCIPreleaseCons(scip, &cons) );
1049  }
1050  else
1051  {
1052  SCIP_Bool infeasible;
1053 
1054  SCIP_CALL( SCIPflushRowExtensions(scip, row) );
1055  SCIP_CALL( SCIPaddRow(scip, row, FALSE, &infeasible) );
1056  SCIP_CALL( SCIPreleaseRow(scip, &row) );
1057 
1058  assert(!infeasible);
1059  }
1060  }
1061 
1062  if( isactive )
1063  {
1064  stacklength = 0;
1065  goto ENDOFLOOP;
1066  }
1067  nunsatarcs -= shift;
1068 
1069  continue;
1070  }
1071  else
1072  {
1073  SCIP_Bool insert = TRUE;
1074 
1075  if( is_pseudoroot )
1076  {
1077  int i = start[v];
1078  const int end = start[v + 1];
1079 
1080  assert(end - i == 2);
1081 
1082  for( ; i != end; i++ )
1083  if( rescap[i] != 0.0 )
1084  break;
1085 
1086  if( i == end )
1087  {
1088  if( augmentingcomponent == -1 )
1089  augmentingcomponent = v;
1090 
1091  if( augmentingcomponent != v )
1092  insert = FALSE;
1093  }
1094  }
1095 
1096  if( insert )
1097  {
1098  /* reinsert active vertex */
1099  gnodeact->dist = currscore;
1100  SCIP_CALL( SCIPpqueueInsert(pqueue, gnodeact) );
1101  }
1102  }
1103 
1104  ENDOFLOOP:
1105 
1106  for( int i = 0; i < ncutverts; i++ )
1107  gmark[cutverts[i] + 1] = FALSE;
1108 
1109  for( int i = 0; i < nnodes + 1; i++ )
1110  {
1111  assert(!gmark[i]);
1112  }
1113 
1114  break;
1115  } /* augmentation loop */
1116  } /* dual ascent loop */
1117 
1118  SCIPdebugMessage("DA: dualglobal: %f \n", dualobj);
1119  *objval = dualobj;
1120 
1121  for( int i = ncsredges; i < nedges; i++ )
1122  {
1123  edgearr[i] = i;
1124  rescap[i] = g->cost[i];
1125  }
1126 
1127  /* re-extend rescap array */
1128  for( int i = 0; i < ncsredges; i++ )
1129  {
1130  if( edgearr[i] != i )
1131  {
1132  SCIP_Real bufferedval = rescap[i];
1133  int a = i;
1134 
1135  rescap[i] = g->cost[i];
1136  while( edgearr[a] != a )
1137  {
1138  const int shift = edgearr[a];
1139  const SCIP_Real min = rescap[shift];
1140 
1141  rescap[shift] = bufferedval;
1142  bufferedval = min;
1143  edgearr[a] = a;
1144  a = shift;
1145  }
1146  }
1147  }
1148 
1149 #ifdef BITFIELDSARRAY
1150  SCIPfreeBufferArray(scip, &bitarr);
1151 #endif
1152 
1153  SCIPfreeMemoryArray(scip, &stackarr);
1154  SCIPfreeMemoryArray(scip, &gmark);
1155  SCIPfreeMemoryArray(scip, &start);
1156  SCIPfreeMemoryArray(scip, &tailarr);
1157  SCIPfreeMemoryArray(scip, &edgearr);
1158  SCIPfreeMemoryArray(scip, &active);
1159 
1160  SCIPpqueueFree(&pqueue);
1161  SCIPfreeMemoryArrayNull(scip, &gnodearr);
1162 
1163  /* call Ascend-And-Prune? */
1164  if( daparams->ascendandprune )
1165  {
1166  SCIP_Bool success;
1167 
1168  SCIP_CALL( SCIPStpHeurAscendPruneRun(scip, NULL, g, rescap, unsatarcs, root, &success, TRUE) );
1169  }
1170 
1171  SCIPfreeBufferArray(scip, &unsatarcs);
1172  SCIPfreeBufferArray(scip, &cutverts);
1173  SCIPfreeBufferArray(scip, &unsattails);
1174 
1175  assert(dualascent_allTermsReachable(scip, g, root, rescap));
1176 
1177  return SCIP_OKAY;
1178 }
1179 
1180 
1181 /**@} */
1182 
1183 /**@name Interface methods
1184  *
1185  * @{
1186  */
1187 
1188 
1189 
1190 /** dual ascent heuristic */
1192  SCIP* scip, /**< SCIP data structure */
1193  const GRAPH* g, /**< graph data structure */
1194  const int* result, /**< solution array or NULL */
1195  const DAPARAMS* daparams, /**< parameter */
1196  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs or NULL */
1197  SCIP_Real* objval /**< pointer to store objective value */
1198 )
1199 {
1200  SCIP_Real* RESTRICT rescap;
1201 
1202  assert(scip && g && daparams && objval);
1203 
1204  if( g->knots == 1 )
1205  {
1206  *objval = 0.0;
1207  return SCIP_OKAY;
1208  }
1209 
1210  if( redcost == NULL )
1211  SCIP_CALL( SCIPallocBufferArray(scip, &rescap, g->edges) );
1212  else
1213  rescap = redcost;
1214 
1215  SCIP_CALL( daExec(scip, g, result, daparams, FALSE, rescap, objval) );
1216 
1217  if( redcost == NULL )
1218  SCIPfreeBufferArray(scip, &rescap);
1219 
1220  return SCIP_OKAY;
1221 }
1222 
1223 
1224 /** updates reduced costs with dual ascent heuristic */
1226  SCIP* scip, /**< SCIP data structure */
1227  const GRAPH* g, /**< graph data structure */
1228  const int* result, /**< solution array or NULL */
1229  const DAPARAMS* daparams, /**< parameter */
1230  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs */
1231  SCIP_Real* objval /**< pointer to store objective value */
1232 )
1233 {
1234  assert(scip && g && daparams && objval && redcost);
1235  assert(GE(*objval, 0.0));
1236 
1237  if( g->knots == 1 )
1238  return SCIP_OKAY;
1239 
1240  assert(!dualascent_allTermsReachable(scip, g, daparams->root, redcost));
1241 
1242  printf("dual bound before %f \n", *objval);
1243 
1244  SCIP_CALL( daExec(scip, g, result, daparams, TRUE, redcost, objval) );
1245 
1246  printf("dual bound after %f \n", *objval);
1247 
1248  assert(dualascent_allTermsReachable(scip, g, daparams->root, redcost));
1249 
1250  return SCIP_OKAY;
1251 }
1252 
1253 
1254 
1255 /** dual ascent heuristic for degree constrained problem */
1257  SCIP* scip, /**< SCIP data structure */
1258  const GRAPH* g, /**< graph data structure */
1259  const int* result, /**< solution array or NULL */
1260  const DAPARAMS* daparams, /**< parameter */
1261  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs or NULL */
1262  SCIP_Real* objval /**< pointer to store objective value */
1263 )
1264 {
1265  SCIP_Real* orgcosts;
1266  SCIP_Real* RESTRICT rescap;
1267  const int nnodes = graph_get_nNodes(g);
1268  const int nedges = graph_get_nEdges(g);
1269  const int root = daparams->root;
1270 
1271  assert(scip && g && daparams && objval);
1272  assert(g->stp_type == STP_DCSTP);
1273  assert(g->maxdeg);
1274 
1275  if( nnodes == 1 )
1276  {
1277  *objval = 0.0;
1278  return SCIP_OKAY;
1279  }
1280 
1281  assert(root == g->source || !daparams->addcuts);
1282  assert(graph_knot_isInRange(g, root) && Is_term(g->term[root]));
1283 
1284  SCIP_CALL( SCIPallocBufferArray(scip, &orgcosts, nedges) );
1285  BMScopyMemoryArray(orgcosts, g->cost, nedges);
1286 
1287  for( int k = 0; k < nnodes; k++ )
1288  {
1289  if( g->maxdeg[k] != 1 )
1290  continue;
1291 
1292  if( !Is_term(g->term[k]) || k == root )
1293  continue;
1294 
1295  for( int e = g->outbeg[k]; e != EAT_LAST; e = g->oeat[e] )
1296  {
1297  g->cost[e] = FARAWAY;
1298  }
1299  }
1300 
1301  if( redcost == NULL )
1302  SCIP_CALL( SCIPallocBufferArray(scip, &rescap, nedges) );
1303  else
1304  rescap = redcost;
1305 
1306  SCIP_CALL( daExec(scip, g, result, daparams, FALSE, rescap, objval) );
1307 
1308  if( redcost == NULL )
1309  SCIPfreeBufferArray(scip, &rescap);
1310 
1311  BMScopyMemoryArray(g->cost, orgcosts, nedges);
1312  SCIPfreeBufferArray(scip, &orgcosts);
1313 
1314  return SCIP_OKAY;
1315 }
1316 
1317 
1318 /** dual ascent heuristic for PCSPG and MWCSP */
1320  SCIP* scip, /**< SCIP data structure */
1321  GRAPH* g, /**< graph data structure */
1322  SCIP_Real* redcost, /**< array to store reduced costs or NULL */
1323  SCIP_Real* objval, /**< pointer to store objective value */
1324  SCIP_Bool addcuts, /**< should dual-ascent add Steiner cuts? */
1325  SCIP_Bool ascendandprune, /**< perform ascend-and-prune and add solution? */
1326  int nruns /**< number of dual ascent runs */
1327  )
1328 {
1329  SCIP_CONSHDLR* conshdlr = NULL;
1330  SCIP_PQUEUE* pqueue;
1331  SCIP_VAR** vars;
1332  GRAPH* transgraph;
1333  SCIP_Real min;
1334  SCIP_Real prio1;
1335  SCIP_Real offset;
1336  SCIP_Real dualobj;
1337  SCIP_Real currscore;
1338  SCIP_Real maxdeviation;
1339  SCIP_Real* rescap;
1340  GNODE* gnodeact;
1341  GNODE** gnodearr;
1342  int s;
1343  int i;
1344  int k;
1345  int v;
1346  int a;
1347  int tail;
1348  int pnode;
1349  int shift;
1350  int root;
1351  int nnodes;
1352  int nterms;
1353  int nedges;
1354  int degsum;
1355  int ncutverts;
1356  int pseudoroot;
1357  int nunsatarcs;
1358  int stacklength;
1359  int norgcutverts;
1360  int* cutverts;
1361  int* stackarr;
1362  STP_Bool* origedge;
1363  int* unsatarcs;
1364  STP_Bool firstrun;
1365  STP_Bool* sat;
1366  STP_Bool* active;
1367  const SCIP_Bool addconss = (SCIPgetStage(scip) < SCIP_STAGE_INITSOLVE);
1368  SCIP_Bool consUseInital = TRUE;
1369  enum DACONS_Type constype = dacons_logicor;
1370 
1371  /* should currently not be activated */
1372  assert(addconss || !addcuts);
1373 
1374  assert(g != NULL);
1375  assert(scip != NULL);
1376  assert(nruns >= 0);
1377  assert(objval != NULL);
1378 
1379  if( g->knots == 1 )
1380  {
1381  *objval = 0.0;
1382  return SCIP_OKAY;
1383  }
1384 
1385 
1386  if( addcuts )
1387  {
1388  SCIP_CALL( daconsGetParams(scip, &constype, &consUseInital) );
1389 
1390  vars = SCIPprobdataGetVars(scip);
1391  assert(vars != NULL);
1392  if( !addconss )
1393  {
1394  conshdlr = SCIPfindConshdlr(scip, "stp");
1395  assert(conshdlr != NULL);
1396  }
1397  }
1398  else
1399  {
1400  vars = NULL;
1401  }
1402 
1403  root = g->source;
1404  degsum = 0;
1405  offset = 0.0;
1406  dualobj = 0.0;
1407 
1408  ncutverts = 0;
1409  norgcutverts = 0;
1410  maxdeviation = DEFAULT_DAMAXDEVIATION;
1411 
1412  SCIP_CALL( graph_transPcGetSap(scip, g, &transgraph, &offset) );
1413 
1414  nnodes = transgraph->knots;
1415  nedges = transgraph->edges;
1416  nterms = transgraph->terms;
1417  pseudoroot = nnodes - 1;
1418 
1419  if( redcost == NULL )
1420  {
1421  SCIP_CALL( SCIPallocBufferArray(scip, &rescap, nedges) );
1422  }
1423  else
1424  {
1425  rescap = redcost;
1426  }
1427 
1428  stacklength = 0;
1429  SCIP_CALL( SCIPallocBufferArray(scip, &stackarr, nnodes) );
1430  SCIP_CALL( SCIPallocBufferArray(scip, &sat, nedges) );
1431  SCIP_CALL( SCIPallocBufferArray(scip, &active, nnodes) );
1432  SCIP_CALL( SCIPallocBufferArray(scip, &cutverts, nnodes) );
1433  SCIP_CALL( SCIPallocBufferArray(scip, &gnodearr, nterms - 1) );
1434  SCIP_CALL( SCIPallocBufferArray(scip, &unsatarcs, nedges) );
1435  SCIP_CALL( SCIPallocBufferArray(scip, &origedge, nedges) );
1436 
1437  for( i = 0; i < nedges; i++ )
1438  if( !Is_term(transgraph->term[transgraph->tail[i]]) && transgraph->head[i] == pseudoroot )
1439  origedge[i] = FALSE;
1440  else if( transgraph->tail[i] == pseudoroot && !Is_term(transgraph->term[transgraph->head[i]]) )
1441  origedge[i] = FALSE;
1442  else
1443  origedge[i] = TRUE;
1444 
1445  for( i = 0; i < nterms - 1; i++ )
1446  {
1447  SCIP_CALL( SCIPallocBuffer(scip, &gnodearr[i]) ); /*lint !e866*/
1448  }
1449 
1450  SCIP_CALL( SCIPpqueueCreate( &pqueue, nnodes, 2.0, GNODECmpByDist, NULL) );
1451 
1452  k = 0;
1453  /* mark terminals as active, add all except root to pqueue */
1454  for( i = 0; i < nnodes; i++ )
1455  {
1456  if( Is_term(transgraph->term[i]) )
1457  {
1458  active[i] = TRUE;
1459  assert(transgraph->grad[i] > 0);
1460  if( i != root )
1461  {
1462  gnodearr[k]->number = i;
1463  gnodearr[k]->dist = transgraph->grad[i];
1464 
1465  for( a = transgraph->inpbeg[i]; a != EAT_LAST; a = transgraph->ieat[a] )
1466  if( SCIPisEQ(scip, transgraph->cost[a], 0.0) )
1467  break;
1468 
1469  if( a != EAT_LAST )
1470  gnodearr[k]->dist += transgraph->grad[transgraph->tail[a]] - 1;
1471 
1472  assert(gnodearr[k]->dist > 0);
1473 
1474  SCIP_CALL( SCIPpqueueInsert(pqueue, gnodearr[k++]) );
1475  }
1476  }
1477  else
1478  {
1479  active[i] = FALSE;
1480  }
1481  transgraph->mark[i] = FALSE;
1482  }
1483 
1484  for( i = 0; i < nedges; i++ )
1485  {
1486  rescap[i] = transgraph->cost[i];
1487  if( SCIPisZero(scip, rescap[i]) )
1488  sat[i] = TRUE;
1489  else
1490  sat[i] = FALSE;
1491  }
1492 
1493  /* dual ascent loop */
1494  while( SCIPpqueueNElems(pqueue) > 0 && !SCIPisStopped(scip) )
1495  {
1496  /* get active vertex of minimum score */
1497  gnodeact = (GNODE*) SCIPpqueueRemove(pqueue);
1498 
1499  v = gnodeact->number;
1500  prio1 = gnodeact->dist;
1501 
1502  firstrun = TRUE;
1503  nunsatarcs = 0;
1504 
1505  /* perform augmentation as long as ... */
1506  for( ; ; )
1507  {
1508  assert(stacklength == 0);
1509  /* 1. step: BFS from v (or connected component) on saturated, incoming arcs */
1510 
1511  if( firstrun )
1512  {
1513  degsum = transgraph->grad[v];
1514  ncutverts = 0;
1515  firstrun = FALSE;
1516  nunsatarcs = 0;
1517  transgraph->mark[v] = TRUE;
1518  cutverts[ncutverts++] = v;
1519  stackarr[stacklength++] = v;
1520  }
1521  /* not in first processing of root component: */
1522  else
1523  {
1524  for( i = norgcutverts; i < ncutverts; i++ )
1525  {
1526  s = cutverts[i];
1527  assert(transgraph->mark[s]);
1528  if( active[s] )
1529  {
1530  active[v] = FALSE;
1531  stacklength = 0;
1532  goto ENDOFLOOP;
1533  }
1534 
1535  stackarr[stacklength++] = s;
1536  }
1537  }
1538 
1539  while( stacklength )
1540  {
1541  pnode = stackarr[--stacklength];
1542 
1543  /* traverse incoming arcs */
1544  for( a = transgraph->inpbeg[pnode]; a != EAT_LAST; a = transgraph->ieat[a] )
1545  {
1546  tail = transgraph->tail[a];
1547  if( sat[a] )
1548  {
1549  if( !transgraph->mark[tail] )
1550  {
1551  /* if an active vertex has been hit, break */
1552  if( active[tail] )
1553  {
1554  active[v] = FALSE;
1555  stacklength = 0;
1556  goto ENDOFLOOP;
1557  }
1558 
1559  degsum += transgraph->grad[tail];
1560  transgraph->mark[tail] = TRUE;
1561  cutverts[ncutverts++] = tail;
1562  stackarr[stacklength++] = tail;
1563  }
1564  }
1565  else if( !transgraph->mark[tail] )
1566  {
1567  unsatarcs[nunsatarcs++] = a;
1568  }
1569  }
1570  }
1571 
1572  currscore = degsum - (ncutverts - 1);
1573 
1574  assert(SCIPisGE(scip, currscore, prio1));
1575 
1576  /* augmentation criteria met? */
1577  if( SCIPisLE(scip, (currscore - prio1) / prio1, maxdeviation) || (SCIPpqueueNElems(pqueue) == 0) )
1578  {
1579  SCIP_Bool in = FALSE;
1580  SCIP_ROW* row;
1581  SCIP_CONS* cons = NULL;
1582 
1583  /* 2. pass: get minimum residual capacity among cut-arcs */
1584 
1585  /* adjust array of unsatisfied arcs */
1586  min = FARAWAY;
1587  shift = 0;
1588 
1589  for( i = 0; i < nunsatarcs; i++ )
1590  {
1591  a = unsatarcs[i];
1592  if( transgraph->mark[transgraph->tail[a]] )
1593  {
1594  shift++;
1595  }
1596  else
1597  {
1598 
1599  assert(!sat[a]);
1600  if( SCIPisLT(scip, rescap[a], min) )
1601  min = rescap[a];
1602  if( shift != 0 )
1603  unsatarcs[i - shift] = a;
1604  }
1605  }
1606 
1607  assert(SCIPisLT(scip, min, FARAWAY));
1608  nunsatarcs -= shift;
1609 
1610  if( nunsatarcs > 0)
1611  assert(!transgraph->mark[transgraph->tail[unsatarcs[nunsatarcs-1]]]);
1612 
1613  norgcutverts = ncutverts;
1614 
1615 
1616  /* 3. pass: perform augmentation */
1617 
1618 
1619  /* create constraint/row */
1620 
1621  if( addcuts )
1622  {
1623  if( addconss )
1624  {
1625  SCIP_CALL( daconsCreateEmpty(scip, constype, consUseInital, &cons) );
1626  }
1627  else
1628  {
1629  SCIP_CALL(SCIPcreateEmptyRowConshdlr(scip, &row, conshdlr, "da", 1.0, SCIPinfinity(scip), FALSE, FALSE, TRUE));
1630  SCIP_CALL(SCIPcacheRowExtensions(scip, row));
1631  }
1632  }
1633 
1634  dualobj += min;
1635  for( i = 0; i < nunsatarcs; i++ )
1636  {
1637  a = unsatarcs[i];
1638  if( a == -1 )
1639  continue;
1640 
1641  if( addcuts && origedge[a] )
1642  {
1643  assert(vars && cons);
1644  assert(addconss);
1645 
1646  if( g->tail[a] == root && g->head[a] == v )
1647  in = TRUE;
1648 
1649  if( constype == dacons_linear )
1650  SCIP_CALL( SCIPaddCoefLinear(scip, cons, vars[a], 1.0) );
1651  else if( constype == dacons_logicor )
1652  SCIP_CALL( SCIPaddCoefLogicor(scip, cons, vars[a]) );
1653  else
1654  SCIP_CALL( SCIPaddCoefSetppc(scip, cons, vars[a]) );
1655 
1656 #ifdef SCIP_DISABLED
1657  if( addconss )
1658  SCIP_CALL( SCIPaddCoefLinear(scip, cons, vars[a], 1.0) );
1659  else
1660  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[a], 1.0) );
1661 #endif
1662  }
1663  rescap[a] -= min;
1664 
1665  assert(SCIPisGE(scip, rescap[a], 0.0));
1666 
1667  if( SCIPisEQ(scip, rescap[a], 0.0) )
1668  {
1669  sat[a] = TRUE;
1670  if( !(transgraph->mark[transgraph->tail[a]]) )
1671  {
1672  tail = transgraph->tail[a];
1673  degsum += transgraph->grad[tail];
1674  transgraph->mark[tail] = TRUE;
1675  cutverts[ncutverts++] = tail;
1676  }
1677  }
1678  }
1679 
1680  if( addcuts )
1681  {
1682  assert(vars != NULL);
1683 
1684  if( !in )
1685  {
1686  for( i = g->outbeg[root]; i != EAT_LAST; i = g->oeat[i] )
1687  {
1688  if( g->head[i] == v )
1689  {
1690  if( constype == dacons_linear )
1691  SCIP_CALL( SCIPaddCoefLinear(scip, cons, vars[i], 1.0) );
1692  else if( constype == dacons_logicor )
1693  SCIP_CALL( SCIPaddCoefLogicor(scip, cons, vars[i]) );
1694  else
1695  SCIP_CALL( SCIPaddCoefSetppc(scip, cons, vars[i]) );
1696 
1697 #ifdef SCIP_DISABLED
1698  if( addconss )
1699  SCIP_CALL( SCIPaddCoefLinear(scip, cons, vars[i], 1.0) );
1700  else
1701  SCIP_CALL( SCIPaddVarToRow(scip, row, vars[i], 1.0) );
1702 #endif
1703  }
1704  }
1705  }
1706 
1707  if( addconss )
1708  {
1709  SCIP_CALL( SCIPaddCons(scip, cons) );
1710  SCIP_CALL( SCIPreleaseCons(scip, &cons) );
1711  }
1712  else
1713  {
1714  SCIP_Bool infeasible;
1715  assert(row != NULL);
1716 
1717  SCIP_CALL( SCIPflushRowExtensions(scip, row) );
1718  SCIP_CALL( SCIPaddRow(scip, row, FALSE, &infeasible) );
1719  SCIP_CALL( SCIPreleaseRow(scip, &row) );
1720 
1721  assert(!infeasible);
1722  }
1723  }
1724 
1725  continue;
1726  }
1727  else
1728  {
1729  /* reinsert active vertex */
1730  gnodeact->dist = currscore;
1731  SCIP_CALL( SCIPpqueueInsert(pqueue, gnodeact) );
1732  }
1733 
1734  ENDOFLOOP:
1735 
1736  for( i = 0; i < ncutverts; i++ )
1737  transgraph->mark[cutverts[i]] = FALSE;
1738 
1739  break;
1740  } /* augmentation loop */
1741  } /* dual ascent loop */
1742 
1743 
1744  *objval = dualobj + offset;
1745  SCIPdebugMessage("DA: dualglobal: %f \n", *objval + SCIPprobdataGetOffset(scip));
1746 
1747  /* call dual Ascend-And-Prune? */
1748  if( ascendandprune )
1749  {
1750  SCIP_Bool success;
1751  SCIP_CALL( SCIPStpHeurAscendPruneRun(scip, NULL, g, rescap, unsatarcs, -1, &success, TRUE));
1752  }
1753 
1754  /* free memory */
1755  SCIPpqueueFree(&pqueue);
1756 
1757  for( i = nterms - 2; i >= 0; i-- )
1758  SCIPfreeBuffer(scip, &gnodearr[i]);
1759 
1760  SCIPfreeBufferArray(scip, &origedge);
1761  SCIPfreeBufferArray(scip, &unsatarcs);
1762  SCIPfreeBufferArray(scip, &cutverts);
1763  SCIPfreeBufferArray(scip, &gnodearr);
1764  SCIPfreeBufferArray(scip, &active);
1765  SCIPfreeBufferArray(scip, &sat);
1766  SCIPfreeBufferArray(scip, &stackarr);
1767 
1768  assert(dualascent_allTermsReachable(scip, transgraph, root, rescap));
1769 
1770  if( redcost == NULL )
1771  SCIPfreeBufferArray(scip, &rescap);
1772 
1773  graph_free(scip, &transgraph, TRUE);
1774 
1775  return SCIP_OKAY;
1776 }
1777 
1778 
1779 /** path based dual ascent heuristic */
1781  SCIP* scip, /**< SCIP data structure */
1782  const GRAPH* transgraph, /**< transformed SAP graph */
1783  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs */
1784  SCIP_Real* objval, /**< pointer to store (dual) objective value */
1785  const int* result /**< solution array or NULL */
1786 )
1787 {
1788  DAPATHS dapaths = { NULL, NULL, NULL, NULL, -1.0, -1, -1, .isUnrootedPcMw = TRUE };
1789 
1790  assert(scip && transgraph && redcost && objval);
1791 
1792  SCIP_CALL( dapathsInit(scip, transgraph, &dapaths) );
1793 
1794  dapathsRunShortestPaths(transgraph, &dapaths, redcost, objval);
1795 
1796  dapathsFreeMembers(scip, &dapaths);
1797 
1798  return SCIP_OKAY;
1799 }
1800 
1801 
1802 /** path based dual ascent heuristic */
1804  SCIP* scip, /**< SCIP data structure */
1805  GRAPH* graph, /**< graph */
1806  SCIP_Real* RESTRICT redcost, /**< array to store reduced costs */
1807  SCIP_Real* objval, /**< pointer to store (dual) objective value */
1808  const int* result /**< solution array or NULL */
1809 )
1810 {
1811  DAPATHS dapaths = { NULL, NULL, NULL, NULL, -1.0, -1, -1, .isUnrootedPcMw = FALSE };
1812 
1813  assert(scip && graph && redcost && objval);
1814  assert(!graph_pc_isPcMw(graph) || graph_pc_isRootedPcMw(graph));
1815 
1816  if( graph_pc_isPcMw(graph) )
1817  graph_pc_2transcheck(scip, graph);
1818 
1819  SCIP_CALL( dapathsInit(scip, graph, &dapaths) );
1820 
1821  if( result )
1822  {
1823  SCIP_CALL( dapathsSortStarts(scip, graph, result, &dapaths) );
1824  }
1825 
1826  dapathsRunShortestPaths(graph, &dapaths, redcost, objval);
1827 
1828  dapathsFreeMembers(scip, &dapaths);
1829 
1830  return SCIP_OKAY;
1831 }
1832 
1833 
1834 /** can all terminal be reached via reduced costs from given root? */
1836  SCIP* scip, /**< SCIP */
1837  const GRAPH* g, /**< graph */
1838  int root, /**< root for reduced costs */
1839  const SCIP_Real* redcost /**< reduced costs */
1840  )
1841 {
1842  int* RESTRICT queue;
1843  STP_Bool* RESTRICT scanned;
1844  int qsize;
1845  const int nnodes = graph_get_nNodes(g);
1846  int termscount;
1847 
1848  assert(scip && redcost);
1849  assert(graph_knot_isInRange(g, root));
1850  assert(Is_term(g->term[root]));
1851 
1852  SCIP_CALL_ABORT( SCIPallocMemoryArray(scip, &queue, nnodes ) );
1853  SCIP_CALL_ABORT( SCIPallocMemoryArray(scip, &scanned, nnodes) );
1854 
1855  BMSclearMemoryArray(scanned, nnodes);
1856 
1857  termscount = 1;
1858  qsize = 0;
1859  scanned[root] = TRUE;
1860  queue[qsize++] = root;
1861 
1862  /* DFS */
1863  while( qsize > 0 )
1864  {
1865  const int k = queue[--qsize];
1866  scanned[k] = TRUE;
1867 
1868  for( int a = g->outbeg[k]; a != EAT_LAST; a = g->oeat[a] )
1869  {
1870  const int head = g->head[a];
1871 
1872  if( SCIPisZero(scip, redcost[a]) )
1873  {
1874  /* vertex not visited yet? */
1875  if( !scanned[head] )
1876  {
1877  scanned[head] = TRUE;
1878  queue[qsize++] = head;
1879 
1880  if( Is_term(g->term[head]) )
1881  termscount++;
1882  }
1883  }
1884  }
1885  }
1886 
1887  SCIPfreeMemoryArray(scip, &scanned);
1888  SCIPfreeMemoryArray(scip, &queue);
1889 
1890  return (termscount == g->terms);
1891 }
1892 
1893 
1894 
1895 /**@} */
SCIP_RETCODE SCIPaddCoefLogicor(SCIP *scip, SCIP_CONS *cons, SCIP_VAR *var)
int graph_get_nEdges(const GRAPH *)
Definition: graph_stats.c:548
int SCIPpqueueNElems(SCIP_PQUEUE *pqueue)
Definition: misc.c:1468
static volatile int nterms
Definition: interrupt.c:38
int * visitlist
Definition: graphdefs.h:314
SCIP_RETCODE SCIPcacheRowExtensions(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1626
int *RESTRICT head
Definition: graphdefs.h:224
SCIP_Bool graph_typeIsDirected(const GRAPH *)
Definition: graph_stats.c:83
SCIP_RETCODE SCIPaddCoefSetppc(SCIP *scip, SCIP_CONS *cons, SCIP_VAR *var)
Definition: cons_setppc.c:9372
SCIP_STAGE SCIPgetStage(SCIP *scip)
Definition: scip_general.c:356
int source
Definition: graphdefs.h:195
SCIP_Bool graph_pc_isRootedPcMw(const GRAPH *)
#define SCIPfreeMemoryArrayNull(scip, ptr)
Definition: scip_mem.h:72
SCIP_RETCODE dualascent_execDegCons(SCIP *scip, const GRAPH *g, const int *result, const DAPARAMS *daparams, SCIP_Real *RESTRICT redcost, SCIP_Real *objval)
Definition: dualascent.c:1256
#define Is_term(a)
Definition: graphdefs.h:102
SCIP_CONSHDLR * SCIPfindConshdlr(SCIP *scip, const char *name)
Definition: scip_cons.c:877
SCIP_RETCODE SCIPflushRowExtensions(SCIP *scip, SCIP_ROW *row)
Definition: scip_lp.c:1649
static SCIP_RETCODE dapathsInit(SCIP *scip, const GRAPH *graph, DAPATHS *dapaths)
Definition: dualascent.c:318
#define SCIPfreeMemoryArray(scip, ptr)
Definition: scip_mem.h:71
#define EPSILON
Definition: lpi_qso.c:93
int terms
Definition: graphdefs.h:192
SCIP_RETCODE SCIPaddVarToRow(SCIP *scip, SCIP_ROW *row, SCIP_VAR *var, SCIP_Real val)
Definition: scip_lp.c:1686
SCIP_VAR ** SCIPprobdataGetVars(SCIP *scip)
static SCIP_RETCODE daExec(SCIP *scip, const GRAPH *g, const int *result, const DAPARAMS *daparams, SCIP_Bool updateRescaps, SCIP_Real *RESTRICT rescap, SCIP_Real *objval)
Definition: dualascent.c:623
SCIP_RETCODE SCIPcreateConsSetcover(SCIP *scip, SCIP_CONS **cons, const char *name, int nvars, SCIP_VAR **vars, SCIP_Bool initial, SCIP_Bool separate, SCIP_Bool enforce, SCIP_Bool check, SCIP_Bool propagate, SCIP_Bool local, SCIP_Bool modifiable, SCIP_Bool dynamic, SCIP_Bool removable, SCIP_Bool stickingatnode)
Definition: cons_setppc.c:9317
#define SCIPallocMemoryArray(scip, ptr, num)
Definition: scip_mem.h:55
SCIP_Bool SCIPisGE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
SCIP_Real damaxdeviation
Definition: dualascent.h:45
int *RESTRICT maxdeg
Definition: graphdefs.h:206
void graph_knot_printInfo(const GRAPH *, int)
Definition: graph_stats.c:151
reduction and dual-cost based primal heuristic for Steiner problems
void graph_free(SCIP *, GRAPH **, SCIP_Bool)
Definition: graph_base.c:796
#define FALSE
Definition: def.h:87
#define DAPATHS_HITLIMIT_PCMW
Definition: dualascent.c:43
Includes dual-ascent for classic Steiner tree and some variants.
int *RESTRICT inpbeg
Definition: graphdefs.h:202
SCIP_Real SCIPinfinity(SCIP *scip)
Problem data for stp problem.
#define TRUE
Definition: def.h:86
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:54
includes various files containing graph methods used for Steiner tree problems
static SCIP_Bool daNodeIsActive(const int *active, int realtail, int dfsbase)
Definition: dualascent.c:472
SCIP_Bool * nodes_abort
Definition: dualascent.c:66
static void dapathsRunShortestPaths(const GRAPH *graph, DAPATHS *dapaths, SCIP_Real *RESTRICT redcost, SCIP_Real *objval)
Definition: dualascent.c:417
SCIP_RETCODE dualascent_update(SCIP *scip, const GRAPH *g, const int *result, const DAPARAMS *daparams, SCIP_Real *RESTRICT redcost, SCIP_Real *objval)
Definition: dualascent.c:1225
#define STP_DCSTP
Definition: graphdefs.h:47
void * SCIPpqueueFirst(SCIP_PQUEUE *pqueue)
Definition: misc.c:1454
SCIP_RETCODE graph_transPcGetSap(SCIP *, GRAPH *, GRAPH **, SCIP_Real *)
Definition: graph_trans.c:931
SCIP_Bool isUnrootedPcMw
Definition: dualascent.c:70
static GRAPHNODE ** active
#define EAT_LAST
Definition: graphdefs.h:58
void SCIPpqueueFree(SCIP_PQUEUE **pqueue)
Definition: misc.c:1263
#define SCIPdebugMessage
Definition: pub_message.h:87
#define FARAWAY
Definition: graphdefs.h:89
SCIP_Real maxdist
Definition: dualascent.c:67
void SCIPsortDownRealInt(SCIP_Real *realarray, int *intarray, int len)
SCIP_Bool SCIPisEQ(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
#define SCIPfreeBufferArray(scip, ptr)
Definition: scip_mem.h:127
int number
Definition: misc_stp.h:55
Constraint handler for the set partitioning / packing / covering constraints .
#define DA_MAXDEVIATION_LOWER
Definition: dualascent.c:40
static SCIP_RETCODE daconsCreateEmpty(SCIP *scip, enum DACONS_Type constype, SCIP_Bool consUseInital, SCIP_CONS **cons)
Definition: dualascent.c:81
SCIP_Real * node_distance
Definition: graphdefs.h:316
SCIP_RETCODE SCIPaddCoefLinear(SCIP *scip, SCIP_CONS *cons, SCIP_VAR *var, SCIP_Real val)
SCIP_RETCODE dualascent_exec(SCIP *scip, const GRAPH *g, const int *result, const DAPARAMS *daparams, SCIP_Real *RESTRICT redcost, SCIP_Real *objval)
Definition: dualascent.c:1191
static SCIP_RETCODE daInit(SCIP *scip, const GRAPH *g, int root, SCIP_Bool is_pseudoroot, int *gmark, int *RESTRICT active, SCIP_PQUEUE *pqueue, GNODE *gnodearr, int *augmentingcomponent)
Definition: dualascent.c:538
int *RESTRICT mark
Definition: graphdefs.h:198
void graph_dijkLimited_reset(const GRAPH *, DIJK *)
Definition: graph_util.c:2105
int *RESTRICT oeat
Definition: graphdefs.h:231
SCIP_RETCODE dualascent_paths(SCIP *scip, GRAPH *graph, SCIP_Real *RESTRICT redcost, SCIP_Real *objval, const int *result)
Definition: dualascent.c:1803
#define GE(a, b)
Definition: portab.h:85
int stp_type
Definition: graphdefs.h:264
SCIP_RETCODE SCIPStpHeurAscendPruneRun(SCIP *scip, SCIP_HEUR *heur, const GRAPH *g, const SCIP_Real *redcosts, int *result, int root, SCIP_Bool *solfound, SCIP_Bool addsol)
SCIP_RETCODE SCIPaddCons(SCIP *scip, SCIP_CONS *cons)
Definition: scip_prob.c:2769
SCIP_Bool SCIPisLT(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
static void dapathsFreeMembers(SCIP *scip, DAPATHS *dapaths)
Definition: dualascent.c:458
void graph_pc_2transcheck(SCIP *, GRAPH *)
SCIP_Real SCIPprobdataGetOffset(SCIP *scip)
Constraint handler for logicor constraints (equivalent to set covering, but algorithms are suited fo...
#define SCIPallocBuffer(scip, ptr)
Definition: scip_mem.h:113
SCIP_Bool is_pseudoroot
Definition: dualascent.h:44
SCIP_Real dist
Definition: misc_stp.h:56
int *RESTRICT grad
Definition: graphdefs.h:201
void graph_getCsr(const GRAPH *, int *RESTRICT, int *RESTRICT, int *RESTRICT, int *)
SCIP_RETCODE SCIPgetBoolParam(SCIP *scip, const char *name, SCIP_Bool *value)
Definition: scip_param.c:241
internal miscellaneous methods
#define NULL
Definition: lpi_spx1.cpp:155
struct dual_ascent_paths DAPATHS
SCIP_RETCODE SCIPgetIntParam(SCIP *scip, const char *name, int *value)
Definition: scip_param.c:260
#define EQ(a, b)
Definition: portab.h:79
int knots
Definition: graphdefs.h:190
#define SCIP_CALL(x)
Definition: def.h:384
#define LT(a, b)
Definition: portab.h:81
SCIP_RETCODE SCIPaddRow(SCIP *scip, SCIP_ROW *row, SCIP_Bool forcecut, SCIP_Bool *infeasible)
Definition: scip_cut.c:241
void * SCIPpqueueRemove(SCIP_PQUEUE *pqueue)
Definition: misc.c:1434
SCIP_Bool dualascent_allTermsReachable(SCIP *scip, const GRAPH *g, int root, const SCIP_Real *redcost)
Definition: dualascent.c:1835
unsigned char STP_Bool
Definition: portab.h:34
int GNODECmpByDist(void *first_arg, void *second_arg)
SCIP_RETCODE graph_dijkLimited_init(SCIP *, const GRAPH *, DIJK **)
Definition: graph_util.c:1989
#define DA_EPS
Definition: dualascent.c:42
#define SCIPallocBufferArray(scip, ptr, num)
Definition: scip_mem.h:115
static SCIP_RETCODE daconsGetParams(SCIP *scip, enum DACONS_Type *constype, SCIP_Bool *consUseInital)
Definition: dualascent.c:112
#define GT(a, b)
Definition: portab.h:84
void graph_dijkLimited_free(SCIP *, DIJK **)
Definition: graph_util.c:2143
#define SCIP_Bool
Definition: def.h:84
int *RESTRICT ieat
Definition: graphdefs.h:230
static SCIP_RETCODE daUpdateRescaps(SCIP *scip, const GRAPH *g, const int *edgemap, int ncsredges, SCIP_Real *rescap)
Definition: dualascent.c:508
int *RESTRICT tail
Definition: graphdefs.h:223
#define MAX(x, y)
Definition: tclique_def.h:83
int *RESTRICT term
Definition: graphdefs.h:196
#define BMScopyMemoryArray(ptr, source, num)
Definition: memory.h:127
Constraint handler for linear constraints in their most general form, .
SCIP_RETCODE SCIPcreateConsLogicor(SCIP *scip, SCIP_CONS **cons, const char *name, int nvars, SCIP_VAR **vars, SCIP_Bool initial, SCIP_Bool separate, SCIP_Bool enforce, SCIP_Bool check, SCIP_Bool propagate, SCIP_Bool local, SCIP_Bool modifiable, SCIP_Bool dynamic, SCIP_Bool removable, SCIP_Bool stickingatnode)
DACONS_Type
Definition: dualascent.c:57
#define CONNECT
Definition: graphdefs.h:87
#define DEFAULT_DAMAXDEVIATION
Definition: dualascent.c:39
SCIP_Bool ascendandprune
Definition: dualascent.h:42
SCIP_RETCODE SCIPpqueueCreate(SCIP_PQUEUE **pqueue, int initsize, SCIP_Real sizefac, SCIP_DECL_SORTPTRCOMP((*ptrcomp)), SCIP_DECL_PQUEUEELEMCHGPOS((*elemchgpos)))
Definition: misc.c:1236
SCIP_RETCODE SCIPcreateConsLinear(SCIP *scip, SCIP_CONS **cons, const char *name, int nvars, SCIP_VAR **vars, SCIP_Real *vals, SCIP_Real lhs, SCIP_Real rhs, SCIP_Bool initial, SCIP_Bool separate, SCIP_Bool enforce, SCIP_Bool check, SCIP_Bool propagate, SCIP_Bool local, SCIP_Bool modifiable, SCIP_Bool dynamic, SCIP_Bool removable, SCIP_Bool stickingatnode)
static void dapathsUpdate(const GRAPH *g, const DAPATHS *dapaths, SCIP_Real *RESTRICT redcost, SCIP_Real *objval)
Definition: dualascent.c:356
SCIP_RETCODE SCIPreleaseRow(SCIP *scip, SCIP_ROW **row)
Definition: scip_lp.c:1553
#define SCIPfreeBuffer(scip, ptr)
Definition: scip_mem.h:125
#define DA_MAXDEVIATION_UPPER
Definition: dualascent.c:41
SCIP_RETCODE SCIPpqueueInsert(SCIP_PQUEUE *pqueue, void *elem)
Definition: misc.c:1335
void graph_pathInLimitedExec(const GRAPH *, const SCIP_Real *, const SCIP_Bool *, int, DIJK *, SCIP_Real *)
Definition: graph_path.c:610
SCIP_Real * cost
Definition: graphdefs.h:221
SCIP_Bool graph_edge_isInRange(const GRAPH *, int)
Definition: graph_stats.c:110
SCIP_RETCODE SCIPreleaseCons(SCIP *scip, SCIP_CONS **cons)
Definition: scip_cons.c:1110
static void daInitRescaps(const GRAPH *g, const int *edgemap, int ncsredges, SCIP_Real *RESTRICT rescap, SCIP_Real *dualobj)
Definition: dualascent.c:491
SCIP_Bool graph_pc_isPcMw(const GRAPH *)
SCIP_VAR * a
Definition: circlepacking.c:57
#define SCIP_Real
Definition: def.h:177
SCIP_Bool SCIPisStopped(SCIP *scip)
Definition: scip_general.c:694
static void dapathsSetRunParams(const GRAPH *graph, DAPATHS *dapaths)
Definition: dualascent.c:262
static void dapathsInitRedCosts(const GRAPH *graph, SCIP_Real *RESTRICT redcost, SCIP_Real *objval)
Definition: dualascent.c:341
int *RESTRICT outbeg
Definition: graphdefs.h:204
int edges
Definition: graphdefs.h:219
SCIP_RETCODE dualascent_execPcMw(SCIP *scip, GRAPH *g, SCIP_Real *redcost, SCIP_Real *objval, SCIP_Bool addcuts, SCIP_Bool ascendandprune, int nruns)
Definition: dualascent.c:1319
int graph_get_nNodes(const GRAPH *)
Definition: graph_stats.c:536
SCIP_Bool graph_knot_isInRange(const GRAPH *, int)
Definition: graph_node.c:52
SCIP_Bool SCIPisZero(SCIP *scip, SCIP_Real val)
SCIP_Bool SCIPisLE(SCIP *scip, SCIP_Real val1, SCIP_Real val2)
#define UNKNOWN
Definition: sepa_mcf.c:4095
#define nnodes
Definition: gastrans.c:65
#define BMSclearMemoryArray(ptr, num)
Definition: memory.h:123
#define LE_FEAS_EPS(a, b, eps)
Definition: portab.h:75
#define SCIP_CALL_ABORT(x)
Definition: def.h:363
SCIP_RETCODE SCIPcreateEmptyRowConshdlr(SCIP *scip, SCIP_ROW **row, SCIP_CONSHDLR *conshdlr, const char *name, SCIP_Real lhs, SCIP_Real rhs, SCIP_Bool local, SCIP_Bool modifiable, SCIP_Bool removable)
Definition: scip_lp.c:1382
static SCIP_RETCODE dapathsSortStarts(SCIP *scip, const GRAPH *graph, const int *result, DAPATHS *dapaths)
Definition: dualascent.c:143
DIJK * dijklimited
Definition: dualascent.c:63
#define DAPATHS_HITLIMIT
Definition: dualascent.c:44
SCIP_RETCODE dualascent_pathsPcMw(SCIP *scip, const GRAPH *transgraph, SCIP_Real *RESTRICT redcost, SCIP_Real *objval, const int *result)
Definition: dualascent.c:1780
SCIP callable library.