Scippy

SCIP

Solving Constraint Integer Programs

reduce_path.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 scipopt.org. */
13 /* */
14 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
15 
16 /**@file reduce_path.c
17  * @brief Path deletion reduction test for Steiner problems
18  * @author Daniel Rehfeldt
19  *
20  * This file implements an improved version of the so called "path substitution test" by Polzin and Vahdati Daneshmand.
21  *
22  * A list of all interface methods can be found in reduce.h.
23  *
24  */
25 
26 /*---+----1----+----2----+----3----+----4----+----5----+----6----+----7----+----8----+----9----+----0----+----1----+----2*/
27 
28 //#define SCIP_DEBUG
29 
30 #include <assert.h>
31 #include "reduce.h"
32 #include "extreduce.h"
33 
34 #define SP_MAXNDISTS 10000
35 #define SP_MAXLENGTH 11
36 #define SP_MAXDEGREE 8
37 #define SP_MAXDEGREE_FAST 5
38 #define SP_MAXDEGREE_PC 10
39 #define SP_MAXNSTARTS 10000
40 #define VNODES_UNSET -1
41 #define SP_MAXNPULLS 4
42 
43 
44 /** path replacement */
45 typedef struct path_replacement
46 {
47  EXTPERMA* extperma; /**< extension data or NULL */
48  SDPROFIT* sdprofit; /**< SD profit */
49  STP_Vectype(SCIP_Real) firstneighborcosts; /**< edge costs from tail of path to non-path neighbors */
50  STP_Vectype(int) firstneighbors; /**< first neighbors (only used in case SDs are being used) */
51  STP_Vectype(int) currneighbors; /**< current neighbors */
52  STP_Vectype(int) pathedges; /**< edges of path */
53  STP_Vectype(int) visitednodes; /**< visited nodes */
54  SCIP_Real* RESTRICT sp_dists; /**< distances to neighbors of path start node */
55  int* RESTRICT sp_starts; /**< CSR like starts for each node in the path, pointing to sp_dists */
56  int* RESTRICT nodes_index; /**< maps each node to index in 0,1,2,..., or to VNODES_UNSET, VNODES_INPATH */
57  SCIP_Bool* RESTRICT nodes_isTerm; /**< terminal node? */
58  SCIP_Bool* RESTRICT nodeindices_isPath; /**< is a node index in the path? */
59  SCIP_Bool* RESTRICT firstneighbors_isKilled; /**< for each first neighbor index: is ruled-out? */
60  SCIP_Real pathcost; /**< cost of path */
61  int nfirstneighbors; /**< number of neighbors of first path node */
62  int pathtail; /**< first node of path */
63  int failneighbor; /**< temporary */
64  int nnodes; /**< number of nodes */
65  int maxdegree;
66  SCIP_Bool probIsPc; /**< prize-collecting problem? */
67  SCIP_Bool useSd; /**< use special distances? */
68 } PR;
69 
70 
71 /*
72  * Local methods
73  */
74 
75 /** gets profit for given node */
76 inline static
78  const SDPROFIT* sdprofit, /**< the SD profit */
79  const int* nodes_index, /**< index array */
80  int node, /**< node to get profit for */
81  int nonsource /**< node that should not be a source */
82 )
83 {
84  const int source1 = sdprofit->nodes_biassource[node];
85 
86  assert(GE(sdprofit->nodes_bias[node], 0.0));
87  assert(LE(sdprofit->nodes_bias[node], FARAWAY));
88  assert(GE(sdprofit->nodes_bias2[node], 0.0));
89  assert(LE(sdprofit->nodes_bias2[node], FARAWAY));
90  assert(GE(sdprofit->nodes_bias[node], sdprofit->nodes_bias2[node]));
91 
92  if( source1 == node )
93  return sdprofit->nodes_bias[node];
94 
95  if( source1 != nonsource && nodes_index[source1] < 0 )
96  {
97  return sdprofit->nodes_bias[node];
98  }
99  else
100  {
101  const int source2 = sdprofit->nodes_biassource2[node];
102 
103  if( source2 != nonsource && nodes_index[source2] < 0 )
104  {
105  return sdprofit->nodes_bias2[node];
106  }
107  }
108 
109  return 0.0;
110 }
111 
112 
113 /** gets head of path */
114 static inline
116  const GRAPH* g, /**< graph data structure */
117  const PR* pr /**< path replacement data structure */
118  )
119 {
120  const int npathedges = StpVecGetSize(pr->pathedges);
121  return g->head[pr->pathedges[npathedges - 1]];
122 }
123 
124 
125 /** deletes given edge */
126 static inline
128  SCIP* scip, /**< SCIP */
129  int edge, /**< edge to delete */
130  GRAPH* g, /**< graph data structure */
131  PR* pr /**< path replacement data structure */
132  )
133 {
134  if( pr->useSd )
135  {
136  EXTPERMA* const extperma = pr->extperma;
137  assert(extperma && extperma->distdata_default);
138 
139  if( extperma->solIsValid )
140  {
141  const int* const result = extperma->result;
142  if( result[edge] == CONNECT || result[flipedge(edge)] == CONNECT )
143  extperma->solIsValid = FALSE;
144  }
145 
146  extreduce_edgeRemove(scip, edge, g, extperma->distdata_default, extperma);
147  }
148  else
149  {
150  graph_edge_delFull(scip, g, edge, TRUE);
151  }
152 }
153 
154 
155 /** removes nodes of degree one and unmarks */
156 static
158  SCIP* scip, /**< SCIP */
159  EXTPERMA* extperma, /**< extension data */
160  GRAPH* g /**< graph data structure (in/out) */
161 )
162 {
163  const int nnodes = graph_get_nNodes(g);
164  SCIP_Bool rerun = TRUE;
165 
166  while( rerun )
167  {
168  rerun = FALSE;
169  for( int i = 0; i < nnodes; ++i )
170  {
171  if( g->grad[i] == 1 && !Is_term(g->term[i]) )
172  {
173  const int sibling = g->head[g->outbeg[i]];
174 
175  extreduce_edgeRemove(scip, g->outbeg[i], g, extperma->distdata_default, extperma);
176  assert(g->grad[i] == 0);
177 
178  g->mark[i] = FALSE;
179 
180  if( Is_term(g->term[sibling]) )
181  continue;
182 
183  if( g->grad[sibling] == 0 )
184  g->mark[sibling] = FALSE;
185  else if( g->grad[sibling] == 1 )
186  rerun = TRUE;
187  }
188  }
189  }
190 }
191 
192 /** tries to rule out neighbors of path head */
193 static inline
195  SCIP* scip, /**< SCIP */
196  const GRAPH* g, /**< graph data structure */
197  PR* pr, /**< path replacement data structure */
198  SCIP_Bool* needFullRuleOut /**< pointer to store whether full rule-out is needed */
199  )
200 {
201  const SCIP_Real* const sp_dists = pr->sp_dists;
202  const SCIP_Real pathcost = pr->pathcost;
203  const int ncurrneighbors = StpVecGetSize(pr->currneighbors);
204  const int nfirst = pr->nfirstneighbors;
205  int failneighbor = VNODES_UNSET;
206  const SCIP_Bool useSd = pr->useSd;
207  DISTDATA* const distdata = useSd ? pr->extperma->distdata_default : NULL;
208 
209  assert(FALSE == *needFullRuleOut);
210  SCIPdebugMessage("starting head rule-out with pathcost=%f \n", pr->pathcost);
211 
212  for( int i = 0; i < ncurrneighbors; i++ )
213  {
214  int j;
215  const int neighbor = pr->currneighbors[i];
216  const int sp_start = pr->sp_starts[pr->nodes_index[neighbor]];
217  int nfails = 0;
218 
219  /* NOTE: we need to be able to rule out all but one of the path-tail neighbors */
220 
221  assert(pr->nodes_index[neighbor] >= 0);
222 
223  for( j = 0; j < nfirst; j++ )
224  {
225  if( pr->firstneighbors_isKilled[j] )
226  continue;
227 
228  SCIPdebugMessage("dist for neighbor=%d, idx=%d: %f \n", neighbor, j, sp_dists[sp_start + j]);
229  if( GT(sp_dists[sp_start + j], pathcost) )
230  {
231  if( useSd )
232  {
233  const int fneighbor = pr->firstneighbors[j];
234  SCIP_Real sd = extreduce_distDataGetSdDoubleEq(scip, g, pathcost, fneighbor, neighbor, distdata);
235 
236  if( LT(sd, pathcost) )
237  continue;
238 
239  if( EQ(sd, pathcost) )
240  {
241  const int startedge = pr->pathedges[0];
242  sd = extreduce_distDataGetSdDoubleForbiddenSingle(scip, g, startedge, neighbor, fneighbor, distdata);
243  if( LE(sd, pathcost) )
244  continue;
245  }
246  }
247 
248  /* second fail? */
249  if( nfails++ > 0 )
250  break;
251  }
252  }
253 
254  if( nfails > 1 )
255  {
256  SCIPdebugMessage("neighbor %d not head-ruled-out \n", neighbor);
257 
258  /* second time that a current neighbor could not be ruled-out? */
259  if( failneighbor != VNODES_UNSET )
260  {
261  pr->failneighbor = VNODES_UNSET;
262  *needFullRuleOut = TRUE;
263  return;
264  }
265  failneighbor = neighbor;
266  }
267  }
268 
269  pr->failneighbor = failneighbor;
270 }
271 
272 
273 /** tries to rule out neighbors of path tail */
274 static inline
276  SCIP* scip, /**< SCIP */
277  const GRAPH* g, /**< graph data structure */
278  PR* pr, /**< path replacement data structure */
279  SCIP_Bool* isRuledOut /**< ruled-out? */
280 )
281 {
282  const SCIP_Real* const sp_dists = pr->sp_dists;
283  const int* const sp_starts = pr->sp_starts;
284  const int pathhead = pathGetHead(g, pr);
285  const int pathhead_idx = pr->nodes_index[pathhead];
286  const SCIP_Bool useSd = pr->useSd;
287  DISTDATA* const distdata = useSd ? pr->extperma->distdata_default : NULL;
288 
289  assert(!*isRuledOut);
290 
291  SCIPdebugMessage("try tail rule-out with pathcost=%f \n", pr->pathcost);
292 
293  for( int i = 0; i < pr->nfirstneighbors; i++ )
294  {
295  if( !pr->firstneighbors_isKilled[i] )
296  {
297  const SCIP_Real extpathcost = pr->pathcost + pr->firstneighborcosts[i];
298  const SCIP_Real altpathcost = sp_dists[sp_starts[pathhead_idx] + i];
299 
300  SCIPdebugMessage("extpathcost for idx=%d: %f \n", i, extpathcost);
301  SCIPdebugMessage("althpathost for idx=%d: %f \n", i, altpathcost);
302 
303  if( GT(altpathcost, extpathcost) )
304  {
305  if( useSd )
306  {
307  const int fneighbor = pr->firstneighbors[i];
308  SCIP_Real sd = extreduce_distDataGetSdDoubleEq(scip, g, extpathcost, fneighbor, pathhead, distdata);
309 
310  if( LT(sd, extpathcost) )
311  continue;
312 
313  if( EQ(sd, extpathcost) )
314  {
315  const int startedge = pr->pathedges[0];
316  sd = extreduce_distDataGetSdDoubleForbiddenSingle(scip, g, startedge, pathhead, fneighbor, distdata);
317  if( LE(sd, extpathcost) )
318  continue;
319  }
320  }
321 
322  SCIPdebugMessage("no rule-out for initial neighbor idx=%d \n", i);
323  return;
324  }
325  }
326  }
327 
328  *isRuledOut = TRUE;
329 }
330 
331 
332 /** tries to rule out neighbors of path tail */
333 static inline
335  SCIP* scip, /**< SCIP */
336  const GRAPH* g, /**< graph data structure */
337  PR* pr, /**< path replacement data structure */
338  SCIP_Bool* isRuledOut /**< ruled-out? */
339 )
340 {
341  const SCIP_Real* const sp_dists = pr->sp_dists;
342  const int* const sp_starts = pr->sp_starts;
343  const int pathhead = pathGetHead(g, pr);
344  const int pathhead_idx = pr->nodes_index[pathhead];
345  int nfails = 0;
346  const SCIP_Bool useSd = pr->useSd;
347  DISTDATA* const distdata = useSd ? pr->extperma->distdata_default : NULL;
348 
349  SCIPdebugMessage("try full tail rule-out with pathcost=%f \n", pr->pathcost);
350 
351  for( int i = 0; i < pr->nfirstneighbors; i++ )
352  {
353  const SCIP_Real altpathcost = sp_dists[sp_starts[pathhead_idx] + i];
354 
355  SCIPdebugMessage("althpathost for idx=%d: %f \n", i, altpathcost);
356 
357  if( GT(altpathcost, pr->pathcost) )
358  {
359  if( useSd )
360  {
361  const int fneighbor = pr->firstneighbors[i];
362  SCIP_Real sd = extreduce_distDataGetSdDoubleEq(scip, g, pr->pathcost, fneighbor, pathhead, distdata);
363 
364  if( LT(sd, pr->pathcost) )
365  {
366  pr->firstneighbors_isKilled[i] = TRUE;
367  continue;
368  }
369 
370  if( EQ(sd, pr->pathcost) )
371  {
372  const int startedge = pr->pathedges[0];
373  sd = extreduce_distDataGetSdDoubleForbiddenSingle(scip, g, startedge, pathhead, fneighbor, distdata);
374  if( LE(sd, pr->pathcost) )
375  {
376  pr->firstneighbors_isKilled[i] = TRUE;
377  continue;
378  }
379  }
380  }
381 
382  SCIPdebugMessage("no full rule-out for initial neighbor idx=%d \n", i);
383  nfails++;
384  }
385  else
386  {
387  pr->firstneighbors_isKilled[i] = TRUE;
388  }
389  }
390 
391  if( nfails > 1 )
392  {
393  SCIPdebugMessage("...no full rule-out \n");
394  assert(!*isRuledOut);
395  }
396  else
397  {
398  *isRuledOut = TRUE;
399  }
400 }
401 
402 
403 /** adds new path node */
404 static inline
406  SCIP* scip, /**< SCIP */
407  int node, /**< node to add */
408  PR* pr /**< path replacement data structure */
409  )
410 {
411  assert(pr->nodes_index[node] == VNODES_UNSET);
412 
413  pr->nodes_index[node] = StpVecGetSize(pr->visitednodes);
414  StpVecPushBack(scip, pr->visitednodes, node);
415  pr->nodeindices_isPath[pr->nodes_index[node]] = TRUE;
416 }
417 
418 /** adds new NON-path node */
419 static inline
421  SCIP* scip, /**< SCIP */
422  int node, /**< node to add */
423  PR* pr /**< path replacement data structure */
424  )
425 {
426  assert(pr->nodes_index[node] == VNODES_UNSET);
427 
428  pr->nodes_index[node] = StpVecGetSize(pr->visitednodes);
429  StpVecPushBack(scip, pr->visitednodes, node);
430  pr->nodeindices_isPath[pr->nodes_index[node]] = FALSE;
431 }
432 
433 
434 /** adds new path head edge to failed node */
435 static inline
437  SCIP* scip, /**< SCIP */
438  const GRAPH* g, /**< graph */
439  PR* pr /**< path replacement data structure */
440  )
441 {
442  const DCSR* const dcsr = g->dcsr_storage;
443  const RANGE* const dcsr_range = dcsr->range;
444  const int* const dcsr_heads = dcsr->head;
445  const int pathhead = pathGetHead(g, pr);
446  const int extnode = pr->failneighbor;
447  int j;
448 
449  assert(pr->nodes_index[extnode] >= 0);
450  assert(!pr->nodeindices_isPath[pr->nodes_index[extnode]]);
451  assert(extnode >= 0);
452 
453  SCIPdebugMessage("adding new path edge to node %d \n", extnode);
454 
455  for( j = dcsr_range[pathhead].start; j != dcsr_range[pathhead].end; j++ )
456  {
457  if( dcsr_heads[j] == extnode )
458  break;
459  }
460  assert(j != dcsr_range[pathhead].end);
461  assert(EQ(dcsr->cost[j], g->cost[dcsr->edgeid[j]]));
462 
463  StpVecPushBack(scip, pr->pathedges, dcsr->edgeid[j]);
464  pr->pathcost += dcsr->cost[j];
465  pr->nodeindices_isPath[pr->nodes_index[extnode]] = TRUE;
466 
467  if( pr->probIsPc )
468  {
469  assert(LE(g->prize[pathhead], dcsr->cost[j]));
470  pr->pathcost -= g->prize[pathhead];
471  }
472 
473  assert(pathGetHead(g, pr) == extnode);
474 }
475 
476 
477 /** adds first two path nodes */
478 static inline
480  SCIP* scip, /**< SCIP */
481  const GRAPH* g, /**< graph */
482  int startedge_tail, /**< tail of start edge */
483  int startdege_head, /**< head of start edge */
484  PR* pr /**< path replacement data structure */
485  )
486 {
487  const DCSR* const dcsr = g->dcsr_storage;
488  const RANGE* const dcsr_range = dcsr->range;
489  const int* const dcsr_heads = dcsr->head;
490  SCIP_Real* RESTRICT sp_dists = pr->sp_dists;
491  int* RESTRICT sp_starts = pr->sp_starts;
492  int starts_final = pr->nfirstneighbors;
493 
494  /* add tail node */
495  addPathNode(scip, startedge_tail, pr);
496  sp_starts[starts_final + 1] = sp_starts[starts_final] + pr->nfirstneighbors;
497  for( int i = sp_starts[starts_final], j = 0; i != sp_starts[starts_final + 1]; i++, j++ )
498  sp_dists[i] = pr->firstneighborcosts[j];
499 
500  /* add head node*/
501  addPathNode(scip, startdege_head, pr);
502  starts_final++;
503  sp_starts[starts_final + 1] = sp_starts[starts_final] + pr->nfirstneighbors;
504  for( int i = sp_starts[starts_final]; i != sp_starts[starts_final + 1]; i++ )
505  sp_dists[i] = FARAWAY;
506 
507  /* update distances by using single edges */
508  for( int j = dcsr_range[startdege_head].start; j != dcsr_range[startdege_head].end; j++ )
509  {
510  const int head = dcsr_heads[j];
511  const int head_index = pr->nodes_index[head];
512 
513  /* unvisited or path tail? */
514  if( head_index < 0 || head == startedge_tail )
515  continue;
516 
517  assert(0 <= head_index && head_index < pr->nfirstneighbors);
518  assert(EQ(sp_dists[sp_starts[starts_final] + head_index], FARAWAY));
519  SCIPdebugMessage("setting distance from first path node(%d) for initial neighbor (idx=%d) to %f \n",
520  startdege_head, head_index, dcsr->cost[j]);
521  sp_dists[sp_starts[starts_final] + head_index] = dcsr->cost[j];
522  }
523 }
524 
525 
526 /** collects neighbors */
527 static inline
529  SCIP* scip, /**< SCIP */
530  const GRAPH* g, /**< graph data structure */
531  PR* pr /**< path replacement data structure */
532  )
533 {
534  const DCSR* const dcsr = g->dcsr_storage;
535  const RANGE* const dcsr_range = dcsr->range;
536  SCIP_Real* RESTRICT sp_dists = pr->sp_dists;
537  const int* const dcsr_heads = dcsr->head;
538  int* RESTRICT nodes_index = pr->nodes_index;
539  int* RESTRICT sp_starts = pr->sp_starts;
540  const int basenode = pathGetHead(g, pr);
541  const int nfirstneighbors = pr->nfirstneighbors;
542 
543  SCIPdebugMessage("extending path to node %d \n", basenode);
544 
545  StpVecClear(pr->currneighbors);
546 
547  for( int i = dcsr_range[basenode].start; i != dcsr_range[basenode].end; i++ )
548  {
549  const int head = dcsr_heads[i];
550  int head_index = nodes_index[head];
551 
552  if( head_index >= 0 && pr->nodeindices_isPath[head_index] )
553  continue;
554 
555  assert(head_index < StpVecGetSize(pr->visitednodes));
556  SCIPdebugMessage("adding neighbor %d head_index=%d\n", head, head_index);
557 
558  StpVecPushBack(scip, pr->currneighbors, head);
559 
560  if( head_index == VNODES_UNSET )
561  {
562  head_index = StpVecGetSize(pr->visitednodes);
563  SCIPdebugMessage("mapping new neighbor %d->%d \n", head, head_index);
564  nodes_index[head] = head_index;
565  StpVecPushBack(scip, pr->visitednodes, head);
566  pr->nodeindices_isPath[head_index] = FALSE;
567  assert(head_index + 1 < SP_MAXNSTARTS);
568 
569  sp_starts[head_index + 1] = sp_starts[head_index] + nfirstneighbors;
570 
571  for( int j = sp_starts[head_index]; j != sp_starts[head_index + 1]; j++ )
572  sp_dists[j] = FARAWAY;
573  }
574  }
575 }
576 
577 
578 /** updates distances from neighbors */
579 static inline
581  SCIP* scip, /**< SCIP */
582  const GRAPH* g, /**< graph data structure */
583  PR* pr /**< path replacement data structure */
584  )
585 {
586  STP_Vectype(int) currneighbors;
587  const SDPROFIT* const sdprofit = pr->sdprofit;
588  const DCSR* const dcsr = g->dcsr_storage;
589  const RANGE* const dcsr_range = dcsr->range;
590  const SCIP_Real* const dcsr_costs = dcsr->cost;
591  SCIP_Real* RESTRICT sp_dists = pr->sp_dists;
592  const int* const dcsr_heads = dcsr->head;
593  const int* const nodes_index = pr->nodes_index;
594  const int* const sp_starts = pr->sp_starts;
595  const int nfirstneighbors = pr->nfirstneighbors;
596  const int nneighbors = StpVecGetSize(pr->currneighbors);
597  const int pathtail = pr->pathtail;
598  const int pathhead = pathGetHead(g, pr);
599  const int nloops = MIN(SP_MAXNPULLS, nneighbors);
600 
601  /* NOTE: to keep the following loop easy */
602  StpVecPushBack(scip, pr->currneighbors, pathhead);
603  currneighbors = pr->currneighbors;
604 
605  for( int loop = 0; loop < nloops; loop++ )
606  {
607  SCIP_Bool hasUpdates = FALSE;
608  /* compute distances to each of the initial neighbors */
609  for( int iter = 0; iter <= nneighbors; iter++ )
610  {
611  const int node = currneighbors[iter];
612  const int node_start = sp_starts[nodes_index[node]];
613 
614  assert(nodes_index[node] >= 0);
615 
616  for( int i = dcsr_range[node].start; i != dcsr_range[node].end; i++ )
617  {
618  const int head = dcsr_heads[i];
619 
620  if( nodes_index[head] != VNODES_UNSET )
621  {
622  SCIP_Real head_profit = 0.0;
623  const int head_start = sp_starts[nodes_index[head]];
624 
625  if( head == pathtail && node == pathhead )
626  continue;
627 
628  if( nodes_index[head] > nfirstneighbors )
629  {
630  assert(head != pathtail);
631  head_profit = getSdProfit(sdprofit, nodes_index, head, node);
632  }
633 
634  for( int k = 0; k < nfirstneighbors; k++ )
635  {
636  SCIP_Real newdist = dcsr_costs[i] + sp_dists[head_start + k];
637  SCIP_Real profitBias = MIN(head_profit, dcsr_costs[i]);
638 
639  profitBias = MIN(profitBias, sp_dists[head_start + k]);
640  //printf("%f \n", profitBias);
641 
642  newdist -= profitBias;
643 
644  if( LT(newdist, sp_dists[node_start + k]) )
645  {
646  SCIPdebugMessage("(l=%d) updating distance for node %d to orgindex %d to %f \n",
647  loop, node, k, newdist);
648  sp_dists[node_start + k] = newdist;
649  hasUpdates = TRUE;
650  }
651  }
652  }
653  }
654  }
655 
656  if( !hasUpdates )
657  break;
658  }
659 
660  StpVecPopBack(currneighbors);
661 }
662 
663 
664 /** preprocess for doing path extension later on */
665 static inline
667  SCIP* scip, /**< SCIP */
668  const GRAPH* g, /**< graph data structure */
669  int startedge, /**< first edge */
670  PR* pr /**< path replacement data structure */
671  )
672 {
673  const int basetail = g->tail[startedge];
674  const int basehead = g->head[startedge];
675  const DCSR* const dcsr = g->dcsr_storage;
676  const RANGE* const dcsr_range = dcsr->range;
677  const int* const dcsr_heads = dcsr->head;
678  const SCIP_Real* const dcsr_costs = dcsr->cost;
679  SCIP_Real* RESTRICT sp_dists = pr->sp_dists;
680  int* RESTRICT sp_starts = pr->sp_starts;
681  const SCIP_Bool isPc = pr->probIsPc;
682 
683  assert(0 == StpVecGetSize(pr->pathedges));
684  assert(0 == StpVecGetSize(pr->visitednodes));
685  assert(0 == StpVecGetSize(pr->currneighbors));
686  assert(0 == StpVecGetSize(pr->firstneighborcosts));
687  assert(0 == StpVecGetSize(pr->firstneighbors));
688 
689  SCIPdebugMessage("---Checking edge %d->%d \n\n", basetail, basehead);
690 
691  pr->pathcost = g->cost[startedge];
692  pr->pathtail = basetail;
693  StpVecPushBack(scip, pr->pathedges, startedge);
694 
695  for( int i = dcsr_range[basetail].start; i != dcsr_range[basetail].end; i++ )
696  {
697  const int head = dcsr_heads[i];
698  assert(pr->nodes_index[head] == VNODES_UNSET);
699  assert(!isPc || !graph_pc_knotIsDummyTerm(g, head));
700 
701  if( head == basehead )
702  continue;
703 
704  SCIPdebugMessage("mapping first neighbor %d->%d \n", head, StpVecGetSize(pr->visitednodes));
705 
706  addNonPathNode(scip, head, pr);
707  StpVecPushBack(scip, pr->currneighbors, head);
708  StpVecPushBack(scip, pr->firstneighborcosts, dcsr_costs[i]);
709 
710  if( pr->useSd )
711  StpVecPushBack(scip, pr->firstneighbors, head);
712  }
713 
714  pr->nfirstneighbors = StpVecGetSize(pr->currneighbors);
715  SCIPdebugMessage("Having %d initial neighbors \n", pr->nfirstneighbors);
716 
717  sp_starts[0] = 0;
718  for( int i = 0; i < pr->nfirstneighbors; i++ )
719  {
720  const int neighbor = pr->currneighbors[i];
721  const SCIP_Real basecost = pr->firstneighborcosts[i];
722  sp_starts[i + 1] = sp_starts[i] + pr->nfirstneighbors;
723 
724  /* set 2-edge distances via path tail */
725  if( isPc )
726  {
727  for( int j = sp_starts[i], k = 0; j != sp_starts[i + 1]; j++, k++ )
728  {
729  assert(LE(g->prize[basetail], MIN(basecost, pr->firstneighborcosts[k])));
730  sp_dists[j] = basecost + pr->firstneighborcosts[k] - g->prize[basetail];
731  }
732 
733  assert(EQ(sp_dists[sp_starts[i] + i], 2.0 * basecost - g->prize[basetail]));
734  }
735  else
736  {
737  for( int j = sp_starts[i], k = 0; j != sp_starts[i + 1]; j++, k++ )
738  sp_dists[j] = basecost + pr->firstneighborcosts[k];
739 
740  assert(EQ(sp_dists[sp_starts[i] + i], 2.0 * basecost));
741  }
742 
743 #ifdef SCIP_DEBUG
744  for( int j = sp_starts[i]; j != sp_starts[i + 1]; j++ )
745  printf("%d->%d dist=%f \n", i, j - sp_starts[i], sp_dists[j]);
746 #endif
747 
748  /* set self-distance */
749  sp_dists[sp_starts[i] + i] = 0.0;
750 
751  /* update distances by using single edges */
752  for( int j = dcsr_range[neighbor].start; j != dcsr_range[neighbor].end; j++ )
753  {
754  const int head = dcsr_heads[j];
755  const int head_index = pr->nodes_index[head];
756 
757  /* unvisited or path tail? */
758  if( head_index < 0 || head == basetail )
759  continue;
760 
761  assert(0 <= head_index && head_index < pr->nfirstneighbors);
762 
763 #ifdef SCIP_DEBUG
764  if( LT(dcsr_costs[j], sp_dists[sp_starts[i] + head_index]) )
765  printf("updating %d->%d distold=%f distnew=%f \n", i, head_index, sp_dists[sp_starts[i] + head_index], dcsr_costs[j]);
766 #endif
767  if( LT(dcsr_costs[j], sp_dists[sp_starts[i] + head_index]) )
768  sp_dists[sp_starts[i] + head_index] = dcsr_costs[j];
769  }
770  }
771 
772  for( int i = 0; i < pr->nfirstneighbors; i++ )
773  pr->firstneighbors_isKilled[i] = FALSE;
774 
775  addInitialPathNodes(scip, g, basetail, basehead, pr);
776 
777  /* NOTE: firstneighborcosts are only used for a subpath, thus we can already subtract the prize here */
778  if( isPc )
779  {
780  for( int i = 0; i < pr->nfirstneighbors; i++ )
781  {
782  pr->firstneighborcosts[i] -= g->prize[basetail];
783  assert(GE(pr->firstneighborcosts[i], 0.0));
784  }
785  }
786 }
787 
788 
789 /** tries to eliminate edges of path starting from edge */
790 static inline
792  SCIP* scip, /**< SCIP */
793  const GRAPH* g, /**< graph data structure */
794  PR* pr, /**< path replacement data structure */
795  SCIP_Bool* isExendible, /**< extendible? */
796  SCIP_Bool* isRedundant /**< redundant? */
797  )
798 {
799  const int npathedges = StpVecGetSize(pr->pathedges);
800  const int pathhead = pathGetHead(g, pr);
801  SCIP_Bool needCombRuleOut = FALSE;
802  SCIP_Bool isSingleRuledOut = FALSE;
803  SCIP_Bool isCombRuledOut = FALSE;
804 
805  assert(*isExendible);
806  assert(!(*isRedundant));
807 
808  if( npathedges > SP_MAXLENGTH || g->grad[pathhead] > pr->maxdegree || pr->nodes_isTerm[pathhead] )
809  {
810  *isExendible = FALSE;
811  return;
812  }
813 
814  pathneighborsCollect(scip, g, pr);
815  pathneighborsUpdateDistances(scip, g, pr);
816 
817  if( StpVecGetSize(pr->currneighbors) == 0 )
818  {
819  // todo really true?
820  *isRedundant = TRUE;
821  return;
822  }
823 
824  ruleOutFromHead(scip, g, pr, &needCombRuleOut);
825 
826  /* NOTE: if we have exactly one neighbor, and a failed neighbor, they are the same and we have to extend */
827  if( StpVecGetSize(pr->currneighbors) == 1 && pr->failneighbor != VNODES_UNSET )
828  {
829  addPathHeadEdge(scip, g, pr);
830  return;
831  }
832 
833  /* we try combination rule-out anyway! */
834  ruleOutFromTailCombs(scip, g, pr, &isCombRuledOut);
835 
836  if( needCombRuleOut && !isCombRuledOut )
837  {
838  *isExendible = FALSE;
839  return;
840  }
841 
842  ruleOutFromTailSingle(scip, g, pr, &isSingleRuledOut);
843 
844  if( !isSingleRuledOut )
845  {
846  *isExendible = FALSE;
847  return;
848  }
849 
850  if( isCombRuledOut )
851  {
852  *isRedundant = TRUE;
853  return;
854  }
855 
856  if( pr->failneighbor == VNODES_UNSET )
857  {
858  *isRedundant = TRUE;
859  return;
860  }
861 
862  /* if neither of the above cases holds, we extend along the failed neighbor */
863  addPathHeadEdge(scip, g, pr);
864 }
865 
866 
867 /** initializes */
868 static
870  SCIP* scip, /**< SCIP data structure */
871  const GRAPH* g, /**< graph data structure */
872  EXTPERMA* extperma, /**< extension data or NULL */
873  PR** pathreplace /**< to initialize */
874  )
875 {
876  PR* pr;
877  const int nnodes = graph_get_nNodes(g);
878 
879  SCIP_CALL( SCIPallocMemory(scip, pathreplace) );
880  pr = *pathreplace;
881 
882  pr->useSd = (extperma != NULL);
883  pr->probIsPc = graph_pc_isPc(g);
884  pr->maxdegree = pr->probIsPc ? SP_MAXDEGREE_PC : SP_MAXDEGREE;
885  if( extperma && extperma->mode == extred_fast )
886  pr->maxdegree = SP_MAXDEGREE_FAST;
887 
888  SCIP_CALL( SCIPallocMemoryArray(scip, &pr->nodes_isTerm, nnodes) );
889  SCIP_CALL( SCIPallocMemoryArray(scip, &pr->nodeindices_isPath, nnodes) );
890  SCIP_CALL( SCIPallocMemoryArray(scip, &pr->nodes_index, nnodes) );
891  SCIP_CALL( SCIPallocMemoryArray(scip, &pr->sp_starts, SP_MAXNSTARTS) );
892  SCIP_CALL( SCIPallocMemoryArray(scip, &pr->sp_dists, SP_MAXNDISTS) );
893  SCIP_CALL( SCIPallocMemoryArray(scip, &pr->firstneighbors_isKilled, pr->maxdegree + 1) );
894 
895  SCIP_CALL( reduce_sdprofitInit(scip, g, &pr->sdprofit) );
896 
897  for( int i = 0; i < nnodes; i++ )
898  pr->nodes_index[i] = VNODES_UNSET;
899 
900  graph_getIsTermArray(g, pr->nodes_isTerm);
901 
902  pr->extperma = extperma;
903  pr->currneighbors = NULL;
904  pr->firstneighborcosts = NULL;
905  pr->pathedges = NULL;
906  pr->firstneighbors = NULL;
907  pr->visitednodes = NULL;
908  pr->nfirstneighbors = -1;
909  pr->nnodes = nnodes;
910  StpVecReserve(scip, pr->currneighbors, pr->maxdegree + 1);
911  StpVecReserve(scip, pr->firstneighborcosts, pr->maxdegree + 1);
912  StpVecReserve(scip, pr->firstneighbors, pr->maxdegree + 1);
913  StpVecReserve(scip, pr->pathedges, SP_MAXLENGTH);
914  StpVecReserve(scip, pr->visitednodes, SP_MAXLENGTH);
915 
916  return SCIP_OKAY;
917 }
918 
919 
920 /** cleans temporary data */
921 static
922 void prClean(
923  PR* pr /**< to be cleaned */
924  )
925 {
926  const int nvisited = StpVecGetSize(pr->visitednodes);
927 
928  for( int i = 0; i < nvisited; i++ )
929  {
930  const int node = pr->visitednodes[i];
931  assert(node >= 0 && node < pr->nnodes && pr->nodes_index[node] != VNODES_UNSET);
932  pr->nodes_index[node] = VNODES_UNSET;
933  }
934 
935  StpVecClear(pr->firstneighbors);
936  StpVecClear(pr->firstneighborcosts);
937  StpVecClear(pr->currneighbors);
938  StpVecClear(pr->pathedges);
939  StpVecClear(pr->visitednodes);
940 
941 #ifndef NDEBUG
942  for( int i = 0; i < pr->nnodes; i++ )
943  assert(pr->nodes_index[i] == VNODES_UNSET);
944 
945  for( int i = 0; i < SP_MAXNSTARTS; i++ )
946  pr->sp_starts[i] = -1;
947 #endif
948 }
949 
950 
951 /** frees */
952 static
953 void prFree(
954  SCIP* scip, /**< SCIP data structure */
955  PR** pathreplace /**< to be freed */
956  )
957 {
958  PR* pr = *pathreplace;
959 
960  StpVecFree(scip, pr->pathedges);
961  StpVecFree(scip, pr->currneighbors);
962  StpVecFree(scip, pr->visitednodes);
963  StpVecFree(scip, pr->firstneighborcosts);
964  StpVecFree(scip, pr->firstneighbors);
965 
966  SCIPfreeMemoryArray(scip, &pr->firstneighbors_isKilled);
967  SCIPfreeMemoryArray(scip, &pr->sp_dists);
968  SCIPfreeMemoryArray(scip, &pr->sp_starts);
969  SCIPfreeMemoryArray(scip, &pr->nodes_index);
970  SCIPfreeMemoryArray(scip, &pr->nodeindices_isPath);
971  SCIPfreeMemoryArray(scip, &pr->nodes_isTerm);
972  reduce_sdprofitFree(scip, &pr->sdprofit);
973 
974  SCIPfreeMemory(scip, pathreplace);
975 }
976 
977 
978 /** is execution of path replacement reduction method promising? */
979 static
981  const GRAPH* g /**< graph structure */
982  )
983 {
984  // todo
985  return TRUE;
986 }
987 
988 
989 /** tries to eliminate edges of path starting from edge */
990 static inline
992  SCIP* scip, /**< SCIP data structure */
993  int startedge, /**< edge to start from (head) */
994  PR* pr, /**< path replacement data structure */
995  GRAPH* g, /**< graph data structure */
996  int* nelims /**< pointer to number of reductions */
997  )
998 {
999  SCIP_Bool pathIsExtendable = TRUE;
1000  SCIP_Bool pathIsRedundant = FALSE;
1001  const int tail = g->tail[startedge];
1002  const int head = g->head[startedge];
1003  const int maxdeg = pr->maxdegree;
1004  assert(StpVecGetSize(pr->pathedges) == 0);
1005 
1006  if( g->grad[tail] > maxdeg || g->grad[head] > maxdeg || pr->nodes_isTerm[tail] || pr->nodes_isTerm[head] )
1007  return SCIP_OKAY;
1008 
1009  if( g->grad[tail] <= 1 )
1010  return SCIP_OKAY;
1011 
1012  pathExendPrepare(scip, g, startedge, pr);
1013 
1014  while( pathIsExtendable )
1015  {
1016  assert(!pathIsRedundant);
1017  pathExend(scip, g, pr, &pathIsExtendable, &pathIsRedundant);
1018 
1019  if( pathIsRedundant )
1020  {
1021  SCIPdebugMessage("deleting edge %d-%d \n", g->tail[startedge], g->head[startedge]);
1022 
1023  deleteEdge(scip, startedge, g, pr);
1024 
1025  (*nelims)++;
1026  break;
1027  }
1028  }
1029 
1030  prClean(pr);
1031 
1032  return SCIP_OKAY;
1033 }
1034 
1035 
1036 /** executes reduction method */
1037 static
1039  SCIP* scip, /**< SCIP data structure */
1040  PR* pr, /**< path replacement data structure */
1041  GRAPH* g, /**< graph data structure */
1042  int* nelims /**< pointer to number of reductions */
1043  )
1044 {
1045  const int nedges = graph_get_nEdges(g);
1046  const SCIP_Bool isPc = pr->probIsPc;
1047  SD* const sddata = pr->useSd ? pr->extperma->distdata_default->sdistdata : NULL;
1048  // todo have edge stack?
1049 
1050  for( int e = 0; e < nedges; e++ )
1051  {
1052  /* deleted edge? */
1053  if( g->oeat[e] == EAT_FREE )
1054  continue;
1055 
1056  /* dummy edge? */
1057  if( isPc && graph_pc_edgeIsExtended(g, e) )
1058  continue;
1059 
1060  if( sddata && reduce_sdgraphEdgeIsInMst(sddata->sdgraph, e) )
1061  continue;
1062 
1063  SCIP_CALL( processPath(scip, e, pr, g, nelims) );
1064 
1065  // todo if edgestack: check whether edge was killed, otherwise go different!
1066  // afterwards deactivate edge!
1067  }
1068 
1069 
1070  return SCIP_OKAY;
1071 }
1072 
1073 /*
1074  * Interface methods
1075  */
1076 
1077 
1078 /** paths elimination while using special distances
1079  * NOTE: SD-repair needs to be set-up! */
1081  SCIP* scip, /**< SCIP data structure */
1082  GRAPH* g, /**< graph data structure */
1083  EXTPERMA* extperma, /**< extension data */
1084  int* nelims /**< pointer to number of reductions */
1085  )
1086 {
1087  PR* pathreplace;
1088  SCIP_Bool hasDcsr = (g->dcsr_storage != NULL);
1089 
1090  assert(scip && g && extperma);
1091 
1092  if( !prIsPromising(g) )
1093  {
1094  return SCIP_OKAY;
1095  }
1096 
1097  if( !hasDcsr )
1098  SCIP_CALL( graph_init_dcsr(scip, g) );
1099 
1100  SCIP_CALL( prInit(scip, g, extperma, &pathreplace) );
1101  SCIP_CALL( pathreplaceExec(scip, pathreplace, g, nelims) );
1102  prFree(scip, &pathreplace);
1103 
1104  if( !hasDcsr )
1105  graph_free_dcsr(scip, g);
1106 
1107  if( *nelims > 0 )
1108  {
1109  deletenodesDeg1(scip, extperma, g);
1110  }
1111 
1112  assert(graph_valid(scip, g));
1113 
1114  return SCIP_OKAY;
1115 }
1116 
1117 
1118 /** paths elimination */
1120  SCIP* scip, /**< SCIP data structure */
1121  GRAPH* g, /**< graph data structure */
1122  int* nelims /**< pointer to number of reductions */
1123  )
1124 {
1125  PR* pathreplace;
1126 
1127  if( !prIsPromising(g) )
1128  {
1129  return SCIP_OKAY;
1130  }
1131 
1132  SCIP_CALL( graph_init_dcsr(scip, g) );
1133  SCIP_CALL( prInit(scip, g, NULL, &pathreplace) );
1134 
1135  SCIP_CALL( pathreplaceExec(scip, pathreplace, g, nelims) );
1136 
1137  prFree(scip, &pathreplace);
1138  graph_free_dcsr(scip, g);
1139 
1140  if( *nelims > 0 )
1141  {
1142  reduce_nodesDeg1(scip, g);
1143  }
1144 
1145  assert(graph_valid(scip, g));
1146 
1147  return SCIP_OKAY;
1148 }
#define STP_Vectype(type)
Definition: stpvector.h:44
static void addPathHeadEdge(SCIP *scip, const GRAPH *g, PR *pr)
Definition: reduce_path.c:436
SDPROFIT * sdprofit
Definition: reduce_path.c:48
int graph_get_nEdges(const GRAPH *)
Definition: graph_stats.c:548
#define StpVecGetSize(vec)
Definition: stpvector.h:139
#define StpVecClear(vec)
Definition: stpvector.h:134
void reduce_nodesDeg1(SCIP *, GRAPH *)
int *RESTRICT head
Definition: graphdefs.h:224
enum EXTRED_MODE mode
#define Is_term(a)
Definition: graphdefs.h:102
EXTPERMA * extperma
Definition: reduce_path.c:47
SCIP_RETCODE reduce_pathreplaceExt(SCIP *scip, GRAPH *g, EXTPERMA *extperma, int *nelims)
Definition: reduce_path.c:1080
#define SCIPfreeMemoryArray(scip, ptr)
Definition: scip_mem.h:71
SCIP_Bool graph_valid(SCIP *, const GRAPH *)
Definition: graph_base.c:1480
SCIP_Real extreduce_distDataGetSdDoubleEq(SCIP *, const GRAPH *, SCIP_Real, int, int, DISTDATA *)
void graph_getIsTermArray(const GRAPH *, SCIP_Bool *)
Definition: graph_base.c:540
#define SCIPallocMemoryArray(scip, ptr, num)
Definition: scip_mem.h:55
#define SP_MAXDEGREE
Definition: reduce_path.c:36
static void ruleOutFromTailCombs(SCIP *scip, const GRAPH *g, PR *pr, SCIP_Bool *isRuledOut)
Definition: reduce_path.c:334
static SCIP_RETCODE prInit(SCIP *scip, const GRAPH *g, EXTPERMA *extperma, PR **pathreplace)
Definition: reduce_path.c:869
void graph_free_dcsr(SCIP *, GRAPH *)
Definition: graph_util.c:1807
#define EAT_FREE
Definition: graphdefs.h:57
#define FALSE
Definition: def.h:87
struct path_replacement PR
#define TRUE
Definition: def.h:86
enum SCIP_Retcode SCIP_RETCODE
Definition: type_retcode.h:54
SCIP_Real * cost
Definition: graphdefs.h:164
#define StpVecPopBack(vec)
Definition: stpvector.h:182
static SCIP_RETCODE processPath(SCIP *scip, int startedge, PR *pr, GRAPH *g, int *nelims)
Definition: reduce_path.c:991
#define StpVecReserve(scip, vec, _size_)
Definition: stpvector.h:186
static void pathExend(SCIP *scip, const GRAPH *g, PR *pr, SCIP_Bool *isExendible, SCIP_Bool *isRedundant)
Definition: reduce_path.c:791
static void addNonPathNode(SCIP *scip, int node, PR *pr)
Definition: reduce_path.c:420
#define SCIPdebugMessage
Definition: pub_message.h:87
#define FARAWAY
Definition: graphdefs.h:89
#define VNODES_UNSET
Definition: reduce_path.c:40
SCIP_RETCODE graph_init_dcsr(SCIP *, GRAPH *)
Definition: graph_util.c:1721
int *RESTRICT mark
Definition: graphdefs.h:198
static void pathExendPrepare(SCIP *scip, const GRAPH *g, int startedge, PR *pr)
Definition: reduce_path.c:666
static SCIP_Real getSdProfit(const SDPROFIT *sdprofit, const int *nodes_index, int node, int nonsource)
Definition: reduce_path.c:77
int *RESTRICT oeat
Definition: graphdefs.h:231
This file implements extended reduction techniques for several Steiner problems.
SCIP_Bool graph_pc_edgeIsExtended(const GRAPH *, int)
#define LE(a, b)
Definition: portab.h:83
#define GE(a, b)
Definition: portab.h:85
SCIP_Real * prize
Definition: graphdefs.h:210
static void deletenodesDeg1(SCIP *scip, EXTPERMA *extperma, GRAPH *g)
Definition: reduce_path.c:157
int *RESTRICT grad
Definition: graphdefs.h:201
#define NULL
Definition: lpi_spx1.cpp:155
static void prFree(SCIP *scip, PR **pathreplace)
Definition: reduce_path.c:953
void extreduce_edgeRemove(SCIP *, int, GRAPH *, DISTDATA *, EXTPERMA *)
#define EQ(a, b)
Definition: portab.h:79
void reduce_sdprofitFree(SCIP *, SDPROFIT **)
#define SCIP_CALL(x)
Definition: def.h:384
#define LT(a, b)
Definition: portab.h:81
void graph_edge_delFull(SCIP *, GRAPH *, int, SCIP_Bool)
Definition: graph_edge.c:418
static void ruleOutFromHead(SCIP *scip, const GRAPH *g, PR *pr, SCIP_Bool *needFullRuleOut)
Definition: reduce_path.c:194
#define GT(a, b)
Definition: portab.h:84
#define SCIP_Bool
Definition: def.h:84
int *RESTRICT tail
Definition: graphdefs.h:223
#define flipedge(edge)
Definition: graphdefs.h:84
static void pathneighborsUpdateDistances(SCIP *scip, const GRAPH *g, PR *pr)
Definition: reduce_path.c:580
int *RESTRICT term
Definition: graphdefs.h:196
SCIP_Real *RESTRICT nodes_bias2
Definition: reducedefs.h:138
SCIP_Bool reduce_sdgraphEdgeIsInMst(const SDGRAPH *, int)
SCIP_Bool graph_pc_isPc(const GRAPH *)
#define CONNECT
Definition: graphdefs.h:87
SCIP_Real extreduce_distDataGetSdDoubleForbiddenSingle(SCIP *, const GRAPH *, int, int, int, DISTDATA *)
static void pathneighborsCollect(SCIP *scip, const GRAPH *g, PR *pr)
Definition: reduce_path.c:528
#define SCIPfreeMemory(scip, ptr)
Definition: scip_mem.h:69
static void addPathNode(SCIP *scip, int node, PR *pr)
Definition: reduce_path.c:405
static int pathGetHead(const GRAPH *g, const PR *pr)
Definition: reduce_path.c:115
SCIP_RETCODE reduce_pathreplace(SCIP *scip, GRAPH *g, int *nelims)
Definition: reduce_path.c:1119
static void deleteEdge(SCIP *scip, int edge, GRAPH *g, PR *pr)
Definition: reduce_path.c:127
static void prClean(PR *pr)
Definition: reduce_path.c:922
SCIP_Bool graph_pc_knotIsDummyTerm(const GRAPH *, int)
SCIP_Real *RESTRICT nodes_bias
Definition: reducedefs.h:137
SCIP_Real * cost
Definition: graphdefs.h:221
#define SP_MAXNPULLS
Definition: reduce_path.c:41
#define SCIP_Real
Definition: def.h:177
#define StpVecFree(scip, vec)
Definition: stpvector.h:153
#define SP_MAXNSTARTS
Definition: reduce_path.c:39
static SCIP_Bool prIsPromising(const GRAPH *g)
Definition: reduce_path.c:980
int *RESTRICT outbeg
Definition: graphdefs.h:204
#define SP_MAXDEGREE_FAST
Definition: reduce_path.c:37
SCIP_RETCODE reduce_sdprofitInit(SCIP *, const GRAPH *, SDPROFIT **)
#define SP_MAXLENGTH
Definition: reduce_path.c:35
#define SP_MAXNDISTS
Definition: reduce_path.c:34
int graph_get_nNodes(const GRAPH *)
Definition: graph_stats.c:536
#define SCIPallocMemory(scip, ptr)
Definition: scip_mem.h:51
static SCIP_RETCODE pathreplaceExec(SCIP *scip, PR *pr, GRAPH *g, int *nelims)
Definition: reduce_path.c:1038
#define nnodes
Definition: gastrans.c:65
#define StpVecPushBack(scip, vec, value)
Definition: stpvector.h:167
static void ruleOutFromTailSingle(SCIP *scip, const GRAPH *g, PR *pr, SCIP_Bool *isRuledOut)
Definition: reduce_path.c:275
DCSR * dcsr_storage
Definition: graphdefs.h:271
includes various reduction methods for Steiner tree problems
static void addInitialPathNodes(SCIP *scip, const GRAPH *g, int startedge_tail, int startdege_head, PR *pr)
Definition: reduce_path.c:479
#define SP_MAXDEGREE_PC
Definition: reduce_path.c:38